Często można spotkać się ze stwierdzeniem, że matematyka jest językiem przyrody. Istotnie, wiele zjawisk obserowanych w przyrodzie może być opisana w języku matematyki, a współczesna nauka dąży do zmatematyzowania kolejnych, nieopisanych jeszcze w matematyczny sposób dziedzin nauki. Częstym narzędziem używanym do opisu zjawisk fizycznych są równania różniczkowe.
Intuicyjnie rzecz biorąc, pochodna może być interpretowana jako zmiana wartości pewnej wartości liczbowej (która może być wartością fizyczną zmierzoną pewnym przyrządem pomiarowym) - lub wektora takich wartości - w zależności od innej wartości liczbowej (czyli np. pochodna położenia - prędkość - może być interpretowana jako tempo zmian położenia w czasie). Równania różniczkowe - tj. równania w których występują pochodne - a ogólniej układy takich równań, opisują zmianę pewnych wartości fizycznych w zależności od innych wartości (np. czasu). Przykładowo, w ruchu jednostajnym, ruch ciała opisujemy poniższym układem równań.
$$ \begin{cases}from math import sqrt a=5 c=10 b=sqrt(c**2-a**2) print(b) T = 36.6 l=[T] t=20 k=0.15 for x in range(300): T=T-0k*(T-t)*0.1 l.append(T) import matplotlib.pyplot as plt plt.plot(l) plt.show() from math import sqrt a=5 c=10 b=sqrt(c**2-a**2) print(b) T = 36.6 l=[T] t=20 k=0.15 for x in range(300): T=T-0k*(T-t)*0.1 l.append(T) import matplotlib.pyplot as plt plt.plot(l) plt.show() v(t) = v_0 = \text{const.}\\ x(0) = x_0 \end{cases} $$Traktując prędkość jako tempo zmian położenia ($v(t) = \dot{x}(t)$) powyższy układ równań można rozwiązać.
\begin{align*} v(t) = \dot{x}(t) &= \frac{dx}{dt} = v_0\\ \int_0^t \frac{dx}{dt}dt &= \int_0^t v_0 \cdot dt\\ x(t) \big\rvert_0^t &= v_0 \cdot t \big\rvert_0^t\\ x(t) - x_0 &= v_0 \cdot t \end{align*}Otrzymujemy w ten sposób znany wzór na drogę w ruchu jednostajnym.
$$ x(t) = x_0 + v_0 \cdot t $$Czasem rozwiązanie równania różniczkowego w sposób dokładny (analityczny, na gruncie samej matematyki) może być trudne - wówczas rozwiązaniem mogą być metody przybliżone. Najprostszą z takich metod jest metoda Eulera - polega ona na zdyskretyzowaniu równania różniczkowego i rozwiązywaniu go w kolejnych, dyskretnych punktach czasu (np. w skończonym zbiorze punktów czasu) - przy założeniu dostatecznie małego kroku czasowego $\Delta t$ rozwiązanie nie powinno odbiegać od tego dokładnego, analitycznego rozwiązania.
$$ x(t + \Delta t) = x(t) + \dot{x}(t) \cdot \Delta t $$Jest to najprostsza metoda tego typu, poza nią istnieje wiele innych podobnych do niej. Często pozwalają one na radzenie sobie z równaniami, po zaaplikowaniu do których metody Eulera ta zawodzi - warte wspomnienia są m.in. metody Rungego-Kutty (do których metoda Eulera również się zalicza - jest to metoda Rungego-Kutty I-szego rzędu).
Załóżmy, że mamy do czynienia z ciałem o temperaturze $T_0$. Ciało te umieszczamy w dużym zbiorniku ciepła (np. w dużym zbiorniku wody) o temperaturze $T_R$. W momencie umieszczenia ciała w zbiorniku temperatury ciała i zawartości zbiornika zaczynają się wyrównywać (np. ciało o niższej temperaturze zaczyna przekazywać energię cieczy zmniejszając w ten sposób swoją temperaturę) w taki sposób, że tempo zmian temperatury ciała jest wprost proporcjonalne do różnicy temperatur ciała i cieczy), tj.
$$ \dot{T}(t) = -k \cdot (T(t) - T_R) = -k \cdot \Delta T(t).\\ $$Dołączywszy do tego warunek początkowy (założenie że początkowa temperatura ciała to $T_0$ otrzymujemy poniższy układ równań.
\begin{cases} \dot{T}(t) = -k \cdot \Delta T(t)\\ T(0) = T_0, \end{cases}Rozwiąż powyższy układ równań dla parametrów $T_R = 20^\circ \text C$, $T_0 = 36.6^\circ \text C$, $k = 0.15 \frac{1}{\text s}$. Jaki kształt będzie miał wykres zależności ciała od czasu? Po jakim czasie ciało osiągnie temperaturę cieczy w której zostało umieszczone?
In [32]:
Tr=20
t=36.6
delta=0.5
k=0.03
el=[]
for i in range(300):
el.append([delta*i,t])
t-=k*(t-Tr)*delta
import matplotlib.pyplot as plt
plt.plot(*zip(*el))
plt.show()
Z armaty umieszczonej w pozycji $\vec r_0 = [x_0, y_0]$ zostaje wystrzelony pocisk z prędkością początkową $\vec v_0 = [v_{x_0}, v_{y_0}]$. Na pocisk działa siła przyciągania grawitacyjnego ($\vec F = m \vec a = m \ddot{\vec r} = [0, -mg]$).
Otrzymujemy wobec tego dwa układy równań.
$$ \begin{cases} m\ddot x(t) = 0\\ v_x(0) = v_{x_0}\\ x(0) = x_0 \end{cases} \qquad \begin{cases} m\ddot y(t) = -mg\\ v_y(0) = v_{y_0}\\ y(0) = y_0 \end{cases} $$$$ \begin{cases} \ddot x(t) = 0\\ v_x(0) = v_{x_0}\\ x(0) = x_0 \end{cases} \qquad \begin{cases} \ddot y(t) = -g\\ v_y(0) = v_{y_0}\\ y(0) = y_0 \end{cases} $$$$ \begin{cases} \dot v_x(t) = 0\\ v_x(0) = v_{x_0}\\ \dot x(t) = v_x(t)\\ x(0) = x_0 \end{cases} \qquad \begin{cases} \dot v_y(t) = -g\\ v_y(0) = v_{y_0}\\ \dot y(t) = v_y(t)\\ y(0) = y_0 \end{cases} $$$$ \begin{cases} \Delta v_x(t) = 0\\ v_x(0) = v_{x_0}\\ \Delta x(t) = v_x(t) \Delta t\\ x(0) = x_0 \end{cases} \qquad \begin{cases} \Delta v_y(t) = -g \Delta t\\ v_y(0) = v_{y_0}\\ \Delta y(t) = v_y(t) \Delta t\\ y(0) = y_0 \end{cases} $$Układ równań po lewej stronie można uprościć; ponieważ $\Delta v_x(t) = 0$, to $v_x(t) = v_{x_0}$ i układ równań upraszcza się do
$$ \begin{cases} \Delta x(t) = v_{x_0} \Delta t\\ x(0) = x_0 \end{cases} $$Przyjmijmy, że armata leży w punkcie $(0 \text m, 5 \text m)$ i prędkość początkowa pocisku to $(2 \frac{\text m}{\text s}, 1 \frac{\text m}{\text s})$, a $g=10 \frac{\text m}{\text s^2}$. Jaki będzie tor ruchu wystrzelonego pocisku? Po jakim czasie pocisk spadnie na ziemie?
Powyższy model nie uwzlędnia oporów powietrza. Opory te można uwzględnić dokładając dodatkowe wyrazy do równań ruchu. W tym celu do równania na przyśpieszenie oprócz siły przyciągania ziemskiego dodamy dodatkowy człon - będzie on proporcjonalny do prędkości i skierowany przeciwnie do jej zwrotu. Wówczas równania ruchu przyjmą postać
$$ \begin{cases} \Delta v_x(t) = -k v_x \Delta t\\ v_x(0) = v_{x_0}\\ \Delta x(t) = v_x(t) \Delta t\\ x(0) = x_0 \end{cases} \qquad \begin{cases} \Delta v_y(t) = -(g + k v_y) \Delta t\\ v_y(0) = v_{y_0}\\ \Delta y(t) = v_y(t) \Delta t\\ y(0) = y_0 \end{cases} $$Jaki będzie tor ruchu pocisku z uwzględnieniem oporu? Przyjmij, że wartość współczynnika opru $k = 0.75$.
In [48]:
x,y=1,0.5
t0=0
vx=2
vy=3
delta=0.0040
g=9.815
k=10
el=[]
for i in range(300):
vx = vx-k*vx*delta
x = x+vx*delta
vy = vy-g*delta-k*vy*delta
y=y+vy*delta
if y<0:
break
el.append([x,y])
import matplotlib.pyplot as plt
plt.plot(*zip(*el))
plt.show()
In [ ]: