www.pudn.com > NumericalComputingwithMatlabCode.zip > quadgui.m, change:2003-03-24,size:4277b

```function [Q,fcount] = quadgui(F,a,b,tol,varargin)
%QUADGUI  Demonstrate numerical evaluation of a definite integral.
%   Q = QUADGUI(F,A,B) shows the steps in approximating the integral
%   of F(x) from A to B by adaptive extrapolated Simpson's quadrature.
%
%   The shaded area shows the integral over the current subinterval.
%   The color switches to green when the desired accuracy is obtained.
%
%
%   F is a string defining a function of a single variable, an inline
%   function, a function handle, or a symbolic expression involving a
%   single variable.  Arguments to QUADGUI beyond the first four,
%   Q = QUADGUI(F,a,b,tol,p1,p2,...), are passed on to the integrand,
%   F(x,p1,p2,..).
%
%   [Q,fcount] = QUADGUI(F,...) also counts the number of evaluations
%   of F(x).
%
%   Examples:
%      betafun = inline('t^p*(1-t)^q','t','p','q')
%

shg
clf reset

% Default tolerance
if nargin < 4 | isempty(tol)
tol = 1.e-4;
end

% Default function and interval.
if nargin < 3
F = @humps;
a = 0;
b = 1;
end

% Make F callable by feval.
if ischar(F) & exist(F)~=2
F = inline(F);
elseif isa(F,'sym')
F = inline(char(F));
end

% Initialization
c = (a + b)/2;
fa = feval(F,a,varargin{:});
fc = feval(F,c,varargin{:});
fb = feval(F,b,varargin{:});

% Scale the plot
h = b - a;
x = [a c b];
y = [fa fc fb];
maxy = max(y);
miny = min(y);
for k = 1:63
v = feval(F,a+k*h/64,varargin{:});
maxy = real(max(maxy,v));
miny = real(min(miny,v));
end
set(gcf,'userdata',0)
plot(x,y,'.','markersize',6);
hold on
p(1) = fill(a,fa,'k');
p(2) = fill(b,fb,'k');
hold off
s = (maxy - miny)/20;
axis([a b miny-s maxy+s])
q(1) = uicontrol('string','step', ...
'units','normal','pos',[.65 .02 .08 .06], ...
'callback','set(gcf,''userdata'',1)');
q(2) = uicontrol('string','auto', ...
'units','normal','pos',[.75 .02 .08 .06], ...
'callback','set(gcf,''userdata'',2)');
q(3) = uicontrol('string','quit', ...
'units','normal','pos',[.85 .02 .08 .06], ...
'callback','set(gcf,''userdata'',3)');

% Recursive call
[Q,k] = quadguistep(F, a, b, tol, fa, fc, fb, varargin{:});
fcount = k + 3;

% Finish
title(sprintf('Q = %8.4f, fcount = %4.0f',Q,fcount))
delete(p);
delete(q(1:2));
set(q(3),'string','close','callback','close(gcf)')

% ---------------------------------------------------------

% Recursive subfunction used by quadtx.

h = b - a;
c = (a + b)/2;
d = (a + c)/2;
e = (c + b)/2;
fd = feval(F,d,varargin{:});
fe = feval(F,e,varargin{:});
Q1 = h/6 * (fa + 4*fc + fb);
Q2 = h/12 * (fa + 4*fd + 2*fc + 4*fe + fb);

u1 = a:h/64:c;
v1 = polyinterp([a d c],[fa fd fc],u1);
u1 = [a u1 c];
v1 = [0 v1 0];
u2 = c:h/64:b;
v2 = polyinterp([c e b],[fc fe fb],u2);
u2 = [c u2 b];
v2 = [0 v2 0];
if (abs(Q2 - Q1) <= tol)
color = [0 2/3 0];
else
color = [.6 .6 .6];
end
p = flipud(get(gca,'child'));
x = [get(p(1),'xdata') d e];
y = [get(p(1),'ydata') fd fe];
set(p(1),'xdata',x,'ydata',y)
set(p(2),'xdata',u1,'ydata',v1,'facecolor',color)
set(p(3),'xdata',u2,'ydata',v2,'facecolor',color)
set(gca,'xtick',sort(x),'xticklabel',[]);
title(num2str(length(x)))
pause(.25)
while get(gcf,'userdata') == 0
pause(.25)
end
if get(gcf,'userdata') == 1
set(gcf,'userdata',0)
end

if (abs(Q2 - Q1) <= tol) | (get(gcf,'userdata') == 3)
Q  = Q2 + (Q2 - Q1)/15;
fcount = 2;
else
[Qa,ka] = quadguistep(F, a, c, tol, fa, fd, fc, varargin{:});
[Qb,kb] = quadguistep(F, c, b, tol, fc, fe, fb, varargin{:});
Q  = Qa + Qb;
fcount = ka + kb + 2;
end
```