最小二乘法 (C语言)

#include

#include

//////////////////////////////////////////////////////////////////////////// //////////////

//矩阵结构体

struct Matrix

{

int m,n;//m为行数,n为列数

double **pm;//指向矩阵二维数组的指针

};

//初始化矩阵mt，并置矩阵的行为m，列为n

void InitMatrix(struct Matrix *mt,int m,int n)

{

int i;

//指定矩阵的行和列

mt->m=m;

mt->n=n;

//为矩阵分配内存

mt->pm=new double *[m];

for(i=0;i

{

mt->pm[i]=new double[n];

}

}

//用0初始化矩阵mt，并置矩阵的行为m，列为n

void InitMatrix0(struct Matrix *mt,int m,int n)

{

int i,j;

InitMatrix(mt,m,n);

for(i=0;i

for(j=0;jpm[i][j]=0.0;

}

//销毁矩阵mt

void DestroyMatrix(struct Matrix *mt)

{

int i;

//释放矩阵内存

for(i=0;im;i++)

{

delete []mt->pm[i];

}

delete []mt->pm;

}

//矩阵相乘mt3=mt1*mt2

//成功返回1，失败返回0

int MatrixMul(Matrix *mt1,Matrix *mt2,Matrix *mt3)

{

int i,j,k;

double s;

//检查行列号是否匹配

if(mt1->n!=mt2->m||mt1->m!=mt3->m||mt2->n!=mt3->n) retu rn 0;

for(i=0;im;i++)

for(j=0;jn;j++)

{

s=0.0;

for(k=0;kn;k++) s=s+mt1->pm[i][k]*mt2->pm[k][j]; mt3->pm[i][j]=s;

}

return 1;

}

//矩阵转置mt2=T(mt1)

//成功返回1，失败返回0

int MatrixTran(Matrix *mt1,Matrix *mt2)

{

int i,j;

//检查行列号是否匹配

if(mt1->m!=mt2->n||mt1->n!=mt2->m) return 0;

for(i=0;im;i++)

for(j=0;jn;j++) mt2->pm[j][i]=mt1->pm[i][j];

return 1;

}

//控制台显示矩阵mt

void ConsoleShowMatrix(Matrix *mt)

{

int i,j;

for(i=0;im;i++)

{

printf("\n");

for(j=0;jn;j++) printf("%2.4f ",mt->pm[i][j]);

}

printf("\n");

}

//////////////////////////////////////////////////////////////////////////// //////////////

//Guass列主元素消元法求解方程Ax=b,a=(A,b)

int Guassl(double **a,double *x,int n)

{

int i,j,k,numl,*h,t;

double *l,max;

l=new double[n];

h=new int[n];

for(i=0;i

for(i=1;i

{

max=fabs(a[h[i-1]][i-1]);

numl=i-1;

//列元的最大值

for(j=i;j

{

if(fabs(a[h[j]][i-1])>max)

{

numl=h[j];

max=fabs(a[h[j]][i-1]);

}

}

if(max<0.00000000001) return 0;

//交换行号

if(numl>i-1)

{

t=h[i];

h[i]=h[numl];

h[numl]=t;

}

for(j=i;j

for(j=i;j

for(k=i;k

}