Matlab求解邊值問題方法+例題

1、1. 把待解的問題轉(zhuǎn)化為標(biāo)準(zhǔn)邊值問題把待解的問題轉(zhuǎn)化為標(biāo)準(zhǔn)邊值問題 2. 因?yàn)檫呏祮栴}可以多解，所以需要為期望解指定一個初始因?yàn)檫呏祮栴}可以多解，所以需要為期望解指定一個初始猜測解。該猜測解網(wǎng)（猜測解。該猜測解網(wǎng)（Mesh）包括區(qū)間）包括區(qū)間a, b內(nèi)的一組網(wǎng)內(nèi)的一組網(wǎng)點(diǎn)（點(diǎn)（Mesh points）和網(wǎng)點(diǎn)上的解）和網(wǎng)點(diǎn)上的解S(x)3. 根據(jù)原微分方程構(gòu)造殘差函數(shù)根據(jù)原微分方程構(gòu)造殘差函數(shù)4. 利用原微分方程和邊界條件，借助迭代不斷產(chǎn)生新的利用原微分方程和邊界條件，借助迭代不斷產(chǎn)生新的S(x)，使殘差不斷減小，從而獲得滿足精度要求的解使殘差不斷減小，從而獲得滿足精度要求的解 ( , )( (

2、 ), ( )0 yf x yg y ay b( )( )( , ( )r xS xf x S xMatlab求解邊值問題方法：求解邊值問題方法：bvp4c函數(shù)函數(shù)solinit=bvpinit(x,v,parameters)生成生成bvp4c調(diào)用指令所必須的調(diào)用指令所必須的“解猜測網(wǎng)解猜測網(wǎng)”sol=bvp4c(odefun,bcfun,solinit,options,p1,p2,)給出微分方程邊值問題的近似解給出微分方程邊值問題的近似解sxint=deval(sol，xint)計算微分方程積分區(qū)間內(nèi)任何一點(diǎn)的解值計算微分方程積分區(qū)間內(nèi)任何一點(diǎn)的解值Matlab求解邊值問題求解邊值問題的基本

3、指令的基本指令solinit=bvpinit(x,v,parameters) x指定邊界區(qū)間指定邊界區(qū)間a,b上的初始網(wǎng)絡(luò)，通常是等距排列的（上的初始網(wǎng)絡(luò)，通常是等距排列的（1M）一維數(shù)組。）一維數(shù)組。注意：使注意：使x(1)=a，x(end)=b；格點(diǎn)要單調(diào)排列。；格點(diǎn)要單調(diào)排列。 v是對解的初始猜測是對解的初始猜測 solinit（可以取別的任意名）是（可以取別的任意名）是“解猜測網(wǎng)（解猜測網(wǎng)（Mesh）”。它是一個結(jié)構(gòu)體，帶如下兩個域：它是一個結(jié)構(gòu)體，帶如下兩個域：solinit.x是表示初始網(wǎng)格有序節(jié)點(diǎn)的（是表示初始網(wǎng)格有序節(jié)點(diǎn)的（1M）一維數(shù)組，并且）一維數(shù)組，并且solinit.x

4、(1)一定是一定是a，solinit.x(end)一定是一定是b。M不宜取得太大，不宜取得太大，10數(shù)量級左右即可。數(shù)量級左右即可。solinit.y是表示網(wǎng)點(diǎn)上微分方程解的猜測值的（是表示網(wǎng)點(diǎn)上微分方程解的猜測值的（NM）二維數(shù)組。）二維數(shù)組。solinit.y(:,i)表示節(jié)點(diǎn)表示節(jié)點(diǎn)solinit.x(i)處的解的猜測值。處的解的猜測值。初始解生成函數(shù)：初始解生成函數(shù)：bvpinit()sol=bvp4c(odefun,bcfun,solinit,options,p1,p2,)輸入?yún)?shù)：輸入?yún)?shù)： odefun是計算導(dǎo)數(shù)的是計算導(dǎo)數(shù)的m函數(shù)文件。該函數(shù)的基本形式為：函數(shù)文件。該函數(shù)的基本

