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$ Dec 29, 2021 at 19:21
  • 1
    $\begingroup$ Any restrictions on $b$? Maybe for small $b$? $\endgroup$
    – Diger
    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
    Jan 1 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
    Jan 5 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
    Jan 11 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$ Jan 17 at 19:35
  • $\begingroup$ @Rainer Glüge the integral does not converge on $x \in (-\infty, \infty)$. $\endgroup$
    – KStarGamer
    Jan 17 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$ Jan 19 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
    Jan 19 at 11:20
  • $\begingroup$ @Rainer Glüge: See my answer at the bottom. $\endgroup$
    – Diger
    Jan 20 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
    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$ Dec 31, 2021 at 16:56
  • 1
    $\begingroup$ This is the wrong way. A numerical calculation is not of interest. $\endgroup$ Dec 31, 2021 at 20:47

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.