Eksempel 15.2 fall_med_luftmotstand.py

# Program fall_med_luftmotstand.py
import numpy as np
import matplotlib.pyplot as plt

# Definerer variabler
tMax = 10 # tid maks [s]
dt = 0.001 # tidsintervall [s]
n = int(tMax/dt) # antall tidsintervaller
t = np.linspace(0,tMax,n) # tid-liste [s]
g = 9.81 # tyngdeakselerasjonen [m/s^2]
r = 0.01 # luftmotstandskonstant [m^-1]

s = np.zeros(n) # posisjon startliste [m]
v = np.zeros(n) # fart startliste [m/s]
a = np.zeros(n) # akselerasjon startliste [m/s^2]

# Setter startverdier for s, v og a
s[0] = 0 # posisjon startverdi [m]
v[0] = 0 # fart startverdi [m/s]
a[0] = -g + r*v[0]**2 # akselerasjon startverdi [m/s^2]

# Steg for steg beregning av a-, v- og s-lister
for i in range(n-1): 
  s[i+1] = s[i] + v[i]*dt # posisjon-liste [m] 
  v[i+1] = v[i] + a[i]*dt # fart-liste [m/s]
  a[i+1] = -g + r*v[i+1]**2 # akselerasjon-liste [m/s^2] 

# Akselerasjonsgraf 
plt.figure(1)
plt.plot(t,a,'-r')
plt.grid()
plt.title ('Akselerasjon')
plt.xlabel('$t$ / s')
plt.ylabel('$a$ / m/s$^2$')

# Fartsgraf
plt.figure(2)
plt.plot(t,v,'-b')
plt.grid()
plt.title ('Fart')
plt.xlabel('$t$ / s')
plt.ylabel('$v$ / m/s')

# Posisjonsgraf
plt.figure(3)
plt.plot(t,s,'-g')
plt.grid()
plt.title ('Posisjon')
plt.xlabel('$t$ / s')
plt.ylabel('$s$ / m')

plt.show()