5、形式為：dydx=odefun(x,y,parameters,p1,p2,)，在此，自變量，在此，自變量x是標(biāo)量，是標(biāo)量，y，dydx是列向量。是列向量。 bcfun是計算邊界條件下殘數(shù)的是計算邊界條件下殘數(shù)的m函數(shù)文件。其基本形式為：函數(shù)文件。其基本形式為：res=bcfun(ya,yb,parameters,p1,p2,)，文件輸入宗量，文件輸入宗量ya,yb是邊界是邊界條件列向量，分別代表?xiàng)l件列向量，分別代表y在在a和和b處的值。處的值。res是邊界條件滿足時的殘數(shù)是邊界條件滿足時的殘數(shù)列向量。列向量。注意：例如注意：例如odefun函數(shù)的輸入宗量中包含若干函數(shù)的輸入宗量中包含若干“未知

6、未知”和和“已知已知”參數(shù)，那么不管在邊界條件計算中是否用到，它們都應(yīng)作為參數(shù)，那么不管在邊界條件計算中是否用到，它們都應(yīng)作為bcfun的輸入宗量。的輸入宗量。 輸入宗量輸入宗量options是用來改變是用來改變bvp4c算法的控制參數(shù)的。在最基本用法算法的控制參數(shù)的。在最基本用法中，它可以缺省，此時一般可以獲得比較滿意的邊值問題解。如需更中，它可以缺省，此時一般可以獲得比較滿意的邊值問題解。如需更改可采用改可采用bvpset函數(shù)，使用方法同函數(shù)，使用方法同odeset函數(shù)。函數(shù)。 輸入宗量輸入宗量p1,p2等表示希望向被解微分方程傳遞的已知參數(shù)。如果無須等表示希望向被解微分方程傳遞的已知參數(shù)

7、。如果無須向微分方程傳遞參數(shù)，它們可以缺省。向微分方程傳遞參數(shù)，它們可以缺省。邊值問題求解指令：邊值問題求解指令：bvp4c()輸出參數(shù)：輸出參數(shù)：輸出變量輸出變量sol是一個結(jié)構(gòu)體是一個結(jié)構(gòu)體usol.x是指令是指令bvp4c所采用的網(wǎng)格節(jié)點(diǎn)；所采用的網(wǎng)格節(jié)點(diǎn)；usol.y是是y(x)在在sol.x網(wǎng)點(diǎn)上的近似解值；網(wǎng)點(diǎn)上的近似解值；usol.yp是是y(x)在在sol.x網(wǎng)點(diǎn)上的近似解值；網(wǎng)點(diǎn)上的近似解值；usol.parameters是微分方程所包含的未知參數(shù)的近似解值。是微分方程所包含的未知參數(shù)的近似解值。 當(dāng)被解微分方程包含未知參數(shù)時，該域存在。當(dāng)被解微分方程包含未知參數(shù)時，該域存在

8、。邊值問題求解指令：邊值問題求解指令：bvp4c()原方程組等價于以下標(biāo)準(zhǔn)形式的原方程組等價于以下標(biāo)準(zhǔn)形式的方程組：方程組：solinit=bvpinit(linspace(0,1,10),1 0);sol=bvp4c(ODEfun,BCfun,solinit);x=0:0.05:0.5;y=deval(sol,x);xP=0:0.1:0.5;yP=0 -0.41286057 -0.729740656. -0.95385538 -1.08743325 -1.13181116;plot(xP,yP,o,x,y(1,:),r-)legend(Analytical Solution,Numerica

9、l Solution)% 定義定義ODEfun函數(shù)函數(shù)function dydx=ODEfun(x,y)dydx=y(2);y(1)+10;% 定義定義BCfun函數(shù)函數(shù)function bc=BCfun(ya,yb)bc=ya(1);yb(1);10(0)(1)0yyyy求解兩點(diǎn)邊值問題：求解兩點(diǎn)邊值問題： 121,yy yy令：令： 122110 yyyy11(0)0,(1)0yy邊界條件為：邊界條件為： 邊值問題的求解邊值問題的求解原方程組等價于以下標(biāo)準(zhǔn)形式的原方程組等價于以下標(biāo)準(zhǔn)形式的方程組：方程組：solinit=bvpinit(linspace(0,1,10),0 1);sol=b

