线性方程组AX=B 的数值计算方法实验一、 实验描述：随着科学技术的发展，线性代数作为高等数学的一个重要组成部分，在科学实践中得到广泛的应用。
本实验的通过C 语言的算法设计以及编程，来实现高斯消元法、三角分解法和解线性方程组的迭代法（雅可比迭代法和高斯-赛德尔迭代法），对指定方程组进行求解。
二、 实验原理：1、高斯消去法：运用高斯消去法解方程组，通常会用到初等变换，以此来得到与原系数矩阵等价的系数矩阵，达到消元的目的。
初等变换有三种：(a)、(交换变换)对调方程组两行；(b)、用非零常数乘以方程组的某一行；(c)、将方程组的某一行乘以一个非零常数，再加到另一行。
通常利用(c)，即用一个方程乘以一个常数，再减去另一个方程来置换另一个方程。
在方程组的增广矩阵中用类似的变换，可以化简系数矩阵，求出其中一个解，然后利用回代法，就可以解出所有的解。
2、选主元：若在解方程组过程中，系数矩阵上的对角元素为零的话，会导致解出的结果不正确。
所以在解方程组过程中要避免此种情况的出现，这就需要选择行的判定条件。
经过行变换，使矩阵对角元素均不为零。
这个过程称为选主元。
选主元分平凡选主元和偏序选主元两种。
平凡选主元：如果()0p pp a ≠，不交换行；如果()0p pp a =，寻找第p 行下满足()0p pp a ≠的第一行，设行数为k ，然后交换第k 行和第p 行。
这样新主元就是非零主元。
偏序选主元：为了减小误差的传播，偏序选主元策略首先检查位于主对角线或主对角线下方第p 列的所有元素，确定行k ，它的元素绝对值最大。
然后如果k p >，则交换第k 行和第p 行。
通常用偏序选主元，可以减小计算误差。
3、三角分解法：由于求解上三角或下三角线性方程组很容易所以在解线性方程组时，可将系数矩阵分解为下三角矩阵和上三角矩阵。
其中下三角矩阵的主对角线为1，上三角矩阵的对角线元素非零。
有如下定理：如果非奇异矩阵A 可表示为下三角矩阵L 和上三角矩阵U 的乘积： A LU = (1) 则A 存在一个三角分解。
而且，L 的对角线元素为1，U 的对角线元素非零。
得到L 和U 后，可通过以下步骤得到X ：(1)、利用前向替换法对方程组LY B =求解Y 。
(2)、利用回代法对方程组UX Y =求解X 。
4、雅可比迭代：考察一般形式的线性方程组：,1,2,3,N ij j i ja xb i N ==∑ (2)设从(2)中分离变出变量i X 将它改写成：111,2,3,Ni i ij jj iij i x b a x i N a =≠⎛⎫ ⎪=-= ⎪ ⎪⎝⎭∑ (3)据此建立雅可比迭代公式：(k 1)(k)111,2,3,N i i ij jj ii j i x b a x i N a +=≠⎛⎫ ⎪=-= ⎪ ⎪⎝⎭∑ (4)5、高斯-赛德尔迭代：通过对雅可比迭代的观察，由于通常1k x +被认为是比k x 更好的x 的近似值，所以在计算1k y +时用1k x +来替换k x 是合理的。
所以对雅可比迭代进行了改进，得到了高斯-赛德尔迭代，其收敛速度更快。
它的迭代公式为： (k 1)(k 1)111,2,3,Ni i ij jj ii j i x b a x i N a ++=≠⎛⎫ ⎪=-= ⎪ ⎪⎝⎭∑ (5)三、 实验内容：1、许多科学应用包含的矩阵带有很多零。
在实际情况中很重要的三角形线性方程组有如下形式: 111211122232223334322111111N N N N N N N N N N NNd x c x b a x d x c x b a x d x c x b a x d x c x b a x d x b --------+=++=++=++=+=(6)构造一个程序求解三角形线性方程组。
可假定不需要行变换，而且可用第k 行消去第k+1行的k x 。
2、使用C 语言编写一个程序求解线性方程组AX B =，其中：1357213500252631A -⎡⎤⎢⎥-⎢⎥=⎢⎥⎢⎥---⎣⎦ 1234B ⎡⎤⎢⎥⎢⎥=⎢⎥⎢⎥⎣⎦3、使用2中的程序求解线性方程组AX B =，其中=1j ij N N A a i -⨯⎡⎤==⎣⎦ ；而且1ij N B b ⨯⎡⎤=⎣⎦，11b N =当2i ≥时，1(i 1)/(i 1)N i b =--。
对3,7,11N =的情况分别求解。
精确解为[]1111X '=。
对得到的结果与精确解的差异进行解释。
4、修改2中的程序，使得它可以通过重复求解N 个线性方程组J J AC E = 其中1,2,,J N = (7)来得到1A - ，则1212[C C C ][E E E ]N N A = (8) 而且112[C C C ]N A -= (9)5、设有如下三角线性方程组，而且系数矩阵具有严格的对角优势：111211122232223334322111111N N N N N N N N N N NNd x c x b a x d x c x b a x d x c x b a x d x c x b a x d x b --------+=++=++=++=+=(10)（1）、根据解线性方程组的迭代法原理，设计一个算法来求解上述方程组，算法必须有效利用系数矩阵的稀疏性。
(2)、用（1）中的程序解下列的三角线性方程组：(a)12123234345484950495043+43+43+434343m m m m m m m m m m m m m m m m +=+=+=+=++=+= (b)12123234345484950495041+42+41+424142m m m m m m m m m m m m m m m m +=+=+=+=++=+=6、利用高斯-赛德尔迭代法求解下列带状方程：1231234123452345646474849504748495048495012252122521225212521225212252125x x x x x x x x x x x x x x x x x x x x x x x x x x x x x -+=-+-+=-+-+=-+-+=-+-+=-+-=-+=四、计算结果及分析：1、此题由于没有提供具体的实例计算，所以,以《数值方法(MATLAB 版)》(第四版)第107页的11、12题为例，进行代码的计算验证。
此题以高斯消去法作为算法的核心，并依此来进行程序的编写。
实例：第11题:求解下列线性方程组。
121232343427239423102412x x x x x x x x x x +=+-=++=-=第12题: 求解下列线性方程组。
1212323434525934219262x x x x x x x x x x +=-+=--+=+=第11题的正确解为:[]1322x '=- 第12题的正确解为:[]2321x '=-精确的结果。
2、此题以带选主元的三角分解法作为算法的核心，以此进行程序的编写。
验算:将x 代入AX B =中，算出[]1.0000 2.0001 3.0000 4.0000B '= 与实际的B 进行比较，此算法在保留四位小数的时候，得出的解比较的正确。
3经过程序计算可看出，在3,7N =时得出的结果比较的精确，与所给的答案一致。
但在11N =时，得出的结果与实际的结果相差很大。
可能由于计算时迭代次数较多，机器误差累计造成。
4、此题由于没有提供具体的实例计算，所以,以《数值方法(MATLAB 版)》 (第四版)第119页的第4(b)题为例，进行代码的计算验证。
J 1-27A =42125-2⎡⎤⎢⎥⎢⎥⎢⎥⎣⎦ 1E 11J ⎡⎤⎢⎥=⎢⎥⎢⎥⎣⎦经验算，可推出[]1012A C C C -=。
5、此题用高斯-赛德尔迭代法作为算法的核心，以此编写程序。
通过对给出的两个实例的计算：6、用5中的程序，通过高斯-赛德尔迭代求解此带状方程解为：五、实验结论：经过C程序计算的得出的数据，经验算后，都比较好的符合了计算的要求，这同时也证明了算法的正确性。
通过本次实验，了解科学计算中的解线性方程组的一些基本方法，通过亲自设计算法，对这些科学的方法有了更加深刻的认识。
对以后的学习有很大帮助。
附件：一、(高斯消去法)：#include<stdio.h>#include<stdlib.h>//头文件声明//float main(){float *a,*b,*c,*d,*x,**A;FILE *fp;//定义指向文件的指针//float temp,sum;int i,j,N;fopen_s(&fp,"D:\\vs程序\\3.4.6-1.12.txt","w" );//确定文件指针的指向// printf("请输入此方程组的元的个数：\n");scanf_s("%d",&N);a=(float *)malloc((N-1)*sizeof(float));//动态定义一维数组，下同//b=(float *)malloc(N*sizeof(float));c=(float *)malloc((N-1)*sizeof(float));d=(float *)malloc(N*sizeof(float));x=(float *)malloc(N*sizeof(float));A=(float **)malloc(N*sizeof(float *));//动态定义二维数组//for (i=0;i<N;i++)A[i]=(float *)malloc(N*sizeof(float));for(i=0;i<N;i++){for(j=0;j<N;j++)A[i][j]=0;//先使系数矩阵为一个0矩阵//}printf("请输入对角系数d[i]:\n");for(i=0;i<N;i++){scanf_s("%f",&d[i]);//输入对角元素的值//}printf("请输入系数a[i]:\n");for(i=0;i<N-1;i++){scanf_s("%f",&a[i]);//输入对角线下的元素系数//}printf("请输入系数c[i]:\n");for(i=0;i<N-1;i++){scanf_s("%f",&c[i]);//输入对角线下的元素系数//}for(i=0;i<N;i++){A[i][i]=d[i];//给系数矩阵的对角元素重新赋值//}for(i=0;i<N-1;i++){A[i+1][i]=a[i];A[i][i+1]=c[i];}//给系数矩阵重新赋值，使之成为三角形线性方程组//printf("请输入系数b[i]:\n");for(i=0;i<N;i++){scanf_s("%f",&b[i]);//输入题中的b矩阵//}for(i=0;i<N-1;i++){temp=(*(A[i+1]+i))/(*(A[i]+i));//对角元素所在列的下一行元素比上此对角元素，求出其比例关系//b[i+1]=b[i+1]-b[i]*temp;//根据上条指令求出的比例关系，进行行变换//for(j=0;j<N;j++)A[i+1][j]=A[i+1][j]-A[i][j]*temp;//进行行变换，已达到消元的目的// }x[N-1]=b[N-1]/A[N-1][N-1];//经过上面的行变换，求出只有一个未知数的那一行的未知数大小//for(i=2;i<=N;i++){sum=0;sum+=A[N-i][N+1-i]*x[N+1-i];x[N-i]=(b[N-i]-sum)/A[N-i][N-i];}//此循环根据书上的回代公式和上面求出的未知数的值，求出输入的方程组的解//if (fp!=NULL){fprintf(fp,"此方程组的解为：");for(i=0;i<N;i++)fprintf(fp,"\nx[%d]=%9.4f",i,x[i]);//将计算结果，输出到文件指针所指向的文件中去//fclose(fp);//关闭文件//}free(a);//释放动态定义指针所指向的空间，便于下次使用。