三角分解法求线性方程.cpp

计算方法的几个经典算法

//#include<stdio.h>
//#include<math.h>

              /*解方程很快的哦*/ 

#include<iostream.h>
 main()
{
   int n,i,j,k,p;
   double x[100],a[100][100],c;
   cout<<"输入方程组的阶n:";
   cin>>n;
   cout<<"输入系数矩阵:"<<"\n";


   for(i=1;i<=n;i++)
	   for(j=1;j<=n+1;j++)
		   cin>>a[i][j]; 


   if(a[1][1]==0)   //S1
   {
	   cout<<"输出方法失败";
		   return 0;
   } 

   for(i=2;i<=n;i++)
	   a[i][1]=a[i][1]/a[1][1];//S2

   for(k=2;k<=n-1;k++)
   {
	   for(j=k;j<=n+1;j++)
	   {
		   c=0;
		   for(p=1;p<=k-1;p++)
			   c=c+a[k][p]*a[p][j];
		   a[k][j]=a[k][j]-c;
	   }
		if(a[k][k]==0)
		{
		     cout<<"输出方法失败";
		        return 0;

		}
		c=0;
	   for(i=k+1;i<=n;i++)
	   {
		   c=0;
		   for(p=1;p<=k-1;p++)
			   c=c+a[i][p]*a[p][k];
		   a[i][k]=(a[i][k]-c)/a[k][k];
	   }
   }
   for(j=n;j<=n+1;j++)
   {
	   c=0;
	   for(p=1;p<=n-1;p++)
		   c=c+a[n][p]*a[p][j];
	   a[n][j]=a[n][j]-c;
   }
   if(a[n][n]==0)
   {
	   cout<<"输出方法失败";
		        return 0;
   }
   x[n]=a[n][n+1]/a[n][n];
   for(i=n-1;i>=1;i--)
   {
	   c=0;
	   for(j=i+1;j<=n;j++)
		   c=c+a[i][j]*x[j];
	   x[i]=(a[i][n+1]-c)/a[i][i];
   }
   for(i=1;i<=n;i++)
   {
	   cout<<"x["<<i<<"]="<<x[i]<<"   ";
	   if((i%5)==0) cout<<endl;
   }
   return 1;
}