# Numerical evaluation of integral gives warning message

1 view (last 30 days)

Show older comments

Dear all,

I am trying to compute the following integral in MATLAB:

with corresponding code is:

p1_minusSign = @(r) exp(-2*pi * integral(@(x) x/(1-x^4/r^4), r, Inf, 'ArrayValued', 1)) * exp(-pi*r^2)*2*pi*r;

p2_minusSign = integral(p1_minusSign, 0, Inf, 'ArrayValued', 1),

Unfortunately, MATLAB gives me a bunch of warning messages of the following form: "Warning: Minimum step size reached near x = 5833.72. There may be a singularity, or the tolerances may be too tight for this problem."

As far as I can tell, I don't see any potential singularity in my expression. In contrast, if we just change the sign in the demoninator of the inside integral, it computes perfectly:

with corresponding code:

p1_plusSign = @(r) exp(-2*pi * integral(@(x) x/(1+x^4/r^4), r, Inf, 'ArrayValued', 1)) * exp(-pi*r^2)*2*pi*r;

p2_plusSign = integral(p1_plusSign, 0, Inf, 'ArrayValued', 1),

Any ideas on what the error might be and how to correct it?

Many thanks in advance for your help.

##### 8 Comments

Torsten
on 27 Jun 2018

### Answers (2)

John D'Errico
on 27 Jun 2018

Edited: John D'Errico
on 27 Jun 2018

Ignoring the fact that you cannot just factor out -1 from that term in the integral, and get what you think you got, we have the added problem that your inner integral is not even finite.

syms x r

K = x/(1 - x^4/r^4)

K =

-x/(x^4/r^4 - 1)

Now, pick some arbitrary value for r. Lets say 1.

int(subs(K,r,1),[1,inf])

ans =

-Inf

But if you change the denominator in the simple way that is NOT mathematically correct, as you wrote, then it is of course finite, as expected.

K2 = x/(1 + x^4/r^4)

K2 =

x/(x^4/r^4 + 1)

int(subs(K2,r,1),[1,inf])

ans =

pi/8

There is a fundamental difference between those two kernels. So, K2 is a mountain that we can indeed climb. (Sorry about that. Could not resist it.) It would probably be a HW assignment in first year calc, to show the former case is not finite?

Anyway, the answer is your problem lacks a finite solution.

##### 9 Comments

John D'Errico
on 28 Jun 2018

Walter Roberson
on 28 Jun 2018

With the +1 the inner integral goes to

(1/4)*r^2*ln((2*r+1)/(2*r^2+2*r+1))

As r goes to infinity the ratio goes to 0, leading to ln(0) which is negative infinity. But you have exp(-2*pi*that) so you are getting exp(infinity)

Walter Roberson
on 28 Jun 2018

Maple says that the inner integral is infinite for non-negative r.

When you substitute that in to the outer expression, you get

int(exp(-Pi*r^2)*r*infinity, r, 0, infinity)

The exp(-Pi*r^2) term is nonnegative, and r is nonnegative, so we can see by inspection that the result will be infinity.

##### 0 Comments

### See Also

### Categories

### Products

### Community Treasure Hunt

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

Start Hunting!