其实列主元高斯消去法无非就是比之前的高斯消去法多了一个判断主元这个步骤,但是里面还是有一些小细节的,比如:你要求一个3*4的增广矩阵,你的主元只需要选两次,第一次是在第一列的0.1.2里面选,第二次就会在第二列的1.2里面选,这里面需要细心一点不然会“连续互换”。其实代码实现也不会太难(本人能力有限,只会用拍照截图来给大家看一道题,大家见谅)。
下面给大家看一道例题帮助大家理解一下列主元高斯消去法:
当然这道题也是非常简单的,但是计算机要解决的话可谓是难上加难,咱们先看下过程分析:
当然过程很简单,下面咱们来看下代码:
#include<stdio.h> //解线性方程组//
#include<math.h>
#define n 3 //n可以自行改变,可以求任意阶的方程组//
int findmain(double a[n][n+1],int l)
{
double max; //定义max为去主元做准备//
int k; //定义主元所在的行数//
k=l;
max=fabs(a[l][l]); //先将第l行l列的元素定位主元,为下面的比大小做准备//
for(int i=l;i<n-1;i++)
{
if(max<fabs(a[i+1][i]))
{
max=fabs(a[i+1][i]);
k=i+1; //谁是主元就让k指向谁//
}
}
return k; //最后返回主元所在的行数//
}
void f(double a[n][n+1])
{
double b[n]={0}; //存结果的数组//
int p=0,m=0; //m是次数,p是找主元的行数//
double count=0,e=0,q=0; //e是交换中介,q是保存数组元素,count是计算x的结果用到的值//
while(1)
{
p=findmain(a,m); //找主元,找n-1次//
for(int l=0;l<n+1;l++) //将主元行移动到首位(这里的首位随着m的变换而变化)//
{
e=a[p][l];
a[p][l]=a[m][l];
a[m][l]=e;
}
for(int i=m+1;i<n;i++)
{
q=a[i][m]; //保存第i行第m个数要不然后面的元素不会被改变了//
for(int j=0;j<n+1;j++)
{
a[i][j]=a[i][j]-(a[m][j]*q/a[m][m]); //线性代数里面正常变换倒三角//
/**********测试*********/
//printf("a[%d][%d]=a[%d][%d]-(a[%d][%d]*a[%d][%d]/a[%d][%d])\n",i,j,i,j,m,j,i,m,m,m);
//printf("%.2lf=%.2lf-(%.2lf*%.2lf/%.2lf)\n",a[i][j],a[i][j],a[m][j],a[i][m],a[m][m]);
}
}
m++; //m保持向后移动,同时也是矩阵变换次数//
printf("您的第%d次变换后的矩阵为:\n",m);
for(int i=0;i<n;i++) //每次都可以看见当前矩阵情况//
{
for(int j=0;j<n+1;j++)
{
printf("%-15lf",a[i][j]);
}
printf("\n");
}
if(m==n-1) //循环结束条件//
break;
}
int flag=0,w; //循环常量,非常重要//
b[n-1]=a[n-1][n]/a[n-1][n-1]; //先存xn,后面就会计算方便一点//
for(int i=n-2;i>=0;i--) //只需要循环n-2次,因为已经又了b[n-1]//
{
flag++; //循环的次数//
w=n-1; //循环常量//
for(int j=flag;j>0;j--) //根据flag的变化而循环//
{
count=b[w]*a[i][w]+count;
/**********测试*************/
/*printf("b[%d]*a[%d][%d]\n",w,i,w);
printf("%lf*%lf=%lf\n",b[w],a[i][w],b[w]*a[i][w]);
printf("%d.count=%lf\n",flag,count);*/
w--; //使a的取值往前走//
}
b[i]=(a[i][n]-count)/a[i][i];
count=0; //count一定要清0不然的话后面没法计算了(上一次的count就会带到下一次运算当中去)//
// printf("%lf=(%lf-%lf)/%lf\n\n",b[i],a[i][n],count,a[i][i]);
}
printf("\n该方程组的解是:("); //正常输出就行//
for(int i=0;i<n;i++)
{
printf("%lf ",b[i]);
}
printf(")\n");
}
int main()
{
double a[n][n+1];
printf("请输入您的%d*%d阶的增广矩阵\n",n,n+1);
for(int i=0;i<n;i++)
{
for(int j=0;j<n+1;j++)
{
scanf("%lf",&a[i][j]);
}
}
printf("您的初始矩阵为:\n");
for(int i=0;i<n;i++)
{
for(int j=0;j<n+1;j++)
{
printf("%-15lf",a[i][j]);
}
printf("\n");
}
f(a);
return 0;
}
当然给大家附上代码运行结果(大家有什么不会可以在讨论区留言,欢迎大家):