学习笔记：常微分方程和常微分方程组的数值解

为什么要求数值解？

有许多常微分方程从理论上讲有解析解，但无法求出，此时只能寻求数值解。且，数值解法必须提供初值。
基本思想：


MATLAB求解标准形式：
M(x，y)*y’=F(x，y）

例1：y’=2x,初值y(0)=10,积分区间[0,10]

:

// An matlab block
f =@(x,y)2*x     %为什么要写要y？因为匿名函数必须同时接受两个输入（x,y）,即使一个输入未使用也是如此
tspan = [0,10]
y0=10;
ode45(f,tspan,y0)

结果：

例2：y’‘=u(1-y**2)y’+y,初值y(0)=1,y’(0)=0,积分区间[0,20]

首先转化为MATLAB标准求解格式：令y=y,z=y’
则转化为(1)y’=z (2)z=u(1-y**2)*z+y

:

// An matlab block
u=1;
f =@(x,y)fun(y,u);
tspan = [0,20];
y0=[1 0];
[x,y]=ode45(f,tspan,y0);
plot(x,y);

// An matlab block
  function ydot = fun(y,u)
  ydot =[y(2);u*(1-y(1)^2)*y(2)+y(1)];
  end

结果：

x,y形成的图像（蓝色）；x,z形成的图像（红色）：

例3：方程组：(1) y’=a-(b+1)y+(y2)z , (2) z’=by-(y2)*z 初值：y(0)=3,z(0)=4,积分区间[0,10]

:

// An matlab block
  a=100;
b=50;
f =@(x,y)fun(y,a,b);
tspan = [0,10];
y0=[3 4];
[x,y]=ode45(f,tspan,y0);
plot(x,y(:,1),"r",x,y(:,2),"g")

// An matlab block
  function ydot = fun(y,a,b)
  ydot =[a-(b+1)*y(1)+y(1)^2*y(2);b*y(1)-y(1)^2*y(2)];
end

结果：