求C语言课程设计：用高斯列主元消元法解线性方程组

2025-03-01 22:21:32

推荐回答（4个）

回答1：

这里向你推荐一下克鲁特算法（其实就是对高斯列主元消元法进行优化，使之更适合于计算机编程），首先将矩阵A进行LU分解（将系数矩阵分解成一个上三角矩阵和一个下三角矩阵），分解的过程中用到了隐式的主元寻找法，同时利用克鲁特算法可以将两个n*n矩阵压缩到一个n*n矩阵中，大大节省了存储空间提高了计算速度。
方程可化为L*U*x=B，令U*x=y --->L*y=B
然后利用回代先求y，再利用y求x
因为该方法在求解过程中不涉及增广矩阵所以矩阵B几乎不参与什么运算，所以它的计算速度应该能够达到高斯列主元消元法的三倍，但原理与其基本一致。
而且我在程序中使用了动态数组方便你今后进行扩展。

以下程序按照《矩阵论第二版》和《C语言数值计算法方法大全》编写，LU分解部分程序主要参考了《C语言数值计算法方法大全》第二章的程序
如果你需要详细的理论讲解我可以将这两本书和源程序发给你，上面的论述相当详细足够你答辩用的了，我的邮箱hu_hu605@163.com

计算结果：
A矩阵:
2 2 5
3 4 7
1 3 3
B矩阵:
5
6
5
解矩阵:
x 1=-7
x 2=0.333333
x 3=3.66667

Press any key to continue

