9
$\begingroup$

I am trying to integrate

$$\int^1_0 t^{2n+2}\exp\left({\frac{\alpha r_0}{t}}\right)dt$$

which is a simple transformation of

$$\int^{\infty}_1 x^{2n}\exp(-\alpha r_0 x)dx$$

using $t = \frac1{x}$ because it is difficult to numerically approximate improper integrals. This does, however, lead to the problem of evaluating the new integrand near zero. It will be very easy to get the proper number of quadrature nodes seeing as the interval is only of length 1 (so the comparable $dt$ can be made very small), but what sort of considerations should I make when integrating near zero?

On some level, I think that simply taking $\int^1_\epsilon t^{2n+2}\exp({\frac{\alpha r_0} {t}})dt$ is a good idea where $\epsilon$ is some small number. However, what number should I choose? Should it be machine epsilon? Is division by machine epsilon a well quantified number? Furthermore, if division my machine epsilon (or close to it) gives an incredibly large number, then taking $\exp(\frac{1}{\epsilon})$ will become even larger.

How should I account for this? Is there a way to have a well defined numerical integral of this function? If not, what is the best way of integrating the function?

$\endgroup$
8
  • 1
    $\begingroup$ Have you considered using Monte Carlo? $\endgroup$ Dec 13, 2011 at 17:15
  • $\begingroup$ I feel like it would not fix the problem. Monte Carlo integration is often reserved for high dimensional integrals. I would run into the exact same problems with Monte Carlo, I would simply have less control over where my function is being evaluated. $\endgroup$
    – drjrm3
    Dec 13, 2011 at 17:16
  • $\begingroup$ You might be right. $\endgroup$ Dec 13, 2011 at 17:21
  • $\begingroup$ I think it would still be good to have an answer (perhaps to a separate, more general question) explaining how one does numerical integration when the function is divergent at one limit, for the general case where it's not possible to do the integral analytically. Then again, that could just as well be found in Numerical Recipes... $\endgroup$
    – David Z
    Dec 13, 2011 at 19:03
  • $\begingroup$ @Faheem: "Monte Carlo is an extremely bad method; it should be used only when all alternative methods are worse." - Alan Sokal $\endgroup$
    – J. M.
    Dec 14, 2011 at 7:31

2 Answers 2

10
$\begingroup$

This can be done by integration by parts: $$ \int^\infty_1 x e^{-ax} = \frac{-1}{a} x e^{-ax}\mid^\infty_1 - \frac{-1}{a} \int^\infty_1 e^{-ax} = \frac{e^{-a}}{a} + \frac{e^{-a}}{a^2} = \frac{a+1}{a^2} e^{-a} $$ and continuing on by induction $$ \int^\infty_1 x^k e^{-ax} = \frac{-1}{a} x^k e^{-ax}\mid^\infty_1 - \frac{-k}{a} \int^\infty_1 x^{k-1} e^{-ax} = \frac{e^{-a}}{a} + \frac{k}{a} \int^\infty_1 x^{k-1} e^{-ax} $$ so that $$ I(k) = \frac{e^{-a}}{a} + \frac{k}{a} I(k-1) $$ and $I(0) = \frac{e^{-a}}{a}$.

$\endgroup$
4
  • $\begingroup$ absolutely no idea how i overlooked this. thank you. $\endgroup$
    – drjrm3
    Dec 13, 2011 at 18:09
  • 1
    $\begingroup$ Clever substitutions and integration by parts should always be one of the first things you do with unruly integrals. $\endgroup$
    – J. M.
    Dec 14, 2011 at 4:30
  • $\begingroup$ It's often a good idea to ask a computer algebra system whenever you have an integral like this. Maple evaluates "$\int_1^\infty x^{2n} \exp(-\alpha x) \mathrm{d}x \text{ assuming }n :: \text{nonnegint}, \alpha > 0$" immediately into $\Gamma(2n+1, \alpha) \alpha^{-2n-1}$; I'm sure Mathematica does the same. (Still a good idea to verify it numerically, of course, which these guys can typically also do.) $\endgroup$
    – Erik P.
    Dec 15, 2011 at 4:06
  • $\begingroup$ Actually Mathematica chooses to represent the answer as ExpIntegralE[-2 n, a r]. If you run FunctionExpand on it, then it gives the same answer as Maple. $\endgroup$
    – Searke
    Jan 3, 2012 at 21:36
1
$\begingroup$

Have a look at QUADPACK. It has routines for integration over (semi-)infinite domains.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.