数值计算：对称正定矩阵的平方根法解法

这也是一个解方程组得。对象是对称正定矩阵。

（该程序在Visual C++.Net 2003 Windows Server 2003下调试成功）

还是没有检查根的存在性。

//对称正定矩阵的平方根法解法

#include "stdafx.h"

#include "conio.h"

#include

#include

#define M 3 //The Size of Matrix

int _tmain(int argc, _TCHAR* argv[])

{


//Definitions


double A[M][M];


double b[M];


double L[M][M];


double y[M];


double x[M];


int i;


int j;


int k;


double Temp;

//Input the Matrix A


std::cout<<"Please Input the Matrix A:\n";


for (i=0;i


{


for (j=0;j


{


std::cin >>A[i][j];


}

}


//End Input


//Input the Matrix b


std::cout <<"Please Input the Matrix b:\n";


for (i=0;i


{


std::cin >>b[i];


}


//End Input





//Calculate Matrix L


for (j=0;j


{


for (i=0;i <=j;i++)


{


Temp=0;


if (i==j)


{


for (k=0;k <=j-1;k++)


{


Temp=Temp+L[j][k]*L[j][k];


}


L[j][j]=sqrt(A[j][j]-Temp);


}


else


{


for (k=0;k <=i-1;k++)


{


Temp=Temp+L[j][k]*L[i][k];


}


L[j][i]=(A[j][i]-Temp)/L[i][i];


}

}


}


//End Calculate

//Calculate Matrix Y


for (i=0;i


{


Temp=0;


for (k=0;k <=i-1;k++)


{


Temp=Temp+L[i][k]*y[k];


}


y[i]=(b[i]-Temp)/L[i][i];


}

//End Calculate

//Calculate Matrix X


for (i=M-1;i>=0;i--)


{


Temp=0;


for (k=i+1;k


{


Temp=Temp+L[k][i]*x[k];


}


x[i]=(y[i]-Temp)/L[i][i];

}


//End Calculate


//Output Results


for (i=0;i


{


std::cout <<"\nX"<<<"="<<<"\T";


}


//End Output


getch();


return 0;

}