MATLAB应用实例分析例分析文档格式.docx
《MATLAB应用实例分析例分析文档格式.docx》由会员分享，可在线阅读，更多相关《MATLAB应用实例分析例分析文档格式.docx（11页珍藏版）》请在冰豆网上搜索。
i=0.5:
3.5;
p=1./（n*log（1+0.01*i））;
T=log（r'
）*p;
plot（r,T）
xlabel（'
r'
）%给x轴加标题
ylabel（'
T'
）%给y轴加标题
q=ones（1,length（i））;
text（7*q-0.2,[T（14,1:
5）+0.5,T（14,6）-0.1,T（14,7）-0.9],num2str（i'
））
40
350.530
25
20T115
1.51022.5533.50123456789r
图1
1
从图1中既可以看到T随r的变化规律，而且还能看到i的不同取值对T—r曲线的影响（图中的六条曲线分别代表i的不同取值）。
二、已知多项式求根
65432已知多项式为，求其根。
h,x,10x，31x,10x,116x，200x,96
对多项式求根问题，我们常用roots（）函数。
>
h=roots（[1-1031-10-116200-96]）%中括号内为多项式系数由高阶到常数。
计算结果显示为（其中i为虚数单位）:
h=
-2.0000
4.0000
3.0000
2.0000+0.0000i
2.0000-0.0000i
1.0000
如果已知多项式的根，求多项式，用poly（）函数。
对上面得到的h的值求多项式，其MATLAB的表达形式及结果如下:
h=[-2.00004.00003.00002.0000+0.0000i2.0000-0.0000i1.0000];
c=poly（h）
计算结果显示为:
c=
1-1031-10-116200-96
三、方程组的求解
8x，x，6x,7.5,123,求解下面的方程组:
3x，5x，7x,4,123,4x，9x，2x,12123,
对于线性方程组求解，常用线性代数的方法，把方程组转化为矩阵进行计算。
1,x,abax,b,x,a\b
a=[816;
357;
492];
%建立系数矩阵
b=[7.5;
4;
12];
%建立常数项矩阵
x=a\b%求方程组的解
x=
1.2931
0.8972
-0.6236
四、数据拟合与二维绘图
在数学建模竞赛中，我们常会遇到这种数据表格问题，如果我们仅凭眼睛观察，很难看到其中的规律，也就更难写出有效的数学表达式从而建立数学模型。
因此可以利用MATLAB的拟合函数，即polyfit（）
函数，并结合MATLAB的绘图功能（利用plot（）函数），得到直观的表示。
2
例:
在化学反应中，为研究某化合物的浓度随时间的变化规律，测得一组数据如下表:
12345678T（分）
y46.48.08.49.289.59.79.86
910111213141516T（分）
y1010.210.3210.4210.510.5510.5810.6
MATLAB的表达形式如下:
t=[1:
16];
%数据输入
y=[46.488.49.289.59.79.861010.210.3210.4210.510.5510.5810.6];
plot（t,y,'
o'
）%画散点图
p=polyfit（t,y,2）%二次多项式拟合
holdon
xi=linspace（0,16,160）;
%在[0,16]等间距取160个点
yi=polyval（p,xi）;
%由拟合得到的多项式及xi，确定yi
plot（xi,yi）%画拟合曲线图
执行程序得到图2;
11
10
9
8
7
6
5
40246810121416
图2
显示的结果为
p=
0.04451.07114.3252-
2p的值表示二阶拟合得到的多项式为:
y=-0.0445t+1.0711t+4.3252
下面是用lsqcurvefit（）函数，即最小二乘拟合方法的Matlab表达:
x0=[0.1,0.1,0.1];
zuixiao=inline（'
x
（1）*t.^2+x
（2）*t+x（3）'
'
x'
t'
）;
x=lsqcurvefit（zuixiao,x0,t,y）%利用最小二乘拟合
其显示的结果为:
-0.04451.07114.3252
可以看出其得到的结果与polyfit函数的结果相同。
这说明在多项式拟合问题上这两个函数的效果是相同的。
3
下面的一个例子将体现lsqcurvefit（）函数的优势。
例2:
在物理学中，为研究某种材料应力与应变的关系，测得一组数据如下表:
925112516252125262531253625应力σ
0.110.160.350.480.610.710.85应变ε
如果假定应力与应变有如下关系（σ为应力值，ε为应变值）:
ε=a+blnσ试计算a、b的值。
x=[925,1125,1625,2125,2625,3125,3625];
y=[0.11,0.16,0.35,0.48,0.61,0.71,0.85];
plot（x,y,'
）
[p,resid1]=polyfit（x,y,2）
xi=linspace（700,3700,3000）;
yi=polyval（p,xi）;
plot（xi,yi）
x0=[0.1,0.1];
fff=inline（'
a
（1）+a
（2）*log（x）'
a'
[a,resid2]=lsqcurvefit（fff,x0,x,y）plot（xi,fff（a,xi）,'
执行程序得到图3，图中蓝色曲线为利用polyfit（）函数得到的曲线，红色曲线为利用lsqcurvefit（）函
数得到的曲线;
0.90.80.70.60.50.40.30.20.10-0.15001000150020002500300035004000
p=
-0.00000.0004-0.2266resid1=
R:
[3x3double]
df:
4
normr:
0.0331
a=
-3.58100.5344
resid2=
0.0064
其中a的值代表利用lsqcurvefit（）函数得到的关系为:
ε=-3.5810+0.5344σ
4
resid1、resid2分别代表运用polyfit（）函数、lsqcurvefit（）函数得到的残差。
可以看出利用lsqcurvefit（）函数残差更小，即得到了更好的拟合效果。
在数学建模的实际问题中，如果问题的机理不明，我们只能采用polyfit（）函数，即多项式拟合的
方法，以获得近似的数据描述函数;
但如果通过分析，可以得到一些机理，那么采用最小二乘的方法将
得到更好的效果，而且得到的拟合函数也更有意义。
五、隐函数的图形绘制
1plot（）只能绘制显函数图形，对于形如的复杂隐函数，很难转化为,lny，ln（,1，y），x,sin（x）,0y
显函数并利用plot（）函数绘制图形，这时就可以用ezplot（）函数直接绘制其曲线。
MATLAB的表达形式如下:
ezplot（'
1/y-log（y）+log（-1+y）+x-sin（x）'
执行程序得到图5
x=sin（3t）cos（t）,y=sin（3t）sin（t）1/y-log（y）+log（-1+y）+x-sin（x）=060.4
40.2
20
-0.2yy0
-0.4-2-0.6-4-0.8
-6-0.8-0.6-0.4-0.200.20.40.60.8-6-4-20246xx
图5图6如果是形如下面的参数方程,同样可以利用ezplot（）函数绘制其曲x,sin3tcost,y,sin3tsint,t,（0,,）线。
MATLAB的表达形式如下:
sin（3*t）*cos（t）'
sin（3*t）*sin（t）'
[0,pi]）
执行程序得到图6。
六、三维图形绘制
假设有一个时间向量t，对该向量进行下列运算则可以构成三个坐标值向量
x,sint,y,cost,z,t
对于上面的方程可以利用ezplot3（）函数或plot3（）函数绘制三维曲线。
这里仅列举ezplot3（）函数的使用。
ezplot3（'
sin（t）'
cos（t）'
[0,6*pi]）
执行程序得到图7:
3绘制下述曲面:
z（r,,）,rcos（3,）,其中0,r,1,0,,,2,
nr=12;
nth=50;
r=linspace（0,1,nr）;
theta=linspace（0,2*pi,nth）;
[R,T]=meshgrid（r,theta）
x=cos（theta'
）*r;
5
y=sin（theta'
surf（x,y,R.^3.*cos（3*T））
执行程序得到图8。
x=sin（t）,y=cos（t）,z=t
15
10z
010.510.500-0.5-0.5-1-1yx
图7图8
除了surf（）函数还有surfc（）、surfl（）、mesh（）、waterfall（）函数也用于曲面的绘制，具体效果如图9所示，可
以针对自己的需要选取适合的曲面绘制函数。
图9
MATLAB作图之三维绘图示例山体绘制
%三维绘图示例山体绘制
%mesh函数演示
x=1.0:
0.1:
2.0;
y=2.0:
3.0;
[X,Y]=meshgrid（x,y）;
z=[5.115.135.145.135.095.044.984.934.894.854.85
5.395.495.515.465.325.144.944.744.594.494.48
5.615.775.815.715.515.234.904.594.364.214.19
5.735.925.975.865.625.274.884.514.234.054.03
5.745.925.975.867.625.274.884.514.214.044.02
5.635.795.846.7410.539.238.914.594.334.184.16
5.425.535.565.497.355.164.934.734.554.454.44
5.145.185.195.1711.125.054.974.904.844.814.80
4.484.804.794.824.874.945.025.105.165.195.20
4.564.454.434.494.644.845.065.285.455.555.56
4.364.194.164.254.474.765.095.415.665.815.83];
figure
（1）
%画网格图
mesh（X,Y,z）;
6
colormap（[010]）;
x轴'
y轴'
zlabel（'
z轴'
%画表面图
figure
（2）
surf（X,Y,z）;
colormap（[100]）;
1212
1010
88z轴z轴66
4433221.81.82.52.51.61.61.41.41.21.22211y轴y轴x轴x轴
七、二项分布的使用
飞机成功起飞的概率问题:
由16架飞机组成的空军飞行中队要求做好立即起飞的准备，其中一架飞机不能立即起飞的概率为
20%，重新起飞需几分钟的时间，因此一架飞机立刻起飞的概率为0.80。
12架飞机能够成功起飞的概率
为多少,
这是一个概率中的二项分布问题，常用binopdf（）函数。
h=binopdf（12,16,0.80）%二项分布函数的概率值
计算结果为
h=0.2001
另一方面，至少有14架飞机立刻成功起飞的概率为:
h=1-binocdf（13,16,0.80）%或h=sum（binopdf（14:
16,16,0.80）），其中binocdf（）为二项分布的累积概率值，计算结果为h=0.3518。
在实际的数学建模竞赛中，仅罗列一个一个的数据是枯燥而又不直观的，很难吸引人们的注意，也不容
易打动评委们的心;
因此，结合数值计算结果，并合理的利用MATLAB的绘图功能会起到事半功倍的效
果。
下面的程序为运行结果的绘图（图10）表示:
n=1:
16;
h=binopdf（n,16,0.80）;
plot（[n;
n],[zeros（1,16）;
h],'
k'
）%二维绘图函数
text（8-.7:
16-.7,h（8:
16）+.005,num2str（h（8:
16）'
3））%在图中进行注释函数axis（[01700.27]）%坐标轴取值范围函数
Numberofaircraftlaunchedontime'
）%给x轴加标题ylabel（'
probability'
）%给y轴加标题
set（gca,'
XTick'
0:
2:
16）
7
XTickLabel'
{'
0架'
2架'
4架'
6架'
8架'
10架'
12架'
14架'
16架'
}）
13
0.2460.2512
110.2110.20.210
90.1580.120.113probability70.16
0.05550.050.028140.01970.005533002468101214160架2架4架6架8架10架12架14架16架Numberofaircraftlaunchedontime
图10图11根据例5中的数据，用线性回归的方法画出95%的置信区间，见图11，程序如下:
x=[1:
[p,resid1]=polyfit（x,y,2）;
[yhat,w]=polyconf（p,x,resid1,0.05）;
plot（x,yhat,'
k-'
x,yhat-w,'
k--'
x,yhat+w,'
x,y,'
ks'
[x;
x],[yhat;
y],'
残差的进一步分析:
首先计算残差值，然后用函数normplot将其画在图12中，以观察残差是否服从正态分布。
程序为:
normplot（y-polyval（p,x））
whitebg（'
white'
）%设置背景为白色
90
NormalProbabilityPlot800.98
700.95
0.9060
0.7550
400.50Probability300.25
200.10100.05
0.020155160165170175180185190195-1-0.500.5Data
图12图13由图12看出，这些残差点与直线非常接近，由此得出结论，这些残差值非常接近于正态分布，选择的模型是合适的。
8
升级会员 