第五版数值分析实验.docx
《第五版数值分析实验.docx》由会员分享，可在线阅读，更多相关《第五版数值分析实验.docx（13页珍藏版）》请在冰豆网上搜索。
第五版数值分析实验
实验一
实验题目：
编写一个拉格朗日插值函数，对不多于9个点的插值节点都可以求出插值函数，任意给定输入x值都可以求出y值。
例如：
（0,1）,（1,1）,（4,2）,（9,3）,（16,4）,（25,5）,（36,6）,（49,7）,（64,8），这实际上是一个平方根函数随意输入x在[0，64]上的值都可以输出一个y值。
进一步可以在这区间上画出图像。
实验原理：
方法一拉格朗日
functiony=lagranger（x0,y0,x）;
%UNTITLEDSummaryofthisfunctiongoeshere
%Detailedexplanationgoeshere
n=length（x0）;
m=length（x）;
fori=1:
m
z=x（i）;
s=0.0;
fork=1:
n
li=1.0;
forj=1:
n
ifj~=k
li=li*（z-x0（j））/（x0（k）-x0（j））;
end
end
s=li*y0（k）+s;
end
y（i）=s;
end
三.函数使用斐波那锲函数
function[f]=fib（n）
%UNTITLED2Summaryofthisfunctiongoeshere
%Detailedexplanationgoeshere
f=[11];
ifn==1&&n>0
f=[1];
elseifn==2
f=[11];
elsefori=3:
n
f（i）=f（i-2）+f（i-1）;
end
end
end
g=f（n）;
end
>>y=fib（10）
y=
11235813213455
>>plot（y,'DisplayName','y','YDataSource','y'）;figure（gcf）
>>
实验二
1.编写一个函数实现被积分
function[y]=shiyan2（x）
%UNTITLED2Summaryofthisfunctiongoeshere
%Detailedexplanationgoeshere
y=log（x）.*sqrt（x）;
end
2.使用MATLAB画出函数曲线
X=0.01:
0.001:
1
Y=shiyan2（x）
Plot（x,y,’-b’）
Gridon
3.编写复合梯形公式计算积分函数
function[y]=echelon（h）
%UNTITLEDSummaryofthisfunctiongoeshere
%Detailedexplanationgoeshere
n=length（h）;
y=[0,0];
fori=1:
n;
temp=0.001;
result=0;
whiletemp+h（i）<1
x=（shiyan2（temp）+shiyan2（temp+h（i）））*h（i）/2;
result=result+x;
temp=temp+h（i）;
end
result=result+（shiyan2（temp）+shiyan2
（1））*h（i）/2;
y（i）=result;
end
end输入一组步长h
调用符合梯形积分函数计算积分结果
计算积分结果的误差
4.辛普森方法
function[y]=quad（h）
%UNTITLED4Summaryofthisfunctiongoeshere
%Detailedexplanationgoeshere
n=length（h）;
y=[0,0];
fori=1:
n;
temp=0.001;
result=0;
whiletemp+h（i）<1;
x=（shiyan2（temp）+4*shiyan2（temp+h（i）/2）+shiyan2（temp+h（i）））*h（i）/6;
result=result+x;
temp=temp+h（i）;
end
result=result+（shiyan2（temp）+shiyan2
（1））*h（i）/6;
y（i）=result;
end
end
输入一组步长h
调用符合梯形积分函数计算积分结果
Y=quade（h）
计算积分结果的误差
实验三
实验目的：
1.编写程序，完成实验题目
实验题目：
用LU分解和列主元消去法解线性方程组
输出Ax=b中系数A=LU分解的矩阵L和U，解向量x及detA；列主元法的行交换次序，解向量x及detA；比较两种方法所得的结果。
方法一LU分解
function[L,U,x]=lux（A,b）
%UNTITLEDSummaryofthisfunctiongoeshere
%Detailedexplanationgoeshere
[n,n]=size（A）;
p=eye（n）;
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],:
）;
p（[km]）=p（[mk]）;
end
fori=k+1:
n
A（i,k）=A（i,k）/A（k,k）;
j=k+1:
n;
A（i,j）=A（i,j）-A（i,k）*A（k,j）;
end
end
end
L=tril（A,-1）+eye（n,n）;
U=triu（A）;
newb=p*b;
y=zeros（n,1）;
fork=1:
n
j=1:
k-1;
y（k）=（newb（k）-L（k,j）*y（j）/L（k,k））;
end
x=zeros（n,1）;
fork=n:
-1:
1
j=k+1:
n;
x（k）=（y（k）-U（k,j）*x（j））/U（k,k）;
end
方法二高斯消去法
functionx=gausslzy（A,b）
%UNTITLED2Summaryofthisfunctiongoeshere
%Detailedexplanationgoeshere
[n,n]=size（A）;
p=eye（n）;
fork=1:
（n-1）
[r,m]=max（abs（A（k:
n,k）））;
m=m+k-1;
ifm>k,t=A（k,:
）;
A（k,:
）=A（m,:
）;
A（m,:
）=t;
end
A（（k+1）:
n,（k+1）:
（n+1））=A（（k+1）:
n,（k+1）:
（n+1））-A（（k+1）:
n,k）/A（k,k）*A（k,（k+1）:
（n+1））;
A（（k+1）:
n,k）=zeros（n-k,1）;
end
x=zeros（n,1）;
x（n）=A（n,n+1）/A（n,n）;
fork=n-1:
-1:
1
x（k,:
）=（A（k,n+1）-A（k,（k+1）:
n）*x（（k+1）:
n））/A（k,k）;
end
升级会员 