分享

矩阵求逆的算法怎么写

 步梦九渊 2011-09-22
//***************************
//求任何一个矩阵的逆矩阵
//***************************
#include <stdio.h>
#include <malloc.h>

void main( void )
{
     float *buffer,*p;   //定义数组首地址指针变量
     short int row,num; //定义矩阵行数row及矩阵元素个数
     short int i,j;
     float determ;      //定义矩阵的行列式

     float comput_D(float *p,short int n);      //求矩阵的行列式
     float Creat_M(float *p, short int m,short int n,short int k); //求代数余子式
     void Print( float *p,short int n);     //打印n×n的矩阵

     printf("\nPlease input the number of rows: ");
     scanf("%d",&row);
   
     num=2 * row * row;
     buffer = (float *)calloc(num, sizeof(float));     //分配内存单元

     p=buffer;
     if(p != NULL)
     {
         for(i=0;i<row;i++)                   //输入各单元值
         {
             printf("Input the number of %d row ",i+1);
             for(j=0;j<row;j++)
             {
                 scanf("%f",p++);
             }
         }  
     }
     else
         printf( "Can't allocate memory\n" );

     printf("\nThe original matrix is:\n");
     Print(buffer,row);     //打印该矩阵

     determ=comput_D(buffer,row);     //求整个矩阵的行列式
     p=buffer + row * row;
     if (determ != 0)
     {
         for (i=0;i<row; i++)       //求逆矩阵
             for (j=0; j<row; j++)
                    *(p+j*row+i)=   Creat_M(buffer,i,j,row)/determ;   
            
         printf("The determinant is %G\n",determ);

         p=buffer + row * row;
         printf("\nThe inverse matrix is:\n");
         Print(p,row);     //打印该矩阵
     }
     else
         printf("The determnant is 0, and there is no inverse matrix !\n");
     free( buffer );
}
//--------------------------------------------------------
//功能:求矩阵 n X n 的行列式
//入口参数:矩阵首地址 p;矩阵行数 n
//返回值:矩阵的行列式值
//--------------------------------------------------------
float comput_D(float *p,short int n)  
{
     short int i,j,m;         //i--row; j--column
     short int lop=0;
     float result=0;
     float mid=1;
   
     if (n!=1)
     {
         lop=(n==2)?1:n;     //控制求和循环次数,若为2阶,则循环1次,否则为n次

         for(m=0;m<lop;m++)
         {
             mid=1;          //顺序求和
             for(i=0,j=m;i<n;i++,j++)
                 mid = mid * ( *(p+i*n+j%n) );
             result+=mid;
         }

         for(m=0;m<lop;m++)
         {                       
             mid=1;          //逆序相减
             for(i=0,j=n-1-m+n; i<n; i++,j--)
                 mid=mid * ( *(p+i*n+j%n));
             result-=mid;
         }
        }
     else result=*p;
     return(result);
}
//----------------------------------------------------
//功能:求k×k矩阵中元素A(mn)的代数余子式
//入口参数:k×k矩阵首地址;元素A的下标m,n; 矩阵行数 k
//返回值: k×k矩阵中元素A(mn)的代数余子式
//----------------------------------------------------
float Creat_M(float *p, short int m,short int n,short int k)
{
     short int len;
     short int i,j;
     float mid_result=0;
     short int quo=1;
     float *p_creat,*p_mid;

     len=(k-1)*(k-1);
     p_creat = (float *)calloc(len, sizeof(float));     //分配内存单元
     p_mid=p_creat;
     for(i=0;i<k;i++)
         for(j=0;j<k;j++)
         {
             if (i!=m && j!=n)
                 *p_mid++ =* (p+i*k+j);            
         }
     //     Print(p_creat,k-1);
     quo = (m + n) %2==0 ? 1:-1;
     mid_result = (float ) quo * comput_D(p_creat,k-1);
     free(p_creat);
     return(mid_result);
}
//-------------------------------------------
//功能:打印n×n的矩阵
//入口参数:n×n矩阵的首地址;该矩阵的行数 n
//返回值: 无
//-------------------------------------------
void Print( float *p,short int n)   
{
     int i,j;
     for (i=0;i<n;i++)
     {
         for (j=0; j<n;j++)
             printf("%10G ",*p++);
         printf("\n");
     }
     printf("--------------\n");
}

    本站是提供个人知识管理的网络存储空间,所有内容均由用户发布,不代表本站观点。请注意甄别内容中的联系方式、诱导购买等信息,谨防诈骗。如发现有害或侵权内容,请点击一键举报。
    转藏 分享 献花(0

    0条评论

    发表

    请遵守用户 评论公约

    类似文章 更多