常微分組邊值問題打靶法matlab解法

《常微分組邊值問題打靶法matlab解法》由會員上傳分享，免費在線閱讀，更多相關(guān)內(nèi)容在行業(yè)資料-天天文庫。

1、常微分方程組兩點邊值問題的數(shù)值解法----張亞苗2011年9月可化為微分方程組方法一：配置法Matlab程序：functionbvcollationclcsolinit=bvpinit(linspace(0,1,20),[100600]);%sol=bvp4c(@twoode,@twobc,solinit);x=linspace(0,1,20);y=deval(sol,x);y'plot(x,y(1,:),x,y(2,:));end%微分方程組functiondydx=twoode(x,y)dydx=[y(2)4*y(1)];end%邊值條件functionres=twobc(ya,yb)r

2、es=[ya(1)-1yb(1)-3];end運行結(jié)果：1.0000-0.42030.9834-0.21170.9777-0.00550.98280.20070.99880.40911.02590.62201.06440.84191.11471.07101.17741.31211.25311.56771.34271.84071.44722.13411.56782.45121.70572.79541.86263.17072.04013.58112.24024.03132.46524.52612.71755.07123.00005.6724方法二：打靶法程序：積分用ode45,搜索用二分法fun

3、ctionshoot_first_tryclcclearalla=0;b=1;a2=-1;b2=0;s=0;tol=1.0*10e-6;h1=Fun(a2)h2=Fun(b2)s=bisect(@Fun,a2,b2,tol)[t1,y1]=ode45(@f,[a,b],[1,s])plot(t1,y1(:,1),t1,y1(:,2))functionxc=bisect(f,a,b,tol)ifsign(f(a))*sign(f(b))>=0error('f(a)*f(b)<0notsatisfied!')endfa=f(a);fb=f(b);k=0;while(b-a)/2>tolc=(a+

4、b)/2;fc=f(c);iffc==0breakendifsign(f(a))*sign(f(b))<0b=c;fb=fc;elsea=c;fa=fc;endendxc=(a+b)/2;functionz=Fun(s)a=0;b=1;yb=3;h=linspace(a,b,20);y0=[1;s];[t,y]=ode45(@f,h,y0);z=y(end,1)-yb;functionydot=f(t,y)%ydot=[0;0];ydot(1)=y(2);ydot(2)=4*y(1);ydot=[ydot(1);ydot(2)];運行結(jié)果：1.0000-0.50000.9969-0.4749

5、0.9940-0.44990.9913-0.42500.9887-0.40010.9799-0.30170.9736-0.20410.9697-0.10690.9683-0.01000.96920.08680.97260.18390.97840.28140.98670.37970.99740.47881.01060.57921.02640.68111.04470.78461.06560.89011.08920.99781.11551.10801.14461.22101.17661.33701.21151.45641.24951.57941.29051.70641.33481.83771.38

6、251.97351.43352.11431.48822.26031.54662.41201.60892.56981.67512.73391.74562.90491.82043.08321.89983.26921.98403.46332.07313.66612.16743.87812.26714.09982.37244.33172.48374.57452.57114.76372.66214.95962.75695.16252.85555.3726