最小二乘法多项式拟合实验报告Word文件下载.docx
《最小二乘法多项式拟合实验报告Word文件下载.docx》由会员分享，可在线阅读，更多相关《最小二乘法多项式拟合实验报告Word文件下载.docx（9页珍藏版）》请在冰豆网上搜索。
其中 s
n
=å
x ,
k+j
j+k i
i=0
n
j
t=å
yixi
所以，在编程计算时，先计算出方程组○1，再用LU分解法计算求出
ai的值，即可得到拟合多项式。
LU分解法的公式为：
[1 0
0⋯0
][𝑢
11𝑢
12𝑢
13⋯
𝑢
1��
][𝑎
11𝑎
12𝑎
13⋯𝑎
1𝑛
]
𝑙
21
10⋯0
0 𝑢
22𝑢
23⋯
𝑢
2��
𝑎
21𝑎
22𝑎
23⋯𝑎
2𝑛
31
𝑙
321⋱⋮
0 0 𝑢
33⋱ ⋮
31𝑎
32𝑎
33⋯𝑎
3𝑛
⋮ ⋮⋮⋱0
0 0 ⋱ ⋱𝑢
𝑛
‒1��
⋮ ⋮ ⋮⋱ ⋮
n1𝑙
n2⋯⋯1
0 0 ⋯0
𝑛
=𝑎
1𝑎
2𝑎
3⋯𝑎
其中L矩阵和U矩阵的计算公式如下：
第一步，当k=1,有：
{u1j
=a1j
（j=2,3,L,n）li1
=ai1=ai1u11 a11
（i=2,3,L,n）
第二步，当k=2,3,⋯,n-1时，有：
å
{u=a-k-1lu
（j=k+1,k+2L,n）
kj kj krrjr=1
lik
k-1
a -å
lu
（k）
ik irrk
u
=
r=1 （k）
kk
（j=k,k+1L,n）
n-1
最后求 ：
unn=ann-å
lnrurn
r=1
三、设计说明
（一）、数据结构
程序采用一维数组的形式来读取文件中给出的已知点处的值和要计算的未知点处的自变量值，最终的拟合计算结果也是采用一维数组的形式输出到文件中。
拟合多项式的系数a和拟合系数方程组的参数t都是采用一维数组来存储的，而拟合系数方程组中的参数s和L、U矩阵都是用二维数组来表示的。
由于要分别计算2、3、4阶拟合结果，所以数组的规模取为5，矩阵的规模取为5*5.
（二）、算法设计及效率分析
在进行LU分解函数中，在计算L矩阵和U矩阵时，因为当k=2,3,⋯,n-1时，
k-1 k-1
计算å
lirurk和å
lkrurj的循环条件不允许k=1时进入，而正好k=1时，计算li1和u1j不
r=1 r=1
需要å
lkrurj，因而对k=1和k=2,3,⋯,n-1,就可以和在一起计算，这样就减少了
程序的长度。
而在分别计算2、3、4阶拟合系数方程组的参数时，没有很好的利用前一阶计算的，而每次都要重新计算；
而且矩阵是一个堆成矩阵，没有好好利用对称矩阵的特性，导致了重复计算，增加了计算量，降低了程序的效率。
而造成这一结果的原因是自己为了编程的简单而忽视了计算量，在以后的编程时要注意改变这一习惯。
（三）、程序结构
程序主要步骤的流程图如下：
以上流程图对应的源程序中的函数分别如下：
//计算拟合系数方程组中的参数s
voidcomputers（doubles[p][p],doublex1[],intm）
//计算拟合系数方程组中的参数t
voidcomputert（doublet[],doublex1[],doubley1[],intm）
//对拟合系数方程组中的参数s组成的矩阵进行LU分解
voidLV（doubleL[p][p],doubleV[p][p],doubles[p][p],intm）
//计算拟合多项式的系数
voidcomputera（doubleL[p][p],doubleV[p][p],doublet[],doublex[],intm）
//由得到的的拟合多项式计算待求点处的函数值
voidcomputerty2（doublea[],doublex2[],doubley2[],intm）
//保存得到的拟合多项式和计算处的参数voidsave（doublea[],intm,doubley2[]）
四、编码实现
#include<
iostream>
#include<
fstream>
string>
math.h>
#definep5 //拟合方程的阶次+1
#defineq5 //已知点的数目，也是带计算点的数目usingnamespacestd;
ofstream outDatay（"
G:
\\连续系统仿真\\拟合实验\\outy.txt"
）;
//用于保存计算结果intmain（）
{
voidcomputers（doubles[p][p],doublex1[],intm）;
//计算拟合方程组的系数s
voidcomputert（doublet[],doublex1[],doubley1[],intm）;
//计算拟
合方程组的系数t
voidLV（doubleL[p][p],doubleV[p][p],doubles[p][p],intm）;
//LU分解
voidcomputera（doubleL[p][p],doubleV[p][p],doublet[],doublex[],intm）;
voidcomputerty2（doublea[],doublex2[],doubley2[],intm）;
voidsave（doublea[],intm,doubley2[]）;
ifstreaminDatax,outDatax,inDatay;
inti;
double x1[q],x2[q],y1[q],y2[q],s[p][p]={0},t[p]={0},a[p],L[p][p],V[p][p];
inDatax.open（"
\\连续系统仿真\\拟合实验\\inx.txt"
//已知点的自变量x值i=0;
while（!
inDatax.eof（））inDatax>
>
x1[i++];
inDatax.close（）;
i=0;
outDatax.open（"
\\连续系统仿真\\拟合实验\\outx.txt"
//要求的插值点的x值while（!
outDatax.eof（））
outDatax>
x2[i++];
outDatax.close（）;
i=0;
inDatay.open（"
\\连续系统仿真\\拟合实验\\iny.txt"
//已知点的因变量y值while（!
inDatay.eof（））
inDatay>
y1[i++];
inDatay.close（）;
computers（s,x1,2）;
computert（t,x1,y1,2）;
LV（L,V,s,2）;
computera（L,V,t,a,2）;
computerty2（a,x2,y2,2）;
save（a,2,y2）;
computers（s,x1,3）;
computert（t,x1,y1,3）;
LV（L,V,s,3）;
computera（L,V,t,a,3）;
computerty2（a,x2,y2,3）;
save（a,3,y2）;
computers（s,x1,4）;
computert（t,x1,y1,4）;
LV（L,V,s,4）;
computera（L,V,t,a,4）;
computerty2（a,x2,y2,4）;
save（a,4,y2）;
return0;
}
inti,j;
for（i=0;
i<
=m;
i++）
{ for（j=0;
j<
j++）s[i][j]=0;
for（i=0;
for（j=0;
j++）
for（intk=0;
k<
q;
k++）
{s[i][j]+=pow（x1[k],i+j）;
inti;
t[i]=0;
for（intj=0;
t[i]+=y1[j]*pow（x1[j],i）;
doublesum;
for（inti=0;
i++）for（intj=0;
j++）
L[i][j]=0;
V[i][j]=0;
for（i=0;
L[i][i]=1;
//计算U矩阵的第i行
{ sum=0;
for（intk=0;
=i-1;
sum+=L[i][k]*V[k][j];
V[i][j]=s[i][j]-sum;
if（i<
m）
//计算L矩阵的第i列
for（j=1;
sum=0;
sum+=L[j][k]*V[k][i];
L[j][i]=（s[j][i]-sum）/V[i][i];
double y[p],sum;
sum=
升级会员 