1、最小二乘法：1）（用1TAA方法计算逆矩阵）#include <stdio.h>#include <math.h>#include <stdlib.h>#include<string.h>#include <iostream.h>#define N 200#define n 9void Getdata(double sun[N])//从txt文档中读取数据（小数）{char data;char sunpot[10]={0000000000};//为防止结果出现‘烫’字int i=0,j=0;double d;FILE *fp=fopen("新建文本文档.txt","r");if(!fp){printf("can't open file\n");}while(!feof(fp)){data=fgetc(fp);if(data!='\n'){sunpot[i]=data;i++;}else if(data=='\n'){sunpot[i]='\0';//给定结束符d=atof(sunpot);//将字符串转换成浮点数sun[j]=d;j++;i=0;//将i复位}}}void Normal(double sun[N],double sun1[N])//将数据进行标准化{double mean,temp=0,variance=0;int i;for(i=0;i<N;i++)temp=temp+sun[i];mean=temp/N;for(i=0;i<N;i++){sun1[i]=sun[i]-mean;}for(i=0;i<N;i++){variance=variance+sun1[i]*sun1[i];}variance=variance/N;variance=sqrt(variance);for(i=0;i<N;i++){sun1[i]=sun[i]/variance;//减去均值除以均方差进行归一化}}void Matrix(double sun[N],double matrixl[n][n],double matrixr[n])//组建方程的左右两个矩阵{double matrix1[N-n][n];//方程左边系数矩阵double matrix_trans[n][N-n];double matrix2[N-n];//方程右边系数矩阵int i,j,k;double temp;for(i=0;i<N-n;i++)//组建N-nxn维方程矩阵{for(j=0;j<n;j++){matrix1[i][j]=sun[n+i-j-1];//此处与matlab不同，matlab是从1开始的}}for(i=0;i<N-n;i++)matrix2[i]=sun[i+n];for(i=0;i<n;i++)//将N-nxn维矩阵转置{for(j=0;j<N-n;j++){matrix_trans[i][j]=matrix1[j][i];}}for(i=0;i<n;i++)//组建nxn维方阵(A'*A){for(k=0;k<n;k++){temp=0;for(j=0;j<N-n;j++){temp+=matrix_trans[i][j]*matrix1[j][k];}matrixl[i][k]=temp;}}for(i=0;i<n;i++){temp=0;for(j=0;j<N-n;j++){temp+=matrix_trans[i][j]*matrix2[j];}matrixr[i]=temp;}}void ExchangeLine(double a[n][n], int i, int j, int p)//交换矩阵第i行和第j行{int k;double temp;for (k=0;k<=p-1;k++){temp=a[j][k];a[j][k]=a[i][k];a[i][k]=temp;}}void Gauss(double a[n][n],int k,int p)//高斯消元,用第k行消去其他行第一个非零元素{int i,j;double m;for(i=k+1;i<p;i++){if(a[i][k]==0)i++;m = a[i][k] / a[k][k];for (j = k; j < p; j ++){a[i][j] = a[i][j] - m * a[k][j];}}}double Determinant(double a[n][n],int p)//求a的行列式{int k,i,j;double result=1.0;for(k=0;k<p-1;k++)//逐一取对角元素以下矩阵进行消元，n-1是因为最后一个对角元素无需消元{i=k;j=k;while(a[i][j]==0)i++;if((i<p)&&i>k){ExchangeLine(a,k,i,p);result=(-1)*result;//换行之后行列式变号一次}else if(i>=p)//即找不到可以交换的行{result=0;break;}Gauss(a,k,p);//对该列进行消元}for(i=0;i<p;i++){result=result*a[i][i];}return result;}void Adjoint_matrix(double a[n][n],double b[n][n],int p,int i,int j)//求a[i][j]元素的伴随矩阵{double temp[n][n]={{0,0,0},{0,0,0},{0,0,0}};//实际上只要n-1xn-1维，但为了调用方便，此处写为nxn 维int k,m,l,r;b[j][i]=1;for (k = 0, l = 0; k <p-1; k++, l ++)//将原矩阵去除i行j列后复制到temp中{if (l == i) l ++;for (m = 0, r = 0; m <p-1; m ++, r ++){if (r == j) r ++;temp[k][m] = a[l][r];}}b[j][i]=b[j][i]*Determinant(temp,n-1);if ((i + j) % 2 == 1)//(-1)^(i+j)b[j][i]=-b[j][i];}void Inverse(double a[n][n],double b,double c[n][n])//求逆最后一步运算{int i,j;for(i=0;i<n;i++){for(j=0;j<n;j++){c[i][j]=a[i][j]/b;}}}void main(){int p=n,i,j;double sunpot[N],sunpot_normal[N];double fi[n];double matrix[n][n];//要换成返回的方阵double matrixr[n];double adjoint[n][n],inverse[n][n];double determinant,temp;Getdata(sunpot);//从txt文档中读取每年的太阳黑子数据Normal(sunpot,sunpot_normal);//数据归一化Matrix(sunpot_normal,matrix,matrixr);//组建方程for(i=0;i<p;i++){for(j=0;j<p;j++){Adjoint_matrix(matrix,adjoint,p,i,j);}}determinant=Determinant(matrix,p);if (determinant== 0)printf("该矩阵不存在逆矩阵！\n");printf("行列式：\n");printf("%10lf\n",determinant);Inverse(adjoint,determinant,inverse);printf("逆矩阵：\n");for(i=0;i<p;i++){for(j=0;j<p;j++){printf("%12lf",inverse[i][j]);}printf("\n");}printf("fi：\n");for(i=0;i<n;i++){temp=0;for(j=0;j<n;j++){temp+=inverse[i][j]*matrixr[j];}fi[i]=temp;printf("%10lf",fi[i]);}printf("\n");}计算结果：fi：1.327258 -0.622683 0.018571 0.129074 -0.141278 0.103665 -0.046903 0.060403 0.154991下面用matlab进行拟合度计算：clearM=load('太阳黑子数.txt');sunpot=(M(:,2));N=200;n=9;[Y,mu,sigma]=zscore(sunpot);B(1:N-n)=Y(n+1:N);B=B';fi=[1.327258 -0.622683 0.018571 0.129074 -0.141278 0.103665 -0.046903 0.060403 0.154991]; %由C语言计算出的fi值fi=[1;-fi'];ts2=idpoly(fi',[]);B=Y(1:N-n);plot(B);compare(B, ts2, 1);2）用构造矩阵1A I I A-→的方法求逆：（列主元高斯消去法）[;][;]部分主要程序：void BuildMatrix(double a[n][M],double b[n][n]){int i,j;for (i=0;i<n;i++){for (j=0;j<n;j++){a[i][j]=b[i][j];}a[i][i+n]=1.0;}}void ExchangeLine(double a[n][M],int i,int j){int k;double temp;for (k=0;k<M;k++){temp=a[i][k];a[i][k]=a[j][k];a[j][k]=temp;}}void Judge(double a[n][M],int i){int j,m=0;//m=0来表示是否需要换行double temp=a[i][i];for(j=i+1;j<n;j++){if(fabs(a[j][i])>fabs(temp)){temp=a[j][i];m=j;}}if(m)ExchangeLine(a,i,m);}void Gauss(double a[n][M],int k) //高斯消元法,用第k行消去其他行{ { int i,j;double m;for(i=k+1;i<n;i++){if(a[i][k]==0)i++;m = a[i][k] / a[k][k];for (j = k; j < M; j ++){a[i][j] = a[i][j] - m * a[k][j];}}}void Normal(double a[n][M]){int i,j;double temp;for(i=0;i<n;i++){temp=a[i][i];if(temp!=1.0){for(j=i;j<M;j++){a[i][j]=a[i][j]/temp;}}}}void Inverse(double a[n][M]){int i,j,k;double temp;for(i=n-1;i>=0;i--){for(j=i-1;j>=0;j--){temp=a[j][i];for(k=i;k<M;k++){a[j][k]=a[j][k]-temp*a[i][k];}}}}void GetInverse(double a[n][M],double b[n][n]) {int i,j;for(i=0;i<n;i++){for(j=0;j<n;j++){b[i][j]=a[i][j+n];}}}void main(){int i,j;double temp;double sunpot[N];double sunpot_normal[N];double matrixr[n];double matrix1[n][n]={0};//方阵double matrix[n][M]={0};//构造矩阵double inverse[n][n]={0};double fi[n];Getdata(sunpot);//从txt文档中读取每年的太阳黑子数据Normal(sunpot,sunpot_normal);//数据归一化Matrix(sunpot_normal,matrix1,matrixr);//组建方程BuildMatrix(matrix,matrix1);for(i=0;i<n-1;i++){Judge(matrix,i);Gauss(matrix,i);}for(i=0;i<n;i++)//判断是否可逆{if(matrix[i][i]==0){printf("矩阵不可逆！\n");exit(-1);}}Normal(matrix);Inverse(matrix);GetInverse(matrix,inverse);printf("逆矩阵：\n");for(i=0;i<n;i++){for(j=0;j<n;j++){printf("%10lf",inverse[i][j]);}printf("\n");}printf("fi：\n");for(i=0;i<n;i++){temp=0;for(j=0;j<n;j++){temp+=inverse[i][j]*matrixr[j];}fi[i]=temp;printf("%10lf",fi[i]);}}。