MATLAB解微分方程.ppt

MATLAB

ODE初值问题的数值解
PDE问题的数值解问题提出
倒葫芦形状容器壁上的刻度问题.对于如图所示圆柱形状容器壁上的容积刻度,可以利用圆柱体体积公式x只要求解上述方程,就可求出体积V与高度x之间的函数关系,从而可标出容器壁上容积的刻度,但问题是函数D(x)无解析表达式,无法求出其解析解.
包含自变量、未知函数及未知函数的导数或微分的方程称为微分方程。
在微分方程中,自变量的个数只有一个,称为常微分方程。
自变量的个数为两个或两个以上的微分方程叫偏微分方程。
微分方程中出现的未知函数最高阶导数的阶数称为微分方程的阶数。
常微分方程：但能求解析解的常微分方程是有限的，大多数的常微分方程是给不出解析解的.8事实上，从实际问题当中抽象出来的微分方程，通常主要依靠数值解法来解决。微分方程数值方法的基本思想
对常微分方程初值问题的数值解法，就是要算出精确解y(x)在区间a,b上的一系列离散节点

处的函数值

数值解法需要把连续性的问题加以离散化，从而求出离散节点的数值解。描述这类算法，要求给出用已知信息
计算的递推公式。

建立这类递推公式的基本方法是在这些节点上用数值积分、数值微分、泰勒展开等离散化方法，对初值问题

数值解和精确解单步法和多步法欧拉（Euler）方法同理，在x=xn处，用差商代替导数：例:用Euler方法求解常微分方程初值问题h=0.2;y(1)=0.2;x=0.2:h:3;
forn=1:14
xn=x(n);yn=y(n);
y(n+1)=yn+h*(yn/xn-2*yn*yn);
end
x0=0.2:h:3;y0=x0./(1+x0.^2);
plot(x0,y0,x,y,x,y,'o') xn y(xn) yn yn-y(xn)
0.0 0 0 0
0.2 0.1923 0.2000 0.0077
0.4 0.3448 0.3840 0.0392
0.6 0.4412 0.5170 0.0758
0.8 0.4878 0.5824 0.0946
1.0 0.5000 0.5924 0.0924
1.2 0.4918 0.5705 0.0787
1.4 0.4730 0.5354 0.0624xn y(xn) yn yn-y(xn)
1.6 0.4494 0.4972 0.0478
1.8 0.4245 0.4605 0.0359
2.0 0.4000 0.4268 0.0268
2.2 0.3767 0.3966 0.0199
2.4 0.3550 0.3698 0.0147
2.6 0.3351 0.3459 0.0108
2.8 0.3167 0.3246 0.0079
3.0 0.3000 0.3057 0.0057二阶Runge-Kutta方法即例蛇形曲线的初值问题x0=0;y0=0;h=.2;
x=.2:h:3;
k1=1;
k2=(y0+h*k1)/x(1)-2*(y0+h*k1)^2;
y(1)=y0+.5*h*(k1+k2);常用的一个公式为functionydot=harmonic(t,y)
ydot=[y(2);-y(1)]functionydot=twobody(t,y)
r=sqrt(y(1)^2+y(2)^2);
ydot=[y(3);y(4);-y(1)/r^3;-y(2)/r^3];
LinearizedDifferentialEquationsJ的特征值是
基于龙格－库塔法，MATLAB求常微分方程数值解的函数，一般调用格式为：
[t,y]=ode23('fname',tspan,y0)
[t,y]=ode45('fname',tspan,y0)
其中fname是定义f(t,y)的函数文件名，该函数文件必须返回一个列向量。tspan形式为[t0,tf],表示求解区间。y0是初始状态列向量。t和y分别给出时间向量和相应的状态向量。ode23:Bogacki,Shampine(1989)和Shampine(1994),”23”表示用两算法：一个2阶，一个3阶F=inline('[y(2);-y(1)]','t','y')
ode23(F,[02*pi],[1;0])ode23(@twobody,[02*pi],[1;0;0;1]);y0=[1;0;0;3];ode23(@twobody,[02*pi],y0);y0=[1;0;0;3];[t,y]=ode23(@twobody,[02*pi],y0);
plot(y(:,1),y(:,2));axisequaly0=[1;0;0;3];[t,y]=ode23(@twobody,[02*pi],y0);
plot(y(:,1))plot(y(:,2))Ap