第四讲 bvp4c函数应用

第四讲bvp4c函数的应用,bvp4c针对两点边界值问题进行数值求解，实现对以下n个一阶常微分方程的数值解：x的取值范围为axb区间内，初始条件,bvp4c调用格式,FunctionName函数的接口,functiondydx=FunctionName(x,y,p1,p2,.)其中x是一标量,y为yj列向量,p1,p2,.等为fj的参数，参数值已知，输出dydx是列向量，对应于fj列向量,BCFunction函数的创建,functionRes=BCFunction(ya,yb,p1,p2,.)ya是表示yj(a)的列向量,yb是表示yj(b)的列向量，即使不要求边界条件，已知参数p1,p2等也必须出现在接口定义语句中，输出Res为列向量。,调用bvp4c函数所涉及到的另外两个函数,变量solinit为一结构？，可由函数bvpinit得到：solinit=bvpinit(x,y)，x向量为初始网格点的估计值，向量y为每个yj的估计值。向量x和y的长度互不相关。bvp4c的输出sol是一结构，为特定数量点对应的解。为使曲线变得更光滑，需要在中间插入一些点，使用：sxint=deval(sol,xint),xint为点向量，函数deval根据这些点向量求解。sol为函数bvp4c的输出。,bvp4c函数应用举例,求解下面这个单变量二阶常微分方程将上式重写为一个一阶方程组：,相关函数的创建,下面创建函数和子函数，实现一阶常微分方程的子函数为OdeBvp，记录边界条件的函数为OdeBC,初始估计值由bvpinit给出，在0和1之间选择5个点，假定y1有一个解-0.05，y2有一个解0.1,设k=2,则程序为:functionbvpExamplek=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,:)functiondydx=OdeBvp(x,y,k)dydx=y(2);x-k*y(1);functionres=OdeBC(ya,yb,k)res=ya(1);yb(1);,求解结果显示,solinit和exmpsol结构,solinit=x:00.250000000000000.500000000000000.750000000000001y:2x5doubleexmpsol=x:1x7doubley:2x7doubleyp:2x7doublesolver:bvp4c0,边值问题举例工程梁,材料力学解法回顾,根据材料力学中关于轴力，剪力和弯矩符号的规定：轴力，拉力为正；剪力，绕隔离体顺时针方向转动者为正；弯矩，使梁的下侧纤维受拉者为正。画出该悬臂梁的剪力图和弯矩图。,求解过程,用积分法，可推出梁的斜率（转角）和挠度表达式为,工程梁的微分控制方程推导,工程梁的无量纲方程,如果定截面梁横向位移为w(x)，长为L，弹性模量是E，横截面的惯性矩是I，单位长度载荷是p0q(x),则无量纲方程如下：,梁的微分控制方程转化为一阶常微分方程组,用bvp4c函数求解受均布载荷的悬臂梁,悬臂梁的bvp4c函数求解,functionCantileverBeamsolinit=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,:)functiondydx=BeamODE(x,y)dydx=y(2);y(3);y(4);1;functionbc=BeamCantileverBC(y0,y1)bc=y0(1);y0(2);y1(3);y1(4);,悬臂梁的bvp4c函数求解结果,课堂练习完成简支梁(一端施加一额外力矩Mr=0.8）的bvp4c函数程序,functionSimpleBeamMr=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,:)functiondydx=BeamODE(x,y,Mr)dydx=y(2);y(3);y(4);1;functionbc=SimpleBeamBC(y0,y1,Mr)bc=y0(1);y0(3);y1(1);y1(3)-Mr;,简支梁(一端施加一额外力矩Mr=0.8）的bvp4c求解结果,思考：用ode45函数和fsolve函数能用来求解梁的问题吗？,对于梁而言，本质是一个边值问题，梁的两端各为两个边界条件，ode45只能求解初值问题，但对于梁的任一端，可以先假设另两个初值已知，这样用ode45函数求可以求解控制微