第3章MATLAB高等数学计算

MATLAB在高等数学计算中的应用非常广泛，它是工程科学的基础，具有抽象性、逻辑性，主要涵盖傅里叶变换、拉普拉斯变换、Z变换、极限、导数、微分、积分与重积分、级数求和与泰勒级数展开、代数方程组、常微分方程、向量代数和插值运算等。为了将抽象的结果进行图形化展示，本章增加了绘图一节，包括二维绘图和三维绘图。
3.1傅里叶变换与反变换
傅里叶变换是进行信号处理、图像处理、音视频处理的数学工具。
3.1.1傅里叶变换
语法格式： 

F＝fourier(f,t,w)%求时域函数f(t)的傅里叶变换

说明： 返回结果F是符号变量w的函数，省略参数w则默认返回结果为w的函数; f为t的函数，省略参数t时，则默认自由变量为x。
3.1.2傅里叶反变换
语法格式： 

f=ifourier (F)%求频域函数F的傅里叶反变换f(t)

f=ifourier (F,w,t)%求频域函数F指定变量w和t算子的傅里叶反变换f(t)

【例31】求f（t）=1t2+1的傅里叶变换及f（t）的反变换。
程序命令： 

syms t，w; 

F=fourier(1/(t^2+1),t,w)%傅里叶变换 

ft=ifourier(F,t)%傅里叶反变换 

f=ifourier(F)%傅里叶反变换默认x为自变量


结果： 

F=pi*exp(-abs(w))

 ft =1/(t^2 + 1)

 f=1/(x^2 + 1) 





3.2拉普拉斯变换与反变换
拉普拉斯变换是工程数学中常用的一种积分变换，它是将时间域变换到复数s域的方法，自动控制理论中的传递函数大部分使用s域表示。
3.2.1拉普拉斯变换
语法格式： 

F=laplace(f,t,s)%求时域函数f的拉普拉斯变换F

说明： 拉普拉斯变换也称拉式变换，返回结果F为s的函数。当参数s省略时，返回结果F默认为s的函数; f为t的函数，当参数t省略时，默认自由变量为t。
3.2.2拉普拉斯反变换
语法格式： 

f＝ilaplace(F,s,t)%求F的拉普拉斯反变换f

说明： 把s域转换成t域的函数。
【例32】求f(t)=cos(at)+sin(at)的拉普拉斯变换和反变换。
程序命令： 

symsa，t，s; 

F1=laplace(sin(a*t)+cos(a*t),t,s) 

f=ilaplace(F1)

fx=ilaplace(sym('1/s'))

结果： 

F1 =a/(a^2 + s^2) + s/(a^2 + s^2)

f =cos(a*t) + sin(a*t)

fx=1

3.3Z变换与反变换
Z变换是将离散系统的时域利用差分方程转化到频域数学模型的一种方法，离散系统传递函数使用Z变换表示。
3.3.1Z变换
Z变换是对连续系统进行离散数学变换，常用于求解线性时不变差分方程的解。
语法格式： 

Z= ztrans(f)%求Z变换

3.3.2Z反变换
将离散系统变换成连续系统的变换称为Z反变换。
语法格式： 

fz＝iztrans(z)%求Z反变换

【例33】求f(x)=xe-10x 的Z变换和f(z)=z(z－1)z2+2z+1的反变换。
程序命令： 

syms x,k,z; 

f=x*exp(-x*10); %定义表达式

F=ztrans(f)%求Z变换

Fz=z*(z-1)/(z^2+2*z+1); %定义Z反变换表达式

F1=iztrans(Fz)

结果： 

F =(z*exp(10))/(z*exp(10) - 1)^2

F1 =3*(-1)^n + 2*(-1)^n*(n - 1)

3.4求极限
极限是一种变化状态的描述。若函数中的某一个变量在变大或变小的过程中，逐渐趋近于某一个确定的数值A，但不能到达A，称数值A为极限值。控制系统稳定性分析中，稳态时间就是指在2%或5%误差状态下，到达控制点的时间。
语法格式： 

limit(f,x,a)%求符号函数f(x)的极限值，即计算当变量x趋近于常数a时f(x)函数的极限值

limit(f,a)%求符号函数f(x)的极限值，由于没有指定符号函数f(x)的自变量，则使用该格式。

%此时，符号函数f(x)的变量使用函数findsym(f)确定默认自变量，变量x趋近于a

limit(f)%求符号函数f(x)的极限值。符号函数f(x)的变量使用函数findsym(f)确定默认变量; 

%没有指定变量的目标值时，系统默认变量趋近于0，即a=0的情况

limit(f,x,a,'right')%求符号函数f的极限值，'right'表示变量x从右边趋近于a

limit(f,x,a,'left')%求符号函数f的极限值，'left'表示变量x从左边趋近于a

limit(f,x,a,'inf')%求符号函数f的极限值，'inf'表示变量x趋近于无穷

【例34】已知表达式F1和F2，求F1和F2的极限。

F1=limx→0x(esinx+1)－2(etanx－1)sin3x
F2=limx→∞x+x－xsin(π/6)


程序命令： 

syms x; %定义x为符号变量

f1=(x*(exp(sin(x))+1)-2*(exp(tan(x))-1))/sin(x)^3; 

f2=(sqrt(x+sqrt(x))-sqrt(x))/ sin(pi/6);

F1=limit(f1,x,0)%求函数的极限

F2=limit(f2.x,inf)

结果： 

F1=-1/2

