#include <conio.h>
#include <math.h>
#include <process.h>
#define N 5//N个点
#define T 3 //T次拟合
#define W 1//权函数
#define PRECISION 0.00001
float pow_n(float a,int n)
{
int i
if(n==0)
return(1)
float res=a
for(i=1i<ni++)
{
res*=a
}
return(res)
}
void mutiple(float a[][N],float b[][T+1],float c[][T+1])
{
float res=0
int i,j,k
for(i=0i<T+1i++)
for(j=0j<T+1j++)
{
res=0
for(k=0k<Nk++)
{
res+=a[i][k]*b[k][j]
c[i][j]=res
}
}
}
void matrix_trans(float a[][T+1],float b[][N])
{
int i,j
for(i=0i<Ni++)
{
for(j=0j<T+1j++)
{
b[j][i]=a[i][j]
}
}
}
void init(float x_y[][2],int n)
{
int i
printf("请输入%d个已知点：\n",N)
for(i=0i<ni++)
{
printf("(x%d y%d):",i,i)
scanf("%f %f",&x_y[i][0],&x_y[i][1])
}
}
void get_A(float matrix_A[][T+1],float x_y[][2],int n)
{
int i,j
for(i=0i<Ni++)
{
for(j=0j<T+1j++)
{
matrix_A[i][j]=W*pow_n(x_y[i][0],j)
}
}
}
void print_array(float array[][T+1],int n)
{
int i,j
for(i=0i<ni++)
{
for(j=0j<T+1j++)
{
printf("%-g",array[i][j])
}
printf("\n")
}
}
void convert(float argu[][T+2],int n)
{
int i,j,k,p,t
float rate,temp
for(i=1i<ni++)
{
for(j=ij<nj++)
{
if(argu[i-1][i-1]==0)
{
for(p=ip<np++)
{
if(argu[p][i-1]!=0)
break
}
if(p==n)
{
printf("方程组无解!\n")
exit(0)
}
for(t=0t<n+1t++)
{
temp=argu[i-1][t]
argu[i-1][t]=argu[p][t]
argu[p][t]=temp
}
}
rate=argu[j][i-1]/argu[i-1][i-1]
for(k=i-1k<n+1k++)
{
argu[j][k]-=argu[i-1][k]*rate
if(fabs(argu[j][k])<=PRECISION)
argu[j][k]=0
}
}
}
}
void compute(float argu[][T+2],int n,float root[])
{
int i,j
float temp
for(i=n-1i>=0i--)
{
temp=argu[i][n]
for(j=n-1j>ij--)
{
temp-=argu[i][j]*root[j]
}
root[i]=temp/argu[i][i]
}
}
void get_y(float trans_A[][N],float x_y[][2],float y[],int n)
{
int i,j
float temp
for(i=0i<ni++)
{
temp=0
for(j=0j<Nj++)
{
temp+=trans_A[i][j]*x_y[j][1]
}
y[i]=temp
}
}
void cons_formula(float coef_A[][T+1],float y[],float coef_form[][T+2])
{
int i,j
for(i=0i<T+1i++)
{
for(j=0j<T+2j++)
{
if(j==T+1)
coef_form[i][j]=y[i]
else
coef_form[i][j]=coef_A[i][j]
}
}
}
void print_root(float a[],int n)
{
int i,j
printf("%d个点的%d次拟合的多项式系数为:\n",N,T)
for(i=0i<ni++)
{
printf("a[%d]=%g,",i+1,a[i])
}
printf("\n")
printf("拟合曲线方程为:\ny(x)=%g",a[0])
for(i=1i<ni++)
{
printf(" + %g",a[i])
for(j=0j<ij++)
{
printf("*X")
}
}
printf("\n")
}
void process()
{
float x_y[N][2],matrix_A[N][T+1],trans_A[T+1][N],coef_A[T+1][T+1],coef_formu[T+1][T+2],y[T+1],a[T+1]
init(x_y,N)
get_A(matrix_A,x_y,N)
printf("矩阵A为：\n")
print_array(matrix_A,N)
matrix_trans(matrix_A,trans_A)
mutiple(trans_A,matrix_A,coef_A)
printf("法矩阵为：\n")
print_array(coef_A,T+1)
get_y(trans_A,x_y,y,T+1)
cons_formula(coef_A,y,coef_formu)
convert(coef_formu,T+1)
compute(coef_formu,T+1,a)
print_root(a,T+1)
}
void main()
{
process()
}
]]>
</Content>
<PostDateTime>2007-4-19 19:23:57</PostDateTime>
</Reply>
<Reply>
<PostUserNickName></PostUserNickName>
<rank>一级(初级)</rank>
<ranknum>user1</ranknum>
<credit>100</credit>
<ReplyID>40389872</ReplyID>
<TopicID>5478010</TopicID>
<PostUserId>1526752</PostUserId>
<PostUserName>jiangxc2004</PostUserName>
<Point>0</Point>
<Content>
<![CDATA[
你可以改一下
不从终端输入,直接在程序中给出参数
请输入5个已知点：
(x0 y0):-2 -0.1
(x1 y1):-1 0.1
(x2 y2):0 0.4
(x3 y3):1 0.9
(x4 y4):2 1.6
矩阵A为：
1 -2 4 -8
1 -1 1 -1
1 0 0 0
1 1 1 1
1 2 4 8
法矩阵为：
5 0 10 0
0 10 0 34
10 0 34 0
0 34 0 130
5个点的3次拟合的多项式系数为:
a[1]=0.408571, a[2]=0.391667, a[3]=0.0857143, a[4]=0.00833333,
拟合曲线方程为:
y(x)=0.408571 + 0.391667*X + 0.0857143*X*X + 0.00833333*X*X*X
]]>
</Content>
<PostDateTime>2007-4-19 19:26:11</PostDateTime>
</Reply>
<Reply>
<PostUserNickName></PostUserNickName>
<rank>一级(初级)</rank>
<ranknum>user1</ranknum>
<credit>100</credit>
<ReplyID>40390406</ReplyID>
<TopicID>5478010</TopicID>
<PostUserId>1526752</PostUserId>
<PostUserName>jiangxc2004</PostUserName>
<Point>0</Point>
<Content>
<![CDATA[
这样就可以直接调用process()函数了!
二次拟合的话就把宏 T 成2
拟合点的数目 N 也可以修改!
也可以去到注释的部分进行返回值的调用!
最小二乘法常用于根据实测数据求线性方程的最近似解。根据如图（图片引用于百度百科）的描述，利用C语言求，使用最小二乘法算法求线性方程的解，程序如下：
#include <stdio.h>#define N 4 //共有4个记录，根据需要增加记录
typedef struct Data{ //定义实验记录结构
int w //实验次数
double x
double y
}DATA
//根据d中的n个DATA记录，计算出线性方程的a,b两值
void getcs(DATA *d,int n,double &a,double &b){
double fi11=0,fi12=0,fi21=0,fi22=0,f1=0,f2=0
int i
for(i=0i<ni++){
fi11+=d[i].w
fi12+=d[i].w*d[i].x
fi21=fi12
fi22+=d[i].w*d[i].x*d[i].x
f1+=d[i].w*d[i].y
f2+=d[i].w*d[i].x*d[i].y
}
//解一元一次方程
b=(f2*fi11/fi21-f1)/(fi22*fi11/fi21-fi12)
a=(f2*fi12/fi22-f1)/(fi21*fi12/fi22-fi11)
}
int main(){
DATA d[N]={ //定义时赋初值，共4个记录
{2,0.1,1.1},
{1,0.2,1.9},
{1,0.3,3.1},
{1,0.4,3.9}
}
double a,b
getcs(d,N,a,b) //计算线性方程参数a,b
printf("线性方程是：Y=%.4lf+%.4lfX\n",a,b)
}
a=(10.0*sumxy-sumx*sumy)/(10.0*sum-sumx*sumx)/*带入a=(NΣxy-ΣxΣy)/(NΣx^2(Σx)^2)NΣx^2(Σx)^2 这边中间到底有没有减号？你注释里面没有减号，但是代码10.0*sum-sumx*sumx中是减法，你看看这边是否有问题
