电子科技大学matlab-上机实验五

Matlab 与数值分析实验报告（五）
沈汉男
一、实验的性质和任务
通过上机实验，使学生对数值积分、微分方程求解方法有一个初步的理解。

二、实验内容
1.选用复合Simpson 公式，计算
并用Matlab 的符号运算工具箱计算其精确值。

比较结果，找出问题原理，提出解决问题的方法。

2.求积分方程
t t e ds s y e e t y --=⎰10
)(12)(的数值解和精确解，分析二者的差异。

3.利用Euler 法对不同的步长求下面初值问题的数值解：
⎩⎨⎧=-=10
)0(20)()('y t y t y 并通过绘图，与方程的解析解进行比较。

三、实验过程
实验1
a.实验过程：
（1）、创建文件simpr1.m 如下：
function s=simpr1(f,a,b,n)h=(b-a)/(2*n);s1=0;s2=0;for k=1:n x=a+h*(2*k-1);s1=s1+feval('f',x);end for k=1:(n-1)x=a+h*2*k;s2=s2+feval('f',x);end s=h*(feval('f',a)+feval('f',b)+4*s1+2*s2)/3
dx
x x ⎰--+1
1)1ln()1ln(
创建文件f.m 如下：
function y=f(x)y=(log(1+x))*(log(1-x));
在主界面中输入
simpr1('f',-1,1,10)
所得结果为：
s =
-Inf
ans =
-Inf
（2）在主界面中输入程序
clear all
syms x
y=int(log(1+x)*log(1-x),'x','-1','1')
eval(y)
所得结果为
y =
2*log(2)^2-pi^2/3-log(16)+4
ans =
-1.1016
b.实验结果分析：
由上面的实验结果可以看出有时候复合Simpson 公式并不能求出所需解实验2：
a.实验过程：原积分方程为
t t e ds s y e e t y --=⎰10
)(12)(对其进行化简和变形，得：
⎪⎩
⎪⎨⎧+==21)0('
e y y y 原问题转化为一个常微分初值问题。

创建文件Euler_1.m:
function E=Euler_1(fun,x0,y0,xN,N) x=zeros(1,N+1);y=zeros(1,N+1);
x(1)=x0;y(1)=y0;
h=(xN-x0)/N;
for n=1:N
x(n+1)=x(n)+h;
y(n+1)=y(n)+h*feval(fun,x(n),y(n)); end
T=[x',y']
plot(0:N,0:N,x,y,'o',x,y,'-')
创建文件f.m：
function z=f(x,y)
z=y;
在主界面中输入
Euler_1('f',-5,(exp(1)+1)/2,5,10),
其结果为：T=
1.0e+003*
-0.00500.0019
-0.00400.0037
-0.00300.0074
-0.00200.0149
-0.00100.0297
00.0595
0.00100.1190
0.00200.2380
0.00300.4759
0.00400.9519
0.0050 1.9038
其解析解为：
t e
y
图像为：
b.实验结论：
Euler公式能够较精确的求出常微分方程的解实验3：
a.实验过程：
创建文件Euler_1.m:
function E=Euler_1(fun,x0,y0,xN,N) x=zeros(1,N+1);y=zeros(1,N+1);
x(1)=x0;y(1)=y0;
h=(xN-x0)/N;
for n=1:N
x(n+1)=x(n)+h;
y(n+1)=y(n)+h*feval(fun,x(n),y(n)); end
T=[x',y']
plot(0:N,0:N,x,y,'o',x,y,'-')
创建文件f.m：
function z=f(x,y)
z=y-20;
在主界面中输入
Euler_1('f',-10,10,10,10)
所的结果为：
T=
-1010
-8-10
-6-70
-4-250
-2-790
0-2410
2-7270
4-21850
6-65590
8-196810
10-590470在主界面中输入
Euler_1('f',-10,10,10,20)所的结果为：T=
-1010
-90
-8-20
-7-60
-6-140
-5-300
-4-620
-3-1260
-2-2540
-1-5100
0-10220
1-20460
2-40940
3-81900
4-163820
5-327660
6-655340
7-1310700
8-2621420
9-5242860
10-10485740在主界面中输入
Euler_1('f',-10,10,10,30)
所的结果为：
T=
1.0e+007*
-0.00000.0000
-0.00000.0000
-0.0000-0.0000
-0.0000-0.0000
-0.0000-0.0000
-0.0000-0.0000
-0.0000-0.0000
-0.0000-0.0000
-0.0000-0.0001
-0.0000-0.0001
-0.0000-0.0002
-0.0000-0.0003
-0.0000-0.0005
-0.0000-0.0008
-0.0000-0.0013
-0.0000-0.0021
0.0000-0.0035
0.0000-0.0059
0.0000-0.0098
0.0000-0.0164
0.0000-0.0273
0.0000-0.0456
0.0000-0.0760
0.0000-0.1266
0.0000-0.2110
0.0000-0.3517
0.0000-0.5862
0.0000-0.9770
0.0000-1.6284
0.0000-2.7140
0.0000-4.5234
原方程解析解为：
t e
=20
y-
其图像为：
b.实验结果分析：
用Euler公式求解常微分方程处置问题时，
步长越短，其数值解越接近其解析解。