How to solve this double integral?

$$f(a,b)=\int _{x=0}^{\infty }\int _{t=-\infty }^{\infty }\exp \left(\frac{-a t^2+i b t}{3 t^2+1}+i t x\right)\frac{x}{3 t^2+1}\mathrm{d}t\mathrm{d}x$$

$$\text{with }a>0,b\in \mathbb{R},i^2=-1$$

Known special solution for ${\bf b=0}$

$$f(a,0)=\frac{\pi}{\sqrt{3}\, {\rm exp}\left(\frac{a}{6}\right)}\left[(a+3) I_0\left(\frac{a}{6}\right)+a I_1\left(\frac{a}{6}\right)\right]$$

where $I_0,I_1$ are Bessel functions of order $0$ and $1$ (proof).

The difference of the special solution to the general solution $f(1,0)-f(1,b)$ is calculated numerically for $a=1$ and $-50<b<50$.

enter image description here

What I tried

I followed the first steps given here. Substitution of $t\rightarrow t \sqrt{3},x\rightarrow x/\sqrt{3}$ removes the factors in the denominator

$$f(a,b)=\sqrt{3}\int _{x=0}^{\infty }\int _{t=-\infty }^{\infty } \exp \left(\frac{1}{t^2+1}\left(\frac{i b t}{\sqrt{3}}-\frac{a t^2}{3}\right)+i t x\right)\frac{x}{t^2+1}\mathrm{d}t\mathrm{d}x$$

$$=\sqrt{3}\int _{x=0}^{\infty }{\rm d}x \frac{x}{\text{exp}(x)}\int _{t=-\infty }^{\infty} \exp \left(\frac{1}{t^2+1}\left(\frac{i b t}{\sqrt{3}}-\frac{a t^2}{3}\right)+i x(t-i) \right)\frac{1}{t^2+1}\mathrm{d}t$$ $$=\frac{\sqrt{3}}{{\rm exp}(a/3)}\int _{x=0}^{\infty }{\rm d}x \frac{x}{\text{exp}(x)}\underbrace{\int _{t=-\infty }^{\infty} \exp \left(\frac{1}{t^2+1}\left(\frac{i b t}{\sqrt{3}}+\frac{a}{3}\right)+i x(t-i) \right)\frac{1}{t^2+1}\mathrm{d}t}_{I(x)}.$$

Now $I(x)$ can be closed in the upper half plane since the contribution along the arc vanishes. Then this $t$-integral encloses the single essential singularity in the upper half plane at $t=i$. Hence we have

$$I(x)=2\pi i \, {\rm Res} \left(\exp \left(\frac{1}{t^2+1}\left(\frac{i b t}{\sqrt{3}}+\frac{a}{3}\right)+i x(t-i) \right)\frac{1}{t^2+1}\right)\Bigg|_{t=i} \, $$