F2=1

3.5求导数
求导数就是求函数的平均变化率，利用MATLAB函数使得求导数变得非常简单。
3.5.1语法格式


diff(s)%没有指定变量和导数阶数，则系统按findsym（）函数确定的默认变量对符号表达式s求

%一阶导数

diff(s,'v')%以v为自变量，对符号表达式s求一阶导数

diff(s,n)%按findsym（）函数确定的默认变量对符号表达式s求n阶导数，n为正整数

diff(s,'v',n)%以v为自变量，对符号表达式s求n阶导数

3.5.2求导数案例
【例35】编写程序求函数f1和f2的导数。

f1=sinx2+3x5+(x+1)3
f1=cosx2+arctan(loge(x))


程序命令： 

syms x; %定义符号变量

f1=sin(x)^2+3*x^5+sqrt((x+1)^3); 

f2=cos(x)^2+atan(log(x))

F1=diff(f1,x)%求函数的极限

F2=diff(f2,x)

结果： 

F1=2*cos(x)*sin(x) + (3*(x + 1)^2)/(2*((x + 1)^3)^(1/2)) + 15*x^4

F2=1/(x*(log(x)^2 + 1))-2*cos(x)*sin(x) 

【例36】已知函数f1,求f1的二阶导数F1及在x=2的值X。

f1=5loge(1+x)
F1=d2fdx
X=d2fdxx=2


程序命令： 

syms x ; 

f1=5/log(1+x); 

F1=diff(f1,x,2)

x=2; 

X=eval(F1)

结果： 

F1=5/(log(x+1)^2*(x+1)^2)+10/(log(x+1)^3*(x+1)^2)

X=1.2983

3.6求积分
MATLAB分别利用int函数和quad函数求积分，int函数可获得解析解，quad函数得到的是小梯形面积求和的解。
3.6.1使用int函数求积分
int函数根据解析的方法求解，对定积分和不定积分均可得到解析解，无任何误差，但速度稍慢。
语法格式： 

int(s)%没有指定积分变量和积分阶数时，系统按findsym（）函数确定的默认变量对被积函数或

%符号表达式s求不定积分

int(s,v)%以v为自变量，对被积函数或符号表达式s求不定积分

int(s,v,a,b)%求以v为自变量的符号表达式s的定积分，a,b分别表示定积分的下限和上限

说明： int(s,v,a,b)求被积函数在区间［a,b］上的定积分。a和b可以是两个具体的数，也可以是符号表达式，还可以是无穷(Inf)。当函数f关于变量x在闭区间［a,b］上可积时，函数返回一个定积分结果； 当a,b中有一个是Inf时，函数返回一个广义积分; 当a,b中有一个为符号表达式时，函数返回一个符号函数。
【例37】求函数∫cos2xsin3xdx的不定积分。
程序命令： 

syms x; 

 f1=cos(2*x)*sin(3*x)

 F1=int(f1)

结果： 

F1 =2*cos(x)^3-cos(x)-(8*cos(x)^5)/5

【例38】求表达式f1和f2的定积分。

f1=∫T/2－T/2(AT2+e－jxt)dt
f2=∫e11x2logexdx


程序命令： 

syms A,t,T,x; 

f1=A*T^2+exp(-j*x*t)

f2=log(x)/x^2

F1=int(f1,t,-T/2, T/2); 

F2=simplify(F1)

F3=int(f2,x,1,exp(1)); 

F4=eval(F3)

结果： 

f1= A*T^2 + exp(-t*x*1i)

f2=log(x)/x^2

F2=A*T^3 + (2*sin((T*s)/2))/x

F4=0.2642

【例39】求f1和f2两个表达式的二重积分。

f1(x)=(x+y)e－xydxdy
f2(x)=loge(x)/y2+y2/x2dxdy,1/2≤x≤2,1≤y≤2


程序命令： 

syms x y; 

f1=(x+1)*exp(-x*y); 

F1=int(int(f1,'x'),'y')

f2=log(x)/y^2+y^2/x^2; 

F2=int(int(f2,x,1/2,2),y,1,2) 

F3=eval(F2)

结果： 

F1=(exp(-x*y)*(x + y))/(x*y)

F2=(5*log(2))/4 + 11/4

F3=3.6164

3.6.2使用quadl函数求积分
计算一元函数的数值积分使用quad函数和quadl函数，它们通过小梯形的面积求和得到积分值，而不是通过解析的方法。当有计算精度限制时，计算速度比int函数快。它们均采用遍历的自适应法计算函数的数值积分，quadl是高阶数值积分，quad是低阶数值积分，它们只能求定积分。
语法格式： 

［Q,Fcnt］=quad(function,a,b)%求function的积分
	
其中，function为被积函数(形式为函数句柄/匿名函数)，a,b分别是积分上下限，［Q,Fcnt］为返回数值积分的结果和函数计算的次数。 
【例310】求表达式f1=∫202x3-x+2dx的定积分。
程序命令： 

F=@(x) 2./(x.^3-x+2); 

［Q ,Fcnt］= quadl(F,0,2)

结果： 

Q=1.7037

Fcnt=48

【例311】已知w=［π/2,π，3π/2］; K=［π/2-1,-2,-3π/2-1］，求表达式Y的定积分。

Y=∫w(1)0x2cos(x)dx-K(1)2+∫w(2)0x2cos(x)dx-K(2)2+

∫w(3)0x2cos(x)dx-K(3)2


