當前位置:編程學習大全網 - 編程語言 - 求C或C++語言編寫的用最小二乘法進行曲線擬合

求C或C++語言編寫的用最小二乘法進行曲線擬合

妳的近似解析表達式為y=at+bt^2+ct^2

是不是想寫成為y=at+bt^2+ct^3

但是實際擬合出來的表達式為y=a[3]+a[2]t+a[1]t^2+a[0]t^3會有個常數項的。

簡單的講,所謂擬合是指已知某函數的若幹離散函數值{f1,f2,…,fn},通過調整該函數中若幹待定系數f(λ1, λ2,…,λ3), 使得該函數與已知點集的差別(最小二乘意義)最小。如果待定函數是線性,就叫線性擬合或者線性回歸(主要在統計中),否則叫作非線性擬合或者非線性回歸。表達式也可以是分段函數,這種情況下叫作樣條擬合。

曲線擬合:

#include <stdio.h>

#include <stdlib.h>

#include <malloc.h>

#include <math.h>

Smooth(double *x,double *y,double *a,int n,int m,double *dt1,double *dt2,double *dt3);

void main()

{

int i ,n ,m ;

double *x,*y,*a,dt1,dt2,dt3,b;

n = 12;// 12個樣點

m = 4; //3次多項式擬合

b = 0; //x的初值為0

/*分別為x,y,a分配存貯空間*/

x = (double *)calloc(n,sizeof(double));

if(x == NULL)

{

printf("內存分配失敗\n");

exit (0);

}

y = (double *)calloc(n,sizeof(double));

if(y == NULL)

{

printf("內存分配失敗\n");

exit (0);

}

a = (double *)calloc(n,sizeof(double));

if(a == NULL)

{

printf("內存分配失敗\n");

exit (0);

}

for(i=1;i<=n;i++)

{

x[i-1]=b+(i-1)*5;

/*每隔5取壹個點,這樣連續取12個點*/

}

y[0]=0;

y[1]=1.27;

y[2]=2.16;

y[3]=2.86;

y[4]=3.44;

y[5]=3.87;

y[6]=4.15;

y[7]=4.37;

y[8]=4.51;

y[9]=4.58;

y[10]=4.02;

y[11]=4.64;

/*x[i-1]點對應的y值是擬合已知值*/

Smooth(x,y,a,n,m,&dt1,&dt2,&dt3); /*調用擬合函數*/

for(i=1;i<=m;i++)

printf("a[%d] = %.10f\n",(i-1),a[i-1]);

printf("擬合多項式與數據點偏差的平方和為:\n");

printf("%.10e\n",dt1);

printf("擬合多項式與數據點偏差的絕對值之和為:\n");

printf("%.10e\n",dt2);

printf("擬合多項式與數據點偏差的絕對值最大值為:\n");

printf("%.10e\n",dt3);

free(x); /*釋放存儲空間*/

free(y); /*釋放存儲空間*/

free(a); /*釋放存儲空間*/

}

Smooth(double *x,double *y,double *a,int n,int m,double *dt1,double *dt2,double *dt3)//(x,y,a,n,m,dt1,dt2,dt3 )

//double *x; /*實型壹維數組,輸入參數,存放節點的xi值*/

//double *y; /*實型壹維數組,輸入參數,存放節點的yi值*/

//double *a; /*雙精度實型壹維數組,長度為m。返回m壹1次擬合多項式的m個系數*/

//int n; /*整型變量,輸入參數,給定數據點的個數*/

//int m; /*整型變量,輸入參數,擬合多項式的項數*/

//double *dt1; /*實型變量,輸出參數,擬合多項式與數據點偏差的平方和*/

//double *dt2; /*實型變量,輸出參數,擬合多項式與數據點偏差的絕對值之和*/

//double *dt3; /*實型變量,輸出參數,擬合多項式與數據點偏差的絕對值最大值*/

{

int i ,j ,k ;

double *s,*t,*b,z,d1,p,c,d2,g,q,dt;

/*分別為s ,t ,b分配存貯空間*/

s = (double *)calloc(n,sizeof(double));

if(s == NULL)

{

printf("內存分配失敗\n");

exit (0);

}

t = (double *)calloc(n,sizeof(double));

if(t == NULL)

{

printf("內存分配失敗\n");

exit (0);

}

b = (double *)calloc(n,sizeof(double));

if(b == NULL)

{

printf("內存分配失敗\n");

exit (0);

}

z = 0;

for(i=1;i<=n;i++)

z=z+x[i-1]/n; /*z為各個x的平均值*/

b[0]=1;

d1=n;

p=0;

c=0;

for(i=1;i<=n;i++)

{

p=p+x[i-1]-z;

c=c+y[i-1];

}

c=c/d1;

p=p/d1;

a[0]=c*b[0];

if(m>1)

{

t[1]=1;

t[0]=-p;

d2=0;

c=0;

g=0;

for(i=1;i<=n;i++)

{

q=x[i-1]-z-p;

d2=d2+q*q;

c=y[i-1]*q+c;

g=(x[i-1]-z)*q*q+g;

}

c=c/d2;

p=g/d2;

q=d2/d1;

d1=d2;

a[1]=c*t[1];

a[0]=c*t[0]+a[0];

}

for(j=3;j<=m;j++)

{

s[j-1]=t[j-2];

s[j-2]=-p*t[j-2]+t[j-3];

if(j>=4)

for(k=j-2;k>=2;k--)

s[k-1]=-p*t[k-1]+t[k-2]-q*b[k-1];

s[0]=-p*t[0]-q*b[0];

d2=0;

c=0;

g=0;

for(i=1;i<=n;i++)

{

q=s[j-1];

for(k=j-1;k>=1;k--)

q=q*(x[i-1]-z)+s[k-1];

d2=d2+q*q;

c=y[i-1]*q+c;

g=(x[i-1]-z)*q*q+g;

}

c=c/d2;

p=g/d2;

q=d2/d1;

d1=d2;

a[j-1]=c*s[j-1];

t[j-1]=s[j-1];

for(k=j-1;k>=1;k--)

{

a[k-1]=c*s[k-1]+a[k-1];

b[k-1]=t[k-1];

t[k-1]=s[k-1];

}

}

*dt1=0;

*dt2=0;

*dt3=0;

for(i=1;i<=n;i++)

{

q=a[m-1];

for(k=m-1;k>=1;k--)

q=q*(x[i-1]-z)+a[k-1];

dt=q-y[i-1];

if(fabs(dt)>*dt3)

*dt3=fabs(dt);

*dt1=*dt1+dt*dt;

*dt2=*dt2+fabs(dt);

}

/*釋放存儲空間*/

free(s);

free(t);

free(b);

return(1);

}

  • 上一篇:使用微信小程序應用場景有哪些?
  • 下一篇:常州大學機器人產業學院好不好
  • copyright 2024編程學習大全網