Matlab求解边值问题方法+例题
重要：MATLAB常微分方程(组)数值解法
Matlab常微分方程求解问题分类
边值问题：
初值问题：
• 定解附加条件在自变量 的一端
• 一般形式为： y' f (x, y)
y(a)
y0
• 初值问题的数值解法一 般采用步进法，如 Runge-Kutta法
➢ 在自变量两端均给定附加 条件
y' f (x, y)
➢ 一般形式：y(a)y1, y(b)y2
1.根据常微分方程要求的求解精度与速度要求
求解初值问题：
y
'
y
2x y
y ( 0 ) 1
(0x1)
比较ode45和ode23的求解精度和速度
ode45和ode23的比较-1
function xODE clear all clc
format long
y0 = 1; [x1,y1] = ode45(@f,[0,1],y0); [x2,y2] = ode23(@f,[0,1],y0); plot(x1,y1,'k-',x2,y2,'b--') xlabel('x') ylabel('y')
rD = k(3)*C(2)-k(5)*C(4);
rE = k(4)*C(3)+k(5)*C(4);
% Mass balances dCdt = [rA; rB; rC; rD; rE];
三个串联的CSTR等温反应器（例4-3）
function IsothermCSTRs clear all clc CA0 = 1.8; % kmol/m^3 CA10 = 0.4; % kmol/m^3 CA20 = 0.2; % kmol/m^3 CA30 = 0.1; % kmol/m^3 k = 0.5; % 1/min tau = 2; stoptime = 2.9; % min [t,y] = ode45(@Equations,[0 stoptime],[CA10 CA20 CA30],[],k,CA0,tau); disp(' Results:') disp(' t CA1 CA2 CA3') disp([t,y]) plot(t,y(:,1),'k--',t,y(:,2),'b:',t,y(:,3),'r-') legend('CA_1','CA_2','CA_3') xlabel('Time (min)') ylabel('Concentration') % -----------------------------------------------------------------function dydt = Equations(t,y,k,CA0,tau) CA1 = y(1); CA2 = y(2); CA3 = y(3); dCA1dt = (CA0-CA1)/tau - k*CA1; dCA2dt = (CA1-CA2)/tau - k*CA2; dCA3dt = (CA2-CA3)/tau - k*CA3; dydt = [dCA1dt; dCA2dt; dCA3dt];
非线性微分方程边值问题打靶算法Matlab程序
非线性微分方程边值问题打靶算法Matlab程序【非线性微分方程边值问题打靶算法】参见:// matlabsky /thread-828-1-1.html【线性微分方程边值问题打靶算法】参见 :// matlabsky /thread-827-1-1.html 【线性微分方程边值问题有限差分算法】参见:// matlabsky /thread-829-1-1.html下面我们讲解下【非线性微分方程边值问题打靶算法】对于边值问题线性边值问题时，p、q和r只是x的函数，但是非线性时它们是x,y和y'的函数我们将问题转换为如下初值问题现在我们需要做的就是找到那个m，使得y(b)=beta，换句话说y现在由两个变量最终控制，即m和x。
对于这个求m的问题，我们可以使用牛顿迭代法得到（1）先给出一个初值m=1（2）计算出对应的y(b)为yb，这个直接使用ode45计算，得到y(end,1)就是yb（3）根据beta和yb的差异更新m（4）将新的m带入（2）重新计算yb并更新m（5）如此迭代直到m稳定或者在允许误差范围内下面详细描述了Matlab代码的运行过程以下内容需要回复才能看到复制内容到剪贴板代码:function [t,x]=nlineshoot(funcn,funcv,tspan,bound,tol,varargin)%对微分方程y''=p*y'+q*y+r,a<t<b,u(a)=alpha,u(b)=beta,其中p,q和r为x,y和y'的函数%只要对应修改p,q,r并将y'用x(2)替换，y用x(1)替换，就可以了%funcn=@(t,x)[x(2);p*x(2)+q*x(1)+r];%funcv=@(t,x)[x(2);p*x(2)+q*x(1)+r;x(4);x(3)*dF/dy+x(4)*dF/y'];%%% Example%F=y''=2y*y',y(0)=-1,y(pi/2)=1%%由方程有p=2y=2*x(1),q=r=0%dF/dy=2*y'=2*x(2) %y'用x(2)替换，y用x(1)替换%dF/dy'=2*y=2*x(1)%%funcn=@(t,x)[x(2);2*x(1)*x(2)];%funcv=@(t,x)[x(2);2*x(1)*x(2);x(4);2*x(2)*x(3)+2*x(1)*x(4)]; %tspan=[0 pi/2];%bound=[-1 1];%tol=1e-8;%[t,y]=nlineshoot(funcn,funcv,tspan,bound,tol);%plot(t,y)%legend('x1','x2')%%by dynamic%see also :// matlabsky%2009.3.12%%lb=bound(1);ub=bound(2);m0=0;m=1;while norm(m-m0)>tolm0=m;[t,x]=ode45(funcv,tspan,[lb;m;0;1],varargin);m=m0-(x(end,1)-ub)/x(end,3);end[t,x]=ode45(funcn,tspan,[lb;m]);下面的图形是代码中附带实例的运行结果线性微分方程边值问题打靶算法Matlab程序【非线性微分方程边值问题打靶算法】参见:// matlabsky /thread-828-1-1.html 【线性微分方程边值问题打靶算法】参见:// matlabsky /thread-827-1-1.html【线性微分方程边值问题有限差分算法】参见:// matlabsky /thread-829-1-1.html注意该算法只能完成二阶常微分方程双边值问题求解，至于其他形式的边值问题必须先转换到二阶形式对于下面的二阶常微分方程利用上面方程的线性结果和两个特殊的初值问题，我们可以构造两个等效的常微分初值方程初值问题，对于初值问题我们就可以直接使用ode**计算器或者龙哥库塔算法求解了。
Matlab求解边值问题方法+例题
y1 y , y2 y1
原方程组等价于以下标准形式的 方程组：
y y 1 2 y2 (1 x 2 ) y1 1
边界条件为：
y1 (0) 1 0 y1 (1) 3 0
边值问题的求解
求解：
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') legend('解曲线','初始网格点解') % 定义ODEfun函数 function dydx=ODEfun(x,y,c) dydx=[y(2);-c*abs(y(1))]; % 定义BCfun函数 function bc=BCfun(ya,yb,c) bc=[ya(1);yb(1)+2];
z(0) 1, z (0) 0, z ( ) 0
本例中，微分方程与参数λ的数值有关。一般而言，对于任意的λ值，该问题无解， 但对于特殊的λ值（特征值），它存在一个解，这也称为微分方程的特征值问题。 对于此问题，可在bvpinit中提供参数的猜测值，然后重复求解BVP得到所需的 参数，返回参数为sol.parameters
边值问题的求解
如果取1，计算结果如何?
infinity=6; solinit=bvpinit(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'); 2 求解： f ff (1 ( f ) ) 0 fprintf('Cebeci & Keller report f''''(0)=0.92768.\n') fprintf('Value computed here is f''''(0)=%7.5f.\n',f(3,1)) 0.5 clf reset plot(eta,f(2,:)); axis([0 infinity 0 1.4]); title... f (0) 0, f (0) 0, ('Falkner-Skan equation, positive wall shear,\beta=0.5.') 边界条件： f ( ) 1, xlabel('\eta'),ylabel('df/d\eta'),shg % 定义ODEfun函数 function dfdeta=ODEfun(eta,f) beta=0.5; dfdeta=[f(2);f(3);-f(1)*f(3)-beta*(1-f(2)^2)]; % 定义BCfun函数 function res=BCfun(f0,finf) res=[f0(1);f0(2);finf(2)-1];
matlab求边值问题例题
matlab求边值问题例题摘要：一、引言二、MATLAB求边值问题的基本步骤1.建立微分方程模型2.离散化3.求解离散化后的线性方程组4.反演法求解边值问题三、实例分析1.二维热传导方程2.一维波动方程四、总结与展望正文：一、引言边值问题广泛应用于物理、工程、数学等领域，解决边值问题的一种有效方法是利用MATLAB进行数值求解。
本文将简要介绍MATLAB求边值问题的基本步骤，并通过实例分析二维热传导方程和一维波动方程的求解过程，以展示MATLAB在边值问题求解中的应用。
二、MATLAB求边值问题的基本步骤1.建立微分方程模型根据实际问题，建立相应的微分方程模型。
例如，对于二维热传导方程，可表示为：u/x = α * u/t2.离散化将空间和时间进行离散化，采用有限差分法将微分方程转化为线性方程组。
3.求解离散化后的线性方程组利用MATLAB内置的求解线性方程组的函数，如“ solve() ”，求解离散化后的线性方程组。
4.反演法求解边值问题对于具有周期性的边值问题，可以采用反演法进行求解。
具体步骤如下：1) 求解对应的周期边界条件下的一维线性方程组；2) 利用MATLAB的插值函数，如“interp1()”或“spline()”，对解进行插值；3) 通过反演公式，求解原始边值问题。
三、实例分析1.二维热传导方程考虑如下二维热传导方程的边值问题：u/x = α * u/t，u(x, 0) = f(x)，u(0, t) = g(t)，u(a, t) = u(b, t)根据基本步骤，首先建立微分方程模型，然后进行离散化，利用MATLAB 求解离散化后的线性方程组。
最后，通过反演法求解边值问题。
2.一维波动方程考虑如下边值问题：i * u/t = u/x，u(x, 0) = f(x)，u(a, t) = u(b, t)同样地，根据基本步骤，建立微分方程模型，离散化，求解离散化后的线性方程组，最后采用反演法求解边值问题。
(偏)微分方程一类边值问题的数值求解(附matlab程序)
（偏）微分方程一类边值问题的数值求解本文介绍了椭圆型（偏）微分方程一类边值问题的数值求解程序（笔者自编）及其使用方法。
程序基于差分原理，将连续的（偏）微分方程在节点上离散化，最终化为线性代数方程组。
对于原理方面的问题很多微分方程数值解的参考书中都有详尽的描述，但一般缺少具体的实现程序。
笔者认为，对这类数值方法是否理解并能运用的关键之一还是在于是否能写出用于解决问题的程序。
理论基础固然必不可少，程序确往往是我们解决问题的敲门砖。
一维椭圆型方程可表示如下：],[,0)(,0)()(b a x x q x p f qu dxdur dx du p dx d Lu ∈>>=++-=其中L 表示微分算子，很明显这是一个线性算子。
这里要求p ，q 大于零是为了保证最终得到的线性方程组有唯一的非零解，但事实上不满足这个条件可能也是有解的，这涉及到微分方程解的存在性也确定性问题，读者若有兴趣可参考相关书籍。
更具体的，我们可以举出一个一维椭圆型方程的例子来： 例（1）：)cos(10)exp(222x u x dx dux dxu d =++- 即：)cos(10),ex p()(,)(,1)(2x f x x q x x r x p ====⎩⎨⎧==3)3(1)1(u u 利用文中提及的程序，可以将以上问题表述为：syms x %定义符号变量xp=@(x) 1; r=@(x) x.^2; q=@(x) exp(x); f=@(x) 10*cos(x);Odbvp(p,r,q,f,1,3,1,3);11.2 1.4 1.6 1.82 2.2 2.4 2.6 2.83-0.50.511.522.532~Odbvp u=f(x)xu哈哈，一条非常漂亮的曲线。
若不满足p,q>0，我们也举一例： 例（2）syms xp=@(x) 10*cos(x); r=@(x) x.^2;q=@(x) 10*sin(x); f=@(x) 10*cos(x);Odbvp(p,r,q,f,1,3,1,3)1 1.2 1.4 1.6 1.82 2.2 2.4 2.6 2.83-500501001502002503003502~Odbvp u=f(x)xu可以发现，程序仍能求解，但结果的光滑性不好。
偏微分方程边值问题的数值解法
求解偏微分方程的边值问题本实验学习使用MATLAB 的图形用户命令pdetool 来求解偏微分方程的边值问题。
这个工具是用有限元方法来求解的，而且采用三角元。
我们用内个例题来说明它的用法。
一、MATLAB 支持的偏微分方程类型考虑平面有界区域D 上的二阶椭圆型PDE 边值问题：()c u u f α-∇∇+= (1.1)其中 (1) , (2) a,f D c x y ⎛⎫∂∂∇=⨯ ⎪∂∂⎝⎭是上的已知函数(3)是标量或22的函数方阵未知函数为(,) (,)u x y x y D ∈。
它的边界条件分为三类：（1）Direchlet 条件：hu f = (1.2)（2）Neumann 条件： ()n c u qu g ∇+= (1.3)（3）混合边界条件：在边界D ∂上部分为Direchlet 条件，另外部分为Neumann 条件。
其中,,,,h r q g c 是定义在边界D ∂的已知函数，另外c 也可以是一个2*2的函数矩阵，n 是沿边界的外法线的单位向量。
在使用pdetool 时要向它提供这些已知参数。
二、例题例题1 用pdetool 求解 22D 1 D: 10u x y u ∂⎧-∆=+≤⎪⎨=⎪⎩ (1.4)解 ：首先在MATLAB 的工作命令行中键入pdetool ，按回牟键确定，于是出现PDE Toolbox 窗口，选Genenic Scalar 模式．( l ）画区域圆单击椭圆工具按钮，大致在（0，0）位置单击鼠标右键，拖拉鼠标到适当位置松开。
为了保证所绘制的圆是标准的单位园，在所绘园上双击，打开 Object Dialog 对话框，精确地输入圆心坐标X-center 为0 、Y-center 为0 及半径Radius 为l ，然后单击OK 按钮，这样单位画已画好．( 2 ）设置边界条件单击工具边界模式按钮 ，图形边界变红，逐段双击边界，打开Boundary condition 对话框．输入边界条件．对于同一类型的边界，可以按Shift 键，将多个边界同时选择，统一设边界条件．本题选择Dirichlet 条件，输入h 为1 , r 为0。
Matlab第一次上机作业
输入： >>tic, n=9;[u,k]=xsj(n), toc,surf(u)
计算Байду номын сангаас果如下：
u=
0 0 0 0 0 0 0 0 0 0 0 0 0.0256 0.0413 0.0508 0.0560 0.0577 0.0560 0.0508 0.0413 0.0256 0 0 0.0413 0.0686 0.0859 0.0955 0.0986 0.0955 0.0859 0.0686 0.0413 0 0.0508 0.0859 0.1088 0.1216 0.1258 0.1216 0.1088 0.0859 0.0508 0 0 0.0560 0.0955 0.1216 0.1364 0.1412 0.1364 0.1216 0.0955 0.0560 0 0 0.0577 0.0986 0.1258 0.1412 0.1462 0.1412 0.1258 0.0986 0.0577 0 0 0.0560 0.0955 0.1216 0.1364 0.1412 0.1364 0.1216 0.0955 0.0560 0 0 0.0508 0.0859 0.1088 0.1216 0.1258 0.1216 0.1088 0.0859 0.0508 0 0 0.0413 0.0686 0.0859 0.0955 0.0986 0.0955 0.0859 0.0686 0.0413 0 0 0.0256 0.0413 0.0508 0.0560 0.0577 0.0560 0.0508 0.0413 0.0256 0 0 0 0 0 0 0 0 0 0 0 0 0
程序一： function [u,k]=xsbj(n) % xsbj:用块 Jacobi 迭代法求解线性方程组 A*u=f % u:方程组的解； k 迭代次数； n:非边界点数 % a：方程组系数矩阵 Aii 的下对角线元素；b:方程系数矩阵 Aii 的主对角线元素 % c：方程组系数矩阵 Aii 的上对角线元素；d:追赶法所求方程的右端向量 % e:允许误差界；er：迭代误差 f=2*1/(n+1)^(2)*ones(n+2,n+2); a=-1*ones(1,n);b=4*ones(1,n);c=-1*ones(1,n);u=zeros(n+2,n+2);e=0.000000001; for k=1:2000 er=0; ub=u; for j=2:n+1 d(1:n)=f(2:n+1,j)+ub(2:n+1,j-1)+ub(2:n+1,j+1); x=zg(a,b,c,d); u(2:n+1,j)=x'; er=er+norm(ub(:,j)-u(:,j),1); end if er/n^2<e,break; end end 程序二： function x=zg(a,b,c,d) % zg:用追赶法求解线性方程组 n=length(b); % LU 分解 u(1)=b(1); for k=2:n if u(k-1)==0,D=0,return; end l(k)=a(k)/u(k-1); u(k)=b(k)-l(k)*c(k-1); end % 追赶法求解之追过程，求解 Ly=d y(1)=d(1); for k=2:n y(k)=d(k)-l(k)*y(k-1); end % 追赶法求解之追过程，求解 Ux=y if u(n)==0,D=0,return;end x(n)=y(n)/u(n); for k=n-1:-1:1 x(k)=(y(k)-c(k)*x(k+1))/u(k); end 输入：
matlab求解初边值问题的偏微分方程
偏微分方程是描述自然界中动态过程的重要数学工具，在工程领域中，求解偏微分方程是很多实际问题的重要一步。
MATLAB作为一种强大的数学计算软件，提供了丰富的工具和函数来求解各种类型的偏微分方程。
本文将介绍使用MATLAB求解初边值问题的偏微分方程的方法和步骤。
一、MATLAB中的偏微分方程求解工具MATLAB提供了几种可以用来求解偏微分方程的工具和函数，主要包括：1. pdepe函数：用于求解偏微分方程初边值问题的函数，可以处理各种类型的偏微分方程，并且可以自定义边界条件和初始条件。
2. pdepeplot函数：用于绘制pdepe函数求解得到的偏微分方程的解的可视化图形，有助于直观地理解方程的解的特性。
3. pdetool工具箱：提供了一个交互式的图形用户界面，可以用来建立偏微分方程模型并进行求解，适用于一些复杂的偏微分方程求解问题。
二、使用pdepe函数求解偏微分方程初边值问题的步骤对于给定的偏微分方程初边值问题，可以按照以下步骤使用pdepe函数进行求解：1. 定义偏微分方程需要将给定的偏微分方程转化为标准形式，即将偏微分方程化为形式为c(x,t,u)∂u/∂t = x(r ∂u/∂x) + ∂(p∂u/∂x) + f(x,t,u)的形式。
2. 编写边界条件和初始条件函数根据实际问题的边界条件和初始条件，编写相应的函数来描述这些条件。
3. 设置空间网格选择合适的空间网格来离散空间变量，可以使用linspace函数来产生均匀分布的网格。
4. 调用pdepe函数求解偏微分方程将定义好的偏微分方程、边界条件和初始条件函数以及空间网格作为参数传递给pdepe函数，调用该函数求解偏微分方程。
5. 可视化结果使用pdepeplot函数绘制偏微分方程的解的可视化图形，以便对解的性质进行分析和理解。
三、实例分析考虑一维热传导方程初边值问题：∂u/∂t = ∂^2u/∂x^2, 0<x<1, 0<t<1u(0,t) = 0, u(1,t) = 0, u(x,0) = sin(πx)使用MATLAB求解该初边值问题的步骤如下：1. 定义偏微分方程将热传导方程化为标准形式，得到c(x,t,u) = 1, r = 1, p = 1, f(x,t,u) = 0。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
y1 y , y2 y1
原方程组等价于以下标准形式的 方程组：
y y 1 2 y2 (1 x 2 ) y1 1
边界条件为：
y1 (0) 1 0 y1 (1) 3 0
边值问题的求解
求解：
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') legend('解曲线','初始网格点解') % 定义ODEfun函数 function dydx=ODEfun(x,y,c) dydx=[y(2);-c*abs(y(1))]; % 定义BCfun函数 function bc=BCfun(ya,yb,c) bc=[ya(1);yb(取1，计算结果如何?
infinity=6; solinit=bvpinit(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'); 2 求解： f ff (1 ( f ) ) 0 fprintf('Cebeci & Keller report f''''(0)=0.92768.\n') fprintf('Value computed here is f''''(0)=%7.5f.\n',f(3,1)) 0.5 clf reset plot(eta,f(2,:)); axis([0 infinity 0 1.4]); title... f (0) 0, f (0) 0, ('Falkner-Skan equation, positive wall shear,\beta=0.5.') 边界条件： f ( ) 1, xlabel('\eta'),ylabel('df/d\eta'),shg % 定义ODEfun函数 function dfdeta=ODEfun(eta,f) beta=0.5; dfdeta=[f(2);f(3);-f(1)*f(3)-beta*(1-f(2)^2)]; % 定义BCfun函数 function res=BCfun(f0,finf) res=[f0(1);f0(2);finf(2)-1];
sxint=deval(sol，xint) 计算微分方程积分区间内任何一点的解值
初始解生成函数：bvpinit()
solinit=bvpinit(x,v,parameters)
x指定边界区间[a,b]上的初始网络，通常是等距排列的（1×M）一维数组。 注意：使x(1)=a，x(end)=b；格点要单调排列。 v是对解的初始猜测 solinit（可以取别的任意名）是“解猜测网（Mesh）”。 它是一个结构体，带如下两个域： solinit.x是表示初始网格有序节点的（1×M）一维数组，并且 solinit.x(1)一定是a，solinit.x(end)一定是b。M不宜取得太大，10 数量级左右即可。 solinit.y是表示网点上微分方程解的猜测值的（N×M）二维数组。 solinit.y(:,i)表示节点solinit.x(i)处的解的猜测值。
Matlab求解边值问题方法：bvp4c函数
1. 把待解的问题转化为标准边值问题
y f ( x, y) g ( y(a), y(b)) 0
2. 因为边值问题可以多解，所以需要为期望解指定一个初始
猜测解。该猜测解网（Mesh）包括区间[a, b]内的一组网 点（Mesh points）和网点上的解S(x) 3. 根据原微分方程构造残差函数 r ( x) S ( x) f ( x, S ( x)) 4. 利用原微分方程和边界条件，借助迭代不断产生新的S(x)，
边值问题求解指令：bvp4c()
输出参数： 输出变量sol是一个结构体 sol.x是指令bvp4c所采用的网格节点； sol.y是y(x)在sol.x网点上的近似解值； sol.yp是y'(x)在sol.x网点上的近似解值； sol.parameters是微分方程所包含的未知参数的近似解值。 当被解微分方程包含未知参数时，该域存在。
z c z 0 c 1
z (0) 0, z (4) 2
令：
z1 z , z2 z1
边值问题的求解
求解：z ( 2q cos 2 x) z 0
q 15
边界条件：
solinit=bvpinit(linspace(0,pi,10),[1;1],lmb); opts=bvpset('Stats','on'); sol=bvp4c(@ODEfun,@BCfun,solinit,opts); lambda=sol.parameters x=[0:pi/60:pi]; y=deval(sol,x); plot(x,y(1,:),'b-',sol.x,sol.y(1,:),'ro') legend('解曲线','初始网格点解') % 定义ODEfun函数 function dydx=ODEfun(x,y,lmb) q=15; dydx=[y(2);-(lmb-2*q*cos(2*x))*y(1)]; % 定义BCfun函数 function bc=BCfun(ya,yb,lmb) bc=[ya(1)-1;ya(2);yb(2)];
边值问题求解指令：bvp4c()
sol=bvp4c(odefun,bcfun,solinit,options,p1,p2,…) 输入参数：
odefun是计算导数的m函数文件。该函数的基本形式为： dydx=odefun(x,y,parameters,p1,p2,…)，在此，自变量x是标量，y， dydx是列向量。 bcfun是计算边界条件下残数的m函数文件。其基本形式为： res=bcfun(ya,yb,parameters,p1,p2,…)，文件输入宗量ya,yb是边界 条件列向量，分别代表y在a和b处的值。res是边界条件满足时的残数 列向量。注意：例如odefun函数的输入宗量中包含若干“未知”和 “已知”参数，那么不管在边界条件计算中是否用到，它们都应作为 bcfun的输入宗量。 输入宗量options是用来改变bvp4c算法的控制参数的。在最基本用法 中，它可以缺省，此时一般可以获得比较满意的边值问题解。如需更 改可采用bvpset函数，使用方法同odeset函数。 输入宗量p1,p2等表示希望向被解微分方程传递的已知参数。如果无须 向微分方程传递参数，它们可以缺省。
使残差不断减小，从而获得满足精度要求的解
Matlab求解边值问题的基本指令
solinit=bvpinit(x,v,parameters) 生成bvp4c调用指令所必须的“解猜测网” sol=bvp4c(odefun,bcfun,solinit,options,p1,p2,…) 给出微分方程边值问题的近似解
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','Numerical Solution') % 定义ODEfun函数 function dydx=ODEfun(x,y) dydx=[y(2);y(1)+10]; % 定义BCfun函数 function bc=BCfun(ya,yb) bc=[ya(1);yb(1)];
z(0) 1, z (0) 0, z ( ) 0
本例中，微分方程与参数λ的数值有关。一般而言，对于任意的λ值，该问题无解， 但对于特殊的λ值（特征值），它存在一个解，这也称为微分方程的特征值问题。 对于此问题，可在bvpinit中提供参数的猜测值，然后重复求解BVP得到所需的 参数，返回参数为sol.parameters
边值问题的求解
y y 10 求解两点边值问题： y (0) y(1) 0
令：
y1 y , y2 y1
原方程组等价于以下标准形式的 方程组：
y y 1 2 y2 y1 10
边界条件为：
y1 (0) 0, y1 (1) 0
边值问题的求解
y (1 x 2 ) y 1 求解： y (0) 1, y (1) 3
令：
solinit=bvpinit(linspace(0,1,10),[0 1]); sol=bvp4c(@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函数 function dydx=ODEfun(x,y) dydx=[y(2);(1+x^2)*y(1)+1]; % 定义BCfun函数 function bc=BCfun(ya,yb) bc=[ya(1)-1;yb(1)-3];