Simpsons rule integration from zero to infinity

I am new to programming and I am trying to calculate with the simpson's 1/3 rule the integral : x*exp(-x^2)(from zero to inf) / exp(-x^2 (from zero to inf). The code that I used is the following:

function fval = myFunInt(x)
fval = x*exp(-x^2);
end
function gval = mygInt(x)
gval = exp(-x^2);
end
h = (b-a)/2;
I_simp = h/3*(myFunInt(a) + 4*myFunInt(a+h) + myFunInt(a+2*h))

disp(I_simp)

h = (b - a);
ga = mygInt(a);
gb = mygInt(b);

%%Simpson's 1/3 Rule
h = (b-a)/2;
I_simp2 = h/3*(mygInt(a) + 4*mygInt(a+h) + mygInt(a+2*h))

disp(I_simp2)
final = I_simp/I_simp2

The result that I get is not a number(NaN). How should I approach this?

Answer 1

[John D'Errico]

Integration of a function using Simpson's rule all the way to infinity will take a LONG time. The last time I checked, infinity was far off.

But if you put in your upper limit there are explicitly +inf, what do you really expect to see? How many steps fall between 0 and inf?

You need to consider if you REALLY need to go all the way to infinity. Perhaps stopping a bit short might be viable. Consider what the value of the kernel here is at some reasonably large number. At what point will you see an underflow? Is there any reason to go beyond the point where the kernel is zero? Think about what you are doing, and these problems become easy enough.