数值分析

本节介绍的是关于函数的函数(function functions),这些函数是用来处理函数而非数值的;

绘图

优化

求解

MatLab中的函数表达

MatLab中用M文件来表示函数,设有如下函数:

他别表示为一称为hump.m的文件中:

function y = humps(x)
y = 1./((x – 0.3).^2 + 0.01) + 1./((x – 0.9).^2 + 0.04) – 6;

这个函数文件可用于数值分析的函数中.
第二种方法就是创造一个行内对象(inline()),方法如下:

f = inline(‘1./((x–0.3).^2 + 0.01) + 1./((x–0.9).^2 + 0.04)–6’);

用了上面的方法创造了函数文件,我们就可以找出函数在2的值:

f(2.0)
ans =
　　 –4.8552

用创造行内对象的方法还可以创造多参数的函数,如下:

f= inline('y*sin(x)+x*cos(y)','x','y')
f(pi,2*pi)
ans =
　　 3.1416

把函数画出来

fplot()可画出在给定范围内的函数值,如下

fplot('humps',[–5 5])
grid on

可通过限制y轴来放大图形

fplot('humps',[–5 5 –10 25])
grid on

你也可直接在fplot()中传递表达式,如:

fplot('2*sin(x+3)',[–1 1])

更可在一附图中画多个函数,如下

fplot('[2*sin(x+3), humps(x)]',[–1 1])

式中,[2*sin(x+3), humps(x)]组成了一个矩阵,每一列都是对应于x的函数

函数的最小值与解

找出一变量的函数的极值

x = fminbnd(’humps’,0.3,1)
x =
　　 0.6370

你可通过向fminbnd()函数传递一个函数optimset()作为参数来把此过程显示为列表形式:

x = fminbnd(’humps’,0.3,1,optimset(’Display’,’iter’))
Func-count 　x 　　　　　f(x) 　　　Procedure
1 　　　　　0.567376 　　12.9098 　　initial
2 　　　　　0.732624 　　13.7746 　　golden
3 　　　　　0.465248 　　25.1714 　　golden
4 　　　　　0.644416　　 11.2693 　　parabolic
5 　　　　　0.6413 　　　11.2583　　 parabolic
6　　　　　 0.637618 　　11.2529 　　parabolic
7　　　　　 0.636985 　　11.2528　　 parabolic
8 　　　　　0.637019　　 11.2528 　　parabolic
9 　　　　　0.637052 　　11.2528 　　parabolic
x =
　　0.6370

多变量函数极值

先创造一个ｍ文件，three_var.m:

function b = three_var(v)
x = v(1);
y = v(2);
z = v(3);
b = x.^2 + 2.5*sin(y) – z^2*x^2*y^2;

现在,以x = –0.6, y = –1.2,z = 0.135为起始点找出函数的极值:

v = [–0.6 –1.2 0.135];
a = fminsearch('three_var',v)
a =
　　 0.0000 –1.5708 0.1803

设置寻找极值的参数

x = fminbnd(fun,x1,x2,options)或
x = fminsearch(fun,x0,options)

其中，options是优化工具箱中(Optimization Toolbox)中的函数所用的一个结构,可如下设置

options = optimset('Display','iter');

options.Display用来设置是否显示中间过程,如为:"iter"则显示,为"off"则不显示,为"final"则只显示最后结果;
options.To1X设置结果的误差范围,默认值是:1.e–4.
options.MaxFunEval设置函数运行次数的上限,默认fminbnd()是500次,fminsearch()是200*length(x0)次

找出函数的解(零点值)

fzero()找出函数的零点值,你可以给出一个起始点,函数会从点开始搜索直到找到一个异号的值,最终给出解;
如果你知道函数会于哪两点异号,你可以给出一个两点的向量,表明起始值和起始搜索步长

a = fzero('humps',–0.2)
a =
　　–0.1316

验证一下，此函数值的确很接近０，

humps(a)
ans =
　　 8.8818e –16

再看看下面的命令，看看你是否能看懂！

humps(1)
ans =
　　16
humps(–1)
ans =
　　 –5.1378
options = optimset('Display','iter');
a = fzero('humps',[–1 1],options)
Func-count 　x 　　　　　　f(x) 　　　　　　Procedure
1 　　　　　–1 　　　　　–5.13779 　　　initial
1 　　　　　　1 　　　　　16 　　　　　　initial
2 　　　　　–0.513876 　–4.02235 　　　interpolation
3 　　　　　0.243062　　　71.6382　　　　bisection
4 　　　　　–0.473635 　–3.83767　　　　interpolation
5 　　　　　–0.115287 　　0.414441　　　 bisection
6 　　　　　–0.150214　　–0.423446　　　interpolation
7 　　　　　–0.132562　　–0.0226907 　　interpolation
8 　　　　　–0.131666　　–0.0011492 　　interpolation
9 　　　　　–0.131618　　1.88371e–07　　interpolation
10 　　　　–0.131618 　　–2.7935e–11　?interpolation
11 　　　　–0.131618 　　8.88178e–16 　　interpolation
12 　　　　–0.131618 　　–9.76996e–15 　interpolation
a =
　　 –0.1316

　　最优化问题常常需要多次的运算才能汇集到一点，因此，起始点的选择显得至关重要，它不仅能提高效率，更可使我们找到的值不是局部的极值，而是全局的最值；
　　下面是可能遇见的问题及解决的方法

数值积分

　　某区域内函数所围的区间可由数值积分来确定，MatLab中的一维函数积分用quad(),quad8().如:

q = quad(’humps’,0,1)
q =
　　 29.8583

　　函数可有第4个参数,指明误差范围,还可有第5个参数,用来画出图形.

例子:计算曲线的长度

如有一曲线如下确定:

x(t)=sin(2t);y(t)=cos(t);z(t)=t;

此函数可如下画出:

t = 0:0.1:3*pi;
plot3(sin(2*t),cos(t),t)

hcurce()用来计算曲线长度:

function f = hcurve(t)
f = sqrt(4*cos(2*t).^2 + sin(t).^2 + 1);

下面就是计算结果:

len = quad(’hcurve’,0,3*pi)
len =
　　 1.7222e+01

二重积分

result = dblquad(’integrnd’,xmin,xmax,ymin,ymax);

默认状态下,dblquad是用的quad()来计算,为更精确,可

result = dblquad(’integrnd’,xmin,xmax,ymin,ymax,[],’quad8’);