Runge-Kutta duseno nihanje
-
- Prispevkov: 19
- Pridružen: 1.12.2008 22:30
Runge-Kutta duseno nihanje
Lepo pozdravljeni!
Z metodo Runge-Kutta 4. stopnje bi rad rešil diferencialno enačbo (v nadaljevanju d.e.) za dušeno nihanje. Vem da morem d.e. zapisati kot sistem dvem d.e. prvega reda. Problem imam s koeficientov k1, k2, k3 in k4. Prosil bi, za pomoč ker ne vem kako naj upoštevam tisti t+1/2 in y+(1/2)k1. Hvala!
d.e.
\(\ddot{x}=-\frac{b}{m}\dot{x}-\frac{k}{m}x\)
kot sistem d.e. prvega reda
\(\frac{dv}{dt} =-\frac{b}{m}\dot{x}-\frac{k}{m}x\)
\(\frac{dx}{dt}=v\)
Runge-Kutta metoda
\(\dot{y}=f(t,\ y)\)
\(y_{n+1}=y_{n}+\frac{1}{6}(k_{1}+2k_{2}+2k_{3}+k_{4})\)
\(k_{1}=h f(t_{n},\ y_{n})\)
\(k_{2}=h f(t_{n}+\frac{1}{2}h,\ y_{n}+\frac{1}{2}k_{1})\)
\(k_{3}=h f(t_{n}+\frac{1}{2}h,\ y_{n}+\frac{1}{2}k_{2})\)
\(k_{4}=h f(t_{n}+h,\ y_{n}+k_{3})\)
Z metodo Runge-Kutta 4. stopnje bi rad rešil diferencialno enačbo (v nadaljevanju d.e.) za dušeno nihanje. Vem da morem d.e. zapisati kot sistem dvem d.e. prvega reda. Problem imam s koeficientov k1, k2, k3 in k4. Prosil bi, za pomoč ker ne vem kako naj upoštevam tisti t+1/2 in y+(1/2)k1. Hvala!
d.e.
\(\ddot{x}=-\frac{b}{m}\dot{x}-\frac{k}{m}x\)
kot sistem d.e. prvega reda
\(\frac{dv}{dt} =-\frac{b}{m}\dot{x}-\frac{k}{m}x\)
\(\frac{dx}{dt}=v\)
Runge-Kutta metoda
\(\dot{y}=f(t,\ y)\)
\(y_{n+1}=y_{n}+\frac{1}{6}(k_{1}+2k_{2}+2k_{3}+k_{4})\)
\(k_{1}=h f(t_{n},\ y_{n})\)
\(k_{2}=h f(t_{n}+\frac{1}{2}h,\ y_{n}+\frac{1}{2}k_{1})\)
\(k_{3}=h f(t_{n}+\frac{1}{2}h,\ y_{n}+\frac{1}{2}k_{2})\)
\(k_{4}=h f(t_{n}+h,\ y_{n}+k_{3})\)
Re: Runge-Kutta duseno nihanje
t+1/2 nima nobenega vpliva ker cas ne nastopa v tvoji diferencialni enacbi. Ostalo je pa takoj jasno, ce upostevas da je zdaj tvoj y vektor.
\(y=(x,v)\)
\(y'=(v,-b/m\cdot v-k/m\cdot x)\)
Zdaj imas zapisan ta y' kot funkcijo x in v (torej kot funkcijo y).
Zdaj bodo k1,k2,k3,... vektorji iz dveh komponent. Od tukaj naprej samo se izvedes postopek.
\(y=(x,v)\)
\(y'=(v,-b/m\cdot v-k/m\cdot x)\)
Zdaj imas zapisan ta y' kot funkcijo x in v (torej kot funkcijo y).
Zdaj bodo k1,k2,k3,... vektorji iz dveh komponent. Od tukaj naprej samo se izvedes postopek.
-
- Prispevkov: 19
- Pridružen: 1.12.2008 22:30
Re: Runge-Kutta duseno nihanje
Najlepša hvala za odgovor. Torej, če prav razumem se k-ji izračunajo tako?
\(k_{1}=\begin{bmatrix} v\\-\frac{b}{m}\ v-\frac{k}{m}\ x \end{bmatrix}\)
\(k_{2}=\begin{bmatrix} v+\frac{1}{2}\ v\\-\frac{b}{m}\ v-\frac{k}{m}\ x +\frac{1}{2}(-\frac{b}{m}\ v-\frac{k}{m}\ x) \end{bmatrix}\)
itd.
\(k_{1}=\begin{bmatrix} v\\-\frac{b}{m}\ v-\frac{k}{m}\ x \end{bmatrix}\)
\(k_{2}=\begin{bmatrix} v+\frac{1}{2}\ v\\-\frac{b}{m}\ v-\frac{k}{m}\ x +\frac{1}{2}(-\frac{b}{m}\ v-\frac{k}{m}\ x) \end{bmatrix}\)
itd.
Re: Runge-Kutta duseno nihanje
Ne, ti si naredil k2=f(y)+k1. k1 mora biti v argumentu. Torej
\(k_2=\begin{bmatrix}
v+h/2 (-b/m v-k/m x)\\
-b/m (v+h/2 (-b/m v-k/m x))-k/m(x+h/2 v)
\end{bmatrix}\)
Se pravi: kot "x" vstavis prvotni x plus x-komponento od k1 (krat tista polovicka koraka). In kot "v" vstavis prvotni v plus v-komponento od k1 (se pravi drugo).
Seveda se nikoli analiticno ne vstavlja, finta je v tem da to delas z numericnimi podatki. Poracunas stevilke za k1 in potem naprej.
\(k_2=\begin{bmatrix}
v+h/2 (-b/m v-k/m x)\\
-b/m (v+h/2 (-b/m v-k/m x))-k/m(x+h/2 v)
\end{bmatrix}\)
Se pravi: kot "x" vstavis prvotni x plus x-komponento od k1 (krat tista polovicka koraka). In kot "v" vstavis prvotni v plus v-komponento od k1 (se pravi drugo).
Seveda se nikoli analiticno ne vstavlja, finta je v tem da to delas z numericnimi podatki. Poracunas stevilke za k1 in potem naprej.
-
- Prispevkov: 19
- Pridružen: 1.12.2008 22:30
Re: Runge-Kutta duseno nihanje
A tako. V praksi bi to tako zgledalo ne?
Koda: Izberi vse
// funkcija za izracun pospeska
double pos(double v, double x){
return (-b/m*v - k/m*x);
}
// Runge-Kutta metoda
void RK4(double dt, double v, double x){
k1v = pos(v,x)*dt;
k1x = v * dt;
k2v = pos(v + k1v/2 , x + k1x/2) * dt;
k2x = (v + k1v/2)*dt;
k3v = pos(v + k2v/2 , x + k2x/2) * dt;
k3x = (v + k2v/2) * dt;
k4v = pos(v + k3v , x + k3x) * dt;
k4x = (v + k3x) * dt;
vnp1 = v + (k1v + 2*k2v + 2*k3v + k4v) / 6.0;
xnp1 = x + (k1x + 2*k2x + 2*k3x + k4x) / 6.0;
}
Re: Runge-Kutta duseno nihanje
Tako nekako ja. Edino pri k4x ti nastopa notri k3x.
V splosnem, ce izraz za odvod x-a ne bi bil tako zelo enostaven, bi imel pa poleg "pos" se eno funkcijo za odvod komponente x. Tako da bi bili izrazi simetricni. Ponavadi numericne kniznjice, ki implementirajo te integratorje (v bistvu katerega koli), zahtevajo funkcijo, ki sprejme vektor in zapise rezultate v drug vektor. Recimo takole
void f(double t, double y[], double dy[]){
dy[0]=y[1];
dy[1]=-b/m*y[1]-k/m*y[0];
}
ce je y={x,v};
Tako da potem integratorja ne briga kako ti to racunas, pa tudi za mnozenje z dt sam poskrbi (ker ponavadi itak adaptivno spreminja korak glede na ocenjeno napako).
V splosnem, ce izraz za odvod x-a ne bi bil tako zelo enostaven, bi imel pa poleg "pos" se eno funkcijo za odvod komponente x. Tako da bi bili izrazi simetricni. Ponavadi numericne kniznjice, ki implementirajo te integratorje (v bistvu katerega koli), zahtevajo funkcijo, ki sprejme vektor in zapise rezultate v drug vektor. Recimo takole
void f(double t, double y[], double dy[]){
dy[0]=y[1];
dy[1]=-b/m*y[1]-k/m*y[0];
}
ce je y={x,v};
Tako da potem integratorja ne briga kako ti to racunas, pa tudi za mnozenje z dt sam poskrbi (ker ponavadi itak adaptivno spreminja korak glede na ocenjeno napako).
-
- Prispevkov: 19
- Pridružen: 1.12.2008 22:30
Re: Runge-Kutta duseno nihanje
OK, super! Najlepša hvala za pomoč!
Re: Runge-Kutta duseno nihanje
Pozdravljeni!
Imam težave pri numeričnem računanju. Naloga pravi da je potrebno poiskati višino, hitrost in pospešek padalca.
\(\ddot{x} = \dot{x}^2 c \rho S/2m - g\) (1) to sem zapisal za višino
zanima me sledeče: reševal sem v matlabu in po gornji enačbi z metodo ode45 za reševanje de, sem višino izračunal pravilno, hitrost pa ki nastopa poleg pa je napačna. Navodilo pravi da je potrebno izračunati z Runge kutta reda 4. Ali je možno, da hitrost ne pride pravilno izračunana zaradi napake metode ?
Poleg tega pa še to: če izračunam hitrost iz (1), bi ta mogla biti identična hitrosti iz enačbe ali se motim ?? ker meni pridejo odstopanja velikosti 1m/s
\(\dot{v} = v^2 c \rho S/2 -g\)
Hvala za odgovore
Imam težave pri numeričnem računanju. Naloga pravi da je potrebno poiskati višino, hitrost in pospešek padalca.
\(\ddot{x} = \dot{x}^2 c \rho S/2m - g\) (1) to sem zapisal za višino
zanima me sledeče: reševal sem v matlabu in po gornji enačbi z metodo ode45 za reševanje de, sem višino izračunal pravilno, hitrost pa ki nastopa poleg pa je napačna. Navodilo pravi da je potrebno izračunati z Runge kutta reda 4. Ali je možno, da hitrost ne pride pravilno izračunana zaradi napake metode ?
Poleg tega pa še to: če izračunam hitrost iz (1), bi ta mogla biti identična hitrosti iz enačbe ali se motim ?? ker meni pridejo odstopanja velikosti 1m/s
\(\dot{v} = v^2 c \rho S/2 -g\)
Hvala za odgovore
Re: Runge-Kutta duseno nihanje
Za takole enacbo bi morala biti kunge-kutta zelo robustna. Napake bi morale biti kvecjemu nekje na 7. decimalki ali kaj takega.
Sklepam, da resujes sistem
\(\dot{x}=v\)
\(\dot{v}=...\)
da dobis zraven poleg hitrosti se visino (ceprav v principu lahko integriras samo drugo enacbo in dobis potek hitrosti, to pa potem integriras se enkrat na katerikoli nacin, da dobis polozaj).
Pri drugi enacbi ti masa manjka, samo to je najbrz samo pozabljeno pri prepisovanju.
Sklepam, da resujes sistem
\(\dot{x}=v\)
\(\dot{v}=...\)
da dobis zraven poleg hitrosti se visino (ceprav v principu lahko integriras samo drugo enacbo in dobis potek hitrosti, to pa potem integriras se enkrat na katerikoli nacin, da dobis polozaj).
Pri drugi enacbi ti masa manjka, samo to je najbrz samo pozabljeno pri prepisovanju.