第7讲：非数与数值微积分ppt课件

1、第七讲第七讲 非数与非数与 数值微积分数值微积分2n本讲教学目标本讲教学目标n 掌握非数的概念和性质掌握非数的概念和性质n 掌握空数组的概念和创建方法掌握空数组的概念和创建方法n 掌握数值微分、梯度和法线掌握数值微分、梯度和法线n 掌握常见的数值积分方法掌握常见的数值积分方法n 了解二重积分的方法了解二重积分的方法3n7.1 “非数和非数和“空数组空数组”n在进行科学计算时经常遇到除零的现象，在进行科学计算时经常遇到除零的现象，这时就会出现错误信息。这时就会出现错误信息。nMATLAB中规定出现上述情况不出致命错中规定出现上述情况不出致命错误信息，仅在输出中进行警告并指出这是误信息，仅在输出中

2、进行警告并指出这是非数。非数。n我们可以巧妙地利用非数对图形进行镂空、我们可以巧妙地利用非数对图形进行镂空、裁减等操作，从而得到很好的可视化效果。裁减等操作，从而得到很好的可视化效果。4n1.1.1 “非数非数”nMATLAB 规定，形如规定，形如n 都会产生非数都会产生非数(Not a Number)，用，用NaN或或nan 表示。表示。n根据根据 IEEE 数学规范，数学规范，NaN具有以下性质：具有以下性质：n只要只要NaN 参加四则运算或函数作用于非数参加四则运算或函数作用于非数其结果也一定是非数，即具有传递性。其结果也一定是非数，即具有传递性。n非数没有大小的概念，因此不能进行逻辑非

3、数没有大小的概念，因此不能进行逻辑运算也不能比较两个非数的大小。运算也不能比较两个非数的大小。0,005n例例1 非数的产生和性质演示。非数的产生和性质演示。n a=0/0,b=0*log(0),c=inf/inf ,d=sqrt(NaN)n Warning: Divide by zero.n In D:ZJJpaperch23.m at line 62n a = NaN; b = NaN; c = NaN; d = NaNn1非数的传递性非数的传递性n 0*a,sin(a) n ans = NaNn ans = NaN6n2非数的不可比较性非数的不可比较性n a=nan % 非数不能比较非数

4、不能比较n ans = 0 n a=nan % 非数不能进行逻辑运算非数不能进行逻辑运算n bc % 两个非数不能比较大小两个非数不能比较大小 n ans = 1n ans = 0n3非数的属性判断非数的属性判断n class(a)% 数据类型归属数据类型归属n isnan(a)% 唯一正确判断非数的指令唯一正确判断非数的指令 n ans = doublen ans = 1 n从结果可知非数也从结果可知非数也是一个数，数据类是一个数，数据类型是双精度的。型是双精度的。7n7.1.2 “空数组空数组n当数组的某些维数为当数组的某些维数为0时，称该数组为空数时，称该数组为空数组。组。n在复杂的计算

5、中有时会产生空数组，也可以在复杂的计算中有时会产生空数组，也可以利用空数组对其他数组进行赋值。利用空数组对其他数组进行赋值。n空数组与全零元素数组不是一回事空数组与全零元素数组不是一回事n空数组确实是数组，用空数组确实是数组，用which、who、whos可以看出。可以看出。n可以用命令可以用命令 isempty 来判断一个数组是否来判断一个数组是否为空数组。为空数组。8n例例2 创建空数组。创建空数组。n A = rand(2,2,2); % 直接赋值法直接赋值法n A(:,:,:) = ;n ndims(A)n size(A)n ans = 3n ans = 0 2 2n可看到可看到 A

6、仍然是仍然是 3 维数组，但其中有一个维数组，但其中有一个维数是维数是 0，则是空数组。，则是空数组。n(2) 用函数创建空数组用函数创建空数组n A = ones(25,0,4); n if An S1=是数组是数组n elsen S0=是空数组是空数组n endnS0 = 是空数组是空数组9n例例3 “空数组的属性。空数组的属性。n class(A) % 显示数据类别显示数据类别n isnumeric(A) % 是数值数组类吗是数值数组类吗n isempty(A) % 是否为空数组？是否为空数组？n size(A) n ans =doublen ans = 1n ans = 1 n ans

