一阶常微分方程的初值

来源:百度知道 编辑:UC知道 时间:2024/05/16 04:12:45
用龙格-库塔方法或EULAR解以下方程组

dy1/dt =y2
dy2/dt =-y1
dy3/dt =-y2

其中y1(0)=-1, y2(0)=0 y3(0)=1
我要的是人工手算或者用C++或C编程算,不要用MATLAB,也不要什么原理

解:设向量り=(y1,y2,y3)' 其中'表示转置,则原方程组可以写成矩阵形式
dり/dt=Aり , 其中矩阵A=
0 1 0
-1 0 0
0 -1 0

从det(入E-A)=入(入2+1)=0 求得特征值0,i,-i
对于特征值0,求解方程组(A-0E)X=0 得到基础解系X1=(0,0,1)'
对于特征值i,求解方程组(A-iE)X=0 得到基础解系X2=(1,i,-1)
对于特征值-i,求解方程组(A+iE)X=0得到基础解系X3=(-1,i,1)

这样就得到了基本解矩阵
り(t)=(X1e^0,X2e^it,X3e^-it)=
0 e^it -e^it
0 ie^it ie^-it
1 -e^it e^-it

于是就得到了原方程组的一个复通解
y1=C2e^it-C3e^-it
y2=C2ie^it+C3ie^-it
y3=C1-C2e^it+C3e^-it

根据实虚通解对应原理,得到原方程组得实通解
y1=C2cost-C3sint
y2=-C2sint+C3cost
y3=C1-C2cost+C3sint

代入初值条件y1(0)=-1,y2(0)=0,y3(0)=1
得到C1=0,C2=-1,C3=0

于是,原方程的解为
y1=-cost
y2=sint
y3=cost

这里, 为常微分方程的右端函数,而 为所求未知函数的初始值。求解常微分方程初值问题用指令ode23 或ode45。使用这两条命令中的任何一条都必须事先编写好函数文件并保存在工作目录下(如取文件名为yprime.m)。命令的具体使用格式为

[x,y] = ode23('yprime', x0, xn, y0)

其中,yprime 为描述常微分方程右端