$$\exp \left(\frac{1}{t^2+1}\left(\frac{i b t}{\sqrt{3}}+\frac{a}{3}\right)+i x(t-i) \right)\frac{1}{t^2+1}$$ can be written as the series $$\sum_{n,m=0}^{\infty}\frac{(ix)^n}{n!}(t-i)^n\frac{1}{m!} \left(\frac{i b t}{\sqrt{3}}+\frac{a}{3}\right)^m\frac{1}{[(t+i)(t-i)]^{m+1}}$$

  • 2
    $\begingroup$ NIntegrate[NIntegrate[Exp[(-t^2 + I*t)/(3*t^2 + 1) + I*t*x]*x/(3*t^2 + 1), {t, -Infinity,Infinity}, WorkingPrecision -> 20], {x, 0, Infinity}] works. The result is 3.49233 in MMA 12.3 $\endgroup$ Commented Dec 29, 2021 at 19:21
  • 1
    $\begingroup$ Any restrictions on $b$? Maybe for small $b$? $\endgroup$
    – Diger
    Commented Dec 29, 2021 at 23:50
  • 1
    $\begingroup$ I can turn it into a single integral: $f(a,b)=2\sqrt{\pi} e^{-b^2/4}-\pi b \, \text{erfc} \, \left(\frac{b}{2}\right)-\int_{-\infty}^{\infty} \frac{dt}{t^2} \left(\frac{e^\frac{i b t-at^2}{3t^2+1}}{3t^2+1} -e^{i b t-t^2}\right)$. Still working at it. $\endgroup$
    – KStarGamer
    Commented Jan 1, 2022 at 2:37
  • 1
    $\begingroup$ I managed to get $f(3a, \sqrt{3} b) = \sqrt{3} \pi e^{-a/2} I_0 \left(\frac{\sqrt{a^2-b^2}}{2}\right) + \sqrt{3\pi} \left(2e^{-b^2/4} -b\sqrt{\pi} \, \text{erfc} \, \left(\frac{b}{2}\right)\right) -\sqrt{3} \int_{-\infty}^{\infty} \left( e^{\frac{i b t-a}{t^2+1}} - e^{ib/t-1/t^2}\right) \, dt$ but haven’t been able to progress further. $\endgroup$
    – KStarGamer
    Commented Jan 5, 2022 at 18:35
  • 1
    $\begingroup$ @Diger I have added a partial solution. I apologise for having not posted it but I have been preoccupied with other things this week. $\endgroup$
    – KStarGamer
    Commented Jan 11, 2022 at 20:11

3 Answers 3


Partial Solution

We begin by transforming our integrand slightly:

$$f(a, b) = \int_{0}^{\infty} x \int_{-\infty}^{\infty} \exp\left(\frac{-a t^2+i b t}{3t^2+1}\right) \frac{e^{i x t}}{3t^2+1} \, dt \, dx$$ Let $a \mapsto 3a, b \mapsto \sqrt{3} \,b$ to give: $$f(3a, \sqrt{3}\,b) = \int_{0}^{\infty} x \int_{-\infty}^{\infty} \exp\left(\frac{-3a t^2+ \sqrt{3}\, i b t}{3t^2+1}\right) \frac{e^{i x t}}{3t^2+1} \, dt \, dx$$ Now perform the substitution $\sqrt{3} \, t = u \implies dt = \frac{1}{\sqrt{3}} du,$ giving us:

$$f(3a, \sqrt{3} \, b) = \frac{1}{\sqrt{3}}\int_{0}^{\infty} x \int_{-\infty}^{\infty} \exp\left(\frac{-a u^2+ \, i b u}{u^2+1}\right) \frac{e^{i x u/\sqrt{3}}}{u^2+1} \, du \, dx$$ Now separate $\exp\left(\frac{-a u^2+ \, i b u}{u^2+1}\right) = e^{-a} \exp\left(\frac{i b u+a}{u^2+1}\right)$ and perform the substitution $x = \sqrt{3} \, v \implies dx = \sqrt{3} \, dv$, giving us:

$$f(3a, \sqrt{3} \, b) = \sqrt{3}\, e^{-a} \int_{0}^{\infty} v \int_{-\infty}^{\infty} \exp\left(\frac{i b u+a}{u^2+1}\right) \frac{e^{i v u}}{u^2+1} \, du \, dv$$ which then the variables can be relabelled with $v \mapsto x$ and $u \mapsto t$.

Consider the following double integral: $$g(3a, \sqrt{3} \, b) = \sqrt{3} e^{-a} \int_{0}^{\infty} x \int_{-\infty}^{\infty} e^{i x t} e^{a-t^2+i b t} \, dt \, dx$$ $$\int_{-\infty}^{\infty} \ e^{i x t} e^{a-t^2+i b t} \, dt = -\frac{i}{2}\sqrt{\pi}e^{a-\frac{1}{4}(b+x)^2}\text{erfi}\left(\frac{b+2it+x}{2}\right) \Bigg]_{t=-\infty}^{t=\infty} = \sqrt{\pi}e^{a-\frac{1}{4}(b+x)^2}$$

