s***e
发帖数: 911 | 1 我推一下. 你用前先debug一下.
一个n阶ODE:
x^{n}=g(x^{n-1},...,x'',x',x)
可以化成n个一阶方程组:
x'=y1
y1'=y2
y2'=y3
......
y^{n-1}=g(y_{n-1},...,y2,y1,x)
上标代表求导.下标就是index. 你的方程是二价方程
d^2 x/(dt)^2 + b^2 * x = c
所以化下来是:
x'=y
y'=c-b^2 * x
下面是一般理论:
Y是一个N维矢量, F是N维向量函数.
Y'=F(Y)
是一个N维一阶微分方程组. F里面不含自变量, 叫自治系统. 非自治系统和你无关,
就不管了. 假设自变量是t, 那么一阶近似是:
Y_{m+1}=Y_{m}+h*F(Y_{m})
四阶情况如下:
K_1=h*F(Y_{m})
K_2=h*F(Y_{m}+(1/2)K_1)
K_3=h*F(K_{m}+(1/2)K_2)
K_4=h*F(Y_{m}+K_3)
那么,
Y_{m+1}=Y_{m}+K_1/6+K_2/3+K_3/3+K_4/6
要注意K_1,K_2,K_3,K_4都是矢量.
你的问题F( |
|