#include
#include
#include
#include
#include
#include
#include
using namespace std;
#define TINY 1.0e-20 //A small number.
#define N 3
void ludcmp(vector > &a, int n, vector &indx, float &d);//对矩阵进行LU分解
void lubksb(vector > &a, int n, vector &indx, vector &b);//对矩阵进行前向和后向回代
void root(vector > &x,vector &col);//解方程结果保存在y中
void iniv(vector > &x,vector line,int n);//对二维动态数组进行初始化
void main()
{
int i,j,n=N;//输入矩阵的维数
float A[N][N]={{2,2,5},{3,4,7},{1,3,3}};//左边A矩阵
float B[N]={5,6,5};//右边B矩阵
vector > x;//建立动态二维数组存放A，保证你的程序进行扩展时只改A，B，N
vector line;
vector y(n);//建立动态数组存放B
iniv(x,line,n);
y.clear();
for(i=0;i {
y.push_back(B[i]);
for(j=0;j {
x[i].push_back(A[i][j]);
}
}
cout<<"A矩阵:"< for(i=0;i {
for(j=0;j {
cout< }
cout< }
cout<<"B矩阵:"< for(i=0;i {
cout< }
root(x,y);//求根
cout<<"解矩阵:"< for(i=0;i {
cout< }
cout<}
void root(vector > &x,vector &col)
{
int n=x.size(),i=0,j=0;
vector index(n);//用于记录寻找主元素过程中对矩阵的初等变换
index.clear();
float m=1.0;//记录变换方式，此程序中无用
ludcmp(x,n,index,m);//进行LU分解
lubksb(x,n,index,col);//根据分解结果进行回带
}
//以下程序按照《矩阵论第二版》和《C语言数值计算法方法大全》编写，LU分解部分程序主要参考了《C语言数值计算法方法大全》第二章的程序
//如果你需要详细的理论讲解我可以将这两本书和源程序发给你，我的邮箱hu_hu605@163.com
void ludcmp(vector > &a, int n, vector &indx, float &d)
{

int i,imax,j,k;
float big=0,dum=0,sum=0,temp=0;
vector vv(n);
vv.clear();
d=1.0;
for (i=0;i {
big=0.0;
for (j=0;j if ((temp=fabs(a[i][j])) > big)
big=temp;
vv[i]=1.0/big;
}
for (j=0;j {
for (i=0;i {
sum=a[i][j];
for (k=0;k sum -= a[i][k]*a[k][j];
a[i][j]=sum;
}
big=0.0;
for (i=j;i {
sum=a[i][j];
for (k=0;k sum -= a[i][k]*a[k][j];
a[i][j]=sum;
if ( (dum=vv[i]*fabs(sum)) >= big)
{
big=dum;
imax=i;
}
}
if (j != imax)
{
for (k=0;k {
dum=a[imax][k];
a[imax][k]=a[j][k];
a[j][k]=dum;
}
d = -(d);
vv[imax]=vv[j];
}
indx[j]=imax;
if (a[j][j] == 0.0)
a[j][j]=TINY;
if (j != n)
{
dum=1.0/(a[j][j]);
for (i=j+1;i a[i][j] *= dum;
}
}
}
void lubksb(vector > &a, int n, vector &indx, vector &b)
{
int i,ii=0,ip,j;
float sum;
for(i=0;i {
ip=indx[i];
sum=b[ip];
b[ip]=b[i];
b[i]=sum;
}
sum=0;
for (i=1;i {
sum=0;
for(j=0;j {
sum+=a[i][j]*b[j];
}
b[i]=b[i]-sum;
}
b[n-1]=b[n-1]/a[n-1][n-1];
for (i=n-2;i>=0;i--)
{
sum=0;
for(j=i+1;j {
sum+=a[i][j]*b[j];
}
b[i]=(b[i]-sum)/a[i][i];
}
}
void iniv(vector > &x,vector line,int n)
{
int i,j;
for(i=0;i {
x.push_back(line);
for(j=0;j {
x[i].clear();
}
}
}

回答2：

我在上面那位师兄提供的代码的基础上，稍稍做了一点修改，运行结果不仅可以显示X的值，还能显示出消元结束后系数矩阵转变成的上三角阵。代码如下：
#include
#include
#include
#include
#define n 4 /*方程组的阶数*/
#define precision 1e-16
static double aa[n][n+1]={{7.2,2.3,-4.4,0.5,15.1},{1.3,6.3,-3.5,\
2.8,1.8},{5.6,0.9,8.1,-1.3,16.6},{1.5,0.4,3.7,5.9,36.9}};
/*增广矩阵的原始数据*/
void main()
{int i,j,det;double a[n+1][n+2],x[n+1];
int GaussElimination_ColumnSelect();
clrscr();
for(i=1;i<=n;i++) for(j=1;j<=n+1;j++)
/*用a[1][1]~a[n][n+1]存放增广矩阵*/
a[i][j]=aa[i-1][j-1];
det=GaussElimination_ColumnSelect(a,x);
/*调用求解方程的函数，获取返回标志值*/
if(det!=0)
for(i=1;i<=n;i++)
printf("\nx[%d]=%f\n",i,x[i]);printf("\n");
getch();
}
int GaussElimination_ColumnSelect(double a[][n+2],double x[n+1])
/*用列主元高斯消去法求解线性方程组的函数*/
{int i,j,k,r;
double c;
for(k=1;k<=n-1;k++) /*消元过程*/
{r=k;
for(i=k;i<=n;i++) /*选列主元*/
if(fabs(a[i][k])>fabs(a[r][k]))r=i;
if(fabs(a[r][k]){printf("\n det A=0.Elimination faild!");exit(0);}
if(r!=k)
{for(j=k;j<=n+1;j++) /*交换k,r两行*/
{c=a[k][j];a[k][j]=a[r][j];a[r][j]=c;}
}
for(i=k+1;i<=n;i++) /*进行消元计算*/
{c=a[i][k]/a[k][k];
for(j=k;j<=n+1;j++)
a[i][j]=a[i][j]-c*a[k][j];
}
}

printf("After modified,the matrix is::\n");
for(i=1;i<=n;i++){
for(j=1;j<=n;j++){
printf("%f\t",a[i][j]);
if(j%n==0) printf("\n");
}
}

if(fabs(a[n][n]){printf("\n det A=0.Algorithm faild!");exit(0);}
for(k=n;k>=1;k--) /*回代过程*/
{x[k]=a[k][n+1];
for(j=k+1;j<=n;j++)
x[k]=x[k]-a[k][j]*x[j];
x[k]=x[k]/a[k][k];
}
return(1);
}

回答3：

#include
#include
#include
#include
#define n 3 /*方程组的阶数*/
#define precision 1e-16
static double aa[n][n+1]={{2,2,5,5},{3,4,7,6},{1,3,3,5}};
/*增广矩阵的原始数据*/
void main()
{int i,j,det;double a[n+1][n+2],x[n+1];
int GaussElimination_ColumnSelect();

for(i=1;i<=n;i++) for(j=1;j<=n+1;j++)
/*用a[1][1]~a[n][n+1]存放增广矩阵*/
a[i][j]=aa[i-1][j-1];
det=GaussElimination_ColumnSelect(a,x);
/*调用求解方程的函数，获取返回标志值*/
if(det!=0)
for(i=1;i<=n;i++)
printf("\nx[%d]=%f\n",i,x[i]);printf("\n");
getch();
}
int GaussElimination_ColumnSelect(double a[][n+2],double x[n+1])
/*用列主元高斯消去法求解线性方程组的函数*/
{int i,j,k,r;double c;
for(k=1;k<=n-1;k++) /*消元过程*/
{r=k;
for(i=k;i<=n;i++) /*选列主元*/
if(fabs(a[i][k])>fabs(a[r][k]))r=i;
if(fabs(a[r][k]){printf("\n det A=0.Elimination faild!");exit(0);}
if(r!=k)
{for(j=k;j<=n+1;j++) /*交换k,r两行*/
{c=a[k][j];a[k][j]=a[r][j];a[r][j]=c;}
}
for(i=k+1;i<=n;i++) /*进行消元计算*/
{c=a[i][k]/a[k][k];
for(j=k+1;j<=n+1;j++)
a[i][j]=a[i][j]-c*a[k][j];
}
}
if(fabs(a[n][n]){printf("\n det A=0.Algorithm faild!");exit(0);}
for(k=n;k>=1;k--) /*回代过程*/
{x[k]=a[k][n+1];
for(j=k+1;j<=n;j++)
x[k]=x[k]-a[k][j]*x[j];
x[k]=x[k]/a[k][k];
}
return(1);
}

具体思想计算方法上有，在这也说不清。

回答4：

最新问答

憋屁对于人体有害吗？

能告诉我怎样快速学好英语的方法吗？【自学】

江铃2006年5月成立的技术开发中心和江铃2007年11月被称国家级研发中心是否相同?具体解释下,谢谢

你好，我从小肠道不好，小时候就拉，现在就便秘，现在总是几天才上一次厕所，最久的有一个星期，我真害怕

帮忙给两个qq视频秀激活码758098509@qq.com

我是山西对口升学的考生，我被9月份的补录给预录取了不想去想退档，退档后有什么风险急急急。请详细回答

华为荣耀6手机如何设置自动开关机

五菱之光有时加速到三四十码左右时，车子喘抖。松开油门再从新加油又没有。。。请问是哪有故障

"喜欢你"用西班牙语怎么说？

急求答案经济法案例分析（抱歉，只剩下5分了）