matlab中函数拟合方法—个人总结.doc

目 录一维插值方案2二维数据内插值(表格查找)3等高线4三维曲面5等高线26三维曲面27matlab绘制温度场(尚未深入研究)13二维曲线（非线性）拟合步骤18三维曲线（非线性）拟合步骤19三维曲线的画法20三维曲面的画法21画三维图3 只有点的数据,没有函数关系式23空间点拟合的基本原理27空间点拟合的最小二乘法28曲面生成后再进行多项式拟合37六点生成曲面38四点生成平面39用三维离散点拟合光滑曲面140用三维离散点拟合光滑曲面240一维插值方案clearyear = 1900:10:2010;product = 75.995 91.972 105.711 123.203 131.669 150.697 179.323 203.212 226.505 249.633 256.344 267.893 p1995 = interp1(year,product,1995)%使用一维数据内插值(该题中只能在1900和2010之间进行插值,大于2010和小于1900都%无效)命令x = 1900:1:2010y = interp1(year,product,x,spine);plot(year,product,o,x,y)插值说明:interp1(x,Y,xi,method) %用指定的算法计算插值：nearest：最近邻点插值，直接完成计算；linear：线性插值（缺省方式），直接完成计算；spine：三次样条函数插值。对于该方法，命令interp1 调用函数spline、ppval、mkpp、umkpp。这些命令生成一系列用于分段多项式操作的函数。命令spline 用它们执行三次样条函数插值；pchip：分段三次Hermite 插值。对于该方法，命令interp1 调用函数pchip，用于对向量x 与y 执行分段三次内插值。该方法保留单调性与数据的外形；cubic：与pchip操作相同；v5cubic：在MATLAB 5.0 中的三次插值。对于超出x 范围的xi 的分量，使用方法nearest、linear、v5cubic的插值算法，相应地将返回NaN。对其他的方法，interp1 将对超出的分量执行外插值算法。yi = interp1(x,Y,xi,method,extrap) %对于超出x 范围的xi 中的分量将执行特殊的外插值法extrap。yi = interp1(x,Y,xi,method,extrapval) %确定超出x 范围的xi 中的分量的外插值extrapval，其值通常取NaN 或0。例1clear;x = 0:10; y = x.*sin(x);xx = 0:.25:10; yy = interp1(x,y,xx)plot(x,y,kd,xx,yy)interp2二维数据内插值(表格查找)X,Y = meshgrid(-3:.25:3);Z = peaks(X,Y);XI,YI = meshgrid(-3:.125:3);ZZ = interp2(X,Y,Z,XI,YI);surfl(X,Y,Z);hold on;surfl(XI,YI,ZZ+15)axis(-3 3 -3 3 -5 20);shading flathold off功能 三维数据插值interp3（查表）x,y,z,v = flow(20);xx,yy,zz = meshgrid(.1:.25:10, -3:.25:3, -3:.25:3);vv = interp3(x,y,z,v,xx,yy,zz);slice(xx,yy,zz,vv,6 9.5,1 2,-2 .2); shading interp;colormap cool等高线clearZ=peaksfor w=1:1:100V=w/10,0,w/10contour(Z,V)%C=contour(Z,V)%Clabel(C)Hold ontitle(等高线及其标注)endend三维曲面x=0:10y=0:.1:1d,B=meshgrid(x,y)z=1./(B.*d.2+1);surf(B,d,z)x=0:0.05:10y=0:0.05:1X,Y=meshgrid(x,y)Z=( X.3+ 3.*Y.2+5*Y); %Z=( X.2+ 3.*Y.3+5*Y);%surf(X,Y,Z)%一张普通的三维曲面,有时需要旋转一下才能看到下图的结果;x=0:0.05:1y=0:0.05:1X,Y=meshgrid(x,y)Z=( X.2-Y.2);% Z=( 4*X.3*Y-4*X.*Y.3);surf(X,Y,Z) %一张普通的三维曲面,有时需要旋转一下才能看到下图的结果;等高线2clearx=-2:0.1:2y=-2:0.1:2X,Y=meshgrid(x,y)Z=(X.2+Y.2).0.5for w=1:1:100V=w/3,w/pi,w/3contour(Z,V)hold onend三维曲面2clearx=-5:0.05:5y=-5:0.05:5X,Y=meshgrid(x,y)Z=1./(X+1).2+(Y+1).2+1)-1.5./(X-1).2+(Y-1).2+1)mesh(X,Y,Z)clear;A=1.486,3.059,0.1;2.121,4.041,0.1;2.570,3.959,0.1;3.439,4.396,0.1;4.505,3.012,0.1;3.402,1.604,0.1;2.570,2.065,0.1;2.150,1.970,0.1;1.794,3.059,0.2;2.121,3.615,0.2;2.570,3.473,0.2;3.421,4.160,0.2;4.271,3.036,0.2;3.411,1.876,0.2;2.561,2.562,0.2;2.179,2.420,0.2;2.757,3.024,0.3;3.439,3.970,0.3;4.084,3.036,0.3;3.402,2.077,0.3;2.879,3.036,0.4;3.421,3.793,0.4;3.953,3.036,0.4;3.402,2.219,0.4;3.000,3.047,0.5;3.430,3.639,0.5;3.822,3.012,0.5;3.411,2.385,0.5;3.103,3.012,0.6;3.430,3.462,0.6;3.710,3.036,0.6;3.402,2.562,0.6;3.224,3.047,0.7;3.411,3.260,0.7;3.542,3.024,0.7;3.393,2.763,0.7;x=A(:,1);y=A(:,2);z=A(:,3);scatter(x,y,5,z)%散点图figureX,Y,Z=griddata(x,y,z,linspace(1.486,4.271),linspace(1.604,4.276),v4);%插值pcolor(X,Y,Z);shading interp%伪彩色图figure,contourf(X,Y,Z) %等高线图clear;A=1.486,3.059,1858;2.121,4.041, 1858;2.570,3.959, 1858;3.439,4.396, 1858;4.505,3.012, 1858;3.402,1.604, 1858;2.570,2.065, 1858;2.150,1.970, 1858;1.794,3.059,2350;2.121,3.615, 2350;2.570,3.473, 2350;3.421,4.160, 2350;4.271,3.036, 2350;3.411,1.876, 2350;2.561,2.562, 2350;2.179,2.420, 2350;2.757,3.024, 2600;3.439,3.970, 2600;4.084,3.036, 2600;3.402,2.077, 2600;2.879,3.036, 2849;3.421,3.793, 2849;3.953,3.036, 2849;3.402,2.219, 2849;3.000,3.047, 3010;3.430,3.639, 3010;3.822,3.012, 3010;3.411,2.385, 3010;3.103,3.012, 3345;3.430,3.462, 3345;3.710,3.036, 3345;3.402,2.562, 3345;3.224,3.047, 3629;3.411,3.260, 3629;3.542,3.024, 3629;3.393,2.763, 3629;x=A(:,1);y=A(:,2);z=A(:,3);scatter(x,y,5,z)%散点图,5是点的大小figure %打开显示图的界面X,Y,Z=griddata(x,y,z,linspace(1.486,4.271),linspace(1.604,4.276),v4);%插值pcolor(X,Y,Z);shading interp%伪彩色图figure;contourf(X,Y,Z) %等高线图figure;mesh(X,Y,Z)A=1.109,1.059,1718;2.021,0.841, 1758;2.870,0.359, 1858;4.039,0.196, 1838;4.505,3.012, 3345;3.402,1.604, 3347;2.570,2.065, 3629;2.150,1.970, 3330;1.794,3.059,2250;2.121,3.615, 3027;2.570,3.473, 2935;3.421,4.160, 1930;4.271,3.036, 2050;3.411,1.876, 3144;2.561,2.562, 3739;2.179,2.420, 1950;2.757,3.024, 3530;3.439,3.970, 2720;4.084,3.036, 2610;3.402,2.077, 3500;2.879,3.036, 3249;3.421,3.793, 2149;3.953,3.036, 2849;3.402,2.219, 2849;3.000,3.047, 3010;3.430,3.639, 3010;3.822,3.012, 2310;3.411,2.385, 3410;3.103,3.012, 3345;3.430,3.462, 3845;3.710,3.036, 2645;3.402,2.562, 2745;3.224,3.047, 3229;3.411,3.260, 3329;3.542,3.024, 3429;3.393,2.763, 3529;x=A(:,1);y=A(:,2);z=A(:,3);scatter(x,y,5,z)%散点图,5是点的大小figure %打开显示图的界面X,Y,Z=griddata(x,y,z,linspace(1.486,4.271),linspace(1.604,4.276),v4);%插值pcolor(X,Y,Z);shading interp%伪彩色图figure;contourf(X,Y,Z) %等高线图figure;mesh(X,Y,Z)A=1.109,1.059,0.4874;2.021,0.841,0.5643;2.870,0.359,0.4628;4.039,0.196,0.4411;4.505,3.012,0.4845;3.402,1.604,0.7857;3.570,3.565,0.7071;2.150,4.870,0.4284;1.794,3.059,1.0000;2.121,3.615,0.8544;2.570,3.473,1.0000;3.421,4.160,0.5447;4.271,3.036,0.5643;3.411,1.876,0.8771;2.561,2.562,1.0000;2.179,2.420,1.0000;2.757,3.024,1.0000;3.439,3.970,0.6008;4.084,3.036,0.6325;3.402,2.077,0.9713;2.879,3.036,1.0000;3.421,3.793,0.6667;3.953,3.036,0.6727;3.402,2.219,1.0000;3.000,3.047,1.0000;3.430,3.639,0.7036;3.822,3.012,0.7180;3.411,4.215,0.5199;1.103,4.612,0.3962;3.430,3.462,0.7857;3.710,3.036,0.7692;3.802,2.462,0.7670;3.424,3.247,0.8771;3.511,3.060,0.8944;4.342,2.724,0.5522;3.803,2.903,0.7352;x=A(:,1);y=A(:,2);z=A(:,3);scatter(x,y,5,z)%散点图,5是点的大小figure %打开显示图的界面X,Y,Z=griddata(x,y,z,linspace(1.486,4.271),linspace(1.604,4.276),v4);%插值pcolor(X,Y,Z);shading interp%伪彩色图figure;contourf(X,Y,Z) %等高线图figure;mesh(X,Y,Z)matlab绘制温度场(尚未深入研究)clearecho ond1=43;d2=7;dx=0.15;dy=0.1;xy=dx/dy;yx=dy/dx;t=zeros(d1,d2);t1=ones(d1,d2);t0=zeros(d1,d2);x=zeros(d1);y=zeros(d2);x(1)=0;x(2)=dx/2;for i=3:d1-1; x(i)=x(i-1)+dx;endx(d1)=(d1-2)*dx;y(1)=0;y(2)=dy/2;for i=3:d2-1; y(i)=y(i-1)+dy;endy(d2)=(d2-2)*dy;t1=20*ones(d1,d2);t=zeros(d1);dt=0.1;ttt=30;%nnn=tt/dt;echo off%for iii=1:nnn% ttt=iii*dt; tf=30;af=6.6;af=1/af;bta=-6.6;v=0.0625; tin=100;tout=100;d=0.05;l=6; for i=1:42; if x(i)v*ttt+28*d t1(i,1)=300-(300-tout)*(x(i)-v*ttt-28*d)/(l-28*d-v*ttt); else zz=-0.123*(x(i)-v*ttt)/d-3.52*exp(-0.123*(x(i)-v*ttt)/d); t1(i,1)=10060*exp(zz); end end for i=1:41; for j=2:6 t1(i,j)=t1(i,1)-10*(j-1); end endfor iii=1:500; t0=t1; cd=1300;a=0.0003; an=1.69,-0.594,0.401,-0.168,0.027,-0.037,0.046,-0.05, 0.039,-0.012; bn=0, 0.333, 0.017,-0.131,0.054,0.0003,0.007,-0.012,0.026,0; fncp=zeros(d1,d2); for i=1:10; fncp=an(i)*cos(i*pi/9*t1)+bn(i)*sin(i*pi/9*t1)+fncp; end fncp=fncp/4.1868; b=1.3e-6;c=1.5e-9; fnk=(a+b*t1+c*t1.2)*360; for i=2:42; for j=2:6; fnae(i,j)=2*yx*fnk(i,j)*fnk(i+1,j)/(fnk(i,j)+fnk(i+1,j); fnaw(i,j)=2*yx*fnk(i,j)*fnk(i-1,j)/(fnk(i,j)+fnk(i-1,j); fnan(i,j)=2*xy*fnk(i,j)*fnk(i,j+1)/(fnk(i,j)+fnk(i,j+1); fnas(i,j)=2*xy*fnk(i,j)*fnk(i,j-1)/(fnk(i,j)+fnk(i,j-1); end end fnap0=cd*fncp*dx*dy/dt; fnbb=fnap0.*t1; %for i=2:41 % for j=1:5 % t1(i,j)=t(i,1)-10*(j-1); % end % end % t0=t1; kk=af+0.5*dx/fnk(1,2); bb=fnbb(2,2)+tf*dy/kk; ap=fnae(2,2)+fnan(2,2)+fnap0(2,2)+dy/kk+2*xy*fnk(2,1)-bta*dx*dy; fff=fnae(2,2)*t1(3,2)+fnan(2,2)*t1(2,3)+2*xy*fnk(2,1)*t1(2,1)+bb; t1(2,2)=fff/ap; for j=3:5; kk=af+0.5*dx/fnk(1,j); bb=fnbb(2,j)+tf*dy/kk; ap=fnae(2,j)+fnas(2,j)+fnan(2,j)+fnap0(2,j)+dy/kk-bta*dx*dy; fff=fnae(2,j)*t1(3,j)+fnan(2,j)*t1(2,j+1)+fnas(2,j)*t1(2,j-1)+bb; t1(2,j)=fff/ap; end kk1=af+0.5*dx/fnk(1,6); kk2=af+0.5*dy/fnk(1,6); bb=fnbb(2,6)+tf*dy/kk1+tf*dx/kk2; ap=fnae(2,6)+fnas(2,6)+fnap0(2,6)+dy/kk1+dx/kk2-bta*dx*dy; fff=fnae(2,6)*t1(3,6)+fnas(2,6)*t1(2,5)+bb; t1(2,6)=fff/ap; for i=3:40; for j=2:6; if j=2; as=2*xy*fnk(i,1); fff=fnae(i,2)*t1(i+1,2)+fnaw(i,2)*t1(i-1,2)+fnan(i,2)*t1(i,3)+as*t1(i,1)+fnbb(i,2); ap=fnae(i,2)+fnaw(i,2)+fnan(i,2)+as+fnap0(i,2)+dy/kk-bta*dx*dy; t1(i,2)=fff/ap; elseif j=6 kk=af+0.5*dy/fnk(i,6); bb=fnbb(i,6)+tf*dx/kk; fff=fnae(i,6)*t1(i+1,6)+fnaw(i,6)*t1(i-1,6)+fnas(i,6)*t1(i,4)+bb; ap=fnae(i,6)+fnaw(i,6)+fnas(i,6)+fnap0(i,6)+dx/kk-bta*dx*dy; t1(i,5)=fff/ap; else fff=fnae(i,j)*t1(i+1,j)+fnaw(i,j)*t1(i-1,j)+fnan(i,j)*t1(i,j+1)+fnas(i,j)*t1(i,j-1)+fnbb(i,j); ap=fnae(i,j)+fnaw(i,j)+fnan(i,j)+fnas(i,j)+fnap0(i,j)-bta*dx*dy; t1(i,j)=fff/ap; end end end for j=3:5; kk=af+0.5*dx/fnk(42,j); bb=fnbb(41,j)+tf*dy/kk; ap=dy/kk+fnaw(41,j)+fnan(41,j)+fnas(41,j)+fnap0(41,j)-bta*dx*dy; fff=fnaw(41,j)*t1(40,j)+fnan(41,j)*t1(41,j+1)+fnas(41,j)*t1(41,j-1)+bb; t1(41,j)=fff/ap; end kk1=af+0.5*dx/fnk(42,6); kk2=af+0.5*dy/fnk(41,7); bb=fnbb(41,6)+tf*dy/kk1+tf*dx/kk2; ap=fnaw(41,6)+fnas(41,6)+fnap0(41,6)+dy/kk1+dx/kk2-bta*dx*dy; fff=fnaw(41,6)*t1(40,6)+fnas(41,6)*t1(41,5)+bb; t1(41,6)=fff/ap; kk=af+0.5*dx/fnk(42,2); bb=fnbb(41,2)+tf*dy/kk; ap=fnaw(41,2)+fnan(41,2)+fnap0(41,2)+dy/kk+2*xy*fnk(41,1)-bta*dx*dy; fff=fnaw(41,2)*t1(40,2)+fnan(41,2)*t1(41,3)+2*xy*fnk(41,1)*t1(41,1)+bb; t1(41,2)=fff/ap; if abs(abs(t1)-abs(t0)20; break; %else %t1=t0+0.5*(t1-t0); endend%if rem(iii,10)=0 figure subplot(3,2,1) t111=t1; %mesh(t111); %view(2); pcolor(t111); subplot(3,2,5) ccc=contour(t111,200,400,600,800,1000); %ccc=contour(t111,5); clabel(ccc); %end % 将平面网格映射到三维网格面上 x1=1:7; y1=1:43; X1,Y1=meshgrid(y1,x1); %绘制加密前三维网格 subplot(3,2,3); mesh(X1,Y1,t111); % 用二元差值公式将网格加密并绘出温度场的三维图 x2=1:.25:7; y2=1:.25:43; X2,Y2=meshgrid(y2,x2); Z2=griddata(X1,Y1,t111,X2,Y2); t111=round(t111); subplot(3,2,4); mesh(X2,Y2,Z2) % 绘制加密后的平面网格 subplot(3,2,2); pcolor(Z2) % 绘制加密后的温度场等高线图 subplot(3,2,6) CCC=contour(Z2,100,200,300,400,500,600,700,800,900,1000); clabel(CCC); % end%下面输入平板中某点坐标以查看此点温度a1=input(请输入你要查看温度的点的X轴坐标值（范围为1到169）);b1=input(请输入你要查看温度的点的Y轴坐标值（范围为1到25）);disp(输入坐标值为),a1,b1a2=round(a1+3)/4);b2=round(b1+3)/4);disp(计算出网格加密前坐标点的温度为)Q1=(t111(b2,a2)display(计算出网格加密后坐标点的温度为)Q2= Z2(b1,a1) 二维曲线（非线性）拟合步骤1. 建立函数function F = myfun(x,xdata)%F = x(1)*xdata.2 + x(2)*sin(xdata) + x(3)*xdata.3;F = x(1)+x(2)*xdata.2 + x(3)*sin(xdata) + x(4)*xdata.3; % 这个函数的可以是任意的,这个要%根据实际情况自己来选;2.然后给出数据xdata和ydataxdata = -1.0:0.5:2.0;ydata = -4.447,-0.452,0.551,0.048,-0.447,0.549,4.552;x0 = 10, 10, 10,1; %初始估计值,可以随意给,有了这个初值后matlab才能进行迭代, 这是matlab的迭代的初始条件, 这个初始条件与结果无关; % 中的初值的数量要和myfun中x(1) x(2) x(3) x(4)的数量一致才行%拟合函数的表达式:x,resnorm = lsqcurvefit(myfun,x0,xdata,ydata) %回车后;x = 0.5491 -2.9977 -0.0000 1.9991x结果即为myfun 中对应的x(1) x(2) x(3) x(4)的值；这样就得到了所求的多项式的系数, resnorm是在x处残差的平方和, 在workpsace中可以查到, resnorm=sum (fun(x,xdata)-ydata).2)，函数建立成功后,在matlab输入的程序如下;clear;%可以有,也可以没有,最好有.xdata = -1.0:0.5:2.0;ydata = -4.447,-0.452,0.551,0.048,-0.447,0.549,4.552;x0 = 10, 10, 10, 1; %初始估计值x,resnorm = lsqcurvefit(myfun,x0,xdata,ydata)怎么样,很简单吧?如果不想编写函数的话也可以将最后一句话改为:x,resnorm = lsqcurvefit(x,xdata)x(1)+x(2)*xdata.2 + x(3)*sin(xdata) + x(4)*xdata.3,x0,xdata,ydata);这样就不用编写函数了; 三维曲线（非线性）拟合步骤 1 设定目标函数. （M函数书写）% 可以是任意的例如：function f=mydata(a,data) %y的值目标函数值 或者是第三维的,a=a(1) ，a(2) 列向量x=data(1,:); %data 是一2维数组，x=x1y=data(2,:); %data 是一2维数组，x=x2 f=a(1)*x+a(2)*x.*y; %这里的a(1)， a(2)为目标函数的系数值。 f的值相当于ydata的值 2 然后给出数据xdata和ydata的数据和拟合函数lsqcurvefit 例如：x1=1.0500 1.0520 1.0530 1.0900 1.0990 1.1020 1.1240 1.1420. 1.1490 1.0500 1.0520 1.0530 1.0900 1.0990 1.1020 1.1240 1.1420 1.1490;x2=3.8500 1.6500 2.7500 5.5000 7.7000 3.3000 4.9500 8.2500 11.5500. 1.6500 2.7500 3.8500 7.7000 3.3000 5.5000 8.2500 11.5500 4.9500;ydata=56.2000 62.8000 62.2000 40.8000 61.4000 57.5000 44.5000 54.8000. 53.9000 64.2000 62.9000 64.1000 63.0000 62.2000 64.2000 63.6000. 52.5000 62.0000;data=x1;x2; %类似于将x1 x2整合成一个2维数组。 a0= -0.0014,0.07;option=optimset(MaxFunEvals,5000);format long;a,resnorm=lsqcurvefit(mydata,a0,data,ydata,option);yy=mydata(a,data);result=ydata yy (yy-ydata)将数据按列的形式排列:第一列是yadata; 第二列是yadata(根据拟合函数得到的计算值);第三列是上述二者的插值; 这个符号实现了数据按列进行排列的功能;这个符号后面有空格result = 56.200000000000003 57.919539361940053 1.719539361940051 62.799999999999997 60.487127691642989 -2.312872308357008 62.200000000000003 59.314824360714432 -2.885175639285571 40.799999999999997 58.216478553229628 17.416478553229631 61.399999999999999 56.130116906144998 -5.269883093855000 57.500000000000000 61.431449491527204 3.931449491527204 44.500000000000000 60.688766221557501 16.188766221557501 54.799999999999997 57.659418794420404 2.859418794420407 53.899999999999999 53.987090284681393 0.087090284681395 64.200000000000003 60.372133152305260 -3.827866847694743 62.899999999999999 59.258494992850515 -3.641505007149483 64.099999999999994 58.085023760117025 -6.014976239882969 63.000000000000000 55.670452618469561 -7.329547381530439 62.200000000000003 61.264213240642817 -0.935786759357185 64.200000000000003 58.857393913448675 -5.342606086551328 63.600000000000001 56.750601335313952 -6.849398664686049 52.500000000000000 53.658187210710317 1.158187210710317 62.000000000000000 62.038605327908876 0.038605327908876% a的值为拟合的目标函数的参数值 利用lsqcurvefit进行拟合的 它完整的语法形式是：% x,resnorm,residual,exitflag,output,lambda,jacobian =lsqcurvefit(fun,x0,xdata,ydata,lb,ub,options)三维曲线的画法三维空间