Gauss-Laguerre quadrature

339 Views Asked by At

I am trying to compute this integral: $$ \int_{0}^{\infty}\prod_{k = 1}^{d}\left(1 - \,\mathrm{e}^{-a_{k}\,t}\right) \,\mathrm{e}^{-t}\,\mathrm{d}t,\quad \mbox{where}\quad a_{k} > 0, \forall\ k. $$

I can compute this just fine for small values $d$, e.g., less than $100$, using a numerical Gauss-Laguerre quadrature. I am having trouble computing this accurately when $d$ gets larger. Any suggestions on how to solve this problem ?.

3

There are 3 best solutions below

7
On

Let we exploit $\int_{0}^{+\infty}e^{-kx}\,dx=\frac{1}{k}$ by expanding the product:

$$ I(\alpha_1,\ldots,\alpha_d)=\int_{0}^{+\infty}e^{-x}\prod_{k=1}^{d}\left(1-e^{-\alpha_k x}\right)\,dx \\= 1-\sum_{k=1}^{d}\frac{1}{\alpha_k+1}+\sum_{1\leq k_1<k_2\leq d}\frac{1}{\alpha_{k_1}+\alpha_{k_2}+1}-\ldots+(-1)^d\frac{1}{\alpha_1+\ldots+\alpha_d+1}\\=\color{blue}{\int_{0}^{1}\prod_{k=1}^{d}\left(1-x^{\alpha_k}\right)\,dx}$$ the last integral does not have a nice general closed form, but it is easy to approximate through Holder's inequality, since Euler's beta function gives:

$$ \int_{0}^{1}(1-x^{\alpha})^p\,dx = \frac{\Gamma(p+1)\,\Gamma\left(1+\frac{1}{a}\right)}{\Gamma\left(p+1+\frac{1}{a}\right)}.$$

4
On

$\newcommand{\angles}[1]{\left\langle\,{#1}\,\right\rangle} \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{\half}{{1 \over 2}} \newcommand{\ic}{\mathrm{i}} \newcommand{\iff}{\Longleftrightarrow} \newcommand{\imp}{\Longrightarrow} \newcommand{\Li}[1]{\,\mathrm{Li}_{#1}} \newcommand{\ol}[1]{\overline{#1}} \newcommand{\pars}[1]{\left(\,{#1}\,\right)} \newcommand{\partiald}[3][]{\frac{\partial^{#1} #2}{\partial #3^{#1}}} \newcommand{\ul}[1]{\underline{#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}$

With $\ds{a_{k} > 0\,,\ \forall\ k \in \braces{1,2,\ldots,d}}$:

\begin{align} &\color{#f00}{\int_{0}^{\infty}\prod_{k = 1}^{d}\pars{1 - \expo{-a_{k}\,t}} \expo{-t}\,\dd t} = \int_{0}^{\infty}\prod_{k = 1}^{d} \pars{a_{k}t\int_{0}^{1}\expo{-a_{k}\,t\,x_{k}}\,\dd x_{k}}\expo{-t}\,\dd t \\[3mm] = &\ \pars{\prod_{k = 1}^{d}a_{k}} \int_{0}^{1}\cdots\int_{0}^{1}\int_{0}^{\infty}t^{d} \exp\pars{-\bracks{\sum_{j = 1}^{d}a_{j}x_{j} + 1}t}\,\dd t \,\prod_{k = 1}^{d}\dd x_{k} \\[3mm] = &\ \Gamma\pars{d + 1}\pars{\prod_{k = 1}^{d}a_{k}}\,\mathrm{f}_{d}\pars{\vec{a}} \tag{1} \\[5mm] &\ \mbox{where}\ \,\mathrm{f}_{d}\pars{\vec{a}} \equiv \int_{0}^{1}\cdots\int_{0}^{1} {\prod_{k = 1}^{d}\,\dd x_{k} \over \pars{\vec{a}\cdot\vec{x} + 1}^{d + 1}} \\[2mm] & \mbox{and}\quad \left\lbrace\begin{array}{rcl} \ds{\vec{a}} & \ds{\equiv} & \ds{\pars{a_{1},a_{2}\,\ldots,a_{d}}} \\[2mm] \ds{\vec{x}} & \ds{\equiv} & \ds{\pars{x_{1},x_{2}\,\ldots,x_{d}}} \end{array}\right. \end{align}


It's still a 'hard' one. However, when $\ds{d}$ is 'very large' the main contribution arises 'from' $\ds{\vec{a}\cdot\vec{x} \gtrsim 0}$ such that \begin{align} \,\mathrm{f}_{d}\pars{\vec{a}} & \sim \int_{0}^{1}\cdots\int_{0}^{1} \exp\pars{-\bracks{d + 1}\vec{a}\cdot\vec{x}} \prod_{k = 1}^{d}\,\dd x_{k} \\[5mm] & = \prod_{k = 1}^{d}\int_{0}^{1} \exp\pars{-\bracks{d + 1}a_{k}x_{k}}\,\dd x_{k} = \prod_{k = 1}^{d}{1 - \expo{-\pars{d + 1}a_{k}} \over \pars{d + 1}a_{k}} \\[5mm] & = {1 \over \pars{d + 1}^{d}} \prod_{k = 1}^{d}{1 - \expo{-\pars{d + 1}a_{k}} \over a_{k}} \end{align} and $\ds{\pars{~\mbox{see expression}\ \pars{1}~}}$ \begin{align} &\color{#f00}{\int_{0}^{\infty}\prod_{k = 1}^{d}\pars{1 - \expo{-a_{k}\,t}} \expo{-t}\,\dd t} \\[3mm] \sim &\ \root{2\pi}d^{d + 1/2}\,\,\expo{-d}\pars{\prod_{k = 1}^{d}a_{k}} {1 \over \pars{d + 1}^{d}} \prod_{k = 1}^{d}{1 - \expo{-\pars{d + 1}a_{k}} \over a_{k}}\tag{2} \\[3mm] & = \root{2\pi}d^{1/2}{1 \over \pars{1 + 1/d}^{d}}\,\expo{-d} \sim \color{#f00}{\root{2\pi}d^{1/2}\expo{-\pars{d + 1}}} \end{align}

Note that $$ \bracks{1 - \expo{-\pars{d + 1}\min\braces{a_{k}}}\,\,}^{d} < \prod_{k = 1}^{d}\bracks{1 - \expo{-\pars{d + 1}a_{k}}} < \bracks{1 - \expo{-\pars{d + 1}\max\braces{a_{k}}}\,\,}^{d} $$

4
On

The large number of factors will make the integrand look like a Gaussian near its peak. The maximum of the integrand is at $t = t^*$ satisfying the equation:

$$\sum_{j=1}^{d} \frac{a_j}{\exp(a_j t^*)-1} = 1$$

The second derivative of the logarithm of the integrand at the maximum is $-\sigma$ with

$$\sigma = \sum_{j=1}^{d} \frac{a_j^2\exp(a_j t^*)}{\left[\exp(a_j t^*)-1\right]^2}$$

The integrand near its peak is then approximately proportional to $\exp\left[-\frac{\sigma}{2}\left(t-t^*\right)^2\right]$. Using the values for $t^*$ and $\sigma$, you can set up an optimal Gauss-Hermite quadrature scheme.