通过MATLAB求二阶全微分方程解析解.docx
《通过MATLAB求二阶全微分方程解析解.docx》由会员分享，可在线阅读，更多相关《通过MATLAB求二阶全微分方程解析解.docx（14页珍藏版）》请在冰豆网上搜索。
通过MATLAB求二阶全微分方程解析解
1.对于二阶全微分方程a
，不同的a,b,c取值会求出不同的解析解，解析解又是由齐次解和特解组成。
其中，齐次解由特征方程决定，而特解的决定因素则比较复杂。
2.对于二阶全微分方程的分析，我们大致分为三种情况：
b^2-4ac>0（两个不同的实根）
b^2-4ac=0（两个相同的重根）
b^2-4ac<0（两个不同的复数根）
对三种情况进行MATLAB编程，分析齐次解和特解后，再改变W的值，观察解析解的变化
3.b^2-4ac>0的情况
STEP1：
求解析解
s1=dsolve（'D2y+3*Dy+2*y=0','y（0）=2,Dy（0）=0','t'）;
s2=dsolve（'D2y+3*Dy+2*y=sin（t）','y（0）=2,Dy（0）=0','t'）;
s3=dsolve（'D2y+3*Dy+2*y=sin（2*t）','y（0）=2,Dy（0）=0','t'）;
s4=dsolve（'D2y+3*Dy+2*y=sin（5*t）','y（0）=2,Dy（0）=0','t'）;
s5=dsolve（'D2y+3*Dy+2*y=sin（13*t）','y（0）=2,Dy（0）=0','t'）;
s6=dsolve（'D2y+3*Dy+2*y=sin（25*t）','y（0）=2,Dy（0）=0','t'）;
STEP2：
绘制图形
（1）求w=1情况下的通解和齐次解
t=1:
0.1:
10;
s1=4*exp（-t）-2*exp（-2*t）%generalsolution
s2=-3/10*cos（t）+1/10*sin（t）-11/5*exp（-2*t）+9/2*exp（-t）%specialsolution
subplot（2,1,1）;
plot（t,s2）;
xlabel（'t'）
ylabel（'y（t）'）
title（'generalsolution'）
subplot（2,1,2）;
plot（t,s1）;
xlabel（'t'）
ylabel（'y（t）'）
title（'specialsolution'）
Figure1-1
（2）求通解随w变化的规律
.w在（0,1）之间的全微分方程通解
clc
clearall
s1=dsolve（'D2y+3*Dy+2*y=0','y（0）=2,Dy（0）=0','t'）;
l2=dsolve（'D2y+3*Dy+2*y=sin（0.05*t）','y（0）=2,Dy（0）=0','t'）;
l3=dsolve（'D2y+3*Dy+2*y=sin（0.15*t）','y（0）=2,Dy（0）=0','t'）;
l4=dsolve（'D2y+3*Dy+2*y=sin（0.25*t）','y（0）=2,Dy（0）=0','t'）;
l5=dsolve（'D2y+3*Dy+2*y=sin（0.5*t）','y（0）=2,Dy（0）=0','t'）;
l6=dsolve（'D2y+3*Dy+2*y=sin（0.75*t）','y（0）=2,Dy（0）=0','t'）;
t=1:
0.1:
10;
s1_n=eval（s1）;
l2_n=eval（l2）;
l3_n=eval（l3）;
l4_n=eval（l4）;
l5_n=eval（l5）;
l6_n=eval（l6）;
holdon
plot（t,s1_n）;
plot（t,l2_n,'m*'）;
plot（t,l3_n,'rx'）;
plot（t,l4_n,'g^'）;
plot（t,l5_n,'bp'）;
plot（t,l6_n,'ko'）;
holdoff
.w在（1,+
）之间的全微分方程通解
t=1:
0.1:
10;
s1=-2*exp（-2*t）+4*exp（-t）;
s2=-3/10*cos（t）+1/10*sin（t）-11/5*exp（-2*t）+9/2*exp（-t）;
s3=-3/20*cos（2*t）-1/20*sin（2*t）-9/4*exp（-2*t）+22/5*exp（-t）;
s4=-15/754*cos（5*t）-23/754*sin（5*t）-63/29*exp（-2*t）+109/26*exp（-t）;
s5=693/170*exp（-t）-39/29410*cos（13*t）-359/173*exp（-2*t）-167/29410*sin（13*t）;
s6=-1283/629*exp（-2*t）-75/393754*cos（25*t）+2529/626*exp（-t）-623/393754*sin（25*t）;
holdon
plot（t,s1）;
plot（t,s2,'m*'）;
plot（t,s3,'rx'）;
plot（t,s4,'g^'）;
plot（t,s5,'bp'）;
plot（t,s6,'ko'）;
holdoff
结论：
在b^2-4ac>0的情况下，特解的形式是C1*sint+C2*cost,齐次解的形式是
C1*EXP（R1*t）+C2*EXP（R2*t）.若w为正值且随w的增大，通解的形式趋近于齐次解。
5.b^2-4ac=0的情况
STEP1：
求解析解
S1=dsolve（'D2y+y=0','y（0）=2,Dy（0）=1','t'）
S2=dsolve（'D2y+y=sin（t）','y（0）=2,Dy（0）=1','t'）
s3=dsolve（'D2y+y=sin（2*t）','y（0）=2,Dy（0）=1','t'）
s4=dsolve（'D2y+y=sin（6*t）','y（0）=2,Dy（0）=1','t'）
s5=dsolve（'D2y+y=sin（10*t）','y（0）=2,Dy（0）=1','t'）
s6=dsolve（'D2y+y=sin（100*t）','y（0）=2,Dy（0）=1','t'）
l0=dsolve（'D2y+y=sin（0.05*t）','y（0）=2,Dy（0）=1','t'）
l1=dsolve（'D2y+y=sin（0.15*t）','y（0）=2,Dy（0）=1','t'）
l2=dsolve（'D2y+y=sin（0.25*t）','y（0）=2,Dy（0）=1','t'）
l3=dsolve（'D2y+y=sin（0.5*t）','y（0）=2,Dy（0）=1','t'）
l4=dsolve（'D2y+y=sin（0.75*t）','y（0）=2,Dy（0）=1','t'）
STEP2：
绘制图形
（1）求w=1情况下的通解和齐次解
t=1:
0.1:
10;
s1=2*exp（-2*t）+5*exp（-2*t）.*t;
s2=54/25*exp（-2*t）+26/5*exp（-2*t）.*t-4/25*cos（t）+3/25*sin（t）;
subplot（2,1,1）;
plot（t,s1）;
xlabel（'t'）
ylabel（'y（t）'）
title（'homogenioussolution'）
subplot（2,1,2）;
plot（t,s2）;
xlabel（'t'）
ylabel（'y（t）'）
title（'generalsolution'）
（2）
（2）求通解随w变化的规律
.w在（0,1）之间的全微分方程通解
t=1:
0.1:
10;
s1=2*exp（-2*t）+5*exp（-2*t）.*t;
l2=5158402/2563201*exp（-2*t）+8025/1601*exp（-2*t）.*t-32000/2563201*cos（1/20*t）+639600/2563201*sin（1/20*t）;
l3=5273762/2588881*exp（-2*t）+8105/1609*exp（-2*t）.*t-96000/2588881*cos（3/20*t）+636400/2588881*sin（3/20*t）;
l4=8706/4225*exp（-2*t）+329/65*exp（-2*t）.*t-256/4225*cos（1/4*t）+1008/4225*sin（1/4*t）;
l5=610/289*exp（-2*t）+87/17*exp（-2*t）.*t-32/289*cos（1/2*t）+60/289*sin（1/2*t）;
l6=11426/5329*exp（-2*t）+377/73*exp（-2*t）.*t-768/5329*cos（3/4*t）+880/5329*sin（3/4*t）;
holdon
plot（t,s1）;
plot（t,l2,'m*'）;
plot（t,l3,'rx'）;
plot（t,l4,'g^'）;
plot（t,l5,'bp'）;
plot（t,l6,'ko'）;
holdoff
.w在（1,+
）之间的全微分方程通解
t=1:
0.1:
10;
s1=2*exp（-2*t）+5*exp（-2*t）.*t;
s2=54/25*exp（-2*t）+26/5*exp（-2*t）.*t-4/25*cos（t）+3/25*sin（t）;
s3=3522/1681*exp（-2*t）+215/41*exp（-2*t）.*t-160/1681*cos（5/2*t）-36/1681*sin（5/2*t）;
s4=350/169*exp（-2*t）+68/13*exp（-2*t）.*t-12/169*cos（3*t）-5/169*sin（3*t）;
s5=1702/841*exp（-2*t）+150/29*exp（-2*t）.*t-20/841*cos（5*t）-21/841*sin（5*t）;
s6=104942/52441*exp（-2*t）+1160/229*exp（-2*t）.*t-60/52441*cos（15*t）-221/52441*sin（15*t）;
holdon
plot（t,s1）;
plot（t,s2,'m*'）;
plot（t,s3,'rx'）;
plot（t,s4,'g^'）;
plot（t,s5,'bp'）;
plot（t,s6,'ko'）;
holdoff
结论：
W属于（0,1）时，随w的增大在齐次解的旁边波动；w属于（1，+
），随w的增大逐渐趋近于齐次解。
4.b^2-4ac<0的情况
1.[b>0]
s2=dsolve（'D2y+Dy+y=sin（t）','y（0）=2,Dy（0）=1','t'）
s3=dsolve（'D2y+Dy+y=sin（2*t）','y（0）=2,Dy（0）=1','t'）
s4=dsolve（'D2y+Dy+y=sin（2.5*t）','y（0）=2,Dy（0）=1','t'）
s5=dsolve（'D2y+Dy+y=sin（3*t）','y（0）=2,Dy（0）=1','t'）
s6=dsolve（'D2y+Dy+y=sin（3.5*t）','y（0）=2,Dy（0）=1','t'）
s7=dsolve（'D2y+Dy+y=sin（5*t）','y（0）=2,Dy（0）=1','t'）
.w在（0,1）之间的全微分方程通解
.w在（1,+
）之间的全微分方程通解
2.[b<0]
.w在（0,1）之间的全微分方程通解
.w在（1,+
）之间的全微分方程通解
2.[b=0]
STEP1：
求解析解
S1=dsolve（'D2y+y=0','y（0）=2,Dy（0）=1','t'）
S2=dsolve（'D2y+y=sin（t）','y（0）=2,Dy（0）=1','t'）
s3=dsolve（'D2y+y=sin（2*t）','y（0）=2,Dy（0）=1','t'）
s4=dsolve（'D2y+y=sin（6*t）','y（0）=2,Dy（0）=1','t'）
s5=dsolve（'D2y+y=sin（10*t）','y（0）=2,Dy（0）=1','t'）
s6=dsolve（'D2y+y=sin（100*t）','y（0）=2,Dy（0）=1','t'）
l0=dsolve（'D2y+y=sin（0.05*t）','y（0）=2,Dy（0）=1','t'）
l1=dsolve（'D2y+y=sin（0.15*t）','y（0）=2,Dy（0）=1','t'）
l2=dsolve（'D2y+y=sin（0.25*t）','y（0）=2,Dy（0）=1','t'）
l3=dsolve（'D2y+y=sin（0.5*t）','y（0）=2,Dy（0）=1','t'）
l4=dsolve（'D2y+y=sin（0.75*t）','y（0）=2,Dy（0）=1','t'）
STEP2：
绘制图形
（3）求w=1情况下的通解和齐次解
t=1:
0.1:
10;
s1=sin（t）+2*cos（t）;
s2=3/2*sin（t）+2*cos（t）-1/2.*cos（t）.*t;
subplot（2,1,1）;
plot（t,s1）;
xlabel（'t'）
ylabel（'y（t）'）
title（'homogenioussolution'）
subplot（2,1,2）;
plot（t,s2）;
xlabel（'t'）
ylabel（'y（t）'）
title（'generalsolution'）
（2）求通解随w变化的规律
.w在（0,1）之间的全微分方程通解
t=1:
0.1:
10;
l0=379/399*sin（t）+2*cos（t）+400/399*sin（1/20*t）;
l1=331/391*sin（t）+2*cos（t）+400/391*sin（3/20*t）;
l2=11/15*sin（t）+2*cos（t）+16/15*sin（1/4*t）;
l3=1/3*sin（t）+2*cos（t）+4/3*sin（1/2*t）;
l4=-5/7*sin（t）+2*cos（t）+16/7*sin（3/4*t）;
s1=sin（t）+2*cos（t）;
holdon
plot（t,s1）;
plot（t,l0,'m*'）;
plot（t,l1,'rx'）;
plot（t,l2,'g^'）;
plot（t,l3,'bp'）;
plot（t,l4,'ko'）;
holdoff
.w在（1,+
）之间的全微分方程通解
t=1:
0.1:
10;
l2=11/15*sin（t）+2*cos（t）+16/15*sin（1/4*t）;
l3=1/3*sin（t）+2*cos（t）+4/3*sin（1/2*t）;
s1=sin（t）+2*cos（t）;
%s2=3/2*sin（t）+2*cos（t）-1/2*cos（t）*t;
s3=5/3*sin（t）+2*cos（t）-1/3*sin（2*t）;
s4=41/35*sin（t）+2*cos（t）-1/35*sin（6*t）;
holdon
plot（t,s1）;
plot（t,s3,'rx'）;
plot（t,s4,'g^'）;
plot（t,l2,'bp'）;
plot（t,l3,'ko'）;
holdoff
结论：
w=1是特殊情况，s2=3/2*sin（t）+2*cos（t）-1/2*cos（t）*t（见figure）；W属于（0,1）时，随w的增大在齐次解的旁边波动；w属于（1，+
），随w的增大逐渐趋近于齐次解。
升级会员 