《正交多项式最小二乘法拟合(共13页).doc》由会员分享,可在线阅读,更多相关《正交多项式最小二乘法拟合(共13页).doc(13页珍藏版)》请在taowenge.com淘文阁网|工程机械CAD图纸|机械工程制图|CAD装配图下载|SolidWorks_CaTia_CAD_UG_PROE_设计图分享下载上搜索。
1、精选优质文档-倾情为你奉上MATLAB 程序设计实践课程考核 一、编程实现以下科学计算算法,并举一例应用之。(参考书籍精通科学计算,王正林等著,电子工业出版社,年)“正交多项式最小二乘法拟合”正交多项式最小二乘法拟合原理正交多项式做最小二乘法拟合:不要求拟合函数y=f(x)经过所有点(xi,yi),而只要求在给定点xi上残差i=f(xi)-yi按照某种标准达到最小,通常采用欧式范数|2作为衡量标准。这就是最小二乘法拟合。(1)根据作为给定节点x0,x1,xm及权函数(x)0,造出带权函数正交的多项式Pn(x)。注意nm,用递推公式表示Pk(x),即这里的Pk(x)是首项系数为1的k次多项式,根
2、据Pk(x)的正交性,得(2)根据公式(1)和(2)逐步求Pk(x)的同时,相应计算系数并逐步把 Pk(x)累加到S(x)中去,最后就可得到所求的拟合函数曲线流程图开始读取点集x,y和n数据size(x)=size(y)N提示x,y维数不匹配Y利用x点集数据构造范德蒙德系数矩阵V调用QR分解函数求的多项式系数输出多项式系数P结束M文件function p = mypolyfit(x,y,n)%定义mypolyfit为最小二乘拟合函数%P = POLYFIT(X,Y,N)以计算以下多项式系数%P(1)*XN + P(2)*X(N-1) +.+ P(N)*X + P(N+1).if isequal
3、(size(x),size(y) error(MATLAB:polyfit:XYSizeMismatch,. X and Y vectors must be the same size.)end%检验X Y维数是否匹配x = x(:);y = y(:); if nargout 2 mu = mean(x); std(x); x = (x - mu(1)/mu(2);end %利用范德蒙德矩阵构造方程组系数矩阵V(:,n+1) = ones(length(x),1,class(x);for j = n:-1:1 V(:,j) = x.*V(:,j+1);end % 对矩阵进行QR分解以求得多项式
4、系数值Q,R = qr(V,0);ws = warning(off,all); p = R(Q*y); warning(ws);if size(R,2) size(R,1) warning(MATLAB:polyfit:PolyNotUnique, . Polynomial is not unique; degree = number of data points.)elseif condest(R) 1.0e10 if nargout 2 warning(MATLAB:polyfit:RepeatedPoints, . Polynomial is badly conditioned. Rem
5、ove repeated data points.) else warning(MATLAB:polyfit:RepeatedPointsOrRescale, . Polynomial is badly conditioned. Remove repeated data pointsn . or try centering and scaling as described in HELP POLYFIT.) endendr = y - V*p;p = p.; % 将多项式系数默认为行向量.5、运行流程图 调用mypolyfit.m进行运算根据所得到的多项式系数构造拟合函数结束输入变量x=0.5
6、000 1.00001.50002.00002.50003.0000;y=1.75 2.45 3.81 4.80 8.00 8.60开始过程:clearx = 0.5000 1.0000 1.5000 2.0000 2.5000 3.0000y=1.75 2.45 3.81 4.80 8.00 8.60x1=0.5:0.05:3.0;p=mypolyfit(x,y,2)y1=p(3)+p(2)*x1+p(1)*x1.2;plot(x,y,*)hold onplot(x1,y1,r)二、 编程计算以下电路问题例8-1-3如图所示电路,已知R=5,L=3,=5,Uc=10,求R,C,和L,S,并画
7、其相量图。理论分析:根据电路分析Z=R+j*(Xl-Xc)Ic=Uc/Z3;Z3=-j*XcIr=Ur/Z2=Uc/Z2;Z2=RI=Ir+IcUl=I*Z1;Z1=j*XLUs=Ul+Ur计算得Ir =2;Ic =2.00iI =2.00 + 2.00iUl =-6.00 + 6.00iUs =4.00 + 6.00iM文件clearR=5;XL=3;XC=5;UC=10;UR=UC;%为给定元件赋值Z1=j*XL;Z2=R;Z3=-j*XC;%定义各电抗disp(电流)IR=UR/Z2%计算IRIC=UC/Z3%计算ICI=IR+IC%计算Idisp(电压)UL=I*Z1%计算ULUS=U
8、C+UL%计算USdisp(IR IC I UL US)disp(幅值);disp(abs(IR,IC,I,UL,US)disp(相角);disp(angle(IR,IC,I,UL,US)*180/pi)ph=compass(IR,IC,I,UL,US);%分别画出IR,IC,I,UL,US相量图set(ph,linewidth,3)运行流程图:开始读取已知数据构造电抗Z1,Z2,Z3计算Ir ,Ic, I ,Ul ,UsIR=UR/Z2;IC=UC/Z3;I=IR+IC;UL=I*Z1;US=UC+UL结束输出Ir ,Ic, I ,Ul ,Us并画出对应的相量图运行图三、 编程解决以下问题求
9、自然三次样条曲线,经过点(-3,2),(-2,0),(1,3),(4,1),而且自由边界条件(-3)=0,(4)=0。算法说明:三次样条也是分片三次插值函数,它是在给定的区间a,b上的一个划分:a=x0x1=x0(j)&(x(k)x0(j+1) l(k)=j; end endendfor k=1:km sum=(3*(x0(l(k)+1)-x(k)2/h(l(k)2-2*(x0(l(k)+1)-x(k)3/h(l(k)3)*y0(l(k); sum=sum+(3*(x(k)-x0(l(k)2/h(l(k)2-2*(x(k)-x0(l(k)3/h(l(k)3)*y0(l(k)+1); sum=sum+h(l(k)*(x0(l(k)+1)-x(k)2/h(l(k)2-(x0(l(k)+1)-x(k)3/h(l(k)3)*m(l(k); s(k)=sum-h(l(k)*(x(k)-x0(l(k)2/h(l(k)2-(x(k)-x0(l(k)3/h(l(k)3)*m(l(k)+1);end举例:x=-3 -2 1 4y=2 0 3 1x0=-3:0.15:4;y0=myspline(x,y,0,0,x0);plot(x0,y0)hold onplot(x,y,or)legend(计算值,实验值)专心-专注-专业