实验五求解常微分方程实验.doc

求解常微分方程初值问题实验报告
专业：数学与应用数学 年级：09级 班级：AS09101 学号：AS0910136 姓名：张书敏

实验目的
1．常微分方程的数值解法具有十分重要的意义，通过这次试验，了解并掌握求解常微分数值解的几种方法：欧拉方法及改进，龙格-库塔格式；
2．通过本次实验，熟悉求解常微分方程 初值问题的有关方法和理论，主要是欧拉法、改进欧拉法、四阶龙格库塔法，学会编制这两种方法的计算程序；
3．体会这两种解法的功能、优缺点及使用的场合。

一、 实验题目
常微分方程求解
6.1 取h=0.02，编制改进的欧拉公式的MATLAB程序，求下列常微分方程的初值问题的数值解：

6.2 取h=0.1，编制四阶经典龙格-库塔公式的matlab程序，求下列常微分方程的初值问题的数值解：

二、 实验原理
1. 对于欧拉法：
利用 Yn+1 = Yn + hf（Xn, Yn）①
Y0 = Y(X0) ② 二式可以完成计算
需要将微分方程表达式和精度计算表达式作为两个函数保存在m文件里并在程序中调用：
①微分方程（wei_fen）
function z=wei_fen(x,y)

z=(2/x)*y+x*x*exp(x);
end
②精确解计算(jing_que)
function z=jing_que(x)
z=x*x*(exp(x)-exp(1))
end
2. 对于改进的欧拉法：
利用 Yn+1 = Yn + 1/2*K1 + 1/2*K2 ① n = 1, 2, 3……
K1 = hf（Xn, Yn） ②
K2 = hf（Xn + h, Yn + K1） ③ 三式可以完成计算
3. 对于龙格—库塔法：
利用 Yn+1 = Yn + 1/6(K1 + 2K2 + 2K3 +K4) ①
K1 = hf（Xn, Yn） ②
K2 = hf（Xn + 1/2*h, Yn + 1/2*K1） ③
K3 = hf（Xn + 1/2*h, Yn + 1/2*K2） ④
K4 = hf（Xn + h, Yn + K3） ⑤ 四式可以完成计算

三、 实验内容
6.1改进欧拉格式

1.1实验步骤：打开matlab软件，新建一个名为maeuler.m的M文件，编写程序，运行程序，记录结果。
%maeuler.m
function [x, y]=maeuler(dyfun,xspan,y0,h)
format short;
x=xspan(1):h:xspan(2); y(1)=y0;
for n=1:(length(x)-1)
k1=feval(dyfun,x(n),y(n));
y(n+1)=y(n)+h*k1;
k2=feval(dyfun,x(n+1),y(n+1));
y(n+1)=y(n)+h*(k1+k2)/2;
end
x=x'; y=y';

1.2实验程序：

>> clear;dyfun=inline('0.9/(1+2*x)+0*y');
>> [x,y]=maeuler(dyfun,[0,1],1,0.02); [x';y']

1.3实验设备： matlab软件。

6.2 四阶经典龙格-库塔公式

1.1实验步骤：打开matlab软件，新建一个名为marunge4.m的M文件，编写程序，运行程序，记录结果。
%marunge4s.m
function [x,y]=marunge4s(dyfun,xspan,y0,h)
x=xspan(1):h:xspan(2);
y=zeros(length(y0),length(x));
y(:,1)=y0(:);
for n=1:(length(x)-1)
k1=feval(dyfun, x(n), y(:,n));
k2=feval(dyfun, x(n)+h/2, y(:,n)+h/2*k1);
k3=feval(dyfun, x(n)+h/2, y(:,n)+h/2*k2);
k4=feval(dyfun, x(n+1), y(:,n)+h*k3);
y(:,n+1)=y(:,n)+h*(k1+2*k2+2*k3+k4)/6;
end
x=x'; y=y';

1.2实验程序：
>> clear;dyfun=inline('-2/(y-x)');
>> [x,y]=marunge4(dyfun,[0,1],1,0.1); [x',y']
1.3实验设备： matlab软件。

四、 实验结果及分析
6.1、改进欧拉格式结果：
ans =

Columns 1 through 6

0 0.0200 0.0400 0.0600 0.0800 0.1000
1.0000 1.0177 1.0346 1.0510 1.0668 1.0821

Columns 7 through 12

0.1200 0.1400 0.1600 0.1800 0.2000 0.2200
1.0968 1.1111 1.1250 1.1384 1.1514 1.1641

Columns 13 through 18

0.2400 0.2600 0.2800 0.3000 0.3200 0.3400
1.1765 1.1885 1.2001 1.2115 1.2227 1.2335

Columns 19 through 24

0.3600 0.3800 0.4000 0.4200 0.4400 0.4600
1.2441 1.2544 1.2645 1.2744 1.2841 1.2936

Columns 25 through 30

0.4800 0.5000 0.5200 0.5400 0.5600 0.5800
1.3029 1.3120 1.3209 1.3296 1.3382 1.3466

Columns 31 through 36

0.6000 0.6200 0.6400 0.6600 0.6800 0.7000
1.3549 1.3630 1.3709 1.3788 1.3864 1.3940

Columns 37 through 42

0.7200 0.7400 0.7600 0.7800 0.8000 0.8200
1.4014 1.4088 1.4160 1.4231 1.4300 1.4369

Columns 43 through 48

0.8400 0.8600 0.8800 0.9000 0.9200 0.9400
1.4437 1.4503 1.4569 1.4634 1.4698 1.4761

Columns 49 through 51

0.9600 0.9800 1.0000
1.4823 1.4884 1.4944

6.2 四阶经典龙格-库塔公式结果：
ans =
Columns 1 through 6
0 0.1000 0.2000 0.3000 0.4000 0.5000
Columns 7 through 12
0.6000 0.7000 0.8000 0.9000 1.0000 1.0000
Columns 13 through 18
0.7576 -0.3995 0.0795 0.3898 4.5805 4.5306
Columns 19 through 22
4.4787 4.4246 4.3682 4.3092