非矩形区域上图形绘制的MATLAB实现*

翟龙余

(江苏财经职业技术学院基础教学部，江苏淮安223003)

MATLAB系统是目前世界上最流行的一种数值计算和图形、图像处理工具软件，目前已广泛应用于高等数学辅助教学中.比如，使用MATLAB软件辅助教学，可以很容易绘制出复杂的立体图形，把曲线、曲面的形成和变化过程准确地模拟出来，这对提高教学效率和培养学生的空间想象能力可起到事半功倍的效果.

如果要在空间直角坐标系下，绘制z=f(x，y)，其中(x，y) ∈［a，b］×［c，d］所代表的三维曲面图形，可按如下步骤，用MATLAB进行绘图［1］.

(1)生成平面网格坐标矩阵：先将［a，b］在x方向划分出m个等分点，再将［c，d］在y方向划分出n个等分点，最后利用meshgrid函数生成网格点的坐标，横坐标放在x矩阵中，纵坐标放在y矩阵中，两个矩阵都是n×m的矩阵；

(2)计算z坐标的矩阵：在上一步的基础上利用MATLAB的数值计算功能计算出函数z=f(x，y)的函数值，生成z坐标的矩阵；

(3)绘制三维曲面：MATLAB提供了mesh函数和surf函数来绘制三维图形，mesh函数用于绘制三维网格图，surf函数用于绘制三维曲面图.

高等数学中的函数z=f(x，y)的自然定义域不一定是矩形区域［a，b］×［c，d］，很多函数的自然定义域是一般形状的平面区域，还有实际应用问题的研究往往不是在函数的自然定义域上研究，而是在一个与具体问题相关的非矩形区域上研究，这种非矩形区域上函数的三维曲面绘制在MATLAB的相关教材中都没有介绍，本文就此结合数学中的典型实例研究了非矩形区域上空间三维曲面图形的绘制.

1 非矩形区域上的MATLAB图形绘制——自变量赋nan法

文献［2］介绍了非矩形区域上MATLAB图形绘制的一种方法：先把区域分割成若干个窄条矩形，利用MATLAB中绘制空间曲面的库函数作出窄条矩形的图形，最后再进行叠加，得到任意区域上的函数的图形.但是，按照该法，如果区域的边界线较为复杂时，所取点有可能落到区域外，如果该区域是函数的自然定义域时，还有可能绘制不出图形.

在MATLAB中，有一个特殊的变量nan(Not a Number)，表示“非数”，MATLAB在处理数据时，会跳过数据中的nan，所有的MATLAB作图函数都会忽略含有nan的数据点，实现一种所谓“镂空”的效果［3］.受此启发，可以通过如下方法作出边界线已知的任意平面区域上的三维图形：

(1)用上述矩形区域绘图法步骤1的方法，生成一个矩形区域的网格坐标，要求该矩形区域包含所绘图形的绘图区域；

(2)剔除绘图区域以外点：先计算出绘图区域边界线上函数值，然后通过MATLAB循环与条件判断语句，一一判断网格点是否在作图区域内，如不在，把区域外的网格点的某个坐标值重新赋值为nan；

(3)按照上述矩形区域作图法步骤2和3的方法进行绘图，就可绘制出所需图形.

上述绘图方法关键之处是把绘图区域外的点的一个坐标的值即函数的一个自变量的值重新赋值为nan，所以称此绘图方法为自变量赋nan法.

以x型区域为例，说明该绘图法的编程过程：

由y=g(x)，y=h(x)，x=a，z=b围成x型区域D(其中g(x)≥h(x))，现要求绘制二元函数z=f(x，y)，(x，y) ∈D的三维曲面图，具体步骤如下：

(1)任取一矩形区域D0，使D0⊇D，如取D0=［a，b］×［gmax，hmin］，其中gmax是函数g(x) 的最大值，hmin是函数h(x)的最小值；

(2)针对D0区域，通过meshgrid命令生成网格点矩阵x和y；

(3)以网格点的横坐标矩阵x，作为自变量，计算出相应的区域D的两条边界线上的函数值矩阵g和h；

(4)通过循环语句，比较矩阵y的元素与g和h相应的元素之间大小，若y(i，j) ＞g(i，j)，或y(i，j)＜h(i，j) 时，y(i，j) 重新赋值为 nan，即y(i，j)=nan.实际上此时纵坐标为y(i，j)所对应的点在区域D以外；

(5)计算函数z=f(x，y)的函数值；

(6)利用 mesh(x，y，z)(或 surf(x，y，z))命令绘制图形，即可得到所需图形.