程序命令： 

 w=［pi/2,pi,pi*1.5］; 

 K=［pi/2-1,-2,-1.5*pi-1］; 

 F1=@(x)(x.^2.*cos(x) -K(1)).^2; 

 F2=@(x)(x.^2.*cos(x) -K(2)).^2; 

 F3=@(x)(x.^2.*cos(x) -K(3)).^2; 

 y=quadl(F1,0,w(1))+ quadl(F2,0,w(2))+ quadl(F3,0,w(3))

结果： 

y=1.378679143103574e+02

3.7求零点与极值
函数的零点和方程的根是等同的，函数极值一般指极大值或极小值，MATLAB提供了相应函数进行计算。
3.7.1求零点
语法格式： 

x=fzero(fun,x0)%求出函数fun在x0最近的零点

x=fzero(fun,x0,options)%由指定的优化参数options进行最小化

x=fzero(problem)%对problem 指定的求根问题求零点

说明： fzero()函数既可以求某个初始值的零点，也可求区间和函数值的零点。
【例312】求函数f(x)=x5-3x4+2x3+x+3的根。
程序命令： 

 f='x^5-3*x^4+2*x^3+x+3'; 

 x=fzero(f,0)

结果： 

x= -0.7693

因为f(x)是一个多项式，所以可以使用roots命令求出相同的实数零点和复共轭零点,即

p=［1 -3 2 0 1 3］; x=roots(p)

x=

1.8846 +0.58974i

1.8846-0.58974i

2.4286e-16+1i

2.4286e-16 -1i

-0.76929 +0i

【例313】求正弦函数在3附近的零点，并求余弦函数在区间［1，2］的零点。
程序命令： 

fun=@sin; 

fun1=@cos；

x=fzero(fun,3)

x1=fzero(fun1,［1 2］)

结果： 

x =3.1416

x1 =1.5708

3.7.2求极值
fminbnd(f,a,b)函数是对f(x)在［a，b］上求得极小值，求-f(x)的极小值时即可得到f(x)的极大值。当不知道极值所在的范围时，可画出该函数的图形，估计极值范围再使用该命令求取。
语法格式： 

［x,min］=fminbnd(f,a,b)%x为取得极小值的点，min为极小值; f表示函数名，a，b表示

%求取极值的范围

【例314】求表达式f=2e-xsin(x)在［0,8］的极值点。
程序命令： 

 syms x; 

 f='2*exp(-x)*sin(x)'; 

 ［x,min1］=fminbnd(f,0,8)

 ［x,max1］=fminbnd('-2*exp(-x)*sin(x)',0,8)

结果： 

x =3.9270

min1 =-0.0279%极小值点

x =0.7854

max1 =-0.6448

max=-max1= 0.6448%极大值点

3.8求方程的解
方程的求解一般包括对线性方程、符号代数方程和常微分方程的求解。对高次方程求解是比较复杂的计算过程，使用MATLAB函数可替代复杂的编程过程，只用一个函数即可完成方程求解。
3.8.1线性方程组求解
(1) 直接使用左除法求解。
【例315】利用左除法求三元一次方程组的解。

x+y+z=1

3x-y+6z=7

y+3z=4


程序命令： 

A=sym('［111; 3-16; 013］'); 

b =sym('［1; 7; 4］'); 

x=A＼b

结果： 

x =-1/3

0

4/3

(2) 使用solve函数求解。
【例316】使用solve函数求解三元一次方程组。

2x+3y-z=2

8x+2y+3z=4

45x+3y+9z=23


程序命令： 

syms x，y，z；%建立符号变量

 ［x,y,z］=solve(2*x+3*y-z-2, 8*x+2*y+3*z-4, 45*x+3*y+9*z-23)

%使用solve函数求解

结果： 

x=151/273

y=8/39

z =-76/273

3.8.2符号代数方程求解
线性方程组的符号解也可用solve函数求解，当方程组不存在符号解又无其他自由参数时，则给出数值解。
语法格式： 

solve(f,'v')%求一个方程的解

solve(f1,f2,…,fn)%求n个方程的解

说明： f既可以是含等号、符号表达式的方程，也可以不含等号或符号表达式，但所求的解均是令f=0的方程。当参数v省略时，默认为方程中的自由变量; 其输出结果为结构数组类型。
【例317】已知方程为ax2+bx+c=0,x2-x-30=0,求符号方程解及数值解。
程序命令： 

 syms a，b，c，x；

 f1=a*x^2+b*x+c; 

 f2=x^2-x-30; 

 Fx=solve(f1,x) %对默认变量x求解

 Fb=solve(f1,b )%对指定变量b求解

 F2=solve(f2,x)

结果： 

Fx=-(b + (b^2 - 4*a*c)^(1/2))/(2*a)

-(b - (b^2 - 4*a*c)^(1/2))/(2*a)

Fb=-(a*x^2 + c)/x

F2=-5

 6

3.8.3常微分方程(组)求解
在MATLAB中，用大写字母D表示导数。例如，Dy表示y′，D2y表示y″，Dy(0)=5表示y′(0)=5。D3y+D2y+Dy-x+5=0表示微分方程y+y″+y′-x+5=0。
语法格式： 

dsolve( )%求解符号常微分方程的解

dsolve(e,c,v)%求解常微分方程e在初值条件c下的特解

说明： 参数v描述方程中的自变量，省略时默认自变量是t。若没有给出初值条件c，则求的是方程的通解。
dsolve在求常微分方程组时的调用格式为

