/*****************************************************************************/
/* 解线性方程组的直接方法（Gauss、Gausspp、Doolittle） */
/*****************************************************************************/
#include"stdio.h"
#include"math.h"
#define N 100
int n,flag;
float A[N][N],C[N][N],y[N],I=1,O=0;
/*****************************************************************************/
/*****************************************************************************/
float Pivoting(int k)
{
/* 选主元函数*/
int i;
float templ;
for(i=k;i if(fabs(A[i][k])>templ)
{
templ=fabs(A[i][k]); /*将系数矩阵第i列中从第i行到第n行中最大值作为新主元*/
flag=i; /* 标志主元templ所在行,作为行交换函数L_exchange()的形参*/
}
return templ;
}
/*****************************************************************************/
/*****************************************************************************/
void L_exchange(int r,int k)
{
/*行交换函数*/
float B[N];
int i;
for(i=1;i B[i]=A[r][i];
for(i=1;i A[r][i]=A[k][i];
for(i=1;i A[k][i]=B[i];
}
/*****************************************************************************/
/*****************************************************************************/
void print()
{
/*输出函数*/
int i,j;
for(i=1;i {
for(j=1;j printf("%9.3f",A[i][j]);
printf("\n");
}
printf("\n");
}
/*****************************************************************************/
/*****************************************************************************/
void Gresult()
{
/*回代求解函数*/
int k,j;
float x[N];
float memery=0;
x[n]=A[n][n+1]/A[n][n]; /* 回代求解过程 */
for( k=n-1;k>=1;k--)
{
for(j=k+1;j memery=memery+A[k][j]*x[j];
x[k]=(A[k][n+1]-memery)/A[k][k];
memery=0; /* 将计数器memery的值置零以便重新计数*/
}
printf("\nHere are the roots for this linear equations:\n");
for(j=n;j>=1;j--) /*输出方程组的根*/
printf("\n\tx%d=%f\n",j,x[j]);
printf("\n");
}
/*****************************************************************************/
/*****************************************************************************/
void Resultant_Gauss()
{
/*Gauss顺序消去法的消元函数*/
int k,i,j,t;
float templ;
for(k=1;k {
for(t=k;t if(A[t][k]==0)
printf("\nError!!The coefficient matrix A is singular!\n");/*输出方程组系数矩阵奇异的信息*/
for(i=k+1;i {
templ=A[i][k]/A[k][k];
for(j=k;j A[i][j]=A[i][j]-A[k][j]*templ;
}
}
}
/*****************************************************************************/
/*****************************************************************************/
void Resultant_GaussPP()
{
int k,i,j;
float ark;
float templ;
for(k=1;k {
ark=Pivoting(k);
if(ark==0) /*判断方程是否为线性方程，即是否合法*/
printf("\nError!!The coefficient matrix A is singular!"); /*输出方程组系数矩阵奇异的信息*/
else if(flag!=k)
L_exchange(flag,k);
for(i=k+1;i {
templ=A[i][k]/A[k][k];
for(j=k;j A[i][j]=A[i][j]-A[k][j]*templ;
}
}
}
/*****************************************************************************/
/*****************************************************************************/
void H_Doolittle()
{
/* Doolittle分解法*/
int i,r,k;
float temp=0,templ=0;
if(A[1][1]==0)
printf("\nError!!The coefficient matrix is singular!\n");/*输出方程组系数矩阵奇异的信息*/
else
{
for(i=2;i A[i][1]=A[i][1]/A[1][1];
r=2;
while(r {
for(i=r;i {
for(k=1;k templ=templ+(A[r][k]*A[k][i]);
A[r][i]=A[r][i]-templ;
templ=0;
}
for(i=r+1;i {
for(k=1;k templ=templ+(A[i][k]*A[k][r]);
A[i][r]=(A[i][r]-templ)/A[r][r];
templ=0;
}
r++;
}
}
}
/*****************************************************************************/
/*****************************************************************************/
void Doolittle_L()
{
/*Doottle下三角分解及结果输出*/
float x[N]; /*定义数组x[N],用于存放方程组的解*/
float templ=0;
int i, r;
x[n]=y[n]/A[n][n];
printf("\nHere are the roots for this linear equations:\n");
for(r=n-1;r>=1;r--)
{
for(i=r+1;i templ+=A[r][i]*x[i];
x[r]=(y[r]-templ)/A[r][r];
templ=0;
}
for(i=1;i {
printf("\tx%d=%9.3f",i,x[i]);
printf("\n");
}
}
/*****************************************************************************/
/*****************************************************************************/
void Doolittle_U()
{
/*Doottle上三角阵U的求解函数*/
float templ=0;
int i,r;
y[1]=A[1][n+1];
for(r=2;r {
for(i=1;i templ=templ+A[r][i]*y[i];
y[r]=A[r][n+1]-templ;
templ=0;
}
}
/*****************************************************************************/
/*****************************************************************************/
void Print_D_L()
{
/*Doottle分解法,下三角矩阵L的输出函数*/
int i,j;
int x=1,y=1;
printf("\n Here is the Doolittle_L:\n");
for(i=1;i {
for(j=1;j {
if(i>j)
printf("%9.3f",A[i][j]);
else if(i==j)
printf("%9.3f",I);
else
printf("%9.3f",O);
}
printf("\n");
}
}
/*****************************************************************************/
/*****************************************************************************/
void Print_D_U()
{
/*Doottle分解法,上三角矩阵U的输出函数*/
int i,j;
float f=1.0;
printf("\nHere is the Doolittle_U:\n");
for(i=1;i {
for(j=1;j {
if(i {printf("%9.3f",A[i][j]);
f=f*A[i][j]; /*计算系数矩阵A的行列式值*/
}
else
printf("%9.3f",O);
}
printf("\n");
}
if(f==0) /* 判断系数矩阵A的行列式值是否为0 */
printf("Error!!! The coefficient matrix is singular!");
/* 若是，则输出系数矩阵A奇异的信息 */
else /* 若不是，则输出系数矩阵A的行列式值 */
{
printf("The value of determinant A is:"); /*打印输出系数矩阵A的行列式值*/
printf("%9.3f",f);
}
}
/*****************************************************************************/
/*****************************************************************************/
void Gauss()
{
Resultant_Gauss();
Gresult();
}
/*****************************************************************************/
/*****************************************************************************/
GaussPP()
{
Resultant_GaussPP();
Gresult();
}
/*****************************************************************************/
/*****************************************************************************/
void Doolittle()
{
H_Doolittle();
Print_D_L();
Print_D_U();
Doolittle_U();
Doolittle_L();
}
/*****************************************************************************/
/*****************************************************************************/
void Choose()
{
/*菜单函数，负责完成对解题方法的选择以及是否退出程序*/
int s=0;
while(s!=4)
{
printf("If you want to run the program to solve this linear equations,\n");
printf("please make sure that 's' is from 1 to 3; if you want to break \n");
printf("out of program,please input 4 !\n");
printf("\n*******************************************\n");
printf("***\t 1 to choose Gauss().\t\t***\n");
printf("***\t 2 to choose Gausspp().\t\t***\n");
printf("***\t 3 to choose Doottle().\t\t***\n");
printf("***\t 4 to break programer!\t\t***\n");
printf("*******************************************\n");
printf("Please enter your choice: "); /*提示用户选择解题方法*/
scanf("%d",&s);
switch(s)
{
case 1: /*输入1,则选择高斯顺序消去法解线性方程组*/
Gauss();
break;
case 2: /*输入2,则选择高斯列主元素消去法解线性方程组*/
GaussPP();
break;
case 3: /*输入3,则选择Doottle分解法解线性方程组*/
Doolittle();
break;
case 4: /*若输入4,则退出程序*/
break;
}
}
}
/**********************************************************************************/
/**********************************************************************************/
void Input() /*输入函数*/
{
int i,j;
printf("\nPlease input dimension of the linear equations:\n");
printf("* n=");
scanf("%d",&n); /*输入线性方程组的维数*/
printf("\nPlease input augmented matrix (A,b) of the linear equations:\n");
printf("(A,b)=\n\n");
for(i=1;i {
for(j=1;j scanf("%f",&A[i][j]);
}
printf("\n");
}
/*****************************************************************************/
/*****************************************************************************/
void main()
{
Input();
Choose();
}
/*****************************************************************************/
/* THE PROGRAMER END */
/*****************************************************************************/
