线性方程组的数值算法C语言实现(附代码)

合集下载

c语言程序线性方程求解有解释

#include <stdio.h> //标准输入输出头文件#include <stdlib.h> //不是C标准库中的头文件包含了的C语言标准库函数的定义定义了五种类型、一些宏和通用工具函数#include <malloc.h> //包含内存分配函数的函数声明定义常量和表现由堆例程使用的类型#include <math.h> //数学函数库//下面代码为定义为双精度浮点数整形int GS(int,double**,double *,double); //求解函数声明double **TwoArrayAlloc(int,int); //创造二维数组函数声明void TwoArrayFree(double **); //析构二维数组函数声明//注：析构函数当对象脱离其作用域时（对象所在的函数已调用完毕）系统自动执行析构函数。

往往用来做“清理善后”的工作void main() //主函数{int i,j,n; //定义性变量double ep,**a,*b; //定义浮点型变量二维指针一维指针ep = 1e-4; //定义精确度之后用来判断printf("你要解几元线性方程组：\n"); //输出解方程的初步要求scanf("%d",&n); //输出函数输入整形数字//注：求解n元线性方程（即有n个未知数），需要有n个方程组，所以系数总数为n*n,每个方程各有一个常数，一共有n个常数//下面代码分别为系数和常数分配总内存a = TwoArrayAlloc(n,n); //为方程组分配系数所占内存:n*n*sizeof(double)b = (double *)calloc(n,sizeof(double)); //为方程组各个方程分配常数内存:n*sizeof(double)if(b == NULL) //判断内存分配情况{printf("内存分配失败\n"); //输出判断结果exit(1); //内存分配失败，退出程序}//下面代码将用户输入的系数和常数分别一一对应的写入上一步所分配的内存中for(i=0;i<n;i++) //输入n元方程组的系数{printf("请输入第%d行相应的系数：\n",i+1); //说明第几行的方程系数for(j=0;j<n;j++) //循环输入对应行的系数{printf("a[%d][%d]: ",i,j); //输入数在数组中的相对位置scanf("%lf",a[i]+j); //输入系数fflush(stdin); //清空缓存}printf("请输入第%d行相应的常数：\n",i+1); //说明第几行方程组的常数printf("b[%d]: ",i); //循环输入对应行的系数scanf("%lf",b+i); //满足上一个条件后数在数组中的对应位置fflush(stdin); //清空缓存}//注:向下循环如果一直满足j<n则一直向下运行下面步骤一直运行j+1的情况直到不能满足j<n//下面代码根据输入的系数和常数，在屏幕上打印出方程组的具体表达式printf("方程组：\n"); //显示方程组提示语for(i=0;i<n;i++) //打印方程组{for(j=0;j<n;j++) //输入对应行的系数{if(a[i][j]>0) //判断方程组的系数为正{if(j>0)printf(" + "); //从第二行开始如果系数为正输出+号if(a[i][j]!=1) //如果系数不等于1printf("%lfX%d",a[i][j],j+1); //输出系数与Xn的乘积else //否则printf("X%d",j+1); //输出Xn}if(a[i][j]<0) //判断系数小于0{if(j>0)printf(" - "); //如果系数为负输出-号if(a[i][j]!=-1) //如果系数不等于-1printf("%lfX%d",fabs(a[i][j]),j+1); //输出系数的绝对值与Xn的乘积else //否则printf("X%d",j+1); //输出Xn}//注:向下循环如果一直满足i<n 则一直向下运行下面步骤一直运行i+1的情况直到不能满足i<n注：x【n】意思是x包含了n个数据的数组}printf("= %lf\n",b[i]); //打印等号以及常数项}//下面代码判断是否可以用高斯消去法求解，具体逻辑请参考高斯消去法求解线性方程if(!GS(n,a,b,ep)) //判断是否可以用高斯消去法求解，{printf("不可以用高斯消去法求解\n"); //如果不能输出提示语exit(0); //不能用高斯消去法求解，退出程序}printf("该方程组的解为：\n"); //输出该方程组的解打印出结果for(i=0;i<n;i++) //输出结果printf("x%d = %.2f\n",i+1,b[i]); //输出结果TwoArrayFree(a); //释放系数所占内存free (b); //释放常数所占内存}//下面代码为高斯消去法的判断依据函数int GS(int n,double **a,double *b,double ep) //求解函数主体为双精度{int i,j,k,l; //定义整型变量double t; //定义浮点型变量for(k=1;k<=n;k++) //判断每个系数的大小{for(l=k;l<=n;l++) //判断每个系数的大小if(fabs(a[l-1][k-1])>ep) //如果每列系数大于epbreak; //跳出else if(l==n) //如果l=nreturn(0); //返回0if(l!=k) //如果1不等于k{for(j=k;j<=n;j++) //交换系数{t = a[k-1][j-1]; //把第k-1行的每个系数赋给ta[k-1][j-1]=a[l-1][j-1]; //把第l-1行的每个系数赋给k-1行a[l-1][j-1]=t; //把t赋给第1-l行的每个系数}t=b[k-1]; //把每1行的常数项赋给tb[k-1]=b[l-1]; //把l-1行的常数项赋给k-1行b[l-1]=t; //把t的值赋给l-1行的常数}t=1/a[k-1][k-1]; //将t置为Xn分之一用来消元for(j=k+1;j<=n;j++) //消元a[k-1][j-1]=t*a[k-1][j-1]; //将每行的系数乘以该行确定的t 保证Xn的系数为1b[k-1]*=t; //把t的值赋给k-1行的常数for(i=k+1;i<=n;i++) //变换系数{for(j=k+1;j<=n;j++) //变换系数a[i-1][j-1]-=a[i-1][k-1]*a[k-1][j-1]; //把k-1行的每个系数分别赋给j-1和i+1行b[i-1]-=a[i-1][k-1]*b[k-1]; //把k-1行的常数项分别赋给i-1行和b }}for(i=n-1;i>=1;i--) //进行消元for(j=i+1;j<=n;j++) //进行消元b[i-1]-=a[i-1][j-1]*b[j-1]; //把j-1行的常数项分别赋给i-1行和breturn(1); //返回1的步骤}//注:向下循环如果一直满足i>=1则一直向下运行下面步骤一直运行i-1的情况直到不能满足i>=1//下面代码为分配内存的函数double **TwoArrayAlloc(int r,int c) //二维函数创建函数主体{double *x,**y; //定义浮点型变量二维指针一维指针int n; //定义变量x=(double *)calloc(r*c,sizeof(double)); //用x指向r*c个大小为double的内存y=(double **)calloc(r,sizeof(double*)); //用y指向r个double*内存，if(!x||!y) //如果不满足x和y{printf("内存分配失败\n"); //输出判断结果exit(1); //内存分配失败，退出程序}for(n=0;n<=r-1;++n) //输入对应取值y[n]=&x[c*n]; //定义为x在c*n处的地址赋给yreturn (y); //返回y}//注:向下循环如果一直满足n<=r 则一直向下运行下面步骤一直运行n+1的情况直到不能满足n<=r//下面代码为释放内存的函数void TwoArrayFree(double **x) //析构函数{free(x[0]); //释放内存free(x); //释放内存。