We choose this function to allow us to avoid the order two pole at $t=0$ of $f(3a, \sqrt{3} \, b)$ when the integrals are swapped and we have a factor of $\frac{1}{t^2}$.

$$\implies g(3a, \sqrt{3} \, b) =\sqrt{3 \pi} \int_{0}^{\infty} x \exp \left(-\frac{(b+x)^2}{4}\right) \, dx = -\sqrt{3\pi}\left(b\sqrt{\pi}\,\text{erf}\,\left(\frac{b+x}{2}\right)+2e^{-\frac{1}{4}(b+x)^2}\right)\Bigg]_{x=0}^{x=\infty}$$ $$\implies g(3a, \sqrt{3}\, b) = \sqrt{3\pi} \left(2 e^{-b^2/4} - b \sqrt{\pi} \, \text{erfc} \, \left(\frac{b}{2}\right)\right)$$

We now introduce $g(3a, \sqrt{3} \, b)$ in order to allow the order of integration to be interchanged:

$$f(3a, \sqrt{3} \, b)=\sqrt3e^{-a}\int_{0}^{\infty}x\int_{-\infty}^{\infty}e^{i xt}\left(\frac{e^{\frac{a+i bt}{1+t^2}}}{1+t^2}-e^{a-t^2+i bt}\right)\, dt \, dx+\sqrt{3\pi}\left(2e^{-b^2/4}-b\sqrt{\pi}\, \, \text{erfc} \,\left(\frac{b}{2}\right)\right)$$ Since $\int_{0}^{\infty} x^{s-1} e^{i x t}\, dx = (-i t)^{-s} \Gamma (s) \implies \int_{0}^{\infty} x e^{i x t} \, dx = -\frac{1}{t^2}$:

$$\implies f(3a, \sqrt{3} \, b) = -\sqrt{3} e^{-a} \int_{-\infty}^{\infty} \frac{1}{t^2}\left( \frac{e^{\frac{i b t+a}{1+t^2}}}{1+t^2} -e^{a-t^2+i b t} \right)\, dt +\sqrt{3\pi}\left(2e^{-b^2/4}-b\sqrt{\pi}\, \, \text{erfc} \,\left(\frac{b}{2}\right)\right)$$

Now introduce the substitution $\frac{1}{t} = u \implies \frac{1}{t^2} \, dt = -du$. Note the limits will still be from $-\infty$ to $\infty$ because if we split the integral at $0$, we take the limit to $0$ from below for the integral that goes from $-\infty$ to $0$, leaving us with: $$f(3a, \sqrt{3} \, b) = \sqrt{3} e^{-a} \int_{-\infty}^{\infty} \left(\frac{e^{\frac{a+i b/u}{1/u^2 + 1}}}{1+1/u^2} - e^{a+ i b/u - 1/u^2} \right) \, du+\sqrt{3\pi}\left(2e^{-b^2/4}-b\sqrt{\pi}\, \, \text{erfc} \,\left(\frac{b}{2}\right)\right)$$ $$=\sqrt{3} e^{-a} \int_{-\infty}^{\infty} \left(\frac{u^2 e^{\frac{a u^2+i b u}{u^2+1}}}{u^2+1} - e^{a+i b/u-1/u^2}\right) \, du+\sqrt{3\pi}\left(2e^{-b^2/4}-b\sqrt{\pi}\, \, \text{erfc} \,\left(\frac{b}{2}\right)\right)$$ Since $\frac{u^2}{1+u^2} = 1-\frac{1}{u^2+1}$ we have: $$f(3a, \sqrt{3} \, b) = \sqrt{3}\int_{-\infty}^{\infty}\frac{e^{\frac{-a+i bt}{1+t^2}}}{1+t^2}\, dt -\sqrt{3}\int_{-\infty}^{\infty}\left(e^{\frac{-a+i bt}{1+t^2}}-e^{-1/t^2+i b/t}\right) \, dt+\sqrt{3\pi}\left(2e^{-b^2/4}-b\sqrt{\pi}\,\text{erfc}\,\left(\frac{b}{2}\right)\right)$$ $$\sqrt{3}\int_{-\infty}^{\infty}\frac{e^{\frac{-a+i bt}{1+t^2}}}{1+t^2}\, dt = \sqrt{3} e^{-a/2}\int_{-\pi/2}^{\pi/2}e^{-a/2\cos(2u)+ib/2\sin(2u)}\,du= \sqrt{3} \pi e^{-a/2} I_0 \left(\frac{\sqrt{a^2-b^2}}{2}\right)$$

