Matlab常微分方程数值解法

Posted jianle23

tags:

篇首语:本文由小常识网(cha138.com)小编为大家整理,主要介绍了Matlab常微分方程数值解法相关的知识,希望对你有一定的参考价值。

实验目的

Matlab实现欧拉法、后退欧拉法、梯形方法和改进欧拉公式

实验要求

1. 给出欧拉法、后退欧拉法、梯形方法和改进欧拉公式算法

2. Matlab实现欧拉法、后退欧拉法、梯形方法和改进欧拉公式

实验内容

技术图片

 实验步骤

  (1)欧拉法算法,

  技术图片

  MATLAB实现,

技术图片
 1 %数值解常微分方程欧拉算法
 2 %例子:dyfun=inline(y-2*x/y);[x,y]=euler2(dyfun,[0,1],1,0.2);
 3 %输入:函数dfun(x,y),求解区间xspan[x0,xN],初值y0,步长h
 4 %输出:节点x,数值解y
 5 function [x,y]=euler2(dyfun,xspan,y0,h)
 6 x=xspan(1):h:xspan(2);y(1)=y0;
 7 for n=1:length(x)-1
 8     y(n+1)=y(n)+h*feval(dyfun,x(n),y(n));
 9 end
10 x=x;y=y;
11 end
euler2

  求解【题目】,

技术图片

 (2)后退欧拉法算法,

技术图片

  MATLAB实现

技术图片
 1 %数值解常微分后退欧拉法算法
 2 %例子:dfun=inline(x+y,x,y);[x,y]=eulerh1(dfun,0,1,0.02,5)
 3 %输入:函数dfun(x,y),初值x0,y0,步长h,维度N
 4 %输出:结点x和数值解y
 5 function [x,y]=eulerh1(dfun,x0,y0,h,N)
 6 x=zeros(1,N+1);
 7 y=zeros(1,N+1);
 8 x(1)=x0;y(1)=y0;
 9 for n=1:N
10     x(n+1)=x(n)+h;
11     z0=y(n)+h*dfun(x(n),y(n));
12     for k=1:3
13         z1=y(n)+h*dfun(x(n+1),z0);
14         if abs(z1-z0)<1e-6
15             break;
16         end
17         z0=z1;
18     end
19     y(n+1)=z1;
20 end
eulerh1

  求解【题目】,

  技术图片技术图片

 (3)梯形方法算法,

 技术图片  MATLAB实现,

技术图片
 1 %数值解常微分梯形欧拉法算法
 2 %例子:dfun=inline(x+y,x,y);[x,y]=eulert(dfun,0,1,0.02,5)
 3 %输入:函数dfun(x,y),初值x0,y0,步长h,维度N
 4 %输出:结点x和数值解y
 5 function [x,y]=eulert(dfun,x0,y0,h,N)
 6 x=zeros(1,N+1);
 7 y=zeros(1,N+1);
 8 x(1)=x0;y(1)=y0;
 9 for n=1:N
10     x(n+1)=x(n)+h;
11     z0=y(n)+h*dfun(x(n),y(n));
12     for k=1:3
13         z1=y(n)+(h/2)*(dfun(x(n),y(n))+dfun(x(n+1),z0));
14         if abs(z1-z0)<1e-6
15             break;
16         end
17         z0=z1;
18     end
19     y(n+1)=z1;
20 end
eulert

   求解【题目】,

   技术图片技术图片

(2)改进欧拉公式算法,

技术图片

   MATLAB实现,

技术图片
 1 %数值解常微分改进欧拉法算法
 2 %例子:dfun=inline(x+y,x,y);[x,y]=eulerh1(dfun,0,1,0.02,5)
 3 %输入:函数dfun(x,y),初值x0,y0,步长h,维度N
 4 %输出:结点x和数值解y
 5 function [x,y]=eulerg2(dfun,x0,y0,h,N)
 6 x=zeros(1,N+1);
 7 y=zeros(1,N+1);
 8 x(1)=x0;y(1)=y0;
 9 for n=1:N
10     x(n+1)=x(n)+h;
11     ybar=y(n)+h*dfun(x(n),y(n));
12     y(n+1)=y(n)+(h/2)*(dfun(x(n),y(n))+dfun(x(n+1),ybar));
13 end
14 end
eulerg2

  求解【题目】,

  技术图片技术图片

求解结果,

  技术图片技术图片技术图片技术图片

技术图片

 小结

  就给定的题目并没有体现出这些算法之间的差异。

 

以上是关于Matlab常微分方程数值解法的主要内容,如果未能解决你的问题,请参考以下文章

MATLAB常微分方程的数值解法

如何用matlab求解常微分方程?matlab解常微分方程之符号解法介绍

MATLAB解决常微分方程

数学建模暑期集训5:matlab求解常微分方程/偏微分方程

99插值法,函数逼近,曲线拟和,数值积分,数值微分,解线性方程组的直接方法,解线性方程组的迭代法,非线性方程求根,常微分方程的数值解法

MATLAB解常微分方程的数值解