cholesky.cpp

本人的计算方法课程作业

//cholesky.cpp--矩阵的CHOLESKY分解
#include<iostream>
#include<cmath>
using namespace std;
const int ROW=20;
void print_array(int x[ROW][ROW]);
void main()
{
	int a[ROW][ROW],g[ROW][ROW];
	int i,j,k,p;
    for(i=1;i<=ROW;i++)
		for(j=1;j<=ROW;j++)              //求A矩阵的各元素值
		{
			if(j<=i) a[i][j]=j;
            else a[i][j]=i;
		}
    cout<<"在进行分解之前，先输出矩阵A:"<<endl;
	print_array(a);                      //输出A矩阵
	for(i=1;i<=ROW;i++)
	{
		for(k=1;k<i;k++)
		{  
			int t1=0;
		    for(p=1;p<k;p++)
	            t1+=g[i][p]*g[k][p];
	        g[i][k]=(a[i][k]-t1)/g[k][k];   //求G矩阵下三角部分的元素值
		}
		int t2=0;
		for(p=1;p<k;p++)
			t2+=g[k][p]*g[k][p];
		g[k][k]=sqrt(a[k][k]-t2);           //求G矩阵对角线上元素的值
		for(i=1;i<k;i++)
			g[i][k]=0;                      //G矩阵上三角部分的元素均为0
	}
	cout<<"通过分解，可得矩阵G为:"<<endl;
	print_array(g);                         //输出矩阵G
    
}
//print_array()用于输出矩阵
void print_array(int x[ROW][ROW])
{
	int i,j;
	for(i=1;i<=ROW;i++)
	{
		for(j=1;j<=ROW;j++)
			cout<<x[i][j]<<" ";
		cout<<endl;
	}
}