This leaves us with the following final expression:

$$f(3a, \sqrt{3} \, b) = \sqrt{3}\pi e^{-a/2} I_0 \left(\frac{\sqrt{a^2-b^2}}{2}\right) +\sqrt{3\pi}\left(2e^{-b^2/4}-b\sqrt{\pi}\,\text{erfc}\,\left(\frac{b}{2}\right)\right) -\sqrt{3}\int_{-\infty}^{\infty}\left(e^{\frac{-a+i bt}{1+t^2}}-e^{-1/t^2+i b/t}\right) \, dt$$

It is not currently clear to me how to evaluate this final integral. It is fairly trivial to determine the latter $e^{-1/t^2+i b/t}$ half of it through a Mellin transform method, however, I cannot do the same for the former half.


From @Yuri Negometyanov’s excellent answer to OP’s related question we have:

$$f(3a,\sqrt{3}\, b) = - \sqrt{3\pi } \sum_{m=1}^{\infty } \frac{(-a)^m \Gamma \left(m-\frac{1}{2}\right) \,{}_1F_2\left(m-\frac{1}{2};\frac{m+1}{2},\frac{m}{2};-\frac{b^2}{16}\right)}{(m-1)! \, m!}+\frac{\sqrt{3} \pi b^2}{4} \,{}_1F_2\left(\frac{1}{2};\frac{3}{2},2;-\frac{b^2}{16}\right)+\sqrt{3}\pi e^{-a/2} \,{}_0\tilde{F}_1\left(;1;\frac{a^2-b^2}{16} \right)-\sqrt{3} \pi b$$

