C语言实现矩阵的秩求解分析

求解一个$m\times n$阶实矩阵$A$的秩。从数学角度来说矩阵的秩就是矩阵中线性无关的行向量的最大个数，并且初等行变换和列变换均不会影响矩阵的秩的大小，为了实现这一过程常使用全选的高斯消去法将矩阵变为上三角矩阵，则其中非零的行向量个数就是矩阵的秩。

#include "stdio.h"
#include "stdlib.h"
#include "math.h"

double mat[100], eps;
int m,n;
int r_mrank(double *mat,int m, int n,double eps){
    int i,j,k,p,is,js,l,v,rank;
    double tmp, pivot;
    double *cpmat;
    if(mat == NULL)
    {
      printf("matrix pointer is Null.\n");
      return(0);
    }
    cpmat = (double*)malloc(m*n*sizeof(double)); 
    for(i=0; i<m*n; i++)
      cpmat[i] = mat[i];
    p = m<n?m:n;
    rank = 0;
    for(k=0; k<p; k++)
    {
	pivot = 0.0;
	for(i=k; i<m; i++)
	  for(j=k; j<n; j++)
	  {
	    tmp = fabs(cpmat[i*n+j]);
	    if(tmp > pivot)
	    {
	      pivot = tmp;
	      is = i;
	      js = j;
	    }
	  }
	if(pivot < eps)	
	  return(rank);
	rank++;	
	if(is != k)
	  for(j=k; j<n; j++)
	  {
	    l = k*n + j;
	    v = is*n + j;
	    tmp = cpmat[l];
	    cpmat[l] = cpmat[v];
	    cpmat[v] = tmp;
	  }
	if(js != k)
	  for(i=k; i<m; i++)
	  {
	    l = i*n + k;
	    v = i*n + js;
	    tmp = cpmat[l];
	    cpmat[l] = cpmat[v];
	    cpmat[v] = tmp;
	  }
	for(i=k+1; i<m; i++)
	    tmp = cpmat[i*n+k]/cpmat[k*n+k];
	    for(j=k+1; j<n; j++)
		cpmat[i*n+j] -= tmp*cpmat[k*n+j];
	}
    }
    free(cpmat);
    return(rank);
}

int main(){
	int count, rank, i ,j;
	scanf("%d %d %lf",&m,&n,&eps);
	printf("Input the matric please: \n");
	for(count=0;count<m*n;count++){
		scanf("%lf", mat+count);
	}
	printf("The matrix you input is \n");
	for(i=0;i<m;i++){
		for(j=0;j<n;j++)
			printf("%lf  ", *(mat+i*n + j));
		printf("\n");
	}
	rank = r_mrank(mat,m,n,eps);
	printf("the rank of the matrix is: %d", rank);
	return 0;
}