最小二乘法递推最小二乘法C语言源程序

Pn_1[j1] = new double [colPn]; siit0 = new double [colPn]; siit1 = new double [colPn]; /*****************************************************************/ for(j1 = 0; j1 < rowPn; j1++) { for(int j2 = 0; j2 < rowPn; j2++) { if( j1 == j2)Pn[j1][j2] = 1e6; else Pn[j1][j2] = 0; } } for(j1 = 0; j1 < rowPn; j1++) { for(int j2 = 0; j2 < rowPn; j2++) { if( j1 == j2)Pn_1[j1][j2] = 1e6; else Pn_1[j1][j2] = 0; } cout << endl; } for(j1 = 0; j1 < rowPn; j1++) { *(siit1 + j1) = 1e-6; } /******************************************************************/ for(int i = 1; i < 28 ;i++) { En_1 = &a[i][0]; for(j1 = 0; j1 < rowPn; j1++) { for(int j2 = 0; j2 < rowPn; j2++) { Pn[j1][j2] = Pn_1[j1][j2]; } } for(j1 = 0; j1 < rowPn; j1++) { *(Pn_En_1 + j1) = 0; for(int j2 = 0; j2 < rowPn; j2++) { *(Pn_En_1 + j1) += Pn[j1][j2]*(*(En_1 + j2)); } }
wk.baidu.com
//ofstream outFile; //cout << "输入的原矩阵显示如下：" << endl; //outFile.open("ym.txt"); for(i1 = 0; i1 < row1; i1++) for(i2 = 0; i2 < col1; i2++) { cout<<x[i1][i2]<<'\t'; if(col1-1==i2)cout << endl; } //outFile.close; } /*校验矩阵*/ void JY(void) { d = new double *[col]; for(int i3 = 0; i3 < col; i3++) d[i3] = new double [col+col]; double *p1 = *d; double *p2 = *c; for(int i4 = 0; i4 < col; i4++) { p1 = &d[i4][0]; p2 = &c[i4][0]; for(int i5 = 0; i5 < col*2; i5++) { if(i5 < col) *p1 = *p2; else { if(i4 == i5-col) *p1 = 1; else *p1 = 0; } p1++; p2++; } } } /*求逆矩阵第一步*/ void cov1(void) { double *p1 = *d; dou
{ pfai = &fai[j][0]; *pfai = -(*(pyk+1)); *(pfai+1) = -(*pyk); *(pfai+2) = *(puk+1); *(pfai+3) = *puk; pyk++; puk++; } }
//iterative least square method #define rowPn 4 #define colPn 9 void iterative_least_square_method(void) { double **Pn,*En_1,*Pn_En_1,*En_Pn_1,**PE,**Pn_1,**Pn_En_0; double En_Pn_En_1 = 0; double *pPn1; double *pPn2; double *Kn; double *siit0,En_siit0,*En_0,BB,*siit1; Pn = new double *[rowPn]; for(int j1 = 0; j1 < rowPn; j1++) Pn[j1] = new double [rowPn]; Pn_En_1 = new double [rowPn]; En_1 = new double [colPn]; En_1 = &a[colPn][0]; Pn_En_1 = new double [colPn]; Pn_En_0 = new double *[rowPn]; for(j1 = 0; j1 < rowPn; j1++) Pn_En_0[j1] = new double [colPn]; En_Pn_1 = new double [colPn]; Kn = new double [colPn]; PE = new double *[colPn]; for (j1 = 0; j1 < colPn; j1++) PE[j1] = new double [colPn]; Pn_1 = new double *[colPn]; for(j1 = 0; j1 < colPn; j1++)
cout << a[i1][i2] <<'\t'; } cout << endl; } } /*原矩阵的转置*/ void inline Tf(void) { b = new double *[col]; for(i = 0;i < col; i++) b[i] = new double[row]; for(i1 = 0; i1 < col; i1++) for(i2 = 0; i2 < row; i2++) b[i1][i2]=a[i2][i1]; } /*原矩阵的转置与原矩阵相乘*/ void mutply(double **x,double **y) { double *p1 = *x; double *p2 = *y; c = new double *[col]; for(i = 0;i < col; i++) c[i] = new double[col]; for(i1 = 0; i1 < col; i1++) for(i2 = 0; i2 < col; i2++) { c[i1][i2] = 0; p1 = &x[i1][0]; p2 = &y[i2][0]; for(int i3 = 0; i3 < row; i3++) { c[i1][i2] += (*p1)*(*p2); p1++; p2++; } } } /*显示*/ void display(double **x,int row1,int col1) {
//least square method #define eps 1e-10 #define row 29 #define col 4 int i1,i2; double **a,**b,**c,**d,**e,**f,*g; int i; /*绝对值函数*/ double abs(double x) { if(x<eps)x=-x; else x=x; return x; } /*输入原矩阵*/ void inline inputf(void) { //cin >> row; //cin >> col; a = new double *[row]; for(i = 0; i < row; i++) a[i] = new double[col]; for(i1 = 0; i1 < row; i1++) { for(i2 = 0; i2 < col; i2++)//cin >> a[i1][i2]; { a[i1][i2] = fai[i1][i2];
} BB = yk[3+i] - En_siit0; for(j1 = 0; j1 < rowPn; j1++) { *(siit1 + j1) = *(siit0 + j1) + (*(Kn + j1))*BB; cout << *(siit1 + j1) << '\t'; } cout<<endl; } }
for(j1 = 0; j1 < rowPn; j1++) { *(En_Pn_1 + j1) = 0;
for(int j2 = 0; j2 < rowPn; j2++) { *(En_Pn_1 + j1) += *(En_1 + j2)*Pn[j2][j1]; } } En_Pn_En_1 = 0; for(j1 = 0; j1 < rowPn; j1++) En_Pn_En_1 += En_Pn_1[j1]*(*(En_1 + j1)); En_Pn_En_1 += 1; for(j1 = 0; j1 < rowPn; j1++) { *(Kn + j1) = *(Pn_En_1 + j1)/En_Pn_En_1; } for(j1 = 0; j1 < rowPn; j1++) { for(int j2 = 0; j2 < colPn; j2++) { PE[j1][j2] = 0; PE[j1][j2] = Kn[j1]*En_Pn_1[j2]; } } for(j1 = 0; j1 < rowPn; j1++) { for(int j2 = 0; j2 < rowPn; j2++) { Pn_1[j1][j2] = Pn[j1][j2] - PE[j1][j2]; } } for(j1 = 0; j1 < rowPn; j1++) { siit0[j1] = siit1[j1]; } En_siit0 = 0; for(j1 = 0; j1 < rowPn; j1++) { En_siit0 += *(En_1 + j1)*(*(siit0 + j1));
最小二乘法,递推最小二乘法 C 语言源程序
//init input #define row_y 29 #define col_y 4 double uk[31]={ 0,1.147,0.201,-0.787,-1.589,-1.052,0.866,1.152,1.573,0.626,0.433, -0.958,0.810,-0.044,0.947,-1.474,-0.719,-0.086,-1.099,1.450,1.151, 0.485,1.633,0.043,1.326,1.706,-0.340,0.890,1.144,1.177,-0.390}; double yk[31]; double **fai; void initinput(void) { //ifstream inFile; int temp = 0; ofstream outFile; outFile.open("yk.txt"); yk[0] = yk[1] = yk[2] = 0; for(int i = 3; i < 32; i++)yk[i]=-1.642*yk[i-1]-0.715*yk[i-2]+0.39*uk[i-1]+0.35*uk[i-2]; for(i = 1; i < 32; i++) { temp++; outFile<<yk[i]<<endl; //if(i%5==0)outFile<<endl; } outFile<<temp<<endl; outFile.close; cout << "Over!" <<endl; //inFile.open("yk.txt"); fai = new double *[row_y]; for(i = 0; i < row_y; i++) fai[i] =new double [col_y]; double *pyk = yk; double *puk = uk; double *pfai = *fai; pyk = &yk[1]; puk = &uk[1]; for(int j = 0; j < row_y; j++)