实例1 绘制马鞍面z=xy，(x，y)∈ {(x，y)||x|+|y|≤10}.

该图在空间关于x=y面对称，其绘图区域在xοy面上既关于x轴对称，又关于y轴对称.下面的图形绘制充分利用了其对称性的特点.绘图程序如下：

%选取矩形区域［0，10］×［0，10］生成网格点，该区域包含了绘图区域的第一象限部分

x1=0：0.1：10；

y1=0：0.1：10；

n=length(x1)；

m=length(y1)；

［x1，y1］=meshgrid(x1，y1)；

%计算出第一象限边界线上的函数值矩阵；

y=10-x1；

%以下循环将第一象限内绘图区域外的点的纵坐标值重新赋值为nan.

for i=1：m；

for j=1：n；

if y1(i，j) ＞ y(i，j)；

y1(i，j)=nan；

end

end

end

%利用对称性分别计算出第一、二、三、四象限区域上的函数值矩阵 z1、z2 、z3 、z4.

z1=x1.*y1；

z2=-x1.*y1；

z3=z1；

z4=z2；

%在同一图形窗口绘制出一、二、三、四象限绘图区域上的图形，最终完成马鞍面图形的绘制.

mesh(x1，y1，z1)

hold on

mesh(-x1，y1，z2)

mesh(-x1，-y1，z3)

mesh(x1，-y1，z4)

得到如图1所示马鞍面，看到图1就知道了它称为马鞍面的原因.

图1 马鞍面

图2 半椭球面

如果函数z=f(x，y)的绘图区域为y型区域，其绘图方法与上类似，通过比较横坐标的值，把绘图区域外的点的横坐标值重新赋值为nan即可.如果绘图区域是更一般区域，可把其分割成若干个x型区域或y型区域，然后在同一图形窗口分别用上述方法绘制出每个区域图形即可.

2 函数z=f(x，y)的非矩形自然定义域上的三维曲面图形绘制——函数值赋nan法

当函数z=f(x，y)的自然定义域为非矩形区域D时，显然此时图形绘制属于非矩形区域上的图形绘制，因而用上述自变量赋nan法即可绘制出其图形.

但对于自然定义域为非矩形区域D的函数z=f(x，y)，多数情况下，如果给定一个矩形区域D0，D0包含函数z=f(x，y)自然定义域D，即D0⊇D.当(x，y)∈D0且(x，y)∉D时，MATLAB软件计算z=f(x，y)的函数值得到结果是虚数或inf，而MATLAB绘图时，含inf坐标的数据点在图中是不会显示的(也不会出现错误信息)，如果此时我们用nan直接替换函数值矩阵中的虚数，再用mesh(x，y，z)(或surf(x，y，z))命令绘制图形，可以更简便地绘制出其自然定义域D上的图形.该方法关键之处是把函数值矩阵中的虚数重新赋值为nan，所以称此绘图方法为函数值赋nan法，其具体步骤如下：

(1)任取一矩形区域D0，使D0包含自然定义域D；

(2)针对D0区域，通过meshgrid命令生成网格点矩阵x和y；

(3)直接计算函数z=f(x，y)的函数值；

(4)通过循环语句，一一判断z(i，j)的虚部是否为，当其虚部不为时，z(i，j)重新赋值为nan，即z(i，j)=nan.实际上该步骤作用就是剔除自然定义域以外的点；

(5)利用 mesh(x，y，z)(或 surf(x，y，z))命令绘制图形，即可得到所需图形.

解 该函数的定义域为

y=-3：0.0025：3；

m=length(x)；% 计算x中元素的个数；

n=length(y)；

［x，y］=meshgrid(x，y)；% 生成矩形区域 的网格点坐标；

z=2.*(1-x.^2./4-y.^2./9).^(1/2)；% 计算函数值

for i=1：n；

for j=1：m；

if imag(z(i，j))～=0；% 查找z中的虚数；

z(i，j)=nan；% 将z中的虚数重新赋值为nan；

end

end

end

mesh(x，y，z)

其图形如图2所示为半椭球面.

注意点：以上两种绘图方法，由于受MATLAB计算精度以及计算机硬件条件的影响，在编程过程中x与y的步长选择要适当，不要太小也不要太大.

［1］刘卫国.MATLAB程序设计与应用［M］.第二版.北京：高等教育出版社.2011.

［2］任明慧，张月莲.MATLAB中矩形区域约束的突破［J］.湖南文理学院学报，2008，20(3)：9-11.

［3］郑阿奇.MATLAB实用教程［M］.第二版.北京：电子工业出版社，2008.