Improper integral Monte Carlo method

Question

The task is to solve an integral using Monte Carlo method: $$ \int_0^\infty \frac{dx}{(x+1) \sqrt x} = \pi $$ But I have not found anywhere how to solve integrals with infinite limits by the Monte Carlo method. Where should I start? What can I do? Thanks in advance!

Answer 1

$\newcommand{\bbx}[1]{\,\bbox[15px,border:1px groove navy]{\displaystyle{#1}}\,} \newcommand{\braces}[1]{\left\lbrace\,{#1}\,\right\rbrace} \newcommand{\bracks}[1]{\left\lbrack\,{#1}\,\right\rbrack} \newcommand{\dd}{\mathrm{d}} \newcommand{\ds}[1]{\displaystyle{#1}} \newcommand{\expo}[1]{\,\mathrm{e}^{#1}\,} \newcommand{\ic}{\mathrm{i}} \newcommand{\mc}[1]{\mathcal{#1}} \newcommand{\mrm}[1]{\mathrm{#1}} \newcommand{\pars}[1]{\left(\,{#1}\,\right)} \newcommand{\partiald}[3][]{\frac{\partial^{#1} #2}{\partial #3^{#1}}} \newcommand{\root}[2][]{\,\sqrt[#1]{\,{#2}\,}\,} \newcommand{\totald}[3][]{\frac{\mathrm{d}^{#1} #2}{\mathrm{d} #3^{#1}}} \newcommand{\verts}[1]{\left\vert\,{#1}\,\right\vert}$ The usual Monte Carlo procedure is given by $\ds{\int_{a}^{b}\mrm{P}\pars{x}\mrm{f}\pars{x}\dd x \approx {1 \over N}\sum_{k = 1}^{N}\mrm{f}\pars{x_{k}}}$ where

• $\ds{\mrm{P}\pars{x}}$ is a PDF in $\ds{\left[a,b\right)}$.
• $\ds{x_{1}, x_{2},\ldots,x_{N}}$ are generated by the distribution $\ds{\mrm{P}\pars{x}}$.
• $\ds{Usually, N\ \mbox{is a}\ "large\ number"}$.

Given a particular integration $\ds{\int_{a}^{b}\phi\pars{x}\,\dd x}$, you write it as $$ \int_{a}^{b}\mrm{P}\pars{x}\bracks{\phi\pars{x} \over \mrm{P}\pars{x}}\,\dd x \approx {1 \over N}\sum_{k = 1}^{N}{\phi\pars{x_{k}} \over \mrm{P}\pars{x_{k}}}\,,\qquad N \gg 1 $$ where $\ds{P}\pars{x}$ is "conveniently chosen". Note that $\ds{\mrm{P}\pars{x} \geq 0\ \mbox{and}\ \int_{a}^{b}\mrm{P}\pars{x}\dd x = 1}$.

For example,

• $\ds{\int_{4}^{9}x^{2}\,\dd x = 5\int_{4}^{9}{1 \over 5}\,x^{2}\,\dd x \approx 5\bracks{{1 \over 10000}\sum_{k = 1}^{10000}x_{k}^{2}}}$ where $\ds{\braces{x_{k}}}$ are generated uniformly in $\ds{\left[4,9\right)}$.
• $\ds{\int_{0}^{\infty}x^{6}\expo{-x}\,\dd x \approx {1 \over 546989}\sum_{k = 1}^{546989}x_{k}^{6}}$ where $\ds{\braces{x_{k}}}$ are "generated" $\ds{\mbox{with}~\expo{-x}}$.

---

Lets go to the present case ( in general, it's convenient to remove the integrable singularities as $\ds{1/\root{x}}$, but lets keep it for the time being ): \begin{align} \int_{0}^{\infty}{\dd x \over \pars{1 + x}\root{x}} & = \int_{0}^{\infty}\overbrace{1 \over \pars{x + 1}^{2}} ^{\ds{\mrm{P}\pars{x}}}\ {1 + x \over \root{x}}\,\dd x \approx {1 \over 10^{6}}\sum_{n = 0}^{10^{6} - 1} {1 + x_{n} \over \root{x_{n}}} \end{align}

---

The following ${\tt javascript}$ code performs the above task:

// gosrabios10sep2020.js
// Run as node gosrabios10sep2020.js in a terminal
"use strict";
const ITERATIONS = 1000000; // One million
let myRand = (function()
{
 let myR = null, temp = null;
return function()
{
do {
myR = Math.random();
temp = 1.0 - myR;
} while (temp <= 0);
return myR/temp;
};
})();
let total = 0, x = null;
for ( let n = 0 ; n < ITERATIONS ; ++n)
{
x = myRand();
total += (x + 1.0)/Math.sqrt(x);
}
console.log(total/ITERATIONS);

A "typical run" yields $\ds{\large{3.143321704930537}}$.

Answer 2

Do a change of variable. If you use $x=\tan y$, the limits correspond for $y$ from $0$ to $\pi/2$. In this case it would be even better to use $x=\tan^2 y$. Then $\sqrt x=\tan y$, and $$x+1=1+\frac{\sin^2 y}{\cos^2 y}=\frac 1{\cos^2 y}$$ All is left is to write $dx$.