太原理工大学数值计算方法实验报告材料Word下载.docx
《太原理工大学数值计算方法实验报告材料Word下载.docx》由会员分享,可在线阅读,更多相关《太原理工大学数值计算方法实验报告材料Word下载.docx(20页珍藏版)》请在冰点文库上搜索。
(2)计算x2。
x2=x1-f(x1)*(x1-x0)/(f(x1)-f(x0))
(3)x0=x1,x1=x2(4)判断是否达到精度,若是输出x1,若否执行
(2)
主要仪器设备
HP计算机
实验记录
1.二分法
//方程求根(二分法).cpp:
定义控制台应用程序的入口点。
//
#include"
stdafx.h"
#include"
iostream"
usingnamespacestd;
classText
{
public:
floatx,y,a,b,c,n=0;
voidGetab()
{
cout<
<
"
请输入计算区间:
(以空格隔开)"
<
endl;
cin>
>
a>
b;
}
floatGetY(floatx)
y=x*x*x+4*x*x-10;
returny;
floatCalculate(floata,floatb)
c=(a+b)/2;
n++;
if(GetY(c)==0||((b-a)/2)<
0.000005)
{
cout<
c<
"
为方程的解"
return0;
}
if(GetY(a)*GetY(c)<
0)
{
returnCalculate(a,c);
if(GetY(c)*GetY(b)<
returnCalculate(c,b);
}
};
intmain()
方程组为:
f(x)=x^3+4x^2-10=0"
floata,b;
Texttext;
text.Getab();
a=text.a;
b=text.b;
text.Calculate(a,b);
return0;
}
2.割线法:
//方程求根(割线法).cpp:
classA
floatx0,x1,y;
y=x*x*x+4*x*x-10;
voidGetNumber()
cout<
请输入两个初始近似值:
x0;
x1;
voidCalculate(floatx0,floatx1)
floatx2;
x2=x1-(GetY(x1)/(GetY(x1)-GetY(x0))*(x1-x0));
if(x2==x1)
cout<
x2<
else
x2<
returnCalculate(x1,x2);
Atext;
text.GetNumber();
a=text.x0;
b=text.x1;
text.Calculate(a,b);
心得体会
使用不同的方法,可以不同程度的求得方程的解,通过二分法计算的程序实现更加了解二分法的特点,二分法过程简单,程序容易实现,但该方法收敛比较慢一般用于求根的初始近似值,不同的方法速度不同。
面对一个复杂的问题,要学会简化处理步骤,分步骤一点一点的循序处理,只有这样,才能高效的解决一个复杂问题。
实验二线性方程组的直接求解
合理选择利用Gauss消元法、主元素消元法、LU分解法、追赶法求解下列方程组:
(1)了解线性方程组常见的直接解法,如Guass消元法、LU分解法、追赶法。
(2)加深对线性方程组求解方法的认识,掌握算法。
1.高斯分解法:
将原方程组化为三角形方阵的方程组:
lik=aik/akk
aij=aij-lik*akjk=1,2,…,n-1
i=k+1,k+2,…,nj=k+1,k+2,…,n+1
由回代过程求得原方程组的解:
xn=ann+1/ann
xk=(akn+1-∑akjxj)/akk(k=n-1,n-2,…,2,1)
2.LU分解法:
将系数矩阵A转化为A=L*U,L为单位下三角矩阵,U为普通上三角矩阵,然后通过解方程组l*y=b,u*x=y,来求解x.
1.高斯消元法:
stdio.h"
math.h"
#include<
stdio.h>
doublea[5][6],a0[5][6];
doublel[5],tmp;
voidExchange(inti)
intj,l,k;
doublemax=a0[i][i],temp;
j=i;
for(k=i;
k<
=3;
k++)
if(a0[k][i]>
max)
max=a0[k][i];
j=k;
}
for(l=i;
l<
=4;
l++)
temp=a0[i][l];
a0[i][l]=a0[j][l];
a0[j][l]=temp;
for(i=1;
i<
i++)
for(j=1;
j<
j++)
a[i][j]=a0[i][j];
voiddisplayA()
{
inti,j;
printf("
\n"
);
for(j=1;
for(i=1;
printf("
%lf"
a[j][i]);
voidmain()
inti,j,k;
{scanf("
%lf"
&
a[i][j]);
a0[i][j]=a[i][j];
displayA();
printf("
列主元素消元法如下"
//消元过程
k=1;
do
Exchange(k);
for(i=k+1;
l[i]=a0[i][k]/a0[k][k];
printf("
l[%i][%i]=%lf"
i,k,l[i]);
for(j=k;
{
a[i][j]=a0[i][j]-l[i]*a0[k][j];
}
displayA();
k++;
if(k==3)break;
for(i=1;
a0[j][i]=a[j][i];
}while
(1);
//回代过程
l[3]=a[3][4]/a[3][3];
for(k=3;
k>
=1;
k--)
{
tmp=0;
for(j=k+1;
j++)tmp+=a[k][j]*l[j];
l[k]=(a[k][4]-tmp)/a[k][k];
x[%i]=%lf\n"
i,l[i]);
}
2.LU分解法:
#include<
math.h>
inti,j,k,r;
doublem=0,p=0;
doublea[3][3];
voidlu(doublea[3][3])
=2;
i++)
if(a[0][0]!
=0)
a[i][0]=a[i][0]/a[0][0];
for(k=1;
k++)
for(j=k;
j++)
for(r=0;
r<
=k-1;
r++)
m=m+a[k][r]*a[r][j];
}a[k][j]=a[k][j]-m;
m=0;
}for(i=k+1;
{{for(r=0;
p=p+a[i][r]*a[r][k];
}a[i][k]=(a[i][k]-p)/a[k][k];
p=0;
}}}
voidmain(){
staticdoublea[3][3]={{1,2,3},{0,1,2},{2,4,1}};
staticdoubleb[3]={14,8,13};
doublec[3];
doubled[3];
doublef[3][3];
doublem=0;
doublen=0;
intr;
lu(a);
输出U的矩阵为\n"
for(i=0;
for(j=i;
f[i][j]=a[i][j];
%f"
f[i][j]);
输出L的矩阵为\n"
for(j=0;
=i;
if(i==j)
a[i][j]=1;
printf("
a[i][j]);
}
else
c[0]=b[0];
for(r=0;
=i-1;
r=r+1)
m=m+a[i][r]*c[r];
c[i]=b[i]-m;
d[2]=c[2]/f[2][2];
i>
=0;
i=i-1)
for(r=2;
r>
i;
r=r-1)
n=n+f[i][r]*d[r];
d[i]=(c[i]-n)/f[i][i];
n=0;
所求方程组解为x1=%f,x2=%f,x3=%f"
d[0],d[1],d[2]);
/*根据LU分解所得两个矩阵及求解步骤计算所求X一组解*/
对于求解线性方程组的各种直接方法来说各有优缺点,在所有的求解方法中都应该注意其解的精度。
注意不同求解方法的不同误差求法。
编写程序的时候需要一步一步慢慢来,逐步增加自己的算法知识水平和解决问题的能力。
实验三线性方程组的迭代求解
使用雅可比迭代法或高斯-赛德尔迭代法对下列方程组进行求解。
雅可比迭代法:
设线性方程组
Ax=b
的系数矩阵A可逆且主对角元素a11,a22,…,ann均不为零,令
D=diag(a11,a22,…,ann)
并将A分解成
A=(A-D)+D
从而线性方程组可写成
Dx=(D-A)x+b
则有迭代公式
x(k+1)=B1x(k)+f1
其中,B1=I-D-1A,f1=D-1b。
inti;
doublex1[20],x2[20],x3[20];
doublex10,x20,x30;
pleaseinputx1,x2,x3:
scanf("
%lf%lf%lf"
x10,&
x20,&
x30);
nx1[n]x2[n]x3[n]\n"
18;
x1[0]=x10;
x2[0]=x20;
x3[0]=x30;
x1[i+1]=0.1*x2[i]+0.2*x3[i]+0.72;
x2[i+1]=0.1*x1[i]+0.2*x3[i]+0.83;
x3[i+1]=0.2*x1[i]+0.2*x2[i]+0.84;
%5d%5lf%5lf%5lf\n"
i,x1[i],x2[i],x3[i]);
在编写算法是不熟悉,查阅了很多资料,经过反复研究和试验后实现了题目的要求,使用雅克比迭代法和高斯-赛德尔都可以得到方程的解,但相比之下,高斯-赛德尔的迭代次数要比雅克比的迭代次数少,能够更快的达到所求的解的精度。
实验四代数插值和最小二乘法拟合
1.学习使用拉格朗日插值法或牛顿插值法求解方法。
2.了解最小二乘法的多项式拟合的具体计算方法并且注意克服正规方程组的病态。
给定数据点(xi,yi)如下:
xi
0.5
0.6
0.7
0.8
0.9
1.0
yi
1
1.75
1.96
2.19
2.44
2.71
3.00
(1)使用拉格朗日插值法或牛顿插值法,求f(0.856)的近似值.
(2)用最小二乘法拟合数据的(n次)多项式,求f(0.856)的近似值.
(3)对比、分析上两结果
设函数在区间[a,b]上n+1互异节点x0,x1,…,xn上的函数值分别为y0,y1,…,yn,求n次插值多项式Pn(x),满足条件
Pn(xj)=yj,j=0,1,…,n
令
Ln(x)=y0l0(x)+y1l1(x)+…+ynln(x)=∑yili(x)
其中l0(x),l1(x),…,ln(x)为以x0,x1,…,xn为节点的n次插值基函数,则Ln(x)是一次数不超过n的多项式,且满足
Ln(xj)=yj,L=0,1,…,n
再由插值多项式的唯一性,得
Pn(x)≡Ln(x)
实验记录(写出实验内容中的程序代码和运行结果)(可分栏或加页)
拉格朗日插值法:
doublem=1.0,a=0.856,l=0;
doublex[6]={0.50,0.60,0.70,0.80,0.90,1.00};
doubley[6]={1.75,1.96,2.19,2.44,2.71,3.00};
=5;
if(i==j)continue;
m=m*((a-x[j])/(x[i]-x[j]));
l+=y[i]*m;
m=1;
结果为%lf"
l);
最小二乘法:
{doublex[7]={0,0.5,0.6,0.7,0.8,0.9,1.0},
y[7]={1,1.75,1.96,2.19,2.44,2.71,3.00},
a0,a1,sum1=0,sum2=0,sum3=0,sum4=0,sum5=0,l,r;
intm=6,i,k;
7;
sum1+=x[i];
sum2+=x[i]*x[i];
sum3+=y[i];
sum4+=x[i]*y[i];
sum5+=y[i]*y[i];
l=sum1/(m+1);
a1=(sum4-l*sum3)/(sum2-l*sum1);
a0=(sum3-sum1*a1)/(m+1);
doubles=sum3*a0+sum4*a1;
r=sum5-s;
y=a0+a1*x\n"
a0=%fa1=%f\t\n"
a0,a1,r);
doubleq=0.856,p;
p=a0+a1*q;
y=%f\n"
p);
拉格朗日插值的优点是插值多项式特别容易建立,缺点是增加节点是原有多项式不能利用,必须重新建立,即所有基函数都要重新计算,这就造成计算量的浪费。
牛顿插值多项式的优点是增加节点时,原先的差商仍旧不变,仍可以使用。
数据拟合的具体作法是:
对给定的数据(xi,yi)(i=0,1,…,m),在取定的函数类中,求p(x)属于此函数类,使误差ri=p(xi)-yi(i=0,1,…,m)的平方和最小,即:
∑ri2=∑(∑p(xi)-yi)2=min
从几何意义上讲,就是寻求与给定点(xi,yi)(i=0,1,…,m)的距离平方和为最小的曲线y=p(x)。