dsolve(e1,e2,…,en,c1,…,cn,v1,…,vn)


该函数求解常微分方程组e1,…,en在初值条件c1,…,cn下的特解。若不给出初值条件，则求方程组的通解，v1,…,vn给出求解变量，若省略自变量，则默认自变量为t。若找不到解析解，则返回其积分形式。
【例318】求微分方程 dydx+2xy=xe-x2的通解。
程序命令： 

syms x; 

f= 'Dy+2*x*y=x*exp(-x^2)'; 

 y=dsolve(f,x)

结果： 

y=C1*exp(-x^2) + (x^2*exp(-x^2))/2

【例319】求微分方程xdydx+y-ex=0在初值条件y(1)=2e下的特解。
程序命令： 

 syms x; 

 eq1='x*Dy+y-exp(x)=0'; 

cond1='y(1)=2*exp(1)'

y=dsolve(eq1,cond1,x)

结果： 

cond1=y(1)=2*exp(1)

y=(exp(1) + exp(x))/x

【例320】求下列微分方程组的通解。

d2xdt+2dxdt=x(t)+2y(t)-e-t

dydt=4x(t)+3y(t)+4e-t


程序命令： 

［x,y］=dsolve('D2x+2*Dx=x+2*y-exp(-t)','Dy=4*x+3*y+4*exp(-t)')

结果： 

x =

exp(t*(6^(1/2) + 1))*(6^(1/2)/5 - 1/5)*(C2 + exp(- 2*t - 6^(1/2)*t)*((11*6^(1/2))/3 - 37/4)) - exp(-t)*(C1 + 6*t) - exp(-t*(6^(1/2) - 1))*(6^(1/2)/5 + 1/5)*(C3 - exp(6^(1/2)*t - 2*t)*((11*6^(1/2))/3 + 37/4))

 y =

 exp(-t)*(C1 + 6*t) + exp(t*(6^(1/2) + 1))*((2*6^(1/2))/5 + 8/5)*(C2 + exp(- 2*t - 6^(1/2)*t)*((11*6^(1/2))/3 - 37/4)) - exp(-t*(6^(1/2) - 1))*((2*6^(1/2))/5 - 8/5)*(C3 - exp(6^(1/2)*t - 2*t)*((11*6^(1/2))/3 + 37/4))

3.9级数
级数指将数列的项依次用加号连接起来的函数，级数操作一般包括级数求和与级数展开，它是数学基础分析的一个分支。
3.9.1级数求和
级数求和运算是数学中常见的一种运算。例如： 

s=a0+a1x+a2x2+a3x3+…+anxn


函数symsum用于通式项表达式s的求和运算。当确定变量个数为n，变量变化范围为a到b时可进行级数求和。
语法格式： 

(1) symsum(n)%求自变量n从0到n-1的前n项的和

(2) symsum(s,n)%s为级数通式项，求自变量n从0到n-1的前n项的和

(3) symsum(s,a,b)%s为级数通式项，默认变量为n，求从a到b的和

(4) symsum(s,n,a,b)%s为级数通式项，求自变量n从a到b的和

当默认变量不变时，前两种用法基本相同，后两种用法基本相同。
【例321】求级数S=∑∞k=11(n+1)2的前n项和、从1到无穷及其前10项和。
程序命令： 

clc; syms n; 

 s=1/(n+1)^2

S1=symsum(s)

 S2=symsum(s,n)

 S10=symsum(s,1,10)%求级数前10项和

 S20=symsum(s,n,1,10)%求级数前10项和

结果： 

s=1/(n + 1)^2

 S1=piecewise(-1 < n,-psi(1,n + 1),n <= -1,psi(1,-n))

 S2=piecewise(-1 < n,-psi(1,n + 1),n <= -1,psi(1,-n))

 S10=85758209/153679680

 S20=85758209/153679680

说明： 
(1) piecewise表示分段函数,n>-1时，结果为-psi(1,n+1); n≤-1时，结果为
psi(1,-n)。
(2) psi(X)为数组X的每个元素计算psi函数。X必须是非负实数。psi函数也称为双γ函数，是γ函数的对数导数。若 psi(k,X)为在X的元素中计算psi函数的第k个导数。psi(0,X)是双γ函数，psi(1,X) 是三γ函数，psi(2,X)是四γ函数，以此类推。
3.9.2一元函数的泰勒级数展开
语法格式： 

taylor(f)%求f关于默认变量的5阶泰勒级数展开

taylor(f,n)%求f关于默认变量的n-1阶泰勒级数展开

taylor(f,n,v)%求f关于变量v的n-1阶泰勒级数展开

taylor(f,n,v,a)%求f在v=a处的n-1阶泰勒级数展开

【例322】已知函数表达式ex=∑∞n=0xnn!,sinx=∑∞n=0(-1)n(2n+1)!x2n+1,cosx=∑∞n=0(-1)n2n!x2n，求三个函数表达式的泰勒级数展开式。
程序命令： 

 syms x; 

 f1=exp(x); f2=sin(x); f3=cos(x); 

 taylorexpx=taylor(f1)

 taylorsinx=taylor(f2)

 taylorcosx=taylor(f3)

结果： 

taylorexpx =x^5/120 + x^4/24 + x^3/6 + x^2/2 + x + 1

taylorsinx =x^5/120 - x^3/6 + x

taylorcosx =x^4/24 - x^2/2 + 1

