求矩阵的逆矩阵

本文转/参考自：http://www.cnblogs.com/rollenholt/articles/2050662.html

#include<iostream>
#include <stdlib.h>

///////////////////////////////////// 
const int M = 3; //方阵的阶为M
const int N = M*2; //扩展矩阵的列数为N
double **arya; //方阵arya
double **aryab; //扩展矩阵aryab
double **inv; //逆矩阵

///////////////////////////////////////////////////
void init_matr(); //初始化方阵和扩展矩阵
void inv_matr(); //计算逆矩阵并存储
void prin(); //输出逆矩阵
///////////////////////////////////////////////////

int main() {
	init_matr(); //初始化方阵和扩展矩阵
	inv_matr(); //计算逆矩阵并存储
	prin(); //输出逆

	return 0;
}

//初始化方阵和扩展矩阵
void init_matr() { //方阵的阶为M
	std::cout << "扩展矩阵的列为 " << N << '\n';
	std::cout << "方阵的阶为 " << M << '\n' << "input the matrix" << '\n';
	arya = (double **)malloc(M*sizeof(double));
	for(int i=0; i<M; i++)
		*(arya+i) = (double *)malloc(M*sizeof(double));
	aryab = (double **)malloc(M*sizeof(double));
	for(i=0; i<M; i++)
		*(aryab+i) = (double *)malloc(N*sizeof(double));

	for(i=0; i<M; i++)                                
		for(int j=0; j<M; j++) {
			std::cin >> *(*(arya+i)+j);
			*(*(aryab+i)+j) = *(*(arya+i)+j);
		}

	std::cout << "扩展矩阵为" << '\n';
	for(i=0; i<M; i++)
		for(int j=M; j<N; j++) {
			if((j-M)==i)
				*(*(aryab+i)+j) = 1;
			else *(*(aryab+i)+j) = 0;
		}
	for(i=0; i<M; i++) {
		for(int j=0; j<N; j++)
			std::cout << *(*(aryab+i)+j) << "  ";
		std::cout << '\n';
	}
}

//////////////////////////////////////////////////////////////////
//求解逆矩阵
void inv_matr() {
	/*****************求逆模块***********************/
    for(int i=0;i<M;i++)
    {
        if(*(*(aryab+i)+i)==0)
        {
            for(int k=i;k<M;k++)
            {
                if(*(*(aryab+k)+k)!=0)
                {
                    for(int j=0;j<N;j++)
                    {
                        double temp;
                        temp=*(*(aryab+i)+j);
                        *(*(aryab+i)+j)=*(*(aryab+k)+j);
                        *(*(aryab+k)+j)=temp;
                    }
                    break;
                }
            }
            if(k==M)
            {
                std::cout<<"该矩阵不可逆！"<< '\n';
            }
        }
        for(int j=N-1;j>=i;j--)
        {
            *(*(aryab+i)+j)/=*(*(aryab+i)+i);
        }
        for(int k=0;k<M;k++)
        {
            if(k!=i)
            {
                double temp=*(*(aryab+k)+i);
                for(int j=0;j<N;j++)
                {
                    *(*(aryab+k)+j)-=temp*(*(*(aryab+i)+j));
                }
            }
        }
    }
    /*****************求逆模块***********************/

	 //存储逆矩阵
	inv = (double **)malloc(M*sizeof(double));
	for(i=0; i<M; i++)
		*(inv+i) = (double *)malloc(M*sizeof(double));
	for(i=0;i<M;i++)
        for(int j=M;j<N;j++)
				*(*(inv+i)+j-M) = *(*(aryab+i)+j);
}

/////////////////////////////////////////////////////////////////////////
//输出逆矩阵
void prin() {
	std::cout << "逆矩阵为" << '\n';
	for(int i=0; i<M; i++) {
		for(int j=0; j<M; j++)
			std::cout << *(*(inv+i)+j) << "  ";
		std::cout << '\n';
	}
}