数值分析实习题作业.docx
- 文档编号:9375844
- 上传时间:2023-02-04
- 格式:DOCX
- 页数:12
- 大小:249.78KB
数值分析实习题作业.docx
《数值分析实习题作业.docx》由会员分享,可在线阅读,更多相关《数值分析实习题作业.docx(12页珍藏版)》请在冰豆网上搜索。
数值分析实习题作业
MATLAB数值分析程序
第一题:
编写列主元消去法程序,计算如下题目,并比较分析两者的区别。
clc;
clear;
disp('请输入A矩阵');
A=[3.016.031.99;1.274.16-1.23;0.987-4.819.34];
disp('请输入b');
b=[111];
A%显示系数矩阵A
b%显示b
n=length(b);%求n
p=1:
n;
lu=A;
y=[];%记录消去法中常数项b的中间变量
fork=1:
n
[c,i]=max(abs(lu(k:
n,k)));%找出最大值放入C,以及最大值的位置放入i
ik=i+k-1;%求最大值对应元素下标
ifik~=k
m=p(k);p(k)=p(ik);p(ik)=m;%记录第几行作为主元
ck=lu(k,:
);lu(k,:
)=lu(ik,:
);lu(ik,:
)=ck;%换主元
end
ifk==n%循环终止条件
break;
end
lu(k+1:
n,k)=lu(k+1:
n,k)/lu(k,k);%求K+1行对应的系数
lu(k+1:
n,k+1:
n)=lu(k+1:
n,k+1:
n)-lu(k+1:
n,k)*lu(k,k+1:
n);%高斯消去
end
l=diag(ones(n,1))+tril(lu,-1);%生成下三角L
u=triu(lu);%生成上三角阵U
y
(1)=b(p
(1));
fori=2:
n
y(i)=b(p(i))-l(i,1:
i-1)*y(1:
i-1)';%计算列主元消去步骤后的常数项
end
x(n)=y(n)/u(n,n);
fori=n-1:
-1:
1
x(i)=(y(i)-u(i,i+1:
n)*x(i+1:
n)')/u(i,i);%列向量存放x
end
x=x';%求出x
x
detA=det(A);%求行列式A的值
detA
condA=cond(A);%求条件数
condA
第
(2)题:
clc;
clear;
disp('请输入A矩阵');
A=[3.006.031.99;1.274.16-1.23;0.990-4.819.34];
disp('请输入b');
b=[111];
A%显示系数矩阵A
b%显示b
n=length(b);%求n
p=1:
n;
lu=A;
y=[];%记录消去法中常数项b的中间变量
fork=1:
n
[c,i]=max(abs(lu(k:
n,k)));%找出最大值放入C,以及最大值的位置放入i
ik=i+k-1;%求最大值对应元素下标
ifik~=k
m=p(k);p(k)=p(ik);p(ik)=m;%记录第几行作为主元
ck=lu(k,:
);lu(k,:
)=lu(ik,:
);lu(ik,:
)=ck;%换主元
end
ifk==n%循环终止条件
break;
end
lu(k+1:
n,k)=lu(k+1:
n,k)/lu(k,k);%求K+1行对应的系数
lu(k+1:
n,k+1:
n)=lu(k+1:
n,k+1:
n)-lu(k+1:
n,k)*lu(k,k+1:
n);%高斯消去
end
l=diag(ones(n,1))+tril(lu,-1);%生成下三角L
u=triu(lu);%生成上三角阵U
y
(1)=b(p
(1));
fori=2:
n
y(i)=b(p(i))-l(i,1:
i-1)*y(1:
i-1)';%计算列主元消去步骤后的常数项
end
x(n)=y(n)/u(n,n);
fori=n-1:
-1:
1
x(i)=(y(i)-u(i,i+1:
n)*x(i+1:
n)')/u(i,i);%列向量存放x
end
x=x';%求出x
x
detA=det(A);%求行列式A的值
detA
condA=cond(A);%求条件数
condA
结果比较:
有些矩阵,系数的微小扰动会使解向量产生巨大差别
第二题:
对龙格函数做拉格朗日多项式插值,并画出龙格现象的函数
第一步,定义原函数
functionf=f(x)
f=1./(1+25*x.^2);
end
第二步:
计算拉格朗日差值函数
functionlangrange=langrange(x,n)
langrange=0;
xx=linspace(-1,1,n+1);
fori=1:
n+1
lix=1;
forj=1:
n+1
ifj~=i
lix=lix.*((x-xx(j))./(xx(i)-xx(j)));
end
end
langrange=f(xx(i)).*lix+langrange;
end
end
第三步:
定义画龙格函数图像
functionrunge1(n)
x=linspace(-1,1,100);
plot(x,f(x),x,langrange(x,n));
end
第四步:
调用上面三个函数
subplot(3,1,1),runge1(5),title('5个节点');
subplot(3,1,2),runge1(10),title('10个节点');
subplot(3,1,3),runge1(20),title('20个节点');
第三题:
计算龙贝格函数积分
function[R,k,T]=romberg(f,a,b,tol)
formatlonge;
k=0;
n=1;
h=b-a;
T=h/2*(f(a)+f(b));
err=1;
whileerr>=tol
k=k+1;
h=h/2;
tmp=0;
fori=1:
n
tmp=tmp+f(a+(2*i-1)*h);
end
T(k+1,1)=T(k)/2+h*tmp;
forj=1:
k
T(k+1,j+1)=T(k+1,j)+(T(k+1,j)-T(k,j))/(4^j-1);
end
n=n*2;
err=abs(T(k+1,k+1)-T(k,k));
end
R=T(k+1,8);
积分子函数:
functionf=f(x)
f=@(x)1./((sin(x))^2+1/4*(cos(x))^2);
end
主函数:
[R,k,T]=romberg(f,0,pi/2,1e-7)
合并三个函数:
functionromberg()
functionf=f(x)
f=@(x)1./((sin(x))^2+1/4*(cos(x))^2);
end
function[R,k,T]=romberg(f,a,b,tol)
formatlonge;
k=0;%µü´ú´ÎÊý
n=1;%Çø¼ä¸öÊý
h=b-a;
T=h/2*(f(a)+f(b));
err=1;
whileerr>=tol
k=k+1;
h=h/2;
tmp=0;
fori=1:
n
tmp=tmp+f(a+(2*i-1)*h);
end
T(k+1,1)=T(k)/2+h*tmp;
forj=1:
k
T(k+1,j+1)=T(k+1,j)+(T(k+1,j)-T(k,j))/(4^j-1);
end
n=n*2;
err=abs(T(k+1,k+1)-T(k,k));
end
R=T(k+1,8);
end
[R,k,T]=romberg(f,0,pi/2,1e-7)
end
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 数值 分析 实习 作业