第四讲 bvp4c函数应用

单击此处编辑母版标题样式,单击此处编辑母版文本样式,第二级,第三级,第四级,第五级,*,第四讲 bvp4c函数的应用,bvp4c针对,两点边界值问题,进行数值求解，实现对以下n个一阶常微分方程的数值解：,x的取值范围为axb区间内，初始条件,bvp4c调用格式,FunctionName函数的接口,function dydx=FunctionName(x,y,p1,p2,.),其中x是一标量,y为 y,j,列向量, p1,p2,.等为f,j,的参数，参数值已知，输出dydx是列向量，对应于f,j,列向量,BCFunction函数的创建,function Res=BCFunction (ya,yb,p1,p2,.),ya是表示 y,j,(a)的列向量,yb是表示 y,j,(b)的列向量，即使不要求边界条件，已知参数p1,p2等也必须出现在接口定义语句中，输出Res为列向量。,调用bvp4c函数所涉及到的另外两个函数,变量solinit为一结构？,，可由函数,bvpinit,得到：solinit=bvpinit(x,y)，x向量为初始网格点的估计值，向量y为每个y,j,的估计值。向量x和y的长度互不相关。,bvp4c的输出sol是一结构，为特定数量点对应的解。为使曲线变得更光滑，需要在中间插入一些点，使用：sxint=,deval,(sol,xint),xint为点向量，函数deval根据这些点向量求解。sol为函数bvp4c的输出。,bvp4c函数应用举例,求解下面这个单变量二阶常微分方程,将上式重写为一个一阶方程组：,相关函数的创建,下面创建函数和子函数，实现一阶常微分方程的子函数为OdeBvp，记录边界条件的函数为OdeBC,初始估计值由bvpinit给出，在0和1之间选择5个点，假定y1有一个解-0.05，y,2,有一个解0.1,设k=2,则程序为:,function bvpExample,k=2;,solinit=bvpinit(linspace(0,1,5),-0.05,0.1),exmpsol=bvp4c(OdeBvp,OdeBC,solinit,k);,x=linspace(0,1,50);,y=deval(exmpsol,x);,plot(x,y(1,:),function dydx=OdeBvp(x,y,k),dydx=y(2);x-k*y(1);,function res=OdeBC(ya,yb,k),res=ya(1);yb(1);,求解结果显示,solinit和exmpsol结构,solinit =,x: 0 0.25000000000000 0.50000000000000 0.75000000000000 1,y: 2x5 double,exmpsol =,x: 1x7 double,y: 2x7 double,yp: 2x7 double,solver: bvp4c,0,边值问题举例工程梁,材料力学解法回顾,根据材料力学中关于轴力，剪力和弯矩符号的规定：轴力，拉力为正；剪力，绕隔离体顺时针方向转动者为正；弯矩，使梁的下侧纤维受拉者为正。画出该悬臂梁的剪力图和弯矩图,。,求解过程,用积分法，可推出梁的斜率（转角）和挠度表达式为,工程梁的微分控制方程推导,F,P,F,N,+,d,F,N,F,N,F,Q,+,d,F,Q,F,Q,M,M+,d,M,dx,d,x,q,工程梁的无量纲方程,如果定截面梁横向位移为w(x)，长为L，弹性模量是E，横截面的惯性矩是I，单位长度载荷是 p,0,q(x),则无量纲方程如下：,梁的微分控制方程转化为一阶常微分方程组,用bvp4c函数求解受均布载荷的悬臂梁,悬臂梁的bvp4c函数求解,function CantileverBeam,solinit=bvpinit(linspace(0,1,10),0.5,0.5,0.5,0.5);,beamsol=bvp4c(BeamODE, BeamCantileverBC,solinit,);,eta=linspace(0,1,50);,y=deval(beamsol,eta);,plot(eta,y(1,:),function dydx = BeamODE(x, y),dydx = y(2);y(3);y(4);1;,function bc=BeamCantileverBC(y0,y1),bc = y0(1);y0(2);y1(3);y1(4);,悬臂梁的bvp4c函数求解结果,课堂练习完成简支梁,(一端施加一额外力矩Mr=0.8）,的bvp4c函数程序,function SimpleBeam,Mr=0.8;,solinit=bvpinit(linspace(0,1,10),0.5,0.5,0.5,0.5);,beamsol=bvp4c(BeamODE, SimpleBeamBC,solinit,Mr,);,eta=linspace(0,1,50);,y=deval(beamsol,eta);,plot(eta,y(1,:),function dydx = BeamODE(x, y,Mr,),dydx = y(2);y(3);y(4);1;,function bc=SimpleBeamBC(y0,y1,Mr,),bc = y0(1);y0(3);y1(1);y1(3,)-Mr,;,简支梁,(一端施加一额外力矩Mr=0.8）,的bvp4c求解结果,思考：用ode45函数和fsolve函数能用来求解梁的问题吗？,对于梁而言，本质是一个边值问题，梁的两端各为两个边界条件，ode45只能求解初值问题，但对于梁的任一端，可以先假设另两个初值已知，这样用ode45函数求可以求解控制微分方程，但结果中包含两个未知量，然后再利用另一端的已知边界条件，这样便建立了一个非线性方程组，可用fsolve函数求得两个未知量。然后再用ode45求解。便可得解。,结论：因此用ode45和fsolve联立是可以求解梁的问题！,