7、 = 25 0 4n从结果看，数组从结果看，数组A 是一是一个双精度的空数组，具个双精度的空数组，具有维数有维数2504。10n7.2 数值微积分数值微积分n7.2.1 数值微分数值微分n在在MATLAB中，没有直接提供求数值导数中，没有直接提供求数值导数的函数，只有计算向前差分的函数的函数，只有计算向前差分的函数diff，其，其格式为：格式为：nDX=diff(X)：计算向量：计算向量X的向前差分的向前差分n DX(i)=X(i+1)-X(i)，i=1,2,n-1。nDX=diff(X,n)：计算：计算X的的n阶向前差分阶向前差分n diff(X,2)=diff(diff(X)。nDX=di

8、ff(A,n,dim)：计算矩阵：计算矩阵A的的n阶差分，阶差分，dim=1时时(缺省缺省)，按列计算；，按列计算；dim=2，按行，按行计算。计算。11n例例4 对我国对我国5298年国民经济数据年国民经济数据GDP发展趋势发展趋势进行分析，求其一、二阶差分，并作图。进行分析，求其一、二阶差分，并作图。n A=load(e:datagdp.txt)n B=A(:,2)n B1=diff(B);n B2=diff(B,2)n subplot(1,3,1),plot(B),title(发展趋势发展趋势)n subplot(1,3,2),plot(B1),title(发展速度发展速度)n subp

9、lot(1,3,3),plot(B2),title(发展加速度发展加速度)n差分方法应用十分广差分方法应用十分广泛，如在经济数据中泛，如在经济数据中一阶差分表示发展速一阶差分表示发展速度，二阶差分表示经度，二阶差分表示经济发展的加速度。济发展的加速度。12n7.2.2 梯度梯度n在在MATLAB中，使用中，使用gradient函数求解近函数求解近似梯度，其格式为：似梯度，其格式为：nfx, fy = gradient(f) 返回矩阵返回矩阵 f 的数的数值梯度，值梯度，fx相当于相当于df/dx，即在，即在x方向列方向列的差分值。的差分值。nfx, fy = gradient(f, h) 求二

10、元函数的求二元函数的梯度，使用数量梯度，使用数量 h 作为各个方向的间隔点。作为各个方向的间隔点。nfx, fy, fz, = gradient(f, h1, h2, ) 做相似工作。做相似工作。13n例例5 对函数对函数 作梯度图。作梯度图。 n v = -2:0.2:2; x,y = meshgrid(v);n z = x .* exp(-(x.2 + y.2);n subplot(1,2,1);surf(x,y,z),n subplot(1,2,2); px,py = gradient(z,.2,.2);n contour(x,y,z), hold on, quiver(x,y,px,p

11、y)22xyxe14n7.2.3 法线法线n计算与绘制法线的命令为：计算与绘制法线的命令为：n NX, NY, NZ = surfnorm( X, Y, Z )n方向法线绘图命令的基本格式为：方向法线绘图命令的基本格式为：nquiver3(Z,U,V,W)nquiver3(X,Y,Z,U,V,W)nquiver3(.,scale)n其中：其中：X,Y,Z 为曲线上的点；为曲线上的点；U,V,W 为每一为每一点上切面的法线方向；点上切面的法线方向；scale 为实数值，表为实数值，表示方向箭头的指向。示方向箭头的指向。15n例例6 求函数求函数 的曲面图，且作出每点的切的曲面图，且作出每点的切面

12、法向量图形。面法向量图形。nX,Y = meshgrid(-2:0.25:2,-1:0.2:1);nZ = X.* exp(-X.2 - Y.2);nU,V,W = surfnorm(X,Y,Z); % 计算法向量分计算法向量分量量nquiver3(X,Y,Z,U,V,W,0.5); % 作法向量的图作法向量的图形形nhold onnsurf(X,Y,Z);colormap hsvnview(-35,45)naxis (-2 2 -1 1 -.6 .6)nhold off22xyxe16n7.2.4 数值积分数值积分n求解定积分的数值方法多种多样，如下都求解定积分的数值方法多种多样，如下都是经

13、常采用的方法。是经常采用的方法。n简单的梯形法简单的梯形法n辛普生辛普生(Simpson) 法法n牛顿柯特斯牛顿柯特斯(Newton-Cotes)法法n其基本思想都是将积分区间其基本思想都是将积分区间a,b分成分成n个子个子区间区间xi,xi+1，i=1,2,n，其中，其中x1=a，xn+1=b。这样求定积分问题就分解为求和。这样求定积分问题就分解为求和问题。问题。17n1. 变步长辛普生法变步长辛普生法n调用格式为：调用格式为：n q = quad(fun,a,b)n q = quad(fun,a,b,tol)n q,n=quad(fun,a,b,tol,trace)n其中其中fun为被积函

14、数；为被积函数；a, b是积分上下限，是积分上下限，tol为给定的误差限，缺省取为给定的误差限，缺省取0.001；trace若取若取0则不展现积分过程，缺省取则不展现积分过程，缺省取0；取非；取非0则展现。返回参数则展现。返回参数q即定积分值，即定积分值，n为被积为被积函数的调用次数。函数的调用次数。18n例例7 计算积分计算积分n Q = quad(1./(1+25*x.2),-1,1) n ans = 0.5494n% 也可以使用内联函数定义被积函数也可以使用内联函数定义被积函数n F = inline(1./(1+25*x.2);n Q = quad(F,-1,1); n ans = 0

15、.5494 1122511dxx19n例例8 对空间曲线求积分。对空间曲线求积分。n设空间曲线的参数方程为：设空间曲线的参数方程为：n x = t.*sin(t), y = t.*cos(t), z=t n t = 0:pi/50:10*pi;n plot3(t.*sin(t),t.*cos(t),t,LineWidth,2)n f = inline(sqrt(4*(t.*sin(t).2 + (t.*cos(t).2 + 1);n len = quad(f,0,pi)n计算结果为：计算结果为：n len = 8.519120n2. 牛顿柯特斯法牛顿柯特斯法n调用格式为：调用格式为：n I,n

16、=quad8(fun,a,b,tol,trace)n其中参数的含义和其中参数的含义和quad函数相似，只是函数相似，只是tol的缺省值取的缺省值取10-6。 n该函数可以更精确地求出定积分的值，且该函数可以更精确地求出定积分的值，且一般情况下函数调用的步数明显小于一般情况下函数调用的步数明显小于quad函数，从而保证能以更高的效率求出所需函数，从而保证能以更高的效率求出所需的定积分值。的定积分值。21n例例9 分别用函数分别用函数quad和和quad8求定积分的近似值，求定积分的近似值，并在相同的精度下比较函数调用次数。并在相同的精度下比较函数调用次数。n format long; %调用函数

17、调用函数quad求定积分求定积分n fx=inline(exp(-x);n I,n=quad(fx,1,2.5,1e-10)n I =n 0.28579444254766n n =n 65n format long; %调用函数调用函数quad8n fx=inline(exp(-x);n I,n=quad8(fx,1,2.5,1e-10)n I =n 0.28579444254754n n =n 3322n3. 被积函数由一个表格定义被积函数由一个表格定义n用用trapz(X,Y)函数求解，其中向量函数求解，其中向量X,Y定义定义函数关系函数关系Y=f(X)。n例例10 用用trapz函数计算定积分函数计算定积分 。n X=1:0.01:2.5;n Y=exp(-X); % 生成函数关系数据向生成函数关系数据向量量n trapz(X,Y)n ans =n 0.285796824163935 . 21dxex23n4. 二重积分二重积分n与