线性方程组的数值算法C语言实现(附代码)

线性方程组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 存在一个三角分解。

c++解线性方程组的几种方法

//解线性方程组
#include<iostream.h>
#include<iomanip.h>
#include<stdlib.h>
//----------------------------------------------全局变量定义区
const int Number=15; //方程最大个数
}
void gauss_row_xiaoqu() //高斯列主元消去法
{
int i,j,k,maxi;double lik;
cout<<"用Gauss列主元消去法结果如下:\n";
for(k=0;k<lenth-1;k++)
{
j=k;
for(maxi=i=k;i<lenth;i++)
double a[Number][Number],b[Number],copy_a[Number][Number],copy_b[Number]; //系数行列式
int A_y[Number]; //a[][]中随着横坐标增加列坐标的排列顺序,如a[0][0],a[1][2],a[2][1]...则A_y[]={0,2,1...};
}
double quanpailie_A() //计算行列式中一种全排列的值
{
int i,j,l;
double sum=0,p;
for(i=0,l=0;i<lenth;i++)
for(j=0;A_y[j]!=i&&j<lenth;j++)
if(A_y[j]>i) l++;

数值分析线性方程组的实验报告包含代码解析

线性代数方程组的直接解法实验目的：线性方程组求解的直接法编程实现，实验内容：线性方程组求解的高斯消去法算法实现线性方程组求解的主元素消去法算法实现线性方程组求解的LU 分解得方法算法实现线性方程组求解追赶法算法实现实验比较：高斯消去、主元素消去、LU 分解都用实例 ⎪⎩⎪⎨⎧-=-+-=+-=++15.3181533126321321321x x x x x x x x x 这个进行比较。

知识理论解线性方程组的方法大致分为直接法和迭代法。

直接法是指假设计算过程中不产生舍入误差，经过有限次的运算可求的方程组精确解的方法。

方程（2-1）⎪⎪⎩⎪⎪⎨⎧=+++=+++=+++nn nn n n n n n n b x a x a x a b x a x a x a b x a x a x a ...... (22112)222212111212111 记为：AX=b ；一、高斯（Gauss ）消元法(1).Gauss 消元法是最基本的一种方法。

先逐次消去变量，将方程组化成同解的上三角方程组。

消元过程：先逐次消去变量，将方程组化成同解的上三角方程组基本思想回代过程：按方程的相反顺序求解三角形方程组，得到原方程组的解。

(2) Gauss 消去法的求解思路为：若先让第一个方程组保持不变，利用它消去其余方程组中的,使之变成一个关于变元的n-1阶方程组。

按照（1）中的思路继续运算得到更为低阶的方程组。

经过n-1步的消元后，得到一个三角方程。

利用求解公式回代得到线性方程组的解。

①消元过程：第一次消元，,0)1(11≠a设 )1(11)1(11a a l i i =记（i=1,2,3……，n ）.将（2-1）中第i 个方程减去第一个方程乘以1i l （i=1,2,3……，n ），完成第一次消元，⎪⎩⎪⎨⎧=++=++=+++)2()2(2)2(2)2(2)2(22)2(22)1(1)1(12)1(121)1(11...... ... n n nn n n n n n b x a x a b x a x a b x a x a x a 得等价方程组（2-2）)2()2(:b x A =简记为其中：)1(11)1()2(j i ij ij a l a a -=；)1(11)1()2(b l b b i i i -=次消元后经过1-k ：⎪⎪⎪⎪⎪⎪⎪⎩⎪⎪⎪⎪⎪⎪⎪⎨⎧=+++=+++=+++=+++++=++++++++++++++++++++... ...... ......... ......)()(1)(1,)()(1)(,11)(1,1)(,1)()(1)(1,)()2(2)2(21)2(1,2)2(22)2(22)1(1)1(11)1(1,1)1(12)1(121)1(11k n n k nn k k k n k k nk k k n k n k k k k k k k k k k kn k kn k k k k k k kk n n k k k k n n k k k k b x a x a x a b x a x a x a b x a x a x a b x a x a x a x a b x a x a x a x a x a 简记为：)1()1(++=k k b x A其中：)()()1(k kjik k ij k ija l a a -=+ )()()1(k k ik k i k ib l b b -=+),,1,(n k j i +=按上述方法完成n-1次消元后得到。

线性方程组求解.cpp

#include <stdio.h>#include <math.h>#include <stdlib.h>#define MAX 20 //定义最大的方程组中方程个数double a[MAX+1][MAX+2] ; //方程组系数存放的增广矩阵存放数组double b[MAX+1][MAX+2] ;int n ; //方程的实际个数int xxx = 0; //是否初始化void init()//输入线性方程组{xxx= 1;printf("\n\t请输入方程组的方程个数: ");scanf ("%d",&n);if( n > MAX )//输入的数超过了定义的最大方程个数，退出{printf("sorry,can not deal !\n\n");printf("最大能处理的方程组为%d阶方程组\n\n,MAX");exit(0);}printf("\n请依次输入方程组的方程的系数(顺序从上到下，从左至右)：\n"); for(int i = 1 ; i < n+1 ; i++ )//n行n+1列的增广矩阵for(int j = 1 ; j < n+2 ; j++ )scanf ("%lf",&a[i][j]);for(int i = 1 ; i < n+1 ; i++ ) //备份for(int j = 1 ; j < n+2 ; j++ )b[i][j]=a[i][j];}void reset()//将数组a还原成初始状态{for(int i = 1 ; i < n+1 ; i++ )for(int j = 1 ; j < n+2 ; j++ )a[i][j]=b[i][j];}void display()//输出增光矩阵{for(int i = 1 ; i < n+1 ; i++ ){for(int j = 1 ; j < n+2 ; j++ )printf("%.3f ",a[i][j]);printf("\n");}}void print()//输出线性方程组{int count=0;//记录某一行中系数为0的个数for(int i = 1 ; i < n+1 ; i++ )for(int j = 1 ; j < n+2 ; j++ ){if( a[i][j] == 0 && j != n+1)count ++ ;else{if( j == 1 || count == j-1)//系数为第一项或者该项之前的所有项的系数都是0printf ("%.3f X%d ",a[i][j],j);//小数后面保留三位elseif( j == n+1 ){printf(" = %.3f\n",a[i][n+1]);count =0;}elsea[i][j]>0 ? printf (" + %.3f X%d ",fabs(a[i][j]),j): printf (" - %.3f X%d ",fabs(a[i][j]),j);}}}void selectMainEle(int k)//选取第k列的主元{double mainEle = fabs(a[k][k]) ,temp=0;int j = k ;//用变量j记住主元初始位置for(int i = k+1 ; i < n+1 ; i++ )//选出绝对值最大的数if( mainEle < fabs(a[i][k])){mainEle = fabs(a[i][k]) ;j = i; //j记录主元（列中的最大数）的行数}if( k != j ) //初始值不是最大值，k,j换行for(int i = 1 ; i < n+2 ; i++ ){temp = a[k][i] ;a[k][i] = a[j][i] ;a[j][i] = temp ;}}void printit()//输出方程组的根{printf("\n\n求解的方程组的根为：\n");for(int i=1;i<n+1;i++){if(i%4==0) printf("\n");//每一行输出三个解 printf("X%d=%.3f ",i,a[i][n+1]);}printf("\n");}void gause()//高斯列主消元法核心{printf("线性方程组为：\n");print();//输出输入的方程组for( int k=1 ; k < n ; k ++)//消元过程{selectMainEle(k);for(int i = k+1 ; i < n+1 ; i++ )a[i][k] /= a[k][k];for(int i = k+1 ; i < n+1 ; i++ )for(int j = k+1 ; j < n+2 ; j++ )a[i][j] -= a[i][k] * a[k][j] ;printf("\n消元过程如下：\n");for(int i=1;i<n+1;i++){for(int j=1;j<n+1;j++)if(j<k+1&&i>j)printf("%.4f ",0.0);elseprintf("%.4f ",a[i][j]);printf("\n");}}a[n][n+1] /= a[n][n] ;//回代过程for(int i=n-1 ; i>0 ; i--){double sum=0.0;for(int j=i+1 ;j<n+1 ;j++)sum += a[i][j] * a[j][n+1];a[i][n+1] = (a[i][n+1] - sum)/a[i][i];}printit();}void gauseyuedang()//高斯-约当列主元消去法，高斯列主元消去法的改进{printf("线性方程组为：\n");print();//输出输入的方程组double temp;for( int k=1 ; k < n+1 ; k++)//消元，只保留主元{selectMainEle(k);temp = a[k][k];//用temp记住主元的值for(int i=k ; i < n+2 ; i++)//将主元变为1a[k][i] /= temp;//开始消元，出去主元第k列元素全部消去变为零for(int i = 1 ; i < n+1 ; i++ )//行控制if( i != k ) //行变换不包括当前主元所在行for(int j = k+1 ; j < n+2 ; j++ )//列控制a[i][j]=a[i][j]-a[k][j]*a[i][k];printf("\n消元过程如下：\n");for(int i=1;i<n+1;i++){for(int j=1;j<n+1;j++)if(j<k+1&&i!=j)printf("%.4f ",0.0);elseprintf("%.4f ",a[i][j]);printf("\n");}}printit();//输出方程组的解}void LU()//LU分解法{double sum ;printf("线性方程组为：\n");print();//输出输入的方程组for(int i=2;i<n+1;i++)a[i][1] /= a[1][1];for(int r=2 ; r<n+1 ; r++){for(int i=r ; i<n+1 ; i++)//U变换{sum = 0;for(int k=1;k<r;k++)sum += a[r][k]*a[k][i];a[r][i] -= sum;}for(int i=r+1 ; i<n+1 ; i++)//L变换{sum = 0;for(int k=1 ; k<r ; k++)sum += a[i][k]*a[k][r];a[i][r] = (a[i][r] - sum)/a[r][r];}}printf("\n\nLU分解中的L为：\n");for(int i=1 ; i < n+1 ; i++){for(int j=1 ;j < n+1 ; j++)if(i < j) //上部为0printf("%.3f ",0.0);elsei==j ? printf("%.3f ",1.0): printf("%.3f ",a[i][j]); //主对角线为1 printf("\n");}printf("\n\nLU分解中的U为：\n");for(int i=1 ; i < n+1 ; i++){for(int j=1 ;j < n+1 ; j++)if(i > j) //上部为0printf("%.3f ",0.0);elseprintf("%.3f ",a[i][j]);printf("\n");}for(int i=2;i<n+1;i++)//第一次回带求解{sum=0;for(int k=1;k<i;k++)sum+=a[i][k]*a[k][n+1];a[i][n+1] -= sum;}printf("\n\n第一次求解的方程组(LY=b)的根为：\n");for(int i=1;i<n+1;i++){if(i%4==0) printf("\n");//每一行输出三个解printf("Y%d=%.3f ",i,a[i][n+1]);}a[n][n+1]/=a[n][n];//Xn的值for(int i=n-1;i>0;i--){sum=0;for(int k=i+1;k<n+1;k++)sum+=a[i][k]*a[k][n+1];a[i][n+1] = (a[i][n+1] - sum)/a[i][i];}printf("\n\n第二次求解的方程组(UX=y)的根，即原方程组的解为：\n"); for(int i=1;i<n+1;i++){if(i%4==0) printf("\n");//每一行输出三个解printf("X%d=%.3f ",i,a[i][n+1]);}printf("\n");}void Jacobi()//雅可比迭代法{printf("线性方程组为：\n");print();//输出输入的方程组int count=0;///迭代次数double latest[n+1],old[n+1],temp[n+1],sum,precision,max,temp2;for(int i=1;i<n+1;i++)if(a[i][i] == 0){printf("sorry，主对角线中有数据为0，不能用Jacobi迭代法处理！\n请返回重新选择");return ;}printf("\n输入求解精度值: ");scanf("%lf",&precision);for(int i=1 ; i < n+1 ; i++) //初始状态解为0{latest[i] = 0.0;//新值，下一次迭代结果old[i] = 0.0;//旧值，上一次的迭代结果}printf("Jacobi迭代法迭代过程如下:\n");do{ max =0;printf("\n\n第%d次迭代的结果: ",count++);for(int i=1 ; i<n+1 ; i++){if(i%4==0) printf("\n");//每一行输出三个解printf("X%d=%.4f ",i,old[i]);}for(int i = 1 ; i < n+1 ; i++)//迭代过程,完成新的一组解{sum =0.0;for(int j = 1 ; j < n+1 ; j++)if(j != i) sum = sum + a[i][j] * old[j];latest[i]=(a[i][n+1] - sum)/a[i][i];temp[i] = latest[i];//Y用来记录本次迭代的值temp2 = fabs(latest[i] - old[i]); //新值减旧值if(max < temp2) max = temp2 ;}printf("最大的偏差为：%.4f\n",max);for(int i=1 ; i<n+1 ; i++)old[i] = temp[i]; //将新值保留到old数组} while( max > precision && count <= 100);//前后两次解的偏差满足精度要求或者迭代次数超过100次printf("\n\n第%d次(最后的)求解的方程组的根为：\n",count);for(int i=1;i<n+1;i++){if(i%4==0) printf("\n");//每一行输出三个解printf("X%d=%.4f ",i,old[i]);}printf("\n");}void GS()//高斯-赛德尔迭代法{printf("线性方程组为：\n");print();//输出输入的方程组double precision,sum,X[n+1],max,temp,old;int count=0;///迭代次数for(int i=1;i<n+1;i++)if(a[i][i] == 0){printf("sorry，主对角线中有数据为0，不能用高斯赛德尔迭代法处理！\n请返回重新选择");return ;}printf("\n输入求解精度值: ");scanf("%lf",&precision);for(int i = 1 ; i < n+1 ; i++)// 解初始化为0X[i] = 0.0 ;printf("高斯赛德尔迭代法迭代过程如下:\n");do{printf("\n\n第%d次迭代的结果: ",count++);for(int i=1 ; i<n+1 ; i++){if(i%4==0) printf("\n");//每一行输出三个解printf("X%d=%.5f ",i,X[i]);}for(int i = 1 ; i < n+1 ; i++)//迭代过程,完成新的一组解{sum = 0.0;max = 0.0;old =X[i];for(int j = 1 ; j < n+1 ; j++)if(j != i) sum += a[i][j] * X[j] ;X[i] = (a[i][n+1] - sum) /a[i][i];temp = fabs(X[i]-old);if(max < temp) max = temp;}printf("最大的偏差为：%.5f\n",max);} while(max > precision && count <= 100);//前后两次解的偏差满足精度要求或者迭代次数超过100次printf("\n\n第%d次(最后的)求解的方程组的根为：\n",count);for(int i=1;i<n+1;i++){if(i%4==0) printf("\n");//每一行输出三个解printf("X%d=%.5f ",i,X[i]);}printf("\n");}int main(){int xz;while(1){printf("\t\t\t*****************************************\n");printf("\t\t\t**** \t 0.更改线性方程组 ****\n");printf("\t\t\t**** \t 1.高斯列主消元法 ****\n");printf("\t\t\t**** \t 2.高斯-约当列主消去法 ****\n");printf("\t\t\t**** \t 3.LU分解法 ****\n");printf("\t\t\t**** \t 4.雅各比迭代法 ****\n");printf("\t\t\t**** \t 5.高斯-赛德尔迭代法 ****\n");printf("\t\t\t**** \t 6.退出 ****\n");printf("\t\t\t*****************************************\n");printf("\t请根据相应的代号，从上述方法中选择其一： ");scanf("%d",&xz);switch(xz){case 0: init(); break;case 1: if(xxx){reset();gause();}else{printf("sorry,还没有输入线性方程组!\n");init();}break;//高斯列主消元法case 2: if(xxx){reset();gauseyuedang();}else{printf("sorry,还没有输入线性方程组!\n");init();}break;//高斯-约当列主元消去法case 3: if(xxx){reset();LU();}else{printf("sorry,还没有输入线性方程组!\n");init();}break;case 4:if(xxx){reset();Jacobi();}else{printf("sorry,还没有输入线性方程组!\n");init();}break;case 5: if(xxx){reset();GS();}else{printf("sorry,还没有输入线性方程组!\n");init();}break;case 6: exit(0);default:printf("\n\t错误输入o(︶︿︶)o 唉！\n\n");break;}system("pause");system("cls");}}。

解线性方程组

//线性方程组列主元高斯消元法求解#include"stdafx.h"#include<memory>#if_DEBUG#include<iostream>//测试时使用#endifusing namespace std;//绝对值函数__inline double _abs(double v){return v < 0 ? -v : v;}#if_DEBUGvoid test(int n, const double *pd){int i, j, k = 0;for (i = 0; i < n; i++){for (j = 0; j <= n; j++){cout << pd[k] << " ";k++;}cout << endl;}}#endif//线性方程组列主元高斯消元法//n 方程元数；pCoef 系数，必须以行主序方式存放的二维数组；//pOut 长度为n 的一维数组(调用者负责维护)，用于输出数据//返回值：0 成功，-1 无解，1 申请内存失败, 2 不定解。

int GaussJordanElimination(int n, const double *pCoef, double *pOut){double *pcf;int rows = n, columns = n + 1;pcf = new double[rows * columns];if (pcf == nullptr) return 1; //巧妇难为无米之炊，内存都申请不到，还能干嘛！memcpy(pcf, pCoef, (rows * columns) * sizeof(double)); //据说这个运行效率很高int r, c, i; //循环变量int a, b;double x, y;//开始消元，将pcf 方阵区处理成到直角三角形（直角在右上角）矩阵for (r = 0; r < rows - 1; r++){//选取主元a = r; x = _abs(pcf[r * columns + r]);for (i = r + 1; i < rows; i++){ //查找主元在哪行if (x < _abs(pcf[i * columns + r])) a = i;}if (a > r){ //主元不是当前行（r）,比较麻烦，需要将第a 行与第r 行兑换//第r 列前面的就不要对换了，因为这些项已经被消元，变成0 了for (c = r; c < columns; c++){x = pcf[r * columns + c];pcf[r * columns + c] = pcf[a * columns + c];pcf[a * columns + c] = x;}}//开始消元a = r * columns; //记住将主元的行地址偏移量，以提高程序运行效率x = -pcf[a + r]; //要多次使用，记下她，以提高程序运行效率if (x == 0) //主元居然为0，纯粹是想坑爹，岂能上当！continue; //继续后面的消元，以便最终判断是无解还是任意解for (i = r + 1; i < rows; i++){ //正在消元b = i * columns;//记住将要消元的行地址偏移量，以提高程序运行效率y = pcf[b + r]; //要多次使用，记下她，以提高程序运行效率if (y != 0){ //y == 0，本行不需要消元y /= x; //要多次使用，记下她，以提高程序运行效率pcf[b + r] = 0; //肯定为0，不用计算。

C语言解线性方程的四种方法

C语言解线性方程的四种方法C语言解线性方程的四种方法发了好几天编了个解线性方程组的小程序，可第一次实战就大败而归。

经过半天的调试，仍找不出纠正的方法。

因为并不是算法的问题，而是因为自己对编译器处理浮点函数的方法不是很理解。

明明D=0的方阵解出来不等于0了，跟踪调试发现，计算过程程序对数据进行了舍去处理，导致最终结果不对。

不过如果没有浮点型的话，这个程序应该算不错了。

复制代码代码如下:#include#include#include#defineNUM100voidprint(void)/使用说明/{clrscr();printf("\n\n\n\n\n\t\t\t\tIntroduction\n");printf("\t--------------------------------------------------------------\n");printf("\tThisprogramwasdesignforcomputelinearequations. \n");printf("\tThewayofuseitisverysimple.\n");printf("\tFirst:Inputthenumberoftheequation;(Input0toexit)\ n");printf("\tSecond:Inputthecoefficientofeveryeqution;\n");printf("\tThird:Inputtheconstantofeveryeqution;\n");printf("\tLast:Chosethewayyouwantusetosolvetheequtions;\ n");printf("\tThat''sall,inputanykeytorunit...\n");printf("\t-------------------------By__TJX------------------------------\n");getch();}voidchose(void)/选择计算方法/{clrscr();fflush(stdin);printf("\n\n\n\n\n\t\tIntroduction\n");printf("\t\tChosetheway,please.\n");printf("\t\ta:Gausseliminant.\n");printf("\t\tb:Gauss_ydeliminant.\n");printf("\t\tc:Iterativeway.\n");printf("\t\td:Cramerway.\n");printf("\t\te:exit.\n");printf("\t\tBy__TJX\n");printf("\t\tPleasechoosenumber:\n");}voidinput(doublea1,doubleb1[],intnum)/数据输入/ {inti,j,t;doublep;charde1,de2;do{printf("Pleaseinputarraya[%d][%d]:\n",num,num);printf("Warn:Thefirstnumberofthearraymustn''tcontainzero!\ n");for(i=1;i<=num;i++){printf("Pleaseinputarraya[%d][]:\n",i);for(j=1;j<=num;j++){t=0;if(i==1&&j==1){do{if(t==0){scanf("%lf",&a1[i][j]);t++;}else{printf("Theinputisinvalid,pleaseinputagain:\n");scanf("% f",&a1[i][j]);}}while(a1[i][j]==0);}elsescanf("%lf",&a1[i][j]);}}printf("\nPleasecheckthevalueofarraya[%d][%d],pressYtoinp utagain.\n",num,num);do{de1=getch();}while(de1!=''y''&&de1!=''Y''&&de1!=''n''&&de1!=''N'');}while(de1==''y''||de1==''Y'');do{printf("Pleaseinputarrayb[%d]:\n",num);p=b1+1;for(i=1;i<=num;i++)scanf("%lf",p++);printf("\nPleasecheckthevalueofarrayb[%d],pressYtoinputag \n",num);do{de2=getch();}while(de2!=''y''&&de2!=''Y''&&de2!=''n''&&de2!=''N'');}while(de2==''y''||de2==''Y'');}intmax(doublet1,doublex1[],intn)/迭代子函数/{inti,temp=0;for(i=1;i<=n;i++)if(fabs(x1[i]-t1[i])>1e-2){temp=1;break;}/printf("%d",temp);/returntemp;}intddcompute(doublea1,doubleb1[],doublex1[],intn)/迭代法计算/{doublet;inti,j,k=0;doublesum1=0.0,sum2=0.0;t=(double)malloc(nsizeof(double));printf("\nPleaseInputTheInitialValueofx:\n");for(i=1;i<=n;i++)scanf("%lf",&x1[i]);do{k++;for(i=1;i<=n;i++)t[i]=x1[i];for(i=1;i<=n;i++){sum1=0.0;sum2=0.0;for(j=1;j<=i-1;j++)sum1=sum1+a1[i][j]x1[j];/printf("sum1=%0.4f",sum1);/for(j=i+1;j<=n;j++)sum2=sum2+a1[i][j]t[j];/printf("sum2= %0.4f",sum2);}/if(a1[i][i]==0||fabs(sum1)>1e+12||fabs(sum2)>1e+12){printf("\nWarning:Theseequtionscan''tbesolvebythisway!\n PressanyKeytocontinue...");getch();free(t);return0;}x1[i]=(b1[i]-sum1-sum2)/a1[i][i];}}while(max(t,x1,n));/for(i=1;i<=n;i++){if(i%3==0)printf("\n");printf("%.4f",x1[i]);}/free(t);return1;}intgscompute(doublea1,doubleb1[],doublex1[],intn)/高斯消元法计算/{inti,j,k;doublem,sum;for(k=1;k<=n-1;k++)for(i=k+1;i<=n;i++){if(a1[k][k]==0){printf("\nTheseequtionscan''tbesolveisthisw \nPressanyKeytocontinue...");getch();return0;}if((m=0-a1[i][k]/a1[k][k])==0){i++;continue;}else{for(j=k+1;j<=n;j++)a1[i][j]=a1[i][j]+a1[k][j]m;b1[i]=b1[i]+b1[k]m;}}/yixiajisuanxzhi/x1[n]=b1[n]/a1[n][n];for(i=n-1;i>=1;i--){sum=0.0;for(j=n;j>=i+1;j--)sum=sum+a1[i][j]x1[j];x1[i]=(b1[i]-sum)/a1[i][i];}return1;}intgs_ydcompute(doublea1,doubleb1[],doublex1[],intn)/高斯_约当法计算/{inti,j,k;doublem,sum;for(k=1;k<=n;k++){i=1;while(i<=n){if(a1[k][k]==0){printf("\nTheseequtionscan''tbesolveisthisw ay.\nPressanyKeytocontinue...");getch();return0;}if(i!=k){if((m=0-a1[i][k]/a1[k][k])==0){i++;continue;}else{for(j=k+1;j<=n;j++)a1[i][j]=a1[i][j]+a1[k][j]m;b1[i]=b1[i]+b1[k]m;}i++;}elsei++;}}/yixiajisuanxzhi/for(i=n;i>=1;i--)x1[i]=b1[i]/a1[i][i];return1;}doublecomputed(doublea,inth,intl,intc1,intn)/计算系数行列式D值/{inti,j,p=1;doublesum=0.0;if(h==n)sum=1.0;else{i=++h;c1[l]=0;for(j=1;j<=n;j++)if(c1[j])if(a[i][j]==0)p++;else{sum=sum+a[i][j]computed(a,i,j,c1,n)pow(-1,1+p);p++;}c1[l]=1;}returnsum;}voidncompute(doublea,doubleb[],doublex[],intn,intc,double h)/克莱姆法计算/{inti,j;doublet[NUM];for(j=1;j<=n;j++){for(i=1;i<=n;i++){t[i]=a[i][j];a[i][j]=b[i];}x[j]=computed(a,0,0,c,n)/h;for(i=1;i<=n;i++)a[i][j]=t[i];}}main(){doublex[NUM];doubleb[NUM];inti,j=2,n=0;intc;doublehe;charm,decision;doublea;a=(double)malloc(NUMsizeof(double));for(i=0;ia[i]=(double)malloc(NUMsizeof(double));print();do{clrscr();do{if(n>=NUM)printf("nistoolarge,pleaseinputagain:\n");elseprintf("Pleaseinputthetotalnumberoftheequationsn(n scanf("%d",&n);}while(n>NUM);if(n==0){for(i=1;ifree(a);exit(1);}input(a,b,n);c=(int)malloc((n+1)sizeof(int));memset(c,1,(n+1)sizeof(int));he=computed(a,0,0,c,n);if(fabs(he)>1e-4)Other:chose();do{m=getche();}while(m!=''a''&&m!=''b''&&m!=''A''&&m!=''B''&&m!=''c'' &&m!=''C''&&m!=''d''&&m!=''D''&&m!=''e''&&m!=''E'');switch(m){case''a'':;case''A'':j=gscompute(a,b,x,n);break;case''b'':;case''B'':j=gs_ydcompute(a,b,x,n);break;case''c'':;case''C'':j=ddcompute(a,b,x,n);break;case''d'':;case''D'':j=1;ncompute(a,b,x,n,c,he);break;case''e'':;case''E'':j=2;break;default:j=2;break;}if(j==1){clrscr();printf("\n\n\n\n");printf("D=%.4f\n",he);for(i=1;i<=n;i++){if(i%5==0)printf("\n");printf("%.4f",x[i]);}}elseif(j==0){printf("\nTheseequtionscan''tbesolveisthisway.\nPleasechos etheotherway.");gotoOther;}else{for(i=1;ifree(a);free(c);exit(1);}}elseprintf("\n\n\tD=%.4f\nThislinearequationshasn''taccurat eanswer!",he);printf("\nDoyouwanttocontinue?(Y/N)\n");do{decision=getchar();}while(decision!=''y''&&decision!=''Y''& &decision!=''n''&&decision!=''N'');}while(decision==''y''||decision==''Y'');for(i=1;ifree(a);free(c);}如对本文有所疑问，请点击进入脚本之家知识社区提问。

1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。

线性方程组AX=B 的数值计算方法实验一、实验描述：随着科学技术的发展，线性代数作为高等数学的一个重要组成部分，在科学实践中得到广泛的应用。

二、实验原理：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);//释放动态定义指针所指向的空间，便于下次使用。

线性方程组的数值算法C语言实现(附代码)

c语言程序 线性方程求解 有解释

线性方程组的数值算法C语言实现(附代码)

最新高斯消去法和列主元高斯消去法解线性方程组的程序(C语言)

c++解线性方程组的几种方法

数值分析线性方程组的实验报告包含代码解析

线性方程组求解.cpp

解线性方程组

C语言解线性方程的四种方法

c语言程序线性方程求解有解释