改进平方根法解矩阵

 

#include<iostream>
#include<string>
#include<stdio.h>
using namespace std;
    float data[3][3]={2,-1,1,-1,-2,3,1,3,1};//矩阵A 的原始数据
    float datat[3][3]={0};//分解过程的矩阵
    float datal[3][3]={0};//分解过后的L矩阵
    float datad[3][3]={0};//分解过后的D 矩阵
    float y[3]={0};  //定义了一个Y的单行的矩阵
    float b[3]={4,5,6};//所给定的初始值的数据
    float x[3]={0};
    //定义了一个x的解数据
    float datalt[3][3]={0};
    //定义了l的逆矩阵的数据
    float  sumdata(int i,int j)
    {
        float sum=0.0;
       for(int k=0;k<=j-1;k++)
       {
             sum+=datat[i][k]*datal[j][k];
       }
       return sum;
    }
    float dsumdata(int i)
    {
           float sum=0.0;
           for(int k=0;k<=i-1;k++)
            sum+=datat[i][k]*datal[i][k];
           return sum;
    }
    float ysum(int i)//这是求解相关Y的求和
    {
     float sum=0.0;
      for(int k=0;k<=i-1;k++)
         sum+=y[k]*datal[i][k];
      return sum;
    }

    float Xsum(int i)
    {
       float sum=0.0;
       for(int k=i+1;k<3;k++)
          sum+=datalt[i][k]*x[k];
       return sum;


    }
int main(int argc,char *argv)
{
    //这是使用改进平方根法的代码
     int lie=3;
     datad[0][0]=data[0][0];
      for(int i=1;i<lie;i++)
      {
          for(int j=0;j<=i-1;j++)
          {
                datat[i][j]=data[i][j]-sumdata(i,j);
                datal[i][j]=datat[i][j]/datad[j][j];

          }
          datad[i][i]=data[i][i]-dsumdata(i);
      }
      cout<<"This is D matrix"<<endl;
      for(int i=0;i<lie;i++)
      {  //定义了一个Y的单行的矩阵
       cout<<datad[i][i]<<endl;
      }
     cout<<"This is L matrix"<<endl;
        for(int i=0;i<lie;i++)
          datal[i][i]=1;
      for(int i=0;i<lie;i++)
      {
          for(int j=0;j<=i;j++)
             cout<<datal[i][j]<<" ";
             cout<<endl;
      }
      for(int i=0;i<lie;i++)
        for(int j=0;j<=i;j++)
      {
          datalt[j][i]=datal[i][j];
      }
      cout<<"This is L-1 matrix"<<endl;
       for(int i=0;i<lie;i++)
      {
          for(int j=i;j<lie;j++)
             cout<<datalt[i][j]<<" ";
             cout<<endl;
      }
      //根据公式求解Y解向量
        for(int i=0;i<lie;i++)
            y[i]=b[i]-ysum(i);
            cout<<"This is Y matrix"<<endl;
        for(int i=0;i<lie;i++)
                cout<<y[i]<<endl;
                //这是求出D 的逆矩阵的数据
        for(int i=0;i<lie;i++)
            datad[i][i]=1.0/datad[i][i];
            cout<<"This is D inverse"<<endl;
         for(int i=0;i<lie;i++)
                cout<<datad[i][i]<<endl;
         for(int i=0;i<lie;i++)
             datad[i][i]=datad[i][i]*y[i];

              cout<<"This is A matrix of D-1*Y"<<endl;
              for(int i=0;i<lie;i++)
                  cout<<datad[i][i]<<"  ";
                  cout<<endl;



             for(int i=lie-1;i>=0;i--)
                  x[i]=(datad[i][i]-Xsum(i))/1;

            cout<<"************This is a solution vector to X**************"<<endl;
            for(int i=0;i<lie;i++)
                  cout<<x[i]<<"  "<<endl;

}