3.10常用绘图功能
MATLAB具有丰富的绘图功能，提供了一系列的绘图函数，不仅包括常用的二维图，还包括三维函数、三维网格、三维网面及三维立体切片图等，使用系统绘图函数不需要过多编程，只需给出一些基本参数就能得到所需图形。此外，MATLAB还提供了直接对图形句柄进行操作的低层绘图操作，包括图形元素，例如坐标轴、曲线、文字等，系统对每个对象分配了句柄，通过句柄即可对该图形元素进行操作，而不影响其他部分。
3.10.1二维绘图
MATLAB的plot()函数是绘制二维图形最基本的函数，它是针对向量或矩阵列来绘制曲线的，绘制以x轴和y轴为线性尺度的直角坐标曲线。
1. 语法格式


plot(x1,y1,option1,x2,y2,option2,…)

说明： x1,y1,x2,y2给出的数据分别为x轴和y轴坐标值，option定义了图形曲线的颜色、字符和线型，它由一对单引号括起来。可以画一条或多条曲线。若x1和y1都是数组，按列取坐标数据绘图。
2. option的含义
option通常由颜色(见表3.1)、字符(见表3.2)和线型(见表3.3)组成。


表3.1颜色表示




选项
含义
选项
含义
选项
含义
'r'
红色
'w'
白色
'k'
黑色
'g'
绿色
'y'
黄色
'm'
锰紫色
'b'
蓝色
'c'
亮青色
——


表3.2字符表示




选项
含义
选项
含义
选项
含义
'.'
画点号
'o'
画圈符
'd'
画菱形符
'*'
画星号
'+'
画十字符
'p'
画五角形符
'x'
画叉号
's'
画方块符
'h'
画六角形符
'^'
画上三角
'>'
画左三角
——

'^'
画下三角
'<'
画右三角
——


表3.3线型表示




选项含义
选项
含义
'-'
画实线
'.-'
点画线
'--'
画虚线
'：'
画点线
【例323】绘制表达式y=2e-0.5t sin(2πt)的曲线。
程序命令： 

t=0： pi/100： 2*pi; y1=2*exp(-0.5*t).*sin(2*pi*t); 

y2=sin(t); plot(t,y1,'b-',t,y2,'r-o')

结果如图3.1所示。


图 3.1例323函数曲线


【例324】绘制表达式x=tsin3t，y=tsintsint的曲线。
程序命令： 

t=0： 0.1： 2*pi; x=t.*sin(3*t); y=t.*sin(t).*sin(t); 

plot(x,y,'r-p'); 

结果如图3.2所示。


图3.2例324函数曲线


3. 图形屏幕控制命令


figure %打开图形窗口

clf %清除当前图形窗口的内容

hold on%保持当前图形窗口的内容

hold off %解除保持当前图形状态

grid on %给图形加上栅格线

grid off %删除栅格线 

box on%在当前坐标系中显示一个边框

box off%去掉边框

close %关闭当前图形窗口

close all %关闭所有图形窗口

【例325】在不同窗口绘制y1=cos(t),y2=sin2(t)的波形图。
程序命令： 

t=0： pi/100： 2*pi; y1=cos(t); y2=sin(t).^2; 

figure(1); plot(t,y1,'g-p'); grid on; figure(2); plot(t,y2,'r-O'); grid on; 

结果如图3.3所示。


图3.3不同窗口绘图


4. 图形标注


title%图题标注

xlabel%x轴说明

ylabel%y轴说明

zlabel%z轴说明

legend%图例标注，legend函数用于绘制曲线所用线型、颜色或数据点标记图例

legend('字符串1','字符串2',…)%按指定字符串顺序标记当前轴的图例 

legend(句柄,'字符串1','字符串2',…)%指定字符串标记句柄图形对象图例 

legend(M)%用字符M矩阵的每一行字符串作为图形对象标签标记图例

legend(句柄,M)%用字符M矩阵的每一行字符串作为指定句柄的图形对象标签标记图例

text%在图形中指定的位置(x,y)处显示字符串string，格式： text (x,y,'string')

annotation %线条、箭头和图框标注，
例如，annotation('arrow',［0.1,0.45］,［0.3,0.5］)绘制箭头线

5. 字体属性
字体属性如表3.4所示。


表3.4字体属性




属性名
注释
属性名
注释
FontName
字体名称
FontWeight
字形
FontSize
字体大小
FontUnits
字体大小单位
FontAngle
字体角度
Rotation
文本旋转角度
BackgroundColor
背景色
HorizontalAlignment
文字水平方向对齐
EdgeColor
边框颜色
VerticalAlignment
文字垂直方向对齐
说明： 
(1) FontName属性定义名称，其取值是系统支持的一种字体名。
(2) FontSize属性设置文本对象的大小，其单位由FontUnits属性决定，默认值为10磅。
(3) FontWeight属性设置字体粗细，取值可以是normal(默认值)、bold、light或demi。
(4) FontAngle属性设置斜体文字模式，取值可以是normal(默认值)、italic或oblique。
(5) Rotation属性设置字体旋转角度，取值是数值量，默认值为0，取正值时表示逆时针方向旋转，取负值时表示顺时针方向旋转。
(6) BackgroundColor和EdgeColor属性设置文本对象的背景颜色和边框线的颜色，可取值为none(默认值)或颜色字母。
(7) HorizontalAlignment属性设置文本与指定点的相对位置，可取值为left(默认值)、center或right。
6. 坐标轴axis的用法
语法格式： 

