求计算方阵行列式的C代码!

bobnie

目前考虑的算法有两种:
1.按行或者按列展开法;
2.对角化法(化为上三角阵即可);
哪位大侠能赐教?如果使用第二种算法,请同时输出上三角矩阵.谢谢......

albertshu

这是以前刚学c时候写的  结构不怎么好  给楼主参考吧

这里采用的计算方法是对角化法 源程序：

#include<stdio.h>
#include<stdlib.h>
#include<math.h>

#define SIZE sizeof(double)
typedef double real;
int Nn;
real *p;

void main()
{
        void inputm(real *,int);
        void printm(real *,int);
        real jisuan(int);
        real *definem(int);

        int n;
        real answer;

        printf("****determinant solver****\n");
        printf("                       by dragon\n");
        do {
                printf("input n please:");
                scanf("%d",&n);
           }  while (n<=0);
        Nn=n;

        p=definem(n);
        inputm(p,n);
        answer=jisuan(n);

        printf("the answer is %.2f\n",answer);

        system("pause");
}

real *definem(int n)
{
        real *pt;
        pt=(real *)malloc(n*n*SIZE);
        return(pt);
}

void inputm(real *pt,int n)
{
        int i,j; double k;

        for(i=1;i<=n;i++)
        {
                for(j=1;j<=n;j++)
                {
                        printf("a%d%d=",i,j);
                        scanf("%lf",&k);
                        *pt=(real)k;
                        pt++;
                        printf("  ");
                }
                printf("\n");
        }
}

real jisuan(int n)
{
        real a(int,int);
        int swap(int);
        int i,j;
        real r,re,k=1;

        if(n==1)  return (*p);
        if(a(n,n)==0)
                {
                        k=-1;
                        if(swap(n)==0) return 0;
                }

        for(i=1;i<=(n-1);i++)
        {
                r=a(i,n)/a(n,n);
                for(j=1;j<=(n-1);j++)
                        *(p+(i-1)*Nn+j-1)=a(i,j)-r*a(n,j);

        }

        k=k*a(n,n);
        re=k*jisuan(n-1);
        return re;
}

int swap(int n)
{
        real a(int,int);
        int i,j;
        real t;
        for(i=0;i<n;i++)
        {
                if(*(p+i*Nn+n-1)!=0)
                {
                        for(j=0;j<=n;j++)
                        {
                                t=a(i+1,j+1);
                                *(p+i*Nn+j)=a(n,j+1);
                                *(p+(n-1)*Nn+j)=t;
                        }
                        return 1;
                }
        }
        return 0;
}

real a(int x,int y)
{
        real t;
        t=*(p+(x-1)*Nn+y-1);
        return(t);
}

dragon_

哈哈，这个源代码我有。数值分析课的各类算法从头写到尾。
最好使用LU分解法，数值稳定性强，最后输出下三角矩阵的对角线元素乘积即可。

bobnie

跪求二楼大侠给出关键算法的注释,小弟数学能力有些差,

bobnie

跪求二楼大侠给出关键算法的注释,小弟数学能力有些差,

dragon_

对矩阵进行LU分解的子程序如下：

1. void LUDecomposition(double matrix[K][K], int rank)
   
2. {
   
3. int i, j, k, t;
   
4. double tmp;
   
5.    
   
6.     for(k = 0; k < rank; k++)
   
7.     {
   
8.         for(i = k; i < rank; i++)
   
9.         {
   
10.             tmp = 0;
    
11.             for(t = 0; t < k; t++)  tmp+= matrix[i][t] * matrix[t][k];
    
12.             matrix[i][k]-= tmp;
    
13.         }
    
14. 
    
15.         if(k < rank - 1)
    
16.         {
    
17.             for(j = k + 1; j < rank; j++)
    
18.             {
    
19.                 tmp = 0;
    
20.                 for(t = 0; t < k; t++)  tmp+= matrix[k][t] * matrix[t][j];
    
21.                 matrix[k][j] = (matrix[k][j] - tmp) / matrix[k][k];
    
22.             }
    
23.         }
    
24.         else{}
    
25.     }
    
26.    
    
27. }

复制代码

matrix为所求矩阵，rank为矩阵的阶。
最后计算matrix主对角线上元素的积就是矩阵的行列式了。算这个应该没难度吧。
就是这么简单，不像2楼的代码那样繁琐。并且速度快，数值稳定性强。但是使用这个算法时要注意，矩阵matrix的所有顺序主子阵必须都是非奇异的。否则还是要使用Gauss消去法。
关于矩阵的LU分解，请参看任意一本数值分析教材。

[ 本帖最后由 dragon_ 于 2008-6-16 01:51 编辑 ]

dragon_

突然发现2楼的代码里有这么两行……
        printf("****determinant solver****\n");
        printf("                       by dragon\n");
by dragon……
囧死…………………
那不是我写的……

bobnie

狂谢大哥,以后多多向你学习......

bobnie

高斯消去法的程序已经调试成功，有需要的请留下Email

dragon_

说实话，估计EM学数学的人不多，学过数值分析的人就更少。
把用雅各比法求实对称矩阵的特征值和特征向量的C代码
和求矩阵的协方差矩阵的C代码
给我一份吧。谢谢。

[ 本帖最后由 dragon_ 于 2011-1-14 19:54 编辑 ]

twl859588

C用来干这个有点晕，好像MATLAB更好用...