数值分析计算方法实验报告.docx
- 文档编号:10294567
- 上传时间:2023-02-09
- 格式:DOCX
- 页数:29
- 大小:324.69KB
数值分析计算方法实验报告.docx
《数值分析计算方法实验报告.docx》由会员分享,可在线阅读,更多相关《数值分析计算方法实验报告.docx(29页珍藏版)》请在冰豆网上搜索。
数值分析计算方法实验报告
课程
名称
计算法
实验项目
名称
线性程组的数值解法
实验项目类型
验证
演示
综合
设计
其他
指导
教师
小兵
成绩
√
1.实验目的:
(1)高斯列主元消去法求解线性程组的过程
(2)熟悉用迭代法求解线性程组的过程
(3)设计出相应的算法,编制相应的函数子程序
2.实验容
分别用高斯列主元消去法,Jacobi迭代法,Gauss--Saidel迭代法,超松弛迭代法求解线性程组
3.实验环境及实验文件存档名
写出实验环境及实验文件存档名
4.实验结果及分析
输出计算结果,结果分析和小结等。
解:
1.高斯列主元消去法:
%用高斯列主元消去法解实验1
%高斯列主元消元法求解线性程组Ax=b
%A为输入矩阵系数,b为程组右端系数
%程组的解保存在x变量中
formatlong;
A=[2,10,0,-3;-3,-4,-12,13;1,2,3,-4;4,14,9,-13]
b=[10,5,-2,7]'
[m,n]=size(A);
ifm~=n
error('矩阵A的行数和列数必须相同');
return;
end
ifm~=size(b)
error('b的大小必须和A的行数或A的列数相同');
return;
end
ifrank(A)~=rank([A,b])
error('A矩阵的秩和增广矩阵的秩不相同,程不存在唯一解');
return;
end
c=n+1;
A(:
c)=b;
fork=1:
n-1
[r,m]=max(abs(A(k:
n,k)));
m=m+k-1;
if(A(m,k)~=0)
if(m~=k)
A([km],:
)=A([mk],:
);
end
A(k+1:
n,k:
c)=A(k+1:
n,k:
c)-(A(k+1:
n,k)/A(k,k))*A(k,k:
c);
end
end
x=zeros(length(b),1);
x(n)=A(n,c)/A(n,n);
fork=n-1:
-1:
1
x(k)=(A(k,c)-A(k,k+1:
n)*x(k+1:
n))/A(k,k);
end
disp('X=');
disp(x);
formatshort;
输出结果:
2.Jacobi迭代法:
%Jacobi迭代法求解实验1
%A为程组的增广矩阵
clc;
A=[2100-310;-3-4-12135;123-4-2;4149-137]
MAXTIME=50;
eps=1e-5;
[n,m]=size(A);
x=zeros(n,1);
y=zeros(n,1);
k=0;
disp('迭代过程X的值情况如下:
')
disp('X=');
while1
disp(x');
fori=1:
1:
n
s=0.0;
forj=1:
1:
n
ifj~=i
s=s+A(i,j)*x(j);
end
y(i)=(A(i,n+1)-s)/A(i,i);
end
end
fori=1:
1:
n
maxeps=max(0,abs(x(i)-y(i)));
end
ifmaxeps<=eps
fori=1:
1:
n
x(i)=y(i);
end
return;
end
fori=1:
1:
n
x(i)=y(i);
y(i)=0.0;
end
k=k+1;
ifk>MAXTIME
error('超过最大迭代次数,退出');
return;
end
end
输出结果:
由于不收敛,故出现上述情况。
3.Gauss--Saidel迭代法
%Gauss_Seidel迭代法求解实验1
%A为程组的增广矩阵
clc;
formatlong;
A=[2100-310;-3-4-12135;123-4-2;4149-137]
[n,m]=size(A);
Maxtime=50;
Eps=10E-5;
x=zeros(1,n);
disp('x=');
fork=1:
Maxtime
disp(x);
fori=1:
n
s=0.0;
forj=1:
n
ifi~=j
s=s+A(i,j)*x(j);
end
end
x(i)=(A(i,n+1)-s)/A(i,i);
end
ifsum((x-floor(x)).^2) break; end; end; X=x; disp('迭代结果: '); X formatshort; 输出结果: 因为不收敛,故出现上述情况。 4.超松弛迭代法: %SOR法求解实验1 %w=1.45 %程组系数矩阵 clc; A=[2,10,0,-3;-3,-4,-12,13;1,2,3,-4;4,14,9,-13] b=[10,5,-2,7]' w=1.45; Maxtime=100; Eps=1E-5; formatlong; n=length(A); k=0; x=ones(n,1); y=x; disp('迭代过程: '); disp('x='); while1 y=x; disp(x'); fori=1: n s=b(i); forj=1: n ifj~=i s=s-A(i,j)*x(j); end end ifabs(A(i,i))<1E-10|k>=Maxtime error('已达最大迭代次数或矩阵系数近似为0,无法进行迭代'); return; end s=s/A(i,i); x(i)=(1-w)*x(i)+w*s; end ifnorm(y-x,inf) break; end k=k+1; end disp('最后迭代结果: '); X=x' formatshort; 输出结果: 由于不收敛,故出现上述情况。 课程 名称 计算法 实验项目 名称 插值法 实验项目类型 验证 演示 综合 设计 其他 指导 教师 小兵 成绩 √ 1.实验目的: (1)学会拉格朗日插值、牛顿插值等基本法 (2)设计出相应的算法,编制相应的函数子程序 (3)会用这些函数解决实际问题 2.实验容 (1)设计拉格朗日插值算法,编制并调试相应的函数子程序 (2)设计牛顿插值算法,编制并调试相应的函数子程序 (3)给定函数四个点的数据如下: X 1.1 2.3 3.9 5.1 Y 3.887 4.276 4.651 2.117 试用拉格朗日插值确定函数在x=2.101,4.234处的函数值。 4)已知 用牛顿插值公式求 的近似值。 3.实验原理 写出本次实验所用算法的算法步骤叙述或画出算法程序框图 4.实验环境及实验文件存档名 写出实验环境及实验文件存档名 5.实验结果及分析 输出计算结果,结果分析和小结等。 (1)拉格朗日法 程序: functionlagrange(x) formatlong; x0=[1.12.33.95.1]; y0=[3.8874.2764.6512.117]; n=length(x0); s=0; forj=0: (n-1) t=1; fori=0: (n-1) ifi~=j t=t*(x-x0(i+1))/(x0(j+1)-x0(i+1)); end end s=s+t*y0(j+1); end s formatshort; 运行结果: (2)牛顿差值法: functionNewton(x) formatlong;%显示15位 x0=[149];%x的值 y0=[123];%y的值 x=5;%插值点 n=max(size(x0)); y=y0 (1);%迭代初始值 disp(y); s=1; dx=y0; fori=1: n-1%构造差商表 dx0=dx; forj=1: n-i dx(j)=(dx0(j+1)-dx0(j))/(x0(i+j)-x0(j)); end df=dx (1); s=s*(x-x0(i)); y=y+s*df;%计算 disp(y); end 运行结果: 实际结果是2.9979.有一定的误差。 总体来说还是很准确的 课程 名称 计算法 实验项目 名称 数值微积分 实验项目类型 验证 演示 综合 设计 其他 指导 教师 小兵 成绩 √ 1.实验目的: (1)学会复化梯形、复化辛浦生求积公式的应用 (3)设计出相应的算法,编制相应的函数子程序 (4)会用这些函数解决实际问题 2.实验容 (1)设计复化梯形公式求积算法,编制并调试相应的函数子程序 (2)设计复化辛浦生求积算法,编制并调试相应的函数子程序 (4)分别用复化梯形公式和复化辛浦生公式计算定积分 取n=2,4,8,16,精确解为0.9460831 3、实验原理 写出本次实验所用算法的算法步骤叙述或画出算法程序框图 4.实验环境及实验文件存档名 写出实验环境及实验文件存档名 5.实验结果及分析 输出计算结果,结果分析和小结等。 本实验是用MATLAB来完成。 梯形法存档名为trap.m辛浦生法存档名为simpson.m 复化梯形公式求积算法: functionT=trap(n) f=sym('sin(x)/x'); a=0; b=1; h=(b-a)/n; T=0; fork=1: (n-1) x0=a+h*k; T=T+limit(f,x0); end T=h*(limit(f,a)+limit(f,b))/2+h*T; T=double(T); 复化辛浦生法: functionS=simpson(n) f=sym('sin(x)/x'); a=0; b=1; h=(b-a)/(2*n); s1=0; s2=0; fork=1: n x0=a+h*(2*k-1); s1=s1+limit(f,x0); end fork=1: (n-1) x0=a+h*2*k; s2=s2+limit(f,x0); end S=h*(limit(f,a)+limit(f,b)+4*s1+2*s2)/3; S=double(S); 运行结果 课程 名称 计算法 实验项目 名称 常微分程的数值解法 实验项目类型 验证 演示 综合 设计 其他 指导 教师 小兵 成绩 √ 1.实验目的: (1)学会四阶龙格-库塔法的使用 (2)设计出相应的算法,编制相应的函数子程序 (3)会用这些函数解决实际问题 2.实验容 (1)分别取h=0.05,N=10;h=0.025,N=20;h=0.01,N=50,用四阶龙格-库塔法求解微分程初值问题: y’=-50y,y(0)=10 (2)某跳伞者在t=0时刻从飞机上跳出,假设初始时刻的垂直速度为0,且跳伞者垂直下落。 已知空气阻力为F=cv2,其中c为常数,v为垂直速度,向下向为正。 写出此跳伞者的速度满足的微分程;若此跳伞者的质量为M=70kg,且已知c=0.27kg/m,利用四阶龙格-库塔公式计算t<=20s的速度(取h=0.1s) 3.实验原理 写出本次实验所用算法的算法步骤叙述或画出算法程序框图 3、实验结果及分析 输出计算结果,结果分析和小结等。 (1) 程序: functionRK(h,n) F='-50*y'; a=0; b=a+h*n; X1=a: h: b; Y1=zeros(1,n+1); Y1 (1)=10; fori=1: n x=X1(i); y=Y1(i); K1=h*eval(F); x=x+h/2; y=y+K1/2; K2=h*eval(F); x=x; y=Y1(i)+K2/2; K3=h*eval(F); x=X1(i)+h; y=Y1(i)+K3; K4=h*eval(F); Y1(i+1)=Y1(i)+(K1+2*K2+2*K3+K4)/6; end X1 disp(Y1); end h=0.05,N=10时,运行结果: h=0.025,n=20时,运行结果: h=0.01,N=50时,输出结果: (2)微分程为: v’=9.8-0.003857v^2,v(0)=0 程序: functionvelosity(h,n) F='9.8-0.003857*y*y'; a=0; b=a+h*n; t=a: h: b; v=zeros(1,n+1); v (1)=0; fori=1: n x=t(i); y=v(i); K1=h*eval(F); x=x+h/2; y=y+K1/2; K2=h*eval(F); x=x; y=v(i)+K2/2; K3=h*eval(F); x=t(i)+h; y=v(i)+K3; K4=h*eval(F); v(i+1)=v(i)+(K1+2*K2+2*K3+K4)/6; end t disp(v); end 运行结果:
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 数值 分析 计算方法 实验 报告