VC最小二乘法求解方程數(shù)多于未知變量的線性方程組.doc

c語言實(shí)現(xiàn),用最小二乘法求解方程數(shù)多于未知變量的線性方程組的最適解(即矛盾方程組) 收藏 一、代碼：/* * ArgMin.h * Author * Fri Aug 5 17:18:47 2005 * Copyright 2005 * Email */#include stdio.h#include math.h/* Interface */solve the linear equations using argmin methodint ArgMin(double *inMtrx,int COLUMN,int rowNum,int colNum,double *solution);/Solve a series of linear equationsint SolveLinearEqts(double *inMtrx,int COLUMN,int rowNum,int colNum,double *solution);/Elimination of unknownsint ReduceUnknowns(double *mtrx_tmp,int COLUMN,int rowNo,int colNo,int rowNum,int colNum);/* Begin Coding */ArgMin using solvelineareqts. int ArgMin(double *mtrx_tmp,int COLUMN,int rowNum,int colNum,double *solution) int k,l,j; double *eqt; int *q; q=&colNum; eqt=(double *)malloc(colNum*colNum*sizeof(double); /printf(%d,*q); for(k=0;k*q-1;k+) for(l=0;l*q-1;l+) eqtk*COLUMN+l=0; for(j=0;jrowNum;j+) eqtk*COLUMN+l+=*(mtrx_tmp+COLUMN*j+k)*(*(mtrx_tmp+COLUMN*j+l); for(k=0;k*q-1;k+) eqtk*COLUMN+*q-1=0; for(j=0;jrowNum;j+) eqtk*COLUMN+*q-1+=*(mtrx_tmp+COLUMN*j+*q-1)*(*(mtrx_tmp+COLUMN*j+k); /* show the contents of eqts */ /*int m,n; printf(eqt:n); for(m=0;m*q-1;m+) for(n=0;n*q;n+) printf(%.2lft,eqtm*COLUMN+n); printf(n); */ if(!SolveLinearEqts(eqt,COLUMN,(*q-1),*q,solution) return 0; return 1;/Solve the solution of a series of linear equationsint SolveLinearEqts(double *inMtrx,int COLUMN,int rowNum,int colNum,double *solution) int i,j; double tmpSum; if(rowNum!=(colNum-1) printf(Cant solve the equations because equation number ); printf(is not the same as unknow parameters!n); return 0; /reduce unknown parameters for(i=0;i=0;i-) tmpSum=0; for(j=i+1;jrowNum;j+) tmpSum+=*(solution+j)*(*(inMtrx+i*COLUMN+j); if(fabs(inMtrxi*COLUMN+i)0.000001) printf(nNeed more Equations to solve elements in matrix An); printf(Tip: You can try setting a different theta value, n or ); printf(checking the data introduced to functionCalculateCoeff!)n); return 0; *(solution+i)=(*(inMtrx+i*COLUMN+colNum-1)-tmpSum)/(*(inMtrx+i*COLUMN+i); return 1;/*Column number of Array is COLUMN_A, but actually the number of figures in each row is colNum! */int ReduceUnknowns(double *mtrx_tmp,int COLUMN,int rowNo,int colNo,int rowNum,int colNum) int i,j,tmpInt; double tmpDbl1,tmpDbl2; double rowMax,colMax; double *p,*mp; /*Adjust the rowNo whose value is zero in matrix_tmprowNocolNo down to a suitable site*/ p=mtrx_tmp+rowNo*COLUMN+colNo; mp=mtrx_tmp+rowNo*COLUMN+colNo; /select main variable: /*rowMax=fabs(*p); for(j=colNo+1;jcolNum;j+) p+; if(rowMaxfabs(*p) rowMax=fabs(*p); if(rowMax=0) return 0; p+=COLUMN-(colNum-colNo-1); colMax=fabs(*(mtrx_tmp+rowNo*COLUMN+colNo)/rowMax; tmpInt=rowNo; for(i=rowNo+1;irowNum;i+) rowMax=fabs(*p); for(j=colNo+1;jcolNum;j+) p+; if(rowMaxfabs(*p) rowMax=fabs(*p); p+=COLUMN-(colNum-colNo-1); if(rowMax=0) return 0; tmpDbl1=fabs(*(mtrx_tmp+i*COLUMN+colNo)/rowMax; if(colMaxtmpDbl1) colMax=tmpDbl1; tmpInt=i; /change the whole row of tmpInt to rowNo for(j=colNo;jcolNum;j+) tmpDbl1=*(mtrx_tmp+rowNo*COLUMN+j); *(mtrx_tmp+rowNo*COLUMN+j)=*(mtrx_tmp+tmpInt*COLUMN+j); *(mtrx_tmp+tmpInt*COLUMN+j)=tmpDbl1; */ /*Transform the matrix_tmp using linear calculation methods*/ tmpDbl1=*p; p+=COLUMN; for(i=rowNo+1;irowNum;i+) if(fabs(*p)0.000001) continue; tmpDbl2=*p; for(j=colNo;jcolNum;j+) *p=*p-*mp*tmpDbl2/tmpDbl1; p+; mp+; p+=COLUMN-(colNum-colNo); mp-=colNum-colNo; return 1;二、測試/*/ (1)int ArgMin(double *inMtrx,int COLUMN,int rowNum,int colNum,double *solution); 函數(shù)功能：用最小二乘法求解方程數(shù)多于未知變量的線性方程組的最適解。 測試代碼： #define Q 5 #define ROWNUM 6 int main() double solutionQ-1,eqtQQ=0; double aROWNUMQ=3,-2,4,6,-11,4,3,2,9,-2,2,6,8,3,4,2,4,5,3,3, 0,0,8,6,-20,3,4,5,6,0; ArgMin(&a00,Q,ROWNUM,Q,&solution0); for(i=0;iQ-1;i+) printf(%.14lfn,solutioni); return 0; 測試結(jié)果： eqts: 42.000000000 38.000000000 61.000000000 84.000000000 -27.000000000 38.000000000 81.000000000 86.000000000 69.000000000 52.000000000 61.000000000 86.000000000 198.000000000 159.000000000 -161.000000000 84.000000000 69.000000000 159.000000000 207.000000000 -183.000000000 3.00000000000000 2.00000000000000 -1.00000000000000 -2.00000000000000 說明：a, 待求解的方程組數(shù)據(jù)： 3,-2,4,6,-11, 4,3,2,9,-2, 2,6,8,3,4, 2,4,5,3,3, 0,0,8,6,-20, 3,4,5,6,0 以上數(shù)據(jù)中共有四個(gè)變量，六個(gè)方程； b, 測試顯示的數(shù)據(jù)中，“eqts”為通過最小二乘法得出的待求解的線性方程組增廣矩陣；下面四行為計(jì)算求得的解，與實(shí)際完全符合。 (2)int SolveLinearEqts(double *inMtrx,int COLUMN,int rowNum,int colNum,double *solution); 函數(shù)功能：用消元法解NN的線性方程組（即未知數(shù)和方程數(shù)相同的方程組）。 測試代碼： double b45=3,-2,4,6,-11,4,3,2,9,-2,2,6,8,3,4,2,4,5,3,3; double c5; SolveLinearEqts(&b00,5,4,5,&c0); for(i=0;i4;i+) printf(%.14lfn,ci); 測試結(jié)果： 3.00000000000000 2.00000000000000 -1.00000000000000 -2.00000000000000 說明：求得的解與實(shí)際完全符合。/*/ (c)int ArgMin(double *inMtrx,int COLUMN,int rowNum,int colNum,double *solution); 功能說明：用最小二乘法求解方程數(shù)多于未知變量的線性方程組的最適解。 參數(shù)說明： (1)inMtrx: 線性方程組以增廣矩陣的形式用此數(shù)組傳入方程； (2)COLUMN: 數(shù)組inMtrx的列數(shù)為COLUMN； (3)rowNum: 數(shù)組inMtrx中數(shù)據(jù)的行數(shù)； (4)colNum: 數(shù)組inMtrx中數(shù)據(jù)的列數(shù)； (5)solution: 求得方程的解從此參數(shù)傳出； (6)返回值：0,表示方程組無唯一解；1,表示函數(shù)正常結(jié)束。 (d)int SolveLinearEqts(double *inMtrx,int COLUMN,int rowNum,int colNum,double *solution); 功能說明：解NN的線性方程組（即未知數(shù)和方程數(shù)相同的方程組）。 參數(shù)說明： (1)inMtrx: 線性方程組以增廣矩陣的形式用此數(shù)組傳入方程； (2)COLUMN: 數(shù)組inMtrx的列數(shù)為COLUMN； (3)rowNum: 數(shù)組inMtrx中數(shù)據(jù)的行數(shù)； (4)colNum: 數(shù)組inMtrx中數(shù)據(jù)的列數(shù)； (5)solution: 求得方程的解從此參數(shù)傳出； (6)返回值：0,表示方程組無唯一解；1,表示函數(shù)正常結(jié)束。 (e)int ReduceUnknowns(double *mtrx_tmp,int COLUMN,int rowNo,int colNo,int rowNum,int colNum); 功能說明：將儲存在mtrx_tmp中、表示待求解線性方程組的增廣矩陣變換成階梯形矩陣。每調(diào)用函數(shù)一次，通過變換矩陣會將其中元素(rowNo,colNo)對應(yīng)的那一列第rowNo行以下（不包括第rowNo行）的元素全變?yōu)?。 參數(shù)說明： (1)mt