完整word版matlab数值分析例题.docx
- 文档编号:8106857
- 上传时间:2023-05-12
- 格式:DOCX
- 页数:12
- 大小:87.54KB
完整word版matlab数值分析例题.docx
《完整word版matlab数值分析例题.docx》由会员分享,可在线阅读,更多相关《完整word版matlab数值分析例题.docx(12页珍藏版)》请在冰点文库上搜索。
完整word版matlab数值分析例题
1、在MATLAB中用Jacobi迭代法讨论线性方程组,
(1)给出Jacobi迭代法的迭代方程,并判定Jacobi迭代法求解此方程组是否收敛。
(2)若收敛,编程求解该线性方程组.
解
(1):
A=[4-11;4—81;-215]%线性方程组系数矩阵
A=
4-11
4-81
—215
>>D=diag(diag(A))
D=
400
0—80
005
〉〉L=—tril(A,-1)%A的下三角矩阵
L=
000
—400
2—10
〉〉U=-triu(A,1)%A的上三角矩阵
U=
01—1
00—1
000
B=inv(D)*(L+U)%B为雅可比迭代矩阵
B=
00.2500—0。
2500
0.500000.1250
0。
4000—0.20000
〉〉r=eigs(B,1)%B的谱半径
r=
0。
3347〈1
Jacobi迭代法收敛。
(2)在matlab上编写程序如下:
A=[4—11;4-81;—215];
〉〉b=[7—2115]';
>〉x0=[000]’;
〉〉[x,k]=jacobi(A,b,x0,1e—7)
x=
2。
0000
4.0000
3。
0000
k=
17
附jacobi迭代法的matlab程序如下:
function[x,k]=jacobi(A,b,x0,eps)
%采用Jacobi迭代法求Ax=b的解
%A为系数矩阵
%b为常数向量
%x0为迭代初始向量
%eps为解的精度控制
max1=300;%默认最多迭代300,超过300次给出警告
D=diag(diag(A));%求A的对角矩阵
L=-tril(A,—1);%求A的下三角阵
U=—triu(A,1);%求A的上三角阵
B=D\(L+U);
f=D\b;
x=B*x0+f;
k=1;%迭代次数
whilenorm(x-x0)>=eps
x0=x;
x=B*x0+f;
k=k+1;
if(k〉=max1)
disp(’迭代超过300次,方程组可能不收敛’);
return;
end
end
2、设有某实验数据如下:
序号
x
y
序号
x
y
1
—3
-3。
99
8
0.5
1。
3776
2
-2。
5
-3。
3011
9
1
1.5403
3
-2
-2。
4161
10
1。
5
1。
5707
4
—1。
5
—1。
4293
11
2
1.5839
5
—1
—0.4597
12
2。
5
1。
6989
6
—0.5
0。
37758
13
3
2.01
7
0
1
(1)在MATLAB中作图观察离散点的结构,用多项式拟合的方法拟合一个合适的多项式函数;
(2)在MATLAB中作出离散点和拟合曲线图。
解
(1):
首先观察离散点的结构,matlab中的程序如下,
x=[—3-2.5-2-1。
5-1—0。
500.511.522。
53];
>〉y=[—3.99-3。
3011-2。
4161-1.4293-0。
45970。
3775811.37761。
54031。
57071。
58391.69892。
01];
>〉plot(x,y,'r*’)
图形如下:
离散点近似如抛物线,所以用二次多项式拟合,所以matlab程序如下:
x=[-3-2.5—2-1.5-1—0。
500.511.522.53];
〉〉y=[—3.99-3。
3011—2.4161—1。
4293—0。
45970。
3775811。
37761.54031.57071。
58391。
69892。
01];
>〉s=polyfit(x,y,2);
〉〉p=poly2str(s,'t’)
p=
-0。
22214t^2+1t+0.74384
(2)做出离散点与拟合曲线的程序如下:
x=[—3—2。
5—2-1。
5—1—0.500.511.522。
53];
〉>f=[—3.99—3。
3011—2。
4161—1.4293-0。
45970.3775811.37761。
54031.57071.58391。
69892。
01];
>>p=polyfit(x,f,2);
〉〉y=polyval(p,x);
>〉plot(x,f,'r+',x,y,'k');
〉>xlabel('x’);
>〉ylabel('y');
>〉title('拟合’);
〉>gtext(datestr(today))
得出的离散点与拟合曲线图像如下:
3、在MATLAB中用复合Simpson公式编程计算
(要求积分精度为
)
解:
matlab程序如下:
[I,step]=jfSimpson('-x—x^2',—1,2,2,1.0e—5)
I=
-4。
5000
step=
4
附复合Simpson程序如下:
function[I,step]=jfSimpson(f,a,b,type,eps)
%type=1辛普森公式
%type=2复合辛普森公式
if(type==2&&nargin==4)
eps=1.0e—4;%缺省精度为0。
0001
end
I=0;
switchtype
case1,
I=((b-a)/6)*(subs(sym(f),findsym(sym(f)),a)+...
4*subs(sym(f),findsym(sym(f)),(a+b)/2)+。
。
.
subs(sym(f),findsym(sym(f)),b));
step=1;
case2,
n=2;
h=(b—a)/2;
I1=0;
I2=(subs(sym(f),findsym(sym(f)),a)+subs(sym(f),findsym(sym(f)),b))/h;
whileabs(I2—I1)>eps
n=n+1;
h=(b—a)/n;
I1=I2;
I2=0;
fori=0:
n—1
x=a+h*i;
x1=x+h;
I2=I2+(h/6)*(subs(sym(f),findsym(sym(f)),x)+。
。
。
4*subs(sym(f),findsym(sym(f)),(x+x1)/2)+。
.。
subs(sym(f),findsym(sym(f)),x1));
end
end
I=I2;
step=n;
end
4、在MATLAB中用四阶Runge-Kutta法编程求解微分方程初值问题
,
并在MATLAB中画图比较方程的解析解与R—K解的结果.
解:
第1步:
首先先把经典RK4子程序编出来如下,用RK4。
m保存。
RK4。
mmatlab程序如下:
function[t,y]=RK4(func,t0,tt,y0,N,varagin)
%Rk方法计算一阶级微分方程组,
%由微分方程知识可以知道,对于高阶微分方程可以化为一阶微分方程组。
%初始时刻为t0,结束时为tt,初始时刻函数值为y0
%N为步数
ifnargin<4
N=100;
end
%如果没有输入N的话,那么默认计算步长为(tt-t0)/100
y(1,:
)=y0(:
)';
h=(tt—t0)/N;%步长
t=t0+[0:
N]’*h;
fork=1:
N
f1=h*feval(func,t(k),y(k,:
));
f1=f1(:
)';
%RK中的k1
f2=h*feval(func,t(k)+h/2,y(k,:
)+f1/2);
f2=f2(:
)';
%RK中的k2
f3=h*feval(func,t(k)+h/2,y(k,:
)+f2/2);
f3=f3(:
)’;
%RK中的k3
f4=h*feval(func,t(k)+h,y(k,:
)+f3);
f4=f4(:
)';
%RK中的k4
y(k+1,:
)=y(k,:
)+(f1+2*(f2+f3)+f4)/6;
%所求结果
End
第2步:
然后把微分方程的表达式用一个调试程序RK_fun表示,用RK_fun.m保存。
RK_fun.mmatlab程序如下:
%RK4方法的测试函数
functionf=RK_fun(t,y)
f=—y+t*t+3;
第3步:
把用RK方法计算的主函数写出,用RK_main.m保存。
RK_main.mmatlab程序如下:
%rk方法的主函数
x0=0;
xt=3;
y0=1;
%符号计算给出分析解
y=dsolve('Dy=—y+t*t+3’,'y(0)=1','t’);
%RK4方法给出计算结果
[x,yrk]=RK4('RK_fun',x0,xt,y0,10);
%对应于真解的离散数据
yreal=subs(y,x);
tol=yreal—yrk;%RK4方法与分析解的误差
plot(x,yreal,x,yrk,’+',x,tol,'*’)
legend('分析解’,’RK4计算结果',’两者误差')
yrk;
[x,yrk,tol];
%yrk为所要计算的结果
第4步matlab中输入y;yrk得:
y=
t^2-4/exp(t)-2*t+5
〉>yrk
yrk=
1.0000
1.5267
1.9647
2.3837
2.8352
3。
3575
3。
9789
4.7203
5。
5972
6.6213
7.8010
最后调用主函数得:
〉〉RK_main
得到分析结果如下图:
5、在MATLAB中利用牛顿切线迭代法计算非线性方程
在区间[1,2]上的一个根。
解,利用matlab写如下程序得:
root=newtonqxf('x*x*x+4*x*x-10',1,2,1.0e—6)
root=
1.3652
附牛顿切线法程序如下:
functionroot=newtonqxf(f,a,b,eps)
%牛顿法求函数f在区间[a,b]上的一个零点
%f为函数名
%a为区间左端点
%b为区间右端点
%eps为根的精度
%root为求出的函数零点
if(nargin==3)
eps=1。
0e-6;
end
f1=subs(sym(f),findsym(sym(f)),a);
f2=subs(sym(f),findsym(sym(f)),b);
if(f1==0)
root=a;
end
if(f2==0)
root=b;
end
if(f1*f2〉0)
disp(’两端点函数值乘积大于0!
’);
return;
else
tol=1;
fun=diff(sym(f));
fa=subs(sym(f),findsym(sym(f)),a);
fb=subs(sym(f),findsym(sym(f)),b);
dfa=subs(sym(fun),findsym(sym(fun)),a);
dfb=subs(sym(fun),findsym(sym(fun)),b);
if(dfa〉dfb)
root=a—fa/dfa;
else
root=b-fb/dfb;
end
while(tol>eps)
r1=root;
fx=subs(sym(f),findsym(sym(f)),r1);
dfx=subs(sym(fun),findsym(sym(fun)),r1);
root=r1—fx/dfx;
tol=abs(root-r1);
end
end
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 完整 word matlab 数值 分析 例题