axis(［xmin，xmax，ymin，ymax］)或axis(［xmin，xmax，ymin，ymax，zmin，zmax］)

说明： 该函数用来标注输出图线的坐标范围。若给出4个参数则标注二维曲线最大值和最小值，给出6个参数则标注三维曲线最大值和最小值。其中： 

axis equal%将两坐标轴设为相等

axis on/off%显示/关闭坐标轴的显示

axis auto%将坐标轴设置为默认值

axis square%产生两轴相等的正方形坐标系

7. 子图分割
语法格式： 

subplot(n,m,p)

其中，n表示行数，m表示列数，p表示绘图序号，顺序是按从左至右、从上至下排列，它把图形窗口分为n*m个子图，在第p个子图处绘制图形。
【例326】利用子图绘制正弦和余弦图形。
程序命令： 

t=0： pi/100： 2*pi; y1=sin(t); y2=cos(t); y3=sin(t).^2; y4=cos(t).^2; 

subplot(2,2,1),plot(t,y1); title('sin(t)'); subplot(2,2,2),plot(t,y2,'g-p'); title('cos(t)')

subplot(2,2,3),plot(t,y3,'r-O'); title('sin^2(t)'); subplot(2,2,4),plot(t,y4,'k-h'); title('cos^2(t) ')

结果如图3.4所示。


图3.4绘制子图



3.10.2三维绘图
1. 绘制三维空间曲线

与plot()函数相类似，可以使用plot3()函数来绘制一条三维空间的曲线。
语法格式： 

plot3(x,y,z,option)%绘制三维曲线


其中，x,y,z以及选项option与plot()函数中的x,y和选项相类似，多了一个z坐标轴。绘图方法可参考plot()函数的使用方法。option指定颜色、线形等。
【例327】已知函数，绘制［0，2π］区间三维函数曲线。

x=(8+3cos(V))cos(U)
y=(8+3cos(V))sin(U)，0<U,V≤2
z=3sin(V)


程序命令： 

r =linspace(0,2*pi,60); ［u,v］=meshgrid(r); 

x=(8+3*cos(v)).*cos(u); y=(8+3*cos(v)).*sin(u); z=3*sin(v); 

plot3(x,y,z); title('三维空间绘图'); xlabel('x轴'); ylabel('y轴'); zlabel('z轴')

结果如图3.5所示。


图3.5三维空间曲面图


2. 网格矩阵的设置
meshgrid函数产生二维阵列和三维阵列。用户需要知道各个四边形顶点的三维坐标值(x,y,z)。
语法格式： 

［X,Y］=meshgrid(x,y)%向量x，y分别指定x轴向和y轴向的数据点。当x为n维向量，y为

%m维向量时，X,Y均为m*n的矩阵。［X,Y］=meshgrid(x)等效于

%［X,Y］=meshgrid(x,x)

语法格式： 

［X,Y,Z］=meshgrid(x,y,z)%产生x轴、y轴和z轴的三维阵列

3. 绘制三维网格曲面图
语法格式： 

mesh(x，y，z，c)


说明： 
(1) 三维网格图是一些四边形相互连接在一起构成的一种曲面图。
(2) x，y，z是维数相同的矩阵，x和y是网格坐标矩阵，z是网格点上的高度矩阵，c用于指定在不同高度下的颜色范围。
(3) c省略时，c=z，即颜色的设定正比于图形的高度。
(4) 当x和y是向量时，要求x的长度必须等于矩阵z的列长度，y的长度必须等于矩阵z的行长度，x和y元素的组合构成网格点的x轴和y轴坐标，z轴坐标则取自z矩阵，然后绘制三维曲线。
【例328】根据函数z=f(x,y)的x轴和y轴坐标找出z的高度，绘制Z=x2+y2的三维网格图形。
程序命令： 

x=-5： 5; y=x; ［X,Y］=meshgrid(x,y)

Z=X.^2+Y.^2 ; mesh(X,Y,Z)

结果如图3.6所示。


图3.6三维网格曲面图


4. 绘制三维曲面图
语法格式： 

surf(x，y，z，c) 

说明： 参数x,y,z,c同mesh函数，它们均使用网格矩阵meshgrid函数产生坐标，自动着色，其三维阴影曲面四边形的表面颜色分布通过shading命令指定。
【例329】绘制马鞍函数z=f(x,y)=x2-y2在［-10，10］区间的曲面图。
程序命令： 

x=-10： 0.1： 10 ; ［xx,yy］=meshgrid(x); zz =xx .^2-yy .^2; 

surf(xx,yy,zz ); title('马鞍面'); xlabel('x轴'); ylabel('y轴') zlabel('z轴'); grid on; 

结果如图3.7所示。


图3.7三维曲面图


5. 绘制特殊的三维立体图
1) 球面图
MATLAB提供了球面和柱面等标准的三维曲面绘制函数，使用户可以很方便地得到标准三维曲面图。 
语法格式： 

sphere(n)%画n等分球面，默认半径为1，n=20,n表示球面绘制的精度

或

［x,y,z］=sphere(n)%获取球面在三维空间的坐标

2) 柱面图
语法格式： 

cylinder(R,n)%R为半径，n为柱面圆周等分数

或

［x,y,z］=cylinder (R,n)%x，y，z代表空间坐标; 若在调用该函数时不带输出参数，则直接绘制

%所需柱面; n决定了柱面的圆滑程度，其默认值为20; 若n值取得比

%较小，则绘制出多面体的表面图


