Evaluating integral numerically using quadgk

Hi

I am trying to evaluate the following integral numerically in MatLAB: http://www.scribd.com/doc/100400549/mwe

However, the sqrt(x^2 + y^2) in the numerator has changed to x. Here is my attempt (based on this thread: http://mathworks.com/matlabcentral/answers/43929-evaluate-advanced-integral-numerically#answer_53917):

clc
clear all

myquad = @(fun,a,b,tol,trace,varargin)quadgk(@(x)fun(x,varargin{:}),a,b,'AbsTol',tol);

sigma = 1;
kappa = 1e4;
B = 1e6;
beta = 1.0;
f = 1e6;

fct  = @(x,y,z,f) (...
    (x)/sqrt(x.^2+y.^2+z.^2)).*(exp((-x.^2-y.^2-z.^2)/(2*sigma^2)))./ ((f - B*beta*sqrt(x.^2+y.^2+z.^2)).^2 + kappa^2/4 );

result = triplequad(@(x,y,z) fct(x,y,z,f), -Inf, Inf, -Inf, Inf, -Inf, Inf, 1e-16, myquad);

However this gives me a warning: "Rank deficient, rank = 0, tol = NaN.". I'm not sure what else to do here. Am I using triplequad correctly?

Best wishes, Niles.

[Niles Martinsen]

I noticed that the warning only enters when the limits are from -infinity to infinity. If I change them to 0..infinity, then the warning disappears. Unfortunately the integrand is odd, so I have to use -inf..inf.

[Teja Muppirala]

If that changes to an x, then the integral is zero by symmetry.

Answer 1

[Mike Hosea]

Make sure that you are using elementwise operations. See the first / in fct? It should be ./ instead. Check all the other multiplications and divisions. If it's only a scalar multiplication or division, then you don't need to change it, but in fact I find it easier just to use all .* and ./ rather than * or /.