**Info**

This question is closed. Reopen it to edit or answer.

# Problem with quad2d integration

1 view (last 30 days)

Show older comments

Hi! I have a question on the quad2d function. When I run the code below ,it returns Q= 7.1754e-045. But it must be approximately 4.7399 e-043 . There is very little difference but it is very important for my project. When I decrease the value of 'AbsTol' ( for example 1e-85), it returns 4.867 e-043 but the code runs very slowly. I need another way to be fast and approximately 4.7399 e-043. Possible ?

fun = @(b,a) (b ./a).^19.* exp(-(439 ./ a).^b) .* ((439 ./ a).^(b-1)) ... .* 1./a .* 1./69.5391;

Q = quad2d(fun,6.00018,75.5392,248.063,1573.73 , 'AbsTol',1e-15,'Singular',false);

##### 1 Comment

Geoff Hayes
on 12 Dec 2014

### Answers (1)

Mike Hosea
on 13 Dec 2014

Is there a typo somewhere above? Here's a script I called "doit.m":

format short g

fun = @(b,a) (b ./a).^19.* exp(-(439 ./ a).^b) .* ((439 ./ a).^(b-1)) .* 1./a .* 1./69.5391;

tic;

Qquad2dNS = quad2d(fun,6.00018,75.5392,248.063,1573.73,'RelTol',100*eps,'AbsTol',1e-50,'Singular',false,'MaxFunEvals',10000)

toc

tic;

Qquad2dS = quad2d(fun,6.00018,75.5392,248.063,1573.73,'RelTol',100*eps,'AbsTol',1e-50,'Singular',true,'MaxFunEvals',10000)

toc

tic;

QIterated = integral2(fun,6.00018,75.5392,248.063,1573.73,'RelTol',100*eps,'AbsTol',1e-50,'method','iterated')

toc

tic;

QTiled = integral2(fun,6.00018,75.5392,248.063,1573.73,'RelTol',100*eps,'AbsTol',1e-50,'method','tiled')

toc

% Break the region up into n^2 pieces and integrate over each piece.

n = 30;

x = linspace(6.00018,75.5392,n);

y = linspace(248.063,1573.73,n);

Q = 0;

for i = 1:numel(x)-1

for j = 1:numel(y)-1

Q = Q + integral2(fun,x(i),x(i+1),y(j),y(j+1),'RelTol',100*eps,'AbsTol',1e-50,'method','iterated');

end

end

QMIterated = Q

Q = 0;

for i = 1:numel(x)-1

for j = 1:numel(y)-1

Q = Q + integral2(fun,x(i),x(i+1),y(j),y(j+1),'RelTol',100*eps,'AbsTol',1e-50,'method','tiled');

end

end

QMTiled = Q

Answers = [QIterated,QTiled,QMIterated,QMTiled,Qquad2dNS,Qquad2dS]

Max = max(Answers);

Min = min(Answers);

ULPs = (Max - Min)/eps(Max)

>> doit

Qquad2dNS =

2.0639e-18

Elapsed time is 0.763301 seconds.

Qquad2dS =

2.0639e-18

Elapsed time is 1.608291 seconds.

QIterated =

2.0639e-18

Elapsed time is 0.519476 seconds.

QTiled =

2.0639e-18

Elapsed time is 1.408459 seconds.

QMIterated =

2.0639e-18

QMTiled =

2.0639e-18

Answers =

2.0639e-18 2.0639e-18 2.0639e-18 2.0639e-18 2.0639e-18 2.0639e-18

ULPs =

46

There are only 46 units of roundoff difference between the max and min for all these computed results, and the relative tolerance of 100*eps is only trying to get within 100 of the true solution. If that's not mathematically correct, then probably the numerical function isn't being evaluated accurately, which wouldn't surprise me given the large exponents.

##### 0 Comments

### See Also

### Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!