Where $\tilde{F}$ is the regularised hypergeometric function (which in this case can be expressed in terms of Bessel functions).

  • $\begingroup$ Interesting problem. The solution should be much simpler if we change in the OP the boundary of $x$ from $[0,\infty]$ to $[-\infty,\infty]$. I am interested in this case. @KStarGamer could you add this case? $\endgroup$ Commented Jan 17, 2022 at 19:35
  • $\begingroup$ @Rainer Glüge the integral does not converge on $x \in (-\infty, \infty)$. $\endgroup$
    – KStarGamer
    Commented Jan 17, 2022 at 22:32
  • $\begingroup$ the result is $-2\pi b$ but can show this only numerically, a=1.42;b=1.79;in= Exp[(-a*t*t+b*I*t)/(3*t*t+1)+I*t*x]*x/(3*t*t+1);liRe= NIntegrate`LevinIntegrandReduce[ in/.{{x->x},{x->-x}}//Mean//Re//ComplexExpand// Simplify,t]; lo=Normal@KeyDrop[liRe@"Rules","Variables"]/. HoldPattern["DifferentialMatrices"->{dm_,___}]:> "DifferentialMatrix"->dm; lRe[x0_?NumericQ]:= Block[{x=x0}, NIntegrate[ifunc[x,t],{t,-Infinity,Infinity}, Method->{"LevinRule",Sequence@@lo}]]; NumberForm[ 2NIntegrate[lRe[x],{x,0,Infinity},PrecisionGoal->4],10] NumberForm[%/(-2*Pi*b),10] (*0.999999984*) $\endgroup$ Commented Jan 19, 2022 at 11:09
  • $\begingroup$ @Rainer Glüge, indeed you are correct- I had made a mistake. The double integral does indeed converge to $-2\pi b$ from $x, \, t \in (-\infty,\infty)$. $\endgroup$
    – KStarGamer
    Commented Jan 19, 2022 at 11:20
  • $\begingroup$ @Rainer Glüge: See my answer at the bottom. $\endgroup$
    – Diger
    Commented Jan 20, 2022 at 21:23

Not a closed form, but maybe some progress. We write

$$f(3a,\sqrt{3}b)=\sqrt{3}\int _{x=0}^{\infty }\int _{t=-\infty }^{\infty } \exp \left(\frac{i b t-a t^2}{t^2+1}+i t x\right)\frac{x}{t^2+1} \, \mathrm{d}t \, \mathrm{d}x \equiv \sqrt{3} \, I(a)$$ and close the $t$-integral in the upper half-plane, since the integral along the semi-circle vanishes. Introducing a dampening factor $e^{-\epsilon x^2}$ for which we are ultimately interested in the case $\epsilon \rightarrow 0$, we can now interchange the integrals, since the $x$-integral converges properly. Using $$\int_0^\infty x\, e^{-\epsilon x^2 + itx} \, {\rm d}x = -\frac{1}{t^2} + O(\epsilon) \, ,$$ the remaining $t$-integral becomes $$I(a)=-\oint _{-\infty }^{\infty } \frac{\exp \left(\frac{i b t-a t^2}{t^2+1}\right)}{t^2(t^2+1)} \, \mathrm{d}t \, .$$ The singularity at $t=0$ is not an issue, since the original integral (before interchanging the integrals) has the only singularity in the upper half plane at $t=i$ and the closed contour can be deformed at will without crossing $t=i$.

Therefore $$I'(a)=\oint _{-\infty }^{\infty } \frac{\exp \left(\frac{i b t-a t^2}{t^2+1}\right)}{(t^2+1)^2} \, \mathrm{d}t = e^{-a} \frac{{\rm d}}{{\rm d}a} \oint _{-\infty }^{\infty } \frac{\exp \left(\frac{i b t+a}{t^2+1}\right)}{t^2+1} \, \mathrm{d}t = e^{-a} \frac{{\rm d}}{{\rm d}a} \int _{-\infty }^{\infty } \frac{\exp \left(\frac{i b t+a}{t^2+1}\right)}{t^2+1} \, \mathrm{d}t$$ where in the last step the vanishing integral over the semi-circle was removed again, giving an ordinary indefinite integral. Substituting $t=\tan(\phi/2)$ converts this into $$I'(a) = \frac{e^{-a}}{2} \frac{{\rm d}}{{\rm d}a} \, e^{a/2} \int _{-\pi }^{\pi} e^{i b/2 \, \sin(\phi) + a/2 \, \cos(\phi)} \, \mathrm{d}\phi \\ \stackrel{z=e^{i\phi}}{=} \frac{e^{-a}}{2} \frac{{\rm d}}{{\rm d}a} \, e^{a/2} \oint_{|z|=1} \frac{{\rm d}z}{iz} \, e^{\frac{a+b}{4}z + \frac{a-b}{4z}} = \pi {e^{-a}} \frac{{\rm d}}{{\rm d}a} \, e^{a/2} \, {\rm Res}\left( \frac{e^{\frac{a+b}{4} z}e^{\frac{a-b}{4z}}}{z} \right)\Bigg|_{z=0} \\ =\pi {e^{-a}} \frac{{\rm d}}{{\rm d}a} \, e^{a/2} \, {\rm Res}\left( \sum_{n,m=0}^\infty \frac{\left(\frac{a+b}{4}\right)^n\left(\frac{a-b}{4}\right)^m}{n!m!} \, z^{n-m-1} \right)\Bigg|_{z=0} \\ \stackrel{n=m}{=} \pi {e^{-a}} \frac{{\rm d}}{{\rm d}a} \, e^{a/2} \sum_{m=0}^\infty \frac{\left(\frac{\sqrt{a^2-b^2}}{4}\right)^{2m}}{m!^2} = \pi {e^{-a}} \frac{{\rm d}}{{\rm d}a} \, e^{a/2} I_0\left( \frac{\sqrt{a^2-b^2}}{2} \right) \, .$$

Hence $$I(a)=\pi \int {\rm d}a \, {e^{-a}} \frac{{\rm d}}{{\rm d}a} \, e^{a/2} I_0\left( \frac{\sqrt{a^2-b^2}}{2} \right) \\ \stackrel{\text{P.I.}}{=} \pi {e^{-a/2}} I_0\left( \frac{\sqrt{a^2-b^2}}{2} \right) + \pi \int {\rm d}a \, {e^{-a/2}} I_0\left( \frac{\sqrt{a^2-b^2}}{2} \right) + {\rm const.} \, .$$

In the special case $b=0$, the integral has closed form giving the already known result. The constant is obtained from some limiting case e.g. $a=b$.

Notice that you can easily calculate $$\frac{f(3a,\sqrt{3}b)-f(3a,-\sqrt{3}b)}{\sqrt{3}} = -2\pi b$$ by exploiting the identity $$\int_{-\infty}^\infty e^{itx} \, {\rm d}x = 2\pi \delta(t) \, .$$


I tried various approaches to get a result computationally - without success. Then I decided to iterate (using two nested loops) both variables $a,b$, the first variable $a$ from $0.1$ to $10.1$ by steps of $0.5$ and the second variable $b$ from $-5$ to $5$ by steps of $0.5$. Whithin these loops, I am calculating the integral numerically and plot the surface of this function depending from both variables $a$ and $b$. To get a rough idea of that integral, here is my first picture:

enter image description here

Computation took hours and I need to optimize the code. After this I will be able to expand the iteration range. Hope this helps a bit.

Here is the code, if anyone would like to run it:

from mpmath import mp
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import axes3d

arr_x = np.arange(0.1,10.1,0.5)
arr_y = np.arange(-5,5,0.5)
X,Y = np.meshgrid(arr_x,arr_y)
Z = np.zeros([arr_x.size, arr_y.size])
I = mp.sqrt(-1)

for idx_a, a in enumerate(arr_x):
   for idx_b, b in enumerate(arr_y):
       f = lambda x, t: mp.exp((-a*t**2+I*b*t)/(3*t**2+1)+I*t*x)*(x/(3*t**2+1))
       m = mp.quad(f, [0, mp.inf], [-mp.inf, mp.inf])
       mreal = m.real
       mimag = m.imag
       Z[idx_a][idx_b] = mreal

fig = plt.figure(figsize=(10,6))
ax1 = fig.add_subplot(111, projection='3d')
mycmap = plt.get_cmap('gist_earth')
ax1.set_title('gist_earth color map')
surf1 = ax1.plot_surface(X, Y, Z, cmap=mycmap)
fig.colorbar(surf1, ax=ax1, shrink=0.5, aspect=5)

Thanks also to the Stack Overflow community who helped me.

Update (2021-12-31):

I ammended the code, performed a new run and updated the figure above that contains the plot.

  • 1
    $\begingroup$ $a$ is positive so arr_x should be changed. $\endgroup$
    – AHusain
    Commented Dec 31, 2021 at 16:35
  • $\begingroup$ Good point - thank you for this hint. I fixed it and initiated a new run (it is a bit computationally intensive and will take a while). $\endgroup$ Commented Dec 31, 2021 at 16:56
  • 1
    $\begingroup$ This is the wrong way. A numerical calculation is not of interest. $\endgroup$ Commented Dec 31, 2021 at 20:47

You must log in to answer this question.

Not the answer you're looking for? Browse other questions tagged .