《线性方程组AX=B的数值计算方法实验》.docVIP

《数值计算方法》实验报告 线性方程组AX=B的数值计算方法实验 学号：2013040207012 姓名：梁哲豪 实验描述 在自然科学和工程技术中很多问题的解决常常归结为解线性代数方程组。例如电学中的网络问题，船体数学放样中建立三次样条函数问题，用最小二乘法求实验数据的曲线拟合问题，解非线性方程组问题，用差分法或者有限元法解常微分方程，偏微分方程边值问题等都导致求解线性方程组，而且后面几种情况常常归结为求解大型线性方程组。线性代数方面的计算方法就是研究求解线性方程组的一些数值解法与研究计算矩阵的特征值及特征向量的数值方法。 关于线性方程组的数值解法一般有两类： 直接法：若在计算过程中没有舍入误差，经过有限步算术运算，可求得方程组的精确解的方法。 迭代法：用某种极限过程去逐步逼近线性方程组精确解的方法。迭代法具有占存储单元少，程序设计简单，原始系数矩阵在迭代过程中不变等优点，但存在收敛性及收敛速度等问题。 上三角线性方程组的求解： 基本算法： 高斯消元法：将原方程组化为三角形方阵的方程组： a （k=1,2,…,n-1; i=k+1,k+2, …,n ;j=k+1,k+2, …,n+1） 由回代过程求得原方程组的解： x x LU分解法： 将系数矩阵A转化为A=L*U，L为单位下三角矩阵，U为普通上三角矩阵，然后通过解方程组l*y=b，u*x=y,来求解x。 二、实验内容 1、许多科学应用包含的矩阵带有很多零。在实际情况中很重要的三角形线性方程组有如下形式： …… 构造一个程序求解三角形线性方程组。可假定不需要变换。而且可用第k行消去第k+1行的。 核心代码： #includeiostream.h #includemath.h #includeiomanip.h #define N 4//矩阵阶数 void ColPivot(double c[N][N+1],double []); //函数声明 void main(){ int i,j; double x[N]; double c[N][N+1]={1,3,5,7,1, 2,-1,3,5,2, 0,0,2,5,3, -2,-6,-3,1,4}; cout----------------------------------------endl; cout系数矩阵为: \n; for(i=0;iN;i++){ for(j=0;jN;j++) coutsetw(10)c[i][j]; coutendl; } cout右侧矩阵 y 为: \n; for(i=0;iN;i++) coutsetw(10)c[i][N]; coutendl; cout----------------------------------------endl; ColPivot(c,x); //调用函数，进行高斯消去法变换 cout变换后得到的三角矩阵: \n; for(i=0;iN;i++){ for(j=0;jN;j++) coutsetw(10)c[i][j]; coutendl; } cout变换后的右侧矩阵 y 为: \n; for(i=0;iN;i++) coutsetw(10)c[i][N]; coutendl; cout----------------------------------------endl; cout方程的解为: \n; for(i=0;iN;i++) cout x[i]= x[i]endl; cout----------------------------------------endl; } void ColPivot(double c[N][N+1],double x[]) { int i,j,k; double p,max; double t[N]; for(i=0;i=N-2;i++) { max=0; k=i; for(j=i+1;jN;j++) if(fabs(c[j][i])max) { k=j; max=fabs(c[j][i]); //选主元 } if(k!=i) for(j=i;j=N;j++) { p=c[i][j]; c[i][j]=c[k][j]; //选出主元后进行交换 c[k][j]=p; } for(j=i+1;jN;j++) { p=c[j][i]/c[i][i]; for(k=i;k=N;k++) c[j][k]-=p*c[i][k]; //高斯消去，进行计算 } } for(i=0;iN;i++) t[i]=c[i][N]; for(i=N-1;i=0;i--) //利用回代法求最终解 { for(j=N-1;ji;j--) t[i]-=c[i][j]*x[j]; x[i]=t[i]/c[i][i]; } } 运行结果：（以具体方