%% Limitations of quad2d
% quad2d begins by mapping the region of integration to a rectangle. Consequently, it may have trouble integrating over a region that does not have four sides or has a side that cannot be mapped smoothly to a straight line. If the integration is unsuccessful, some helpful tactics are leaving Singular set to its default value of true, changing between Cartesian and polar coordinates, or breaking the region of integration into pieces and adding the results of integration over the pieces.
% For instance:

fun = @(x,y)abs(x.^2 + y.^2 - 0.25);
c = @(x)-sqrt(1 - x.^2);
d = @(x)sqrt(1 - x.^2);
quad2d(fun,-1,1,c,d,'AbsTol',1e-8,...
    'FailurePlot',true,'Singular',false);

% The failure plot shows two areas of difficulty, near the points (-1,0) and (1,0) and near the circle x^2 + y^2 = 0.25.
%
% Changing the value of Singular to true will cope with the geometric singularities at (-1,0) and (1,0). The larger shaded areas may need refinement but are probably not areas of difficulty.

Q = quad2d(fun,-1,1,c,d,'AbsTol',1e-8, ... 
     'FailurePlot',true,'Singular',true);

%
% From here you can take advantage of symmetry:

Q = 4*quad2d(fun,0,1,0,d,'Abstol',1e-8,...
     'Singular',true,'FailurePlot',true)

%
% However, the code is still working very hard near the singularity. It may not be able to provide higher accuracy:

Q = 4*quad2d(fun,0,1,0,d,'Abstol',1e-10,...
     'Singular',true,'FailurePlot',true);

%
% At higher accuracy, a change in coordinates may work better.

polarfun = @(theta,r) fun(r.*cos(theta),r.*sin(theta)).*r;
Q = 4*quad2d(polarfun,0,pi/2,0,1,'AbsTol',1e-10);

%
% It is best to put the singularity on the boundary by splitting the region of integration into two parts:

Q1 = 4*quad2d(polarfun,0,pi/2,0,0.5,'AbsTol',5e-11);
Q2 = 4*quad2d(polarfun,0,pi/2,0.5,1,'AbsTol',5e-11);
Q = Q1 + Q2;

% Copyright 2012 The MathWorks, Inc.