微分方程数值解实验

微分方程数值解课程设计报告

班级：______________ ： _________

学号：___________

成绩：

2017年 6月 21 日

目录

一、摘要 (1)

二、常微分方程数值解 (2)

2.1 4阶Runge-Kutta法和Adams4阶外插法的基本思路 (2)

2.2 算法流程图 (3)

2.3 用matlab编写源程序 (5)

2.4 常微分方程数值解法应用举例 (7)

三、常系数扩散方程的经典差分格式 (9)

3.1 有限差分法的基本思路 (9)

3.2 算法流程图 (10)

3.3 用matlab编写源程序 (10)

3.4 有限差分法应用举例 (12)

四、椭圆型方程的五点差分格式 (14)

4.1 五点差分法的基本思路 (14)

4.2 算法流程图 (15)

4.3 用matlab编写源程序 (15)

4.4 五点差分法应用举例 (16)

五、自我总结 (19)

六、参考文献 (20)

一、摘要

自然界与工程技术中的很多现象，可以归结为微分方程定解问题。其中，常微分方程求解是微分方程的重要基础容。但是，对于许多的微分方程，往往很难得到甚至不存在精确的解析表达式，这时候，数值解提供了一个很好的解决思路。，针对于此，本文对常微分方程数值解法进行了简单研究，主要讨论了一些常用的数值解法，如欧拉法、改进的欧拉法、Runge—Kutta方法、Adams法以及椭圆型方程、抛物型方程的有限差分方法等，通过具体的算例，结合MATLAB求解画图，初步给出了一般常微分方程数值解法的求解过程。同时，通过对各种方法的误差分析，让大家对各种方法的特点和适用围有一个直观的感受。

关键词：微分方程数值解、MATLAB

二、常微分方程数值解

2.1 基本思路

常微分方程数值解法(numerical methods forordinary differential equations)计算数学的一个分支.是解常微分方程各类定解问题的数值方法.现有的解析方法只能用于求解一些特殊类型的定解问题，实用上许多很有价值的常微分方程的解不能用初等函数来表示，常常需要求其数值解.所谓数值解，是指在求解区间一系列离散点处给出真解的近似值.这就促成了数值方法的产生与发展.

常微分方程初值问题的数值解法是求方程（1）的解在点列1(0,1,)n n n x x h n -=+=上的近似值n y ，这里n h 是1n x -到n x 的步长，一般略去下标记为h 。

0(,)

()dy

f x y dx

y x y

⎧=⎪⎨⎪=⎩ （1） 经典的R K -方法是一个四阶的方法，它的计算公式是：

1

12341213

24

3(22)6

(,)(,)

22(,)22

(,)n n n n n n n n n n h y y K K K K K f x y h h K f x y K h h K f x y K K f x h y hK +⎧

=++++⎪⎪

=⎪⎪⎪

=++⎨⎪

⎪

=++⎪⎪

=++⎪⎩ （2） R K -方法的优点是：单步法、精度高，计算过程便于改变步长，缺点是计

算量较大，每前进一步需要计算四次函数值f。在用龙格库塔方法时,要注意n的

选择要合适，n太大，会使计算量加大，n太小，h较大，可能会使误差增大。因此选择合适的n很重要。我们要在考虑精度的基础上，选择合适的n。

2.2 算法步骤

2.2.1、四阶龙格－库塔（R-K）方法流程图：

2.2.2、Adams4阶外插法流程图：

2.2.3、实例求解流程：

2.3 用matlab编写源程序

Matlab程序源代码：

--------------定义Rk4.m文件----------------------------------------

function dy = Rk4 (x,y)

dy=zeros(3,1);

dy(1)=10*(-y(1)+y(2));

dy(2)=28*y(1)-y(2)-y(1)*y(3);

dy(3)=y(1)*y(2)-8*y(3)/3;

end

-------------------------------------------------------------------------------------------

--------------定义Adams4.m文件---------------------------------------

function [x,y,z]=adams4(x1,y1,z1,x2,y2,z2,x3,y3,z3,h)

%Adams外插法

kfy=0;ksy=0;kty=0;

kfz=0;ksz=0;ktz=0;

kfx=10*(y3-x3); %eval_r(abx);

kfy=x3*(28-z3)-y3;%eval_r(aby);

kfz=x3*y3-8/3*z3; %eval_r(abz);

ksx=10*(y2-x2);

ksy=x2*(28-z2)-y2;;

ksz=x2*y2-8/3*z2;

ktx=10*(y1-x1 );

kty=x1*(28-z1)-y1;

ktz=x1*y1-8/3*z1;

x=x3+h/12*(23*kfx-16*ksx+5*ktx);

y=y3+h/12*(23*kfy-16*ksy+5*kty);

z=z3+h/12*(23*kfz-16*ksz+5*ktz);

end

-------------------------------------------------------------------------------------------

--------------定义exe11.m文件------------------------------------------

[t,y]=ode45(Rk4,[0,30],[12,2,9])

suptitle('Runge-Kutta4阶法') %总标题

subplot(2,2,1);

plot(t,y(:,1));grid on;

legend('x关于t 的变化关系图',1);

xlabel('t','FontSize',14);

ylabel('x','FontSize',14);