最小二乘法

示例源代码 test.cpp:  

#include <iostream>
using namespace std;

float gs(float a[10][10],float b[10],int n )
{
    int i,j,k,l;
    float s;
    k=1;
    
    while(k!=n+1)
    {
        if(a[k][k]!=0)
        {
             for(i=k+1;i<=n+1;i++)
             {
                  a[i][k]=a[i][k]/a[k][k];
                  b[i]=b[i]-a[i][k]*b[k];
   
                  for(j=k+1;j<=n+1;j++)
                      a[i][j]=a[i][j]-a[i][k]*a[k][j];
             }
        }
        k=k+1;
    }

    for(k=n+1;k>=1;k--)
    {
         s=0;
         for(l=k+1;l<=n+1;l++)
              s=s+a[k][l]*b[l];
              b[k]=(b[k]-s)/a[k][k];
    }

    return 0;
}

void generateC(float c[10][10],float x[10],int m,int n)//生成中间矩阵C
{
     for(int i=1;i<=m;i++)
     {
          c[i][1]=1;
 
          for(int j=2;j<=n+1;j++)
              c[i][j]=x[i]*c[i][j-1];
     }
}

void juzhenT(float ct[10][10],float c[10][10],int m,int n)//生成矩阵C的转置矩阵
{
     for(int i=1;i<=n+1;i++)
          for(int j=1;j<=m;j++)
               ct[i][j]=c[j][i];
}

int main(void)
{
    float x[10],y[10],c[10][10],ct[10][10],A[10][10],B[10];
    int m,n,i,j;
    cout<<"请输入样本数据的个数和未知数的最高次方数"<<endl;
    cin>>m>>n;
    cout<<"x[i]"<<endl;
    for(i=1;i<=m;i++)
        cin>>x[i];
    cout<<"y[i]"<<endl;
    for(i=1;i<=m;i++)
        cin>>y[i];

    generateC(c,x,m,n);
    juzhenT(ct,c,m,n);
    for(i=1;i<=n+1;i++)
        for(int k=1;k<=n+1;k++)
        {
            A[i][k]=0;
            for(j=1;j<=m;j++)
            {
                A[i][k]=A[i][k]+ct[i][j]*c[j][k];
            }
        }

    for(i=1;i<=n+1;i++)
    {
        B[i]=0;
        for(j=1;j<=m;j++)
            B[i]=B[i]+ct[i][j]*y[j];
    }

    gs(A,B,n );
    for(i=1;i<=n+1;i++)
        cout<<B[i]<<endl;
    cout<<"这组数据的近似函数关系为：y="<<B[1]<<"+"<<B[2]<<"x"<<endl;
    return 0;
}

Makefile:

test:
	g++ test.cpp -o test
clean:
	rm -rf test

编译并运行过程：

$ make
g++ test.cpp -o test
$ ./test 
请输入样本数据的个数和未知数的最高次方数
3
2
x[i]
4
9
16
y[i]
36
88
122
-22.2287
16.4048
-0.461907
这组数据的近似函数关系为：y=-22.2287+16.4048x