# Decimal and very small values returning zeros

9 views (last 30 days)

Show older comments

I am running a integral calculation which is returning all zeros.

So i am running an integral over a meshgrid for space and time which i believe is working correctly.

My time scale is 0 to 1e-11 but when i use these values I get all zeros. If I reduce the it to 0 to 1e-9 I get numerical answers.

Is Matlab having rounding errors? If so how do I tell matlab to keep the sig figs?

t=0:2e-13:9e-11;

x=0:1:5;

[time space]=meshgrid(t,x);

fun=@(c) (pi.*time).^(-1/2).*exp(-(space-c).^2./(time)).*exp(-abs(space)/1);

answer=integral(fun,-inf,inf,'ArrayValued',true)

##### 1 Comment

James Tursa
on 19 Jul 2022

### Accepted Answer

Bruno Luong
on 19 Jul 2022

Edited: Bruno Luong
on 19 Jul 2022

Is is more subtle than I though.

You should not use 'ArrayValued'==true to evaluate a series of integral. The use case is to integrate a vector function depends on the domain, therefore MATLAB stop when the integrate vector converges.

- You should do just regular integral inside a loop, vectorize it is a very BAD idea;
- As your function is has the support so narrow (order of time which is tiny) and you integrate over the entire real axis, matlab indeed canot detect and miss where your function is ~= 0, and integral returns 0 for most of the time. You should help integral by narrowing the integration interval

t=0:2e-13:9e-11;

x=0:1:5;

[time, space]=meshgrid(t,x);

answer = zeros(size(time));

for i=1:numel(time)

ti = time(i);

si = space(i);

fun=@(c) (pi.*ti).^(-1/2).*exp(-(si-c).^2./(ti)).*exp(-abs(si)/1);

d = abs(log(realmin)); % 708.3964, theoretically exp(-x) for x >= d is zero numerically

lo = si-sqrt(d*ti);

hi = si+sqrt(d*ti); % fun(x) for x < lo or x > hi must be tiny

answer(i)=integral(fun,lo,hi);

end

answer

max(answer, [], 'all')

min(answer, [], 'all') % not 0 anymore

### More Answers (2)

Jan
on 19 Jul 2022

Edited: Jan
on 19 Jul 2022

No, Matlab uses the standard IEEE754 conventions and has no rounding errors.

The question is funny: A huge number of scientists use Matlab for decades to get reliable results. You can be sure, that such problem would have been found before.

How do you check, if the output is zero? Remember that the standard output to the command window has a limited size of decimal figures. So maybe all you need is:

format long g

See: format

x = [25 56.31156 255.52675 9876899999];

format short

x % Some zeros are shown, which are not zeros numerically:

format long g

x

##### 0 Comments

John D'Errico
on 19 Jul 2022

Edited: John D'Errico
on 19 Jul 2022

Um, why are you doing this numerically at all??????? With the variable of integration as c, your kernel is:

(pi*t)^(-1/2).*exp(-(x-c).^2./t).*exp(-abs(x))

with limits from -inf to inf. c enters in there only in one place. And this will be a known integral. I'll do it effectively by hand though.

syms t real positive

syms x real

syms c

K = (pi*t)^(-1/2).*exp(-(x-c).^2./t).*exp(-abs(x))

First, I'll transform the problem to make it a little more clear.

u = (x-c)/sqrt(t)

Then we will have

du = -dc/sqrt(t)

Since the limits of integration were -inf to inf, the only impact is to swap the limits, since we have x-c, but then we have an extra factor of -1 in there in the du term. So the net result is we will still integrate over the limits [-inf,inf].

Ku = 1/sqrt(pi) * exp(-u^2)*exp(-abs(x))

But what is the integral of exp(-u^2) over -inf,inf]? (Yes, I know I can do this on paper, but this is a MATLAB forum.)

syms u

int(exp(-u^2),u,[-inf,inf])

Do you see that cancels the 1/sqrt(pi)? And since x is not a function of c, it just comes out of the integral.

The result is, EVERYTHING COLLAPSES. The value of the integral you want to compute is just

exp(-abs(x))

It is not even apparently a function of t. Using a numerical integration to solve the problem is just wrong, since you then have serious numerical problems solving it. Even if the kernel were a little more complicated, I might at worst suggest the use of a Gauss-Hermite integration. As a check, see that Bruno computed a nmerical solution.

format long g

x = (0:5)';

[x,exp(-abs(x))]

Indeed, that is the same as what Bruno found. There is no dependence on t at all.

### See Also

### Community Treasure Hunt

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

Start Hunting!