Asymptotic bound for $\int_0^\infty \int_0^\infty (x+y)^m e^{-\frac{x^2}{2i} - \frac{y^2}{2j}} dx\, dy\;$ for $i$ and $j$ large

Whilst trying to count certain types of bipartite graphs, I'm lead to try to bound the following quantity $$I:=\int_0^\infty \int_0^\infty (x+y)^m e^{-\frac{x^2}{2i} - \frac{y^2}{2j}} dx\,dy$$ where $$i,j$$ and $$m$$ are integers, and I'm interested in the asymptotics for large $$i$$ and $$j$$ and potentially $$m$$ (although it would suffice to have a good upper bound when $$i \approx j$$ and $$m=o(i)$$).

One can derive an exact expression for the integral by multiplying out the terms and using known identities for the quantities $$\int_0^\infty x^k e^{-\frac{x^2}{2i}} dx$$, however the asymptotics of this sum is unclear to me.

It would seem more natural to use a type of `saddle-point' method here, approximate the logarithm of the function around its maximum at $$(x_0,y_0) = \left(i \sqrt{\frac{m}{i+j}},j \sqrt{\frac{m}{i+j}} \right)$$ using the first two terms of the Taylor series, and so evaluate the integral in this region as a standard Gaussian, and then show that the contribution from outside this region is negligible.

This would lead to the following bound, which I would guess is in fact the correct asymptotic order $$I \approx \exp\left(m\log\sqrt{(i+j)(m)}-\frac{m}{2}\right)\pi\sqrt{2ij}.$$ However, I can't get the regions in which the approximation is correct and the region in which the integral is negligible to overlap.

I suspect that this integral will have been considered somewhere in the literature, or at the very least will be susceptible to standard techniques in a field I'm not familiar with.

• Have a look at my new answer. – Claude Leibovici Sep 19 '20 at 14:05

A upper bound

(With the help of Maple)

With the substitution $$u = x+y, v = y$$, we have \begin{align} I &= \int_0^\infty \int_0^u u^m \mathrm{e}^{-(u-v)^2/(2i) - v^2/(2j)} \mathrm{d} v \mathrm{d}u\\ &= \int_0^\infty \sqrt{\frac{\pi ij}{2i+2j}}\, u^m \mathrm{e}^{-\frac{u^2}{2i+2j}} \left[\mathrm{erf}\Big(\tfrac{u}{i}\sqrt{\tfrac{ij}{2i+2j}}\Big) + \mathrm{erf}\Big(\tfrac{u}{j}\sqrt{\tfrac{ij}{2i+2j}}\Big) \right] \mathrm{d}u\\ &= \int_0^\infty \sqrt{\frac{\pi ij}{2i+2j}}\, u^m \mathrm{e}^{-\frac{u^2}{2i+2j}} \mathrm{erf}\Big(\tfrac{u}{i}\sqrt{\tfrac{ij}{2i+2j}}\Big)\mathrm{d}u \\ &\qquad + \int_0^\infty \sqrt{\frac{\pi ij}{2i+2j}}\, u^m \mathrm{e}^{-\frac{u^2}{2i+2j}} \mathrm{erf}\Big(\tfrac{u}{j}\sqrt{\tfrac{ij}{2i+2j}}\Big) \mathrm{d}u\\ &= \sqrt{\pi}(\tfrac{2i+2j}{ij})^{m/2}i^{m+1}\int_0^\infty w^m \mathrm{erf}(w)\mathrm{e}^{-w^2i/j} \mathrm{d} w\\ &\qquad + \sqrt{\pi}(\tfrac{2i+2j}{ij})^{m/2}j^{m+1}\int_0^\infty w^m \mathrm{erf}(w)\mathrm{e}^{-w^2j/i}\mathrm{d} w\\ &= \sqrt{\pi}(\tfrac{2i+2j}{ij})^{m/2}i^{m+1} \Big(\int_0^\infty w^m \mathrm{e}^{-w^2i/j} \mathrm{d} w - \int_0^\infty w^m (1 - \mathrm{erf}(w))\mathrm{e}^{-w^2i/j} \mathrm{d} w\Big)\\ &\qquad + \sqrt{\pi}(\tfrac{2i+2j}{ij})^{m/2}j^{m+1} \Big(\int_0^\infty w^m \mathrm{e}^{-w^2j/i}\mathrm{d} w - \int_0^\infty w^m (1-\mathrm{erf}(w))\mathrm{e}^{-w^2j/i}\mathrm{d} w\Big)\\ &= 2\sqrt{\pi}2^{m/2-1}(i+j)^{m/2}\sqrt{ij}\, \Gamma(\tfrac{m+1}{2})\\ &\qquad - \sqrt{\pi}(\tfrac{2i+2j}{ij})^{m/2}i^{m+1} \int_0^\infty w^m (1 - \mathrm{erf}(w))\mathrm{e}^{-w^2i/j} \mathrm{d} w\\ &\qquad - \sqrt{\pi}(\tfrac{2i+2j}{ij})^{m/2}j^{m+1} \int_0^\infty w^m (1-\mathrm{erf}(w))\mathrm{e}^{-w^2j/i}\mathrm{d} w\\ &\le 2\sqrt{\pi}2^{m/2-1}(i+j)^{m/2}\sqrt{ij}\, \Gamma(\tfrac{m+1}{2})\\ &\qquad - \sqrt{\pi}(\tfrac{2i+2j}{ij})^{m/2}i^{m+1} \int_0^\infty w^m \Big(\sqrt{\frac{2\mathrm{e}}{\pi}}\frac{\sqrt{\beta-1}}{\beta}\mathrm{e}^{-\beta w^2}\Big)\mathrm{e}^{-w^2i/j} \mathrm{d} w\\ &\qquad - \sqrt{\pi}(\tfrac{2i+2j}{ij})^{m/2}j^{m+1} \int_0^\infty w^m \Big(\sqrt{\frac{2\mathrm{e}}{\pi}}\frac{\sqrt{\beta-1}}{\beta}\mathrm{e}^{-\beta w^2}\Big)\mathrm{e}^{-w^2j/i} \mathrm{d} w\\ &= 2\sqrt{\pi}2^{m/2-1}(i+j)^{m/2}\sqrt{ij}\, \Gamma(\tfrac{m+1}{2})\\ &\qquad - \sqrt{\pi}(\tfrac{2i+2j}{ij})^{m/2}i^{m+1} \sqrt{\frac{\mathrm{e}}{2\pi}}\frac{\sqrt{\beta-1}}{\beta} (\beta +\tfrac{i}{j})^{-(m+1)/2}\Gamma(\frac{m+1}{2})\\ &\qquad - \sqrt{\pi}(\tfrac{2i+2j}{ij})^{m/2}j^{m+1} \sqrt{\frac{\mathrm{e}}{2\pi}}\frac{\sqrt{\beta-1}}{\beta} (\beta +\tfrac{j}{i})^{-(m+1)/2}\Gamma(\frac{m+1}{2}) \end{align} where $$\mathrm{erf}(w) = \frac{2}{\sqrt{\pi}}\int_0^w \mathrm{e}^{-t^2}\mathrm{d} t$$ is the error function, and we have used $$1 - \mathrm{erf}(w) \ge \sqrt{\frac{2\mathrm{e}}{\pi}}\frac{\sqrt{\beta-1}}{\beta}\mathrm{e}^{-\beta w^2}$$ (for $$w\ge 0$$, $$\beta > 1$$; see https://en.wikipedia.org/wiki/Error_function). We may choose $$\beta = \frac{5}{4}$$.

• That's very useful. I still have hope that the approximate bound can be determined by a more elementary argument, but as far as I can tell the asymptotics from that first term are precisely what I need. – Joshua Erde Sep 18 '20 at 19:14
• @JoshuaErde I did some numerical experiment. For example, for $m = 10, i = 1000, j = 900, b = 5/4$, denote $I_u = 2\sqrt{\pi}2^{m/2-1}(i+j)^{m/2}\sqrt{ij}\, \Gamma(\tfrac{m+1}{2})$, then $\frac{I_u}{I} \approx 1.0072$; We need to prove (or disprove) that $I_u$ (the first term in my last upper bound expression, ignore the 2nd and 3rd terms) is asymptotic upper bound, as $i, j\to \infty$. – River Li Sep 19 '20 at 0:23
• Oh, perhaps I've misunderstood. I had thought your argument gives $I_u$ as an upper bound for the interval, as the other two terms are negative. Or are you just saying that it's not clear that asymptotically $I_u \approx I$? The upper bound is sufficient for my application, but I agree it would be nice to know it's asymptotically tight. – Joshua Erde Sep 19 '20 at 7:45
• @JoshuaErde Yes, it may hold that $\lim_{i, j\to \infty} \frac{I_u}{I} = 1$. – River Li Sep 19 '20 at 11:54

There is an exact solution for the integral (big surprise for me !).

I try to write the expression for $$-\frac{\sqrt{2} (m+1)}{j}\,I_m=T_1+T_2+T_3+T_4$$ $$T_1=-\frac{\sqrt{\frac{\pi }{2}} (m+1) i^{\frac{m+1}{2}} \Gamma \left(\frac{m+1}{2}\right) \left(\frac{2 j}{i}+2\right)^{m/2}}{\sqrt{j}}$$ $$T_2=2^{\frac{m+1}{2}} j^{m/2} \Gamma \left(\frac{m}{2}+1\right) \, _2F_1\left(\frac{1}{2},1;\frac{m+3}{2};-\frac{j}{i}\right)$$ $$T_3=\frac{2^{\frac{m+1}{2}} (m+1) (i+j) i^{m/2} \Gamma \left(\frac{m}{2}+1\right) \, _2F_1\left(1,\frac{1-m}{2};-\frac{1}{2};-\frac{j}{i}\right)}{j m}$$ $$T_4=-\frac{2^{\frac{m+1}{2}} (m+1) i^{m/2} \Gamma \left(\frac{m}{2}+1\right) (i-j (m-3)) \, _2F_1\left(1,\frac{1-m}{2};\frac{1}{2};-\frac{j}{i}\right)}{j m}$$

I have the feeling that I have been unable to simplify properly.

Edit

Using $$(x+y)^m=\sum_{k=0}^m \binom{m}{k}\, x^{m-k}\,y^k$$ $$\int_0^\infty\int_0^\infty x^{m-k}\, y^k\,e^{-\frac{x^2}{2 i}-\frac{y^2}{2 j}}\,dx\,dy=2^{\frac{m-2}{2}} i^{\frac{m+1-k}{2} }j^{\frac{k+1}{2}} \Gamma \left(\frac{k+1}{2}\right) \Gamma \left(\frac{m+1-k}{2} \right)$$ and then the hypergeometric functions by the summation.

If, as I suggested in a comment, we let $$i=p^2$$ and $$j=a^2p^2$$ $$i^{\frac{m+1-k}{2} }j^{\frac{k+1}{2}}=a^{k+1} p^{m+2}$$ which could be more comfortable.

• Thanks. Do you have a reference for this result, and what is the function $F_1$? – Joshua Erde Sep 18 '20 at 11:43
• @JoshuaErde. These are gaussian hypergeometric functions. The result was not obtained by hand with pen and paper. A CAS did it – Claude Leibovici Sep 18 '20 at 11:45
• Very nice work, Claude. How did you evaluate the hypergeometric functions? Your formulas for $I_1,I_2,...$ are all in terms of $i,j,\pi$ and powers of these values. I am very curious. – K.defaoite Sep 18 '20 at 12:32
• @K.defaoite. Very easy ! Ask a CAS. Cheers :) – Claude Leibovici Sep 18 '20 at 12:56

Wondering if I made or not a mistake somewhere in my previous answer, I restarted using the binomial expansion of $$(x+y)^m$$ and ended with something apparently simpler (but also apparently different). The end result write $$I=\frac{(2i)^{\frac m2}}{4j} \Big[\cdots\Big]$$ with $$\Big[\cdots\Big]=\Gamma \left(\frac{m-2}{2}\right) \left(i (i-j(m-3))-(i+j)^2 \, _2F_1\left(1,\frac{1-m}{2};-\frac{1}{2};-\frac{j}{i}\right)\right)+$$ $$2j \sqrt{\pi ij}\, \left(\frac{i+j}{i}\right)^{\frac m2} \Gamma \left(\frac{m+1}{2}\right)$$

If the first term can be neglected (or similar to the second one), then $$I\sim \sqrt \pi \,2^{\frac{m-2}{2}} \Gamma \left(\frac{m+1}{2}\right) (i+j)^{\frac m2}\sqrt{ij}$$ which looks like what you wrote.

Edit

In order to check, I made $$j=i$$ which makes $$\frac I{2^{\frac{m-4}{2}} i^{\frac{m+2}{2}}}=$$ $$\sqrt{\pi }\, 2^{\frac{m}{2}+1} \Gamma \left(\frac{m+1}{2}\right)-\left(4 \, _2F_1\left(1,\frac{1-m}{2};-\frac{1}{2};-1\right)+m-4\right) \Gamma \left(\frac{m-2}{2}\right)$$

I rewrote it as $$I=\sqrt{\pi }\, 2^{m-1}\, i^{\frac{m+2}{2}} \Gamma \left(\frac{m+1}{2}\right)\, (1-K)$$ with $$K=\frac {\left(4 \, _2F_1\left(1,\frac{1-m}{2};-\frac{1}{2};-1\right)+m-4\right) \Gamma \left(\frac{m}{2}-1\right) } {\sqrt{\pi }\, 2^{\frac{m+2}{2}} \Gamma \left(\frac{m+1}{2}\right) }$$

As shown below, when $$m$$ increases, factor $$K$$ tends very fact to $$-1$$ $$\left( \begin{array}{cc} m & K \\ 3 & -0.883883 \\ 4 & -0.924413 \\ 5 & -0.950175 \\ 6 & -0.966854 \\ 7 & -0.977796 \\ 8 & -0.985044 \\ 9 & -0.989880 \\ 10 & -0.993128 \\ 15 & -0.998968 \\ 20 & -0.999839 \end{array} \right)$$

In other words, at least for $$j=i$$, for large $$m$$, an asymptotics is $$I \sim \sqrt{\pi }\, 2^m\, i^{\frac{m+2}{2}} \,\Gamma \left(\frac{m+1}{2}\right)$$

Using Stirling approximation $$\log(I) =m\log \left(\frac{2 i m}{e}\right)+\log \left(\sqrt{2} \pi i\right)-\frac{1}{12 m}+O\left(\frac{1}{m^3}\right)$$ which seems to be very close to what you wrote.

• Thanks. This seems to agree with River Li's answer in terms of the leading term. – Joshua Erde Sep 21 '20 at 9:45

$$\color{green}{\textbf{Final version of 26.09.20.}}$$

$$\color{brown}{\mathbf{Integration.}}$$

The area of the integration is the first quadrant.

Using the polar coordinates after scaling substitutions, one can get $$I=\int\limits_0^\infty\int\limits_0^\infty(x+y)^me^{^{\Large-\frac{x^2}{2i}-\frac{y^2}{2j}}}\,\text dx\,\text dy =\left|\genfrac{}{}{0}{}{\frac {x\mathstrut}{\sqrt{i\small\mathstrut}}\mapsto x}{\frac y{\sqrt j}\mapsto y}\right|\\ =\sqrt{ij}\int\limits_0^\infty\int\limits_0^\infty(x\sqrt{i\small\mathstrut}+y\sqrt j)^me^{^{\Large-\frac{x^2}{2}-\frac{y^2}{2}}}\,\text dx\,\text dy =\left|\genfrac{}{}{0}{}{x\mapsto \rho\cos\varphi}{y\mapsto \rho\sin\varphi}\right|\\ =\sqrt{ij}\int\limits_0^{\large^\pi/_2}\int\limits_0^\infty\left(\sqrt{i\small\mathstrut}\,\cos\varphi+\sqrt j\sin\varphi\right)^me^{^{\Large-\frac{\rho^2}{2}}}\,\rho^{m+1}\,\text d\rho\,\text d\varphi,$$ $$I=\sqrt{ij(i+j)}I_1 I_2,\tag1$$ where $$I_1 = \int\limits_0^\infty e^{^{\Large-\frac{\rho^2}{2}}}\,\rho^{m+1}\,\text d\rho = 2^{^{\frac{m-1}2}}\Gamma\left(\frac{m+1}2\right),\tag2$$ $$I_2 = \int\limits_0^{\large^\pi/_2}\left(\sqrt{\frac{i\small\mathstrut}{i+j}}\,\cos\varphi+\sqrt{\frac{j\small\mathstrut}{i+j}}\sin\varphi\right)^m\,\text d\varphi = \int\limits_0^{\large^\pi/_2}\cos^m\left(\varphi-\arctan\dfrac ji\right) \,\text d\varphi,$$ $$I_2 = G\left(m,\dfrac\pi2-\arctan\dfrac ji\right) - G\left(m,-\arctan\dfrac ji\right),$$ where $$G(m,z)$$ is the antiderivative (see Appendium below).

The plot of the antiderivatives $$G(2,z) = \frac z2\left(1+\text{sinc }2z\right),$$ $$G(4,z) = \frac z8\left(3+4\text{ sinc }2z+\text{sinc }4z\right),$$ $$G(6,z) = \frac z{32}\left(10+15\text{ sinc }2z+6\text{ sinc }4z+\text{sinc }6z \right),$$ see below.

The plot of $$\;G(20,z):$$

The plot of $$G(m,z),\;m=3,5,7,9:$$

The plot of $$\;G(21,z):$$

The plot of $$\;I_2(21,z),\; z=\arctan \dfrac ji:$$

The highest values of $$G(m,z)$$ can be expressed in the hypergeometric functions.

$$\color{brown}{\mathbf{Conclusions.}}$$

From the consideration above and additional investigations, the next conclusions can be done.

• If $$\;i\to\infty,j\to\infty,\dfrac ji =\text{const},m=\text{const},\;$$ then $$\;I\to \text{const}\cdot\sqrt{ij(i+j)}.$$
• If $$\;m\to\infty,i=\text{const},j=\text{const},\;$$ then $$\;I\to\text{const}\cdot\Gamma\left(\frac{m+1}2\right) (\sqrt2)^m I_2.\;$$
• Integral $$\;|I_2|\lesssim\dfrac\pi{2\sqrt[3]{m+1}}.\;$$ It decreases additionally if $$\;i\gg j\;$$ or $$\;j\gg i.$$

$$\color{brown}{\mathbf{Appendium.\ Antiderivatives\ of\ cos^m(z).}}$$

$$G(m,z) = \int \cos^{m} z\,\text dz,\tag{A1}$$ wherein \begin{align} & G(0,z) = z,\quad G(1,z) = \sin z,\tag{A2}\\[4pt] &2^{2k}\cos^{2k} z = (e^{iz}+e^{-iz})^{2k} = \sum\limits_{j=0}^{k-1}\dbinom{2k}{j}\left(e^{(2k-2j)iz}+e^{(2j-2k)iz}\right)+\dbinom{2k}k\\[4pt] & = 2\sum\limits_{j=0}^{k-1}\dbinom{2k}j\cos(2k-2j)z+\dbinom{2k}k,\\[4pt] &2^{2k+1}\cos^{2k+1} z = (e^{iz}+e^{-iz})^{2k+1} = \sum\limits_{j=0}^k\dbinom{2k+1}{j}\left(e^{(2k+1-2j)iz}+e^{(2j-2k-1)iz}\right)\\[4pt] & = 2\sum\limits_{j=0}^{k}\dbinom{2k+1}{j}\cos(2k+1-2j)z,\\[4pt] &\cos^m z = 2^{-m}\left(2\sum\limits_{j=0}^{\genfrac\lfloor\rfloor{}{}{\large m-1}2}\dbinom{m}{j}\cos(m-2j)z+\left(\genfrac\lceil\rceil{}{}{m-1}2-\genfrac\lfloor\rfloor{}{}{m-1}2\right)\right),\tag{A3}]\\[4pt] & G(m,z) = \dfrac z{2^m}\left(2\sum\limits_{j=0}^{\genfrac\lfloor\rfloor{}{}{\large m-1}2}\dbinom{m}{j}\text{ sinc }((m-2j)z)+\genfrac\lceil\rceil{}{}{m-1}2-\genfrac\lfloor\rfloor{}{}{m-1}2\right).\tag{A4}\\ \end{align}

The table of $$G(m,z)$$ for $$m=2\dots8$$ is presented below. $$\begin{vmatrix} BG(m,z) & m=2 & m=3 & m=4 & m=5 & m=6 & m=7 & m=8 \\ B & 2 & 4 & 8 & 16 & 32 & 64 & 128 \\ z & 1 & 0 & 3 & 0 & 10 & 0 & 35 \\ \sin z & 0 & 3 & 0 & 10 & 0 & 35 & 0 \\ \dfrac{\sin2z}2 & 1 & 0 & 4 & 0 & 15 & 0 & 56 \\ \dfrac{\sin3z}3 & & 1 & 0 & 5 & 0 & 21 & 0 \\ \dfrac{\sin4z}4 & & & 1 & 0 & 6 & 0 & 28 \\ \dfrac{\sin5z}5 & & & & 1 & 0 & 7 & 0 \\ \dfrac{\sin6z}6 & & & & & 1 & 0 & 8 \\ \dfrac{\sin7z}7 & & & & & & 1 & 0 \\ \dfrac{\sin8z}8 & & & & & & & 1 \tag{A5} \end{vmatrix}$$

Not an answer, but this is too long for a comment. Perhaps one can expand $$(x+y)^m$$ into a sum (which will probably work out nicely if $$m\in\Bbb{N}$$, otherwise perhaps not) and use an identity given by Mathematica: $$\int_0^\infty x^n\exp\left(\frac{-ax^2}{2}\right)\exp(-bx)\mathrm{d}x$$ $$=2^{\frac{n-1}{2}} a^{-\left(\frac{n+2}{2}\right)}\left( -b\sqrt{2} \ \Gamma \left(\frac{n+2}{2}\right) \ _{1} F_{1}\left(\left[\frac{n+2}{2} ,\frac{3}{2}\right] ;\frac{b^{2}}{2a}\right) +\sqrt{a} \ \Gamma \left(\frac{n+1}{2}\right) \ _{1} F_{1}\left(\left[\frac{n+1}{2} ,\frac{1}{2}\right] ;\frac{b^{2}}{2a}\right)\right)$$ For $$\operatorname{Re}(a)>0$$ and $$\operatorname{Re}(n)>-1$$. Here $${}_1F_1$$ is the Kummer confluent hypergeometric function of the first kind, defined by the power series $${}_1 F_1\left([\alpha,\beta];z\right)=\sum_{k=0}^\infty \frac{\Gamma(\alpha+k)\Gamma(\beta)}{\Gamma(\beta+k)\Gamma(\alpha)}z^k$$ For $$\alpha,\beta,z \in\Bbb{C}$$.

• Yes, this is the exact expression I talked about in the question, however it's not at all clear how it behaves asymptotically. I would guess that for specific $i$ and $j$ there is some $k$ such that the term with the power $x^k y^{m-k}$ is dominating the whole sum, but I'm hoping there's a more straightforward way to approximate the integral. – Joshua Erde Sep 18 '20 at 12:44

$$I=$$

For $$i\inℕ$$ and $$j\inℕ$$ always true. But $$y\inℝ$$ and $$x\inℝ$$ without restrictions.

If the restriction are presented first to the integral then the solutions looks this way: $$I(x/;ℝ,y/;ℝ,i/;ℕ,j/;ℕ,m/;ℕ)=$$

This is done with Mathematica and Wolfram Language. I used the built-ins $$Assuming$$ and $$Integrate$$. This should be working with CAS that have not the high-level mathmatical knwoledge implemented.