Quad 2d giving inconsistant answers

Hello,

I'm trying to use quad2d to integrate a function whose value is 0 outside of the (-5,5)(-5,5) box. However, I get different answers if I tell quad to integrate that box specifically or a larger box around it. Any idea what's happening? I'm using R2009a. (Sorry for the inefficent Matlab coding, it's partly because I hacked the code to work with quad2d)

quad2d(@(x,y)(of_Nd_ExpHimmelblau(x,y)),-5,5,-5,5,'AbsTol',10^-14,'RelTol',10^-14)

ans =

38.008891501144596

quad2d(@(x,y)(of_Nd_ExpHimmelblau(x,y)),-5,7,-5,11,'AbsTol',10^-14,'RelTol',10^-14)

ans =

38.008544462736467

function outMat=of_Nd_ExpHimmelblau(x1,y1)

[N,D]=size(x1);

x=zeros(N*D,2);

x(:,1)=x1(:);

x(:,2)=y1(:);

beta=0.01;

[numSamples,numDim]=size(x);

xmin=-5;

xmax=5;

out=zeros(numSamples,1);

for jj=1:numSamples

    if sum(x(jj,:)<xmin)>0 || sum(x(jj,:)>xmax)>0

        out(jj)=0;

    else

        for ii=1:(numDim-1)

            xx=x(jj,ii);

            yy=x(jj,ii+1);

            out(jj)=out(jj)+exp(-((xx^2+yy-11)^2+(xx+yy^2-7)^2)*beta);

        end

    end

end

outMat=zeros(N,D);

outMat(:)=out(:);

Answer 1

[Mike Hosea]

Did you turn the warnings off? Your second integration failed:

>> quad2d(@(x,y)(of_Nd_ExpHimmelblau(x,y)),-5,5,-5,5,'AbsTol',10^-14,'RelTol',10^-14)
Warning: RelTol was increased to 100*eps('double') = 2.22045e-14. 
> In quad2d at 173 
ans =
  38.008891501144596
>> quad2d(@(x,y)(of_Nd_ExpHimmelblau(x,y)),-5,7,-5,11,'AbsTol',10^-14,'RelTol',10^-14)
Warning: RelTol was increased to 100*eps('double') = 2.22045e-14. 
> In quad2d at 173 
Warning: Reached the maximum number of function evaluations (2000). The result fails the global error test. 
> In quad2d at 241 
ans =
  38.008544462736474

In general you should use the limits (function handles if need be) to avoid integrating over discontinuities. Setting the integrand to zero outside of the region of interest was a bad habit encouraged of necessity for DBLQUAD, but it is a horrible thing to do, sending accuracy down and run-times up by 10x. If you have sharp edges or discontinuities, especially if you require high accuracy, you really need to break the itegration up into pieces, putting the nasty bits on the boundaries. Of course in this case it is not really necessary since the integrand is zero on the other side of the discontinuity. The reason QUAD2D has more trouble with discontinuities than DBLQUAD, generally speaking, is because it refines the mesh in 2D instead of in one dimension at a time. A curvilinear discontinuity will attract an unnecessarily fine mesh along its length as well as "across" it. -- Mike

Answer 2

[Brendan]

Hi,

Thanks, I had no idea that my warning is turned off. I'm not really sure how I did that given I haven't messed with the 'warning' command since starting the program up. That certainly explains the inconsistency.

Thanks