10、vp4c(ODEfun,BCfun,solinit);x=0:0.1:1;y=deval(sol,x);xP=0:0.1:1.0;yP=1 1.0743 1.1695 1.2869 1.4284. 1.5965 1.7947 2.0274 2.3004 2.6214 3;plot(xP,yP,o,x,y(1,:),r-)legend(Analytical Solution,Numerical Solution,.location,Northwest)legend boxoff% 定義定義ODEfun函數(shù)函數(shù)function dydx=ODEfun(x,y)dydx=y(2);(1+x2)*y(

11、1)+1;% 定義定義BCfun函數(shù)函數(shù)function bc=BCfun(ya,yb)bc=ya(1)-1;yb(1)-3;求解：求解： 121,yy yy令：令： 邊界條件為：邊界條件為： 2(1)1(0)1, (1)3yxyyy12221(1)1 yyyxy03) 1 (01)0(11yy邊值問題的求解邊值問題的求解c=1;solinit=bvpinit(linspace(0,4,10),1 1);sol=bvp4c(ODEfun,BCfun,solinit,c);x=0:0.1:4;y=deval(sol,x);plot(x,y(1,:),b-,sol.x,sol.y(1,:),ro)

12、legend(解曲線解曲線,初始網(wǎng)格點(diǎn)解初始網(wǎng)格點(diǎn)解)% 定義定義ODEfun函數(shù)函數(shù)function dydx=ODEfun(x,y,c)dydx=y(2);-c*abs(y(1);% 定義定義BCfun函數(shù)函數(shù)function bc=BCfun(ya,yb,c)bc=ya(1);yb(1)+2;求解：求解： 121,zz zz令：令： 0 zc z2)4(, 0)0(zz邊值問題的求解邊值問題的求解1csolinit=bvpinit(linspace(0,pi,10),1;1,lmb);opts=bvpset(Stats,on);sol=bvp4c(ODEfun,BCfun,solinit

13、,opts);lambda=sol.parametersx=0:pi/60:pi;y=deval(sol,x);plot(x,y(1,:),b-,sol.x,sol.y(1,:),ro)legend(解曲線解曲線,初始網(wǎng)格點(diǎn)解初始網(wǎng)格點(diǎn)解)% 定義定義ODEfun函數(shù)函數(shù)function dydx=ODEfun(x,y,lmb)q=15;dydx=y(2);-(lmb-2*q*cos(2*x)*y(1);% 定義定義BCfun函數(shù)函數(shù)function bc=BCfun(ya,yb,lmb)bc=ya(1)-1;ya(2);yb(2);求解：求解： 邊界條件：邊界條件： 本例中，微分方程與參數(shù)本

14、例中，微分方程與參數(shù)的數(shù)值有關(guān)。一般而言，對于任意的的數(shù)值有關(guān)。一般而言，對于任意的值，該問題無解，值，該問題無解，但對于特殊的但對于特殊的值（特征值），它存在一個解，這也稱為微分方程的特征值問題。值（特征值），它存在一個解，這也稱為微分方程的特征值問題。對于此問題，可在對于此問題，可在bvpinit中提供參數(shù)的猜測值，然后重復(fù)求解中提供參數(shù)的猜測值，然后重復(fù)求解BVP得到所需的得到所需的參數(shù)，返回參數(shù)為參數(shù)，返回參數(shù)為sol.parameters(2 cos2 )0zqx z15q0)(, 0)0(, 1)0(zzz邊值問題的求解邊值問題的求解infinity=6;solinit=bvpin

15、it(linspace(0,infinity,5),0 0 1);options=bvpset(stats,on);sol=bvp4c(ODEfun,BCfun,solinit,options);eta=sol.x;f=sol.y;fprintf(n);fprintf(Cebeci & Keller report f(0)=0.92768.n)fprintf(Value computed here is f(0)=%7.5f.n,f(3,1)clf resetplot(eta,f(2,:);axis(0 infinity 0 1.4);title.(Falkner-Skan equation, positive wall shear,beta=0.5.)xlabel(eta),ylabel(df/deta),shg% 定義定義ODEfun函數(shù)