文档介绍:数值计算方法
实
验
报
告
班级:计算机31班
姓名: 秦鹏程
学号: 03055014
2005年12月11日
本次上机实验共包括四组题目,分别是关于解线性代数方程组,三次自然样条插值,龙贝格(Romberg)积分法和牛顿迭代法,++语言实现的,并且都能编译通过,运行得到正确答案。
第一组题目: 线性代数方程组的解法
程序算法思想及原理
解线性代数方程组的直接法主要有:高斯消去法和三角分解法,其运算量远比克莱姆法则或约当消去法小。为保证计算过程顺利、稳定,一般需要结合选主元的技术以避免出现小主元,造成过大误差。求解过程总会有舍入误差,误差可能影响的大小取决于方程组固有的条件数cond(A)=||A-1|| ||A||,条件数大时称为病态方程组。实际问题可根据常见病态推测方程是否病态。对病态方程组必须采取慎重措施,才能得出比较准确的解。
本次上机实验中,我总共编程实现了3种算法,分别是高斯消去法,约当消去法,选列主元的高斯消去法,选全主元的高斯消去法,杜里特尔分解法,选主元的杜里特尔分解法。这5种算法都是解线性代数方程组的常用算法,又有自己的局限性。
程序实现细节
高斯消去法
高斯消去法其实就是一种消元法,只是步骤规范,便于程序实现。
本程序中,定义了一个数组a[4][5],用于保存方程组的系数,数组x[4]用于保存最后结果,程序在一开始输出了原来的系数矩阵,便于同消去法得到的矩阵做比较。程序的主要部分就是一个三重循环
for (k=0;k<=n-2;k++)
for (i=k+1;i<=n-1;i++)
{
c=a[i][k]/a[k][k];
a[i][k]=0;
for(j=k+1;j<=n;j++)
a[i][j]=a[i][j]-c*a[k][j];
}
经过上述三重循环,就可以得到一个三角矩阵,以后通过回代便可以得到结果。回代过程为如下:
for(k=n-2;k>=0;k--)
{
xx=0;
for(j=k+1;j<=n-1;j++)
{
xx=xx+a[k][j]*x[j];
}
x[k]=(a[k][n]-xx)/a[k][k];
}
程序为了使中间过程更加透明,除了在一开始就输出原来的系数矩阵以外,还输出了变形后的矩阵,并且在得到结果后,还将得到的x得解带回到原来的方程中(方程1和方程4),从代入后计算的结果可只,得到的解是正确的。
程序运行的结果如下:
选列主元的高斯消去法
高斯消去法消去过程中,第k步求n-k个倍数a[i-1][k-1]/a[k-1][k-1]用到的除数a[k-1][k-1],称为主元。若主元接近0计算机将“溢出”而停止计算或产生较大误差。造成这种结果的原因就是小主元的出现,为了避免出现小主元,可在第k步的第k列元素中选绝对值最大的作为主元,然后交换这个最大的主元所在行和第k行,继续进行消去过程。注意到,交换行相当于改变方程顺序并步影响到原方程的解。
程序分4个部分,第1个部分是找到当前列的最大元素
for (k=0;k<=n-2;k++)
{
max=fabsl(a[k][k]);
for (i=k+1;i<=n-1;i++)//该for语句用于找到当前列的最大元素
{
if(fabsl( a[i][k]) > max)
{
max=fabsl(a[i][k]);
maxi=i;
}
}
第2个部分用于找到后进行交换
if(max!=fabsl(a[k][k])) //找到后进行交换
{
for(j=0;j<=n;j++)
{
tt=a[k][j];
a[k][j]=a[maxi][j];
a[maxi][j]=tt;
}
}
第3个部分才是真正的消去过程
for(i=k+1;i<=n-1;i++)//这个地方应该还有一次for循环,不应影响以前的行变换
{
c=a[i][k]/a[k][k];
a[i][k]=0;
for(j=k+1;j<=n;j++)
a[i][j]=a[i][j]-c*a[k][j];
}
第4个部分是回代
x[n-1]=a[n-1][n]/a[n-1][n-1];//回代
for(k=n-2;k>=0;k--)
{
xx=0;
for(j=k+1;j<=n-1;j++)
{
xx=xx+a[k][j]*x[j];
}
x[k]=(a[k][n]-xx)/a[k][k];
}