3) 利用多峰函数绘图 
多峰函数为

f(x,y)=3(1－x)2e－x2－(y+1)2－10x5－x3－y5e－x2－y2－13e－(x+1)2－y2


语法格式： 

peaks(n)%输出n*n大小的矩阵峰值函数图形

或

［x,y,z］=peaks(n)%x，y，z代表空间坐标


【例330】使用子图分割分别绘制球面图、柱面图、多峰函数图和函数z=f(x,y)=sinx2+y2x2+y2在［-10，10］区间的曲面图。
程序命令： 

x=-10： 0.5： 10 ; ［xx,yy］=meshgrid(x) 

zz =sin(sqrt(xx.^2+yy.^2))./sqrt(xx.^2+yy.^2) 

subplot(2,2,1); surf(xx,yy,zz ); title('函数图'); xlabel('x轴'); ylabel('y轴'); zlabel('z轴'); 

subplot(2,2,2); sphere(20); title('函数图'); xlabel('x轴'); ylabel('y轴'); zlabel('z轴'); 

subplot(2,2,3); cylinder(20); title('函数图'); xlabel('x轴'); ylabel('y轴'); zlabel('z轴'); 

subplot(2,2,4); peaks; title('函数图'); xlabel('x轴'); ylabel('y轴'); zlabel('z轴'); 

结果如图3.8所示。


图3.8例330函数图像


6. 立体切片图
语法格式： 

slice(X,Y,Z,V,x0,y0,z0)%绘制向量x1,y1,z1中的点沿x轴、y轴、z轴方向的切片图，V的

%大小决定了每一点的颜色

slice(V,x0,y0,z0)%按数组x1,y1,z1的网格单位进行切片

说明： 
(1) 该函数是显示三维函数V=V(X,Y,Z)确定的超立体形在x轴、y轴与z轴方向上的若干点切片图。在V=V(X,Y,Z)中的变量X取一定值X1，则函数V=V(X1,Y,Z)变成带颜色的立体曲面切片图，各点的坐标由向量x0,y0,z0指定。参量X,Y,Z为三维数组，用于指定立方体V的坐标。参量X、Y与Z必须有单调的正交间隔，如同用命令meshgrid生成的一样，在每一点上的颜色由对超立体V的三维内插值确定。

(2) slice(V,x0,y0,z0)省略了X,Y,Z，默认取值为X =1： m，Y=1： n，Z=1： p。m,n,p为V的三维数组(阶数)。
【例331】绘制三维函数V=f(x,y,z)=xe-x2-y2-z2的立体切片图，要求切片的坐标为(［-0.5,0.8,1.5］,0.5,［］)。
程序命令： 

［x,y,z］=meshgrid(-2： 0.2： 2,-2： 0.25： 2,-2： 0.2： 2)

V=x.*exp(-x.^2-y.^2-z.^2); 

slice(x,y,z,V,［-0.5,0.8,1.5］,0.5,［］); 

绘图的结果如图3.9所示。


图3.9三维切片图


3.11函数插值
已知函数表达式的情况下，插值就是在已知的数据点之间利用某种算法寻找估计值的过程，即根据一元线性函数表达式f(x)中的两点，找出f(x)在中间点的数值。插值运算可大大减少编程语句，使得程序简洁清晰。
3.11.1一维插值
MATLAB提供的一维插值函数为interp1( )，定义如下： 
若已知在点集x上的函数值y，构造一个解析函数曲线，通过曲线上的点求出它们之间的值，这一过程称为一维插值。
语法格式： 

yi=interp1(x,y,xi); %x,y为已知数据值，xi为插值数据点; 

y1=interp1(x,y,xi,'method'); %x,y为已知数据值，xi为插值点，method为设定插值方法

说明： method常用的设置参数有linear，nearest和spline，分别表示线性插值、最邻近插值和三次样条插值。linear也称为分段线性插值(默认值)，spline函数插值所形成的曲线最平滑、效果最好。
(1) nearest(最邻近插值)： 该方法将内插点设置成最接近于已知数据点的值，其特点是插值速度最快，但平滑性较差。
(2) linear(线性插值)： 该方法连接已有数据点作线性逼近。它是interp1（）函数的默认方法，其特点是需要占用更多的内存，速度比nearest方法稍慢，但是平滑性优于nearest方法。
(3) spline(三次样条插值)： 该方法利用一系列样条函数获得内插数据点，从而确定已有数据点之间的函数。其特点是处理速度慢，但占用内存少，可以产生最光滑的插值结果。
【例332】先绘制0~2π的正弦曲线，按照线性插值、最邻近插值和三次样条插值三种方法，每隔0.5进行插值，绘制插值后曲线并进行对比。
程序命令： 

clc; x=0:2*pi; y=sin(x); 

 xx=0:0.5:2*pi; 

 subplot(2,2,1); plot(x,y); title('原函数图'); 

 y1=interp1(x,y,xx,’linear’); 

 subplot(2,2,2); plot(x,y,'o',xx,y1,'r'); title('线性插值')

 y2=interp1(x,y,xx,'nearest'); 

 subplot(2,2,3); plot(x,y,'o',xx,y2,'r'); title('最邻近插值')

 y3=interp1(x,y,xx,'spline'); 

 subplot(2,2,4); plot(x,y,'o',xx,y3,'r'); title('三次样条插值')

三种插值与原函数插值结果如图3.10所示。


图3.10三种插值与原函数图的比较


结论： 从线性插值、最邻近插值和三次样条插值三种方法可以看出，三次样条插值的曲线效果最平滑。
【例333】设某一天24h内，从零点开始每间隔2h测得的环境温度数据分别为
12,9,9,10,18,24,28,25,20,16,12,11,11,请推测上午9时的温度。
程序命令： 

x=0:2:24; 

y=［12,9 ,9,10,18,24,28,25,20 ,16,12,11,11］; 

x1=9; 

y1=interp1(x,y,x1,'spline'); 

plot(x,y,'b-p',x,y,x1,y1,'r-h'); 

title('插值点绘图'); 

text (x1-3,y1,'插值点'); 

结果曲线如图3.11所示。


图3.11插值绘图


【例334】假设2000年至2020年的产量每间隔2年数据分别为90，105，123，131，150，179，203，226，249，256，267，试估计2015年的产量并绘图。
程序命令： 

clear; 

year=2000:2:2020；

product=［ 90 105 123 131 150 179 203226 249 256 267 ］; 

x =2000:1:2020; 

y=interp1(year,product,x); 

p2015=interp1(year,product,2015)

plot(year,product,'b-O',x,y)

title('2000年-2020年的产量')

插值结果： 

p2015=237.5000

默认插值曲线如图3.12所示。


图3.12默认插值曲线


【例335】对离散分布在y=exsinx函数曲线上的数据，分别进行三次样条插值和线性插值计算，并绘制曲线。
程序命令： 

clear; 

 x=［0 2 4 5 8 12 12.8 17.2 19.9 20］; 

 y=exp(x).*sin(x); 

 xx=0:.25:20; 

 yy=interp1(x,y,xx,’spline’); 

 plot(x,y,'o',xx,yy); hold on

 yy1=interp1(x,y,xx,’linear’); 

 plot(x,y,'o',xx,yy1); hold on

结果： 
插值后曲线如图3.13所示。


图3.13三次样条插值与线性插值绘图


3.11.2二维插值
二维插值函数为二元函数。
语法格式： 

ZZ=interp2(X,Y,Z,X1,Y1)%X和Y分别是m维和n维向量，表示节点; Z为n*m的矩阵，表示

%节点值; X1(行向量)和Y1(列向量) 是插值点的一维数组，为插值

%范围，若插值在范围外的点，则返回 NaN(数值为空)

ZZ=interp2(Z,X1,Y1)%表示X1=1：n,Y1=1：m，其中［m,n］=size(Z)

ZZ=interp2(X,Y,Z,X1,Y1,method)%用指定的方法method计算二维插值，method可以取值为

%
linear(双线性插值算法,默认值)、
nearest（最邻近插值）、

%
spline(三次样条插值)或
cubic（双三次插值）

说明： interp2（）函数能够较好地进行二维插值运算，但是它只能处理以网格形式给出的数据。
【例336】已知工人的平均工资从1980年到2020年逐年提升，计算在 2000 年工作了12年的员工的平均工资。
程序命令：

years=1980:10:2020; 

times= 10:10:30; 

salary= ［1500 1990 2000 3010 3500 4000 4100 4200 4500 5600 7000 8000 9500 10000 12000］; 

S= interp2(service,years,salary,12,2000)

结果： 

S=4120 

【例337】对函数z=f(x,y)=sinx2+y2x2+y2进行不同插值拟合曲面，并比较拟合结果。
程序命令： 

［x,y］=meshgrid(-8:0.8:8); 

 z=sin(sqrt(x.^2+y.^2))./sqrt(x.^2+y.^2)；

 subplot(2,2,1); surf(x,y,z); title('原图'); 

［x1,y1］=meshgrid(-8:0.5:8); 

 z1=interp2(x,y,z,x1,y1); 

 subplot(2,2,2); surf(x1,y1,z1); title('linear'); 

 z2=interp2(x,y,z,x1,y1,'cubic'); 

 subplot(2,2,3); surf(x1,y1,z2); title('cubic'); 

 z3=interp2(x,y,z,x1,y1,'spline'); 

 subplot(2,2,4); surf(x1,y1,z3); title('spline'); 

程序运行结果如图3.14所示。


图3.14四种插值曲面拟合结果


3.11.3三维插值
三维插值函数 interp3（）和n维网格插值函数interpn（）的调用格式与函数interp1（）和interp2（）一致，需要使用三维网格生成函数实现，即［X,Y,Z］=meshgrid(X1,Y1,Z1)，其中X1,Y1,Y1为三维所需要的分割形式，以向量形式给出三维数组，目的是返回X,Y,Z的网格数据。
语法格式： 

interp3(X,Y,Z,V,X1,Y1,Z1,method )

说明： V表示函数，使用方法与interp2()函数一致。
【例338】已知三维函数V(x,y,z)=x2+y2+z2，通过函数生成网格型样本点，根据样本点进行拟合,并使用slice函数绘制拟合图。
程序命令： 

clc; 

［x,y,z］=meshgrid(-1:0.1:1); 

［x1,y1,z1］=meshgrid(-1:0.1:1); 

V=x.^2+y.^2+z.^2; %拟合函数

V1=interp3(x,y,z,V,x1,y1,z1,'spline'); %三维拟合

x0=［-0.5,0.5］; y0=［0.2,-0.2］; z0 =［-1,-0.5,0.5］; %图形位置

slice(x1,y1,z1,V1,x0,y0,z0); title('三维拟合'); 

程序运行结果如图3.15所示。


图3.15三维插值图