21.510.50-0.5-1-1.5-2-2-1.5-1-0.500.511.52图4.6-3
(2)
[x0,y0]=ginput(2); disp([x0,y0])
-0.7926 -0.7843 0.7926 0.7843
(3)
fun='[sin(x(1)-x(2)),cos(x(1)+x(2))]'; %<12> [xy,f,exit]=fsolve(fun,[x0(2),y0(2)]) %<13> Optimization terminated successfully:
First-order optimality less than OPTIONS.TolFun, and no negative/zero curvature detected xy =
0.7854 0.7854 f =
1.0e-006 *
-0.0984 0.2011 exit = 1
〖说明〗
[fun.m]
function ff=fun(x) ff(1)=sin(x(1)-x(2)); ff(2)=cos(x(1)+x(2));
4.7 函数极值点
4.7.1 一元函数的极小值点 4.7.2 多元函数的极小值点
【例4.7.2-1】求f(x,y)?100(y?x)?(1?x)的极小值点。它即是著名的Rosenbrock's \测试函数。该测试函数有一片浅谷,许多算法难以越过此谷。(演示本例搜索过程的文件名为exm04072_1_1.m 。) (1)
ff=inline('100*(x(2)-x(1)^2)^2+(1-x(1))^2','x');
222(2)
x0=[-1.2,1];[sx,sfval,sexit,soutput]=fminsearch(ff,x0) sx =
1.0000 1.0000
11
sfval =
8.1777e-010 sexit = 1
soutput =
iterations: 85 funcCount: 159
algorithm: 'Nelder-Mead simplex direct search'
(3)
[ux,sfval,uexit,uoutput,grid,hess]=fminunc(ff,x0)
Warning: Gradient must be provided for trust-region method; using line-search method instead.
> In D:\\MATLAB6P1\\toolbox\\optim\\fminunc.m at line 211
Optimization terminated successfully:
Current search direction is a descent direction, and magnitude of directional derivative in search direction less than 2*options.TolFun ux =
1.0000 1.0000 sfval =
1.9116e-011 uexit = 1
uoutput =
iterations: 26 funcCount: 162 stepsize: 1.2992
firstorderopt: 5.0020e-004
algorithm: 'medium-scale: Quasi-Newton line search' grid =
1.0e-003 * -0.5002 -0.1888 hess =
820.4028 -409.5496 -409.5496 204.7720
4.8 数值积分
4.8.1 一元函数的数值积分 4.8.1.1 闭型数值积分
?x? 。 【例 4.8.1.1-1】求I??edx,其精确值为0.74684204012(1)
syms x;IS=int('exp(-x*x)','x',0,1)
vpa(IS) IS =
1/2*erf(1)*pi^(1/2) ans =
.74682413281242702539946743613185
(2)
fun=inline('exp(-x.*x)','x');
Isim=quad(fun,0,1),IL=quadl(fun,0,1)
12
Isim =
0.7468 IL =
0.7468
(3)
Ig=gauss10(fun,0,1) Ig =
0.7463
(4)
xx=0:0.1:1.5;ff=exp(-xx.^2); pp=spline(xx,ff); int_pp=fnint(pp);
Ssp=ppval(int_pp,[0,1])*[-1;1] Ssp =
0.7468
(5)
图4.8-1
4.8.1.2 开型数值积分
[gauss10.m]
function g = gauss10(fun,a,b) %GAUSS10(fun,a,b) % fun
%====================================================== x = [0.1488743390;0.4333953941;0.6974095683;... 0.8650633667;0.9739065285];
w = [0.2955242247;0.2692667193;0.2190863625;... 0.1494513492;0.0666713443];
t = .5*(b+a)+.5*(b-a)*[-flipud(x);x]; W = [flipud(w);w];
g = sum(W.*feval(fun,t))*(b-a)/2;
【例 4.8.1.2-1】当f(x)?cos(x)时,比较解析积分和近似积分。 (1)
syms x;F=int('cos(x)','x',-1,1),vpa(F) F =
2*sin(1) ans =
1.6829419696157930133050046432606 (2)
aF=cos(1/sqrt(3))+cos(-1/sqrt(3)) aF =
1.6758
【例4.8.1.2-2】求I(1)
??101?lndx,准确结果是?.88622692? 。 x2 13
syms x;IS=vpa(int('sqrt(log(1/x))','x',0,1)) Warning: Explicit integral could not be found.
> In D:\\MATLAB6P5\\toolbox\\symbolic\\@sym\\int.m at line 58 In D:\\MATLAB6P5\\toolbox\\symbolic\\@char\\int.m at line 9 IS =
.88622692545275801364908374167057
(2)用quad指令求积
ff=inline('sqrt(log(1./x))','x');Isim=quad(ff,0,1) Warning: Divide by zero.
> In D:\\MATLAB6P5\\toolbox\\matlab\\funfun\\inlineeval.m at line 13 In D:\\MATLAB6P5\\toolbox\\matlab\\funfun\\@inline\\feval.m at line 34 In D:\\MATLAB6P5\\toolbox\\matlab\\funfun\\quad.m at line 62 Isim =
0.8862
(3)
Ig=gauss10(ff,0,1) Ig =
0.8861
4.8.2 多重数值积分
4.8.2.1 积分限为常数的二重积分指令
【例 4.8.2.1-1】计算Sx01(1)
syms x y
ssx01=vpa(int(int(x^y,x,0,1),y,1,2)) ssx12=vpa(int(int(x^y,y,0,1),x,1,2))
Warning: Explicit integral could not be found.
> In D:\\MATLAB6P5\\toolbox\\symbolic\\@sym\\int.m at line 58 ssx01 =
.40546510810816438197801311546435 ssx12 =
1.2292741343616127837489278679215
21y1?dy和S?2xydx?dy。 ???xdx?x12?0???1?01??????(2)
zz=inline('x.^y','x','y');
nsx01=dblquad(zz,0,1,1,2) nsx12=dblquad(zz,1,2,0,1) nsx01 = 0.4055 nsx12 =
1.2293
4.8.2.2 内积分限为函数的二重积分
[double_int.m]
function SS=double_int(fun,innlow,innhi,outlow,outhi) %double_int % %fun %innlow,innhi
14
%outlow,outhi
y1=outlow;y2=outhi;x1=innlow;x2=innhi;f_p=fun; SS=quad(@G_yi,y1,y2,[],[],x1,x2,f_p);
[G_yi.m]
function f=G_yi(y,x1,x2,f_p) %G_yi %y %x1,x2 % %f_p
y=y(:);n=length(y);
if isnumeric(x1)==1;xx1=x1*ones(size(y));else xx1=feval(x1,y);end if isnumeric(x2)==1;xx2=x2*ones(size(y));else xx2=feval(x2,y);end for i=1:n
f(i)=quad(f_p,xx1(i),xx2(i),[],[],y(i)); end f=f(:);
【例 4.8.2.2-1】计算I42????(x2?y2)dx?dy。 1???y?(1)
[x_low.m]
function f=x_low(y) f=sqrt(y);
(2)
(3)
ff=inline('x.^2+y.^2','x','y'); SS=double_int(ff,@x_low,2,1,4)
Warning: Minimum step size reached; singularity possible. > In D:\\MATLAB6P5\\toolbox\\matlab\\funfun\\quad.m at line 88 In D:\\MATLAB6p5\\work\\G_yi.m at line 11
In D:\\MATLAB6P5\\toolbox\\matlab\\funfun\\@inline\\feval.m at line 20 In D:\\MATLAB6P5\\toolbox\\matlab\\funfun\\quad.m at line 62 In D:\\MATLAB6p5\\work\\double_int.m at line 8 SS =
9.5810
(4)
Ssym=vpa(int(int('x^2+y^2','x','sqrt(y)',2),'y',1,4)) Ssym =
9.5809523809523809523809523809524
4.8.3 卷积
4.8.3.1 “完整”离散序列的数值卷积 4.8.3.2 “截尾”离散序列的数值卷积
15