I came across a problem on this site (which I now cannot find) asking about the distribution of sums of i.i.d. inverse gamma random variables for paramters $\alpha=\beta=1$. Defining $Z=X+Y$ with $X,Y\sim\operatorname{InvGamma}(1,1)$ we can apply the usual formula to obtain the density
$$
f_Z(z)=\int_0^z(z-x)^{-2}\exp\left(-\tfrac{1}{z-x}\right)x^{-2}\exp\left(-\tfrac{1}{x}\right)\,\mathrm dx,\qquad z>0.
$$
This integral seems deceptively hard which is mainly due to the exponential terms. Now we do have the useful fact that
$$
\int(z-x)^{-2}\exp\left(-\tfrac{1}{z-x}\right)\,\mathrm dx=\exp\left(-\tfrac{1}{z-x}\right),
$$
however, I don't think integrating by parts gets us very far. Another observation is that the integrand is symmetric about $x=z/2$ so that
$$
f_Z(z)=2\int_0^{z/2}(z-x)^{-2}\exp\left(-\tfrac{1}{z-x}\right)x^{-2}\exp\left(-\tfrac{1}{x}\right)\,\mathrm dx.
$$
But again I'm not sure this accomplishes much. How can we deal with the exponential terms to put this integral in a more useful form for evaluating? Taylor series expansions are certainly out. Maybe a clever substitution? The integrand is quite interesting and reminds me of a bump function. Here it is for $z=1$:
Sum of inverse gamma random variables with unity parameters? (challenging integral)
definite integralsexponential functiongamma distributionsmooth-functionsspecial functions
Related Solutions
Let $\mathcal{I}$ denote the value of the definite integral
$$\mathcal{I}:=\int_{0}^{1}\mathrm{d}u\,\sqrt{\frac{u^{2}-2-2\sqrt{u^{4}-u^{2}+1}}{4u^{6}-8u^{4}+8u^{2}-4}}\approx1.5436866339.$$
Note: the denominator of the radicand of the outer square root has the factorization
$$\begin{align} 4u^{6}-8u^{4}+8u^{2}-4 &=4t^{3}-8t^{2}+8t-4;~~~\small{\left[u^{2}=t\right]}\\ &=4\left(t^{3}-2t^{2}+2t-1\right)\\ &=4\left(t-1\right)\left(t^{2}-t+1\right).\\ \end{align}$$
Using the substitution $u^{2}=t$, the integral $\mathcal{I}$ can be rewritten as
$$\begin{align} \mathcal{I} &=\int_{0}^{1}\mathrm{d}u\,\sqrt{\frac{u^{2}-2-2\sqrt{u^{4}-u^{2}+1}}{4u^{6}-8u^{4}+8u^{2}-4}}\\ &=\int_{0}^{1}\mathrm{d}t\,\frac{1}{2\sqrt{t}}\sqrt{\frac{t-2-2\sqrt{t^{2}-t+1}}{4t^{3}-8t^{2}+8t-4}};~~~\small{\left[u=\sqrt{t}\right]}\\ &=\int_{0}^{1}\mathrm{d}t\,\frac{1}{2\sqrt{t}}\sqrt{\frac{t-2-2\sqrt{t^{2}-t+1}}{4\left(t-1\right)\left(t^{2}-t+1\right)}}\\ &=\frac14\int_{0}^{1}\mathrm{d}t\,\sqrt{\frac{2-t+2\sqrt{t^{2}-t+1}}{t\left(1-t\right)\left(t^{2}-t+1\right)}}.\\ \end{align}$$
Now the inner radical is just a square root of a quadratic function, which suggests that the integral might be further simplified using an appropriate Euler substitution.
Consider a substitution given implicitly by the relation
$$\sqrt{t^{2}-t+1}=t+x.$$
Solving for $t$, we obtain $t=\frac{1-x^{2}}{1+2x}$. The integral $\mathcal{I}$ is then transformed to
$$\begin{align} \mathcal{I} &=\frac14\int_{0}^{1}\mathrm{d}t\,\sqrt{\frac{2-t+2\sqrt{t^{2}-t+1}}{t\left(1-t\right)\left(t^{2}-t+1\right)}}\\ &=\frac14\int_{1}^{0}\mathrm{d}x\,\frac{(-2)\left(1+x+x^{2}\right)}{\left(1+2x\right)^{2}}\sqrt{\frac{3\left(1+x\right)}{\left(1-x\right)}\cdot\frac{\left(1+2x\right)}{x\left(2+x\right)}\cdot\frac{\left(1+2x\right)^{2}}{\left(1+x+x^{2}\right)^{2}}};~~~\small{\left[t=\frac{1-x^{2}}{1+2x}\right]}\\ &=\frac12\int_{0}^{1}\mathrm{d}x\,\sqrt{\frac{3\left(1+x\right)}{x\left(1-x\right)\left(2+x\right)\left(1+2x\right)}}\\ &=\frac{\sqrt{3}}{2}\int_{0}^{1}\mathrm{d}x\,\frac{\left(1+x\right)}{\sqrt{x\left(1-x\right)\left(1+x\right)\left(1+2x\right)\left(2+x\right)}}.\\ \end{align}$$
Next, look what happens when we transform the integral using the linear fractional transformation $x=\frac{1-y}{1+y}$:
$$\begin{align} \mathcal{I} &=\frac{\sqrt{3}}{2}\int_{0}^{1}\mathrm{d}x\,\frac{\left(1+x\right)}{\sqrt{x\left(1-x\right)\left(1+x\right)\left(1+2x\right)\left(2+x\right)}}\\ &=\frac{\sqrt{3}}{2}\int_{0}^{1}\mathrm{d}x\,\frac{\left(1+x\right)}{\sqrt{\left(1+x\right)^{6}\left(\frac{1-x}{1+x}\right)\left(\frac{x}{1+x}\right)\left(\frac{1}{1+x}\right)\left(\frac{1+2x}{1+x}\right)\left(\frac{2+x}{1+x}\right)}}\\ &=\frac{\sqrt{3}}{2}\int_{0}^{1}\mathrm{d}x\,\frac{\left(\frac{1}{1+x}\right)^{2}}{\sqrt{\left(\frac{1-x}{1+x}\right)\left(\frac{x}{1+x}\right)\left(\frac{1}{1+x}\right)\left(\frac{1+2x}{1+x}\right)\left(\frac{2+x}{1+x}\right)}}\\ &=\frac{\sqrt{3}}{2}\int_{1}^{0}\mathrm{d}y\,\frac{(-2)}{\left(1+y\right)^{2}}\cdot\frac{\left(\frac{1+y}{2}\right)^{2}}{\sqrt{y\left(\frac{1-y}{2}\right)\left(\frac{1+y}{2}\right)\left(\frac{3-y}{2}\right)\left(\frac{3+y}{2}\right)}};~~~\small{\left[x=\frac{1-y}{1+y}\right]}\\ &=\sqrt{3}\int_{0}^{1}\mathrm{d}y\,\frac{1}{\sqrt{y\left(1-y\right)\left(1+y\right)\left(3-y\right)\left(3+y\right)}}\\ &=\sqrt{3}\int_{0}^{1}\mathrm{d}y\,\frac{1}{\sqrt{y\left(1-y^{2}\right)\left(9-y^{2}\right)}}.\\ \end{align}$$
Recalling Euler's integral formula for the Gauss hypergeometric function: for real argument and parameters,
$$\int_{0}^{1}\mathrm{d}t\,\frac{t^{b-1}\left(1-t\right)^{c-b-1}}{\left(1-zt\right)^{a}}=\operatorname{B}{\left(b,c-b\right)}\,{_2F_1}{\left(a,b;c;z\right)};~~~\small{\left(a,b,c,z\right)\in\mathbb{R}^{4}\land0<b<c\land z<1},$$
(where $\operatorname{B}$ here denotes the usual beta function), we arrive at the following representation for $\mathcal{I}$ as a particular value of ${_2F_1}$:
$$\begin{align} \mathcal{I} &=\sqrt{3}\int_{0}^{1}\mathrm{d}y\,\frac{1}{\sqrt{y\left(1-y^{2}\right)\left(9-y^{2}\right)}}\\ &=\frac{1}{\sqrt{3}}\int_{0}^{1}\mathrm{d}y\,\frac{1}{\sqrt{y}\sqrt{1-y^{2}}\sqrt{1-\frac19y^{2}}}\\ &=\frac{1}{\sqrt{3}}\int_{0}^{1}\mathrm{d}t\,\frac{1}{2\sqrt{t}}\cdot\frac{1}{\sqrt[4]{t}\sqrt{1-t}\sqrt{1-\frac19t}};~~~\small{\left[y=\sqrt{t}\right]}\\ &=\frac{1}{2\sqrt{3}}\int_{0}^{1}\mathrm{d}t\,\frac{t^{b-1}\left(1-t\right)^{c-b-1}}{\left(1-zt\right)^{a}};~~~\small{\left[a:=\frac12,b:=\frac14,c:=\frac34,z:=\frac19\right]}\\ &=\frac{1}{2\sqrt{3}}\operatorname{B}{\left(b,c-b\right)}\,{_2F_1}{\left(a,b;c;z\right)}\\ &=\frac{1}{2\sqrt{3}}\operatorname{B}{\left(\frac14,\frac12\right)}\,{_2F_1}{\left(\frac12,\frac14;\frac34;\frac19\right)}.\\ \end{align}$$
Given $\left(a,b,z\right)\in\mathbb{R}_{>0}\times\mathbb{R}_{>0}\times\left(0,1\right)$, the Gauss hypergeometric function obeys the following two functional relations:
$${_2F_1}{\left(a,b;2b;z\right)}=\left(\frac{1+\sqrt{1-z}}{2}\right)^{-2a}\,{_2F_1}{\left(a,a-b+\frac12;b+\frac12;\left(\frac{1-\sqrt{1-z}}{1+\sqrt{1-z}}\right)^{2}\right)};~~~\small{b<a+\frac12},$$
and
$$\begin{align} {_2F_1}{\left(a,b;\frac12;z\right)} &=\frac{\Gamma{\left(a+\frac12\right)}\,\Gamma{\left(b+\frac12\right)}}{2\,\Gamma{\left(\frac12\right)}\,\Gamma{\left(a+b+\frac12\right)}}\bigg{[}{_2F_1}{\left(2a,2b;a+b+\frac12;\frac{1-\sqrt{z}}{2}\right)}\\ &~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~+{_2F_1}{\left(2a,2b;a+b+\frac12;\frac{1+\sqrt{z}}{2}\right)}\bigg{]}.\\ \end{align}$$
The following pair of identities are then immediate corollaries of the pair above by setting $b=a$: for $0<a\land0<z<1$,
$${_2F_1}{\left(a,\frac12;a+\frac12;z^{2}\right)}=\left(1+z\right)^{-2a}\,{_2F_1}{\left(a,a;2a;\frac{4z}{\left(1+z\right)^{2}}\right)},$$
and
$$\begin{align} {_2F_1}{\left(a,a;\frac12;z\right)} &=\frac{\left[\Gamma{\left(a+\frac12\right)}\right]^{2}}{2\,\Gamma{\left(\frac12\right)}\,\Gamma{\left(2a+\frac12\right)}}\bigg{[}{_2F_1}{\left(2a,2a;2a+\frac12;\frac{1-\sqrt{z}}{2}\right)}\\ &~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~+{_2F_1}{\left(2a,2a;2a+\frac12;\frac{1+\sqrt{z}}{2}\right)}\bigg{]}.\\ \end{align}$$
Continuing with our main evaluation of the integral $\mathcal{I}$, the quadratic transformations of ${_2F_1}$ given above allow us to reduce our integral to standard complete elliptic integrals:
$$\begin{align} \mathcal{I} &=\frac{1}{2\sqrt{3}}\operatorname{B}{\left(\frac14,\frac12\right)}\,{_2F_1}{\left(\frac12,\frac14;\frac34;\frac19\right)}\\ &=\frac{1}{2\sqrt{3}}\operatorname{B}{\left(\frac14,\frac12\right)}\,{_2F_1}{\left(\frac14,\frac12;\frac34;\frac19\right)}\\ &=\frac{1}{2\sqrt{3}}\operatorname{B}{\left(\frac14,\frac12\right)}\,\frac{\sqrt{3}}{2}\,{_2F_1}{\left(\frac14,\frac14;\frac12;\frac34\right)}\\ &=\frac{1}{2\sqrt{3}}\operatorname{B}{\left(\frac14,\frac12\right)}\,\frac{\sqrt{3}}{2}\cdot\frac{\left[\Gamma{\left(\frac34\right)}\right]^{2}}{2\,\Gamma{\left(\frac12\right)}}\bigg{[}{_2F_1}{\left(\frac12,\frac12;1;\frac{2-\sqrt{3}}{4}\right)}+{_2F_1}{\left(\frac12,\frac12;1;\frac{2+\sqrt{3}}{4}\right)}\bigg{]}\\ &=\frac{\Gamma{\left(\frac14\right)}\,\Gamma{\left(\frac34\right)}}{8}\bigg{[}{_2F_1}{\left(\frac12,\frac12;1;\frac{2-\sqrt{3}}{4}\right)}+{_2F_1}{\left(\frac12,\frac12;1;\frac{2+\sqrt{3}}{4}\right)}\bigg{]}\\ &=\frac{\sqrt{2}\,\pi}{8}\bigg{[}{_2F_1}{\left(\frac12,\frac12;1;\frac{2-\sqrt{3}}{4}\right)}+{_2F_1}{\left(\frac12,\frac12;1;\frac{2+\sqrt{3}}{4}\right)}\bigg{]}\\ &=\frac{1}{2\sqrt{2}}\bigg{[}K{\left(\frac{\sqrt{2-\sqrt{3}}}{2}\right)}+K{\left(\frac{\sqrt{2+\sqrt{3}}}{2}\right)}\bigg{]}\\ &=\frac{1}{2\sqrt{2}}\bigg{[}K{\left(\frac{\sqrt{2-\sqrt{3}}}{2}\right)}+K^{\prime}{\left(\frac{\sqrt{2-\sqrt{3}}}{2}\right)}\bigg{]}\\ &=\frac{1}{2\sqrt{2}}\bigg{[}K{\left(\sin{\frac{\pi}{12}}\right)}+K^{\prime}{\left(\sin{\frac{\pi}{12}}\right)}\bigg{]}.\\ \end{align}$$
where here $K{(k)}$ is the complete elliptic integral of the first kind defined as a function of elliptic modulus $k$ by
$$K{(k)}:=\int_{0}^{1}\mathrm{d}x\,\frac{1}{\sqrt{(1-x^{2})(1-k^{2}x^{2})}};~~~\small{-1<k<1},$$
and $K^{\prime}{(k)}$ is the complementary complete elliptic integral of the first kind and is defined in terms of $K$ by
$$K^{\prime}{(k)}:=K{\left(\sqrt{1-k^{2}}\right)}.$$
We can complete our calculation by recognizing that the modulus $k=\sin{\frac{\pi}{12}}$ is in fact the third elliptic integral singular value, $k_{3}$.
Finally, we obtain:
$$\begin{align} \mathcal{I} &=\frac{1}{2\sqrt{2}}\bigg{[}K{\left(\sin{\frac{\pi}{12}}\right)}+K^{\prime}{\left(\sin{\frac{\pi}{12}}\right)}\bigg{]}\\ &=\frac{1}{2\sqrt{2}}\bigg{[}K{\left(k_{3}\right)}+K^{\prime}{\left(k_{3}\right)}\bigg{]}\\ &=\frac{1}{2\sqrt{2}}\bigg{[}1+\frac{K^{\prime}{\left(k_{3}\right)}}{K{\left(k_{3}\right)}}\bigg{]}K{\left(k_{3}\right)}\\ &=\frac{1+\sqrt{3}}{2\sqrt{2}}\,K{\left(k_{3}\right)}\\ &=\frac{1+\sqrt{3}}{2\sqrt{2}}\cdot\frac{\sqrt[4]{3}}{6}\operatorname{B}{\left(\frac12,\frac16\right)}\\ &=\frac{1+\sqrt{3}}{2^{5/2}\,3^{3/4}}\cdot\frac{\Gamma{\left(\frac12\right)}\,\Gamma{\left(\frac16\right)}}{\Gamma{\left(\frac23\right)}}\\ &=\frac{1+\sqrt{3}}{2^{11/6}\,3^{3/4}}\cdot\frac{\pi\,\Gamma{\left(\frac13\right)}}{\left[\Gamma{\left(\frac23\right)}\right]^{2}}\\ &=\frac{\sqrt{3+2\sqrt{3}}}{2^{10/3}}\cdot\frac{\left[\Gamma{\left(\frac13\right)}\right]^{3}}{\pi}.\blacksquare\\ \end{align}$$
The integrand has interesting behavior around $y=x$, so let's split the integral there into two pieces. Let's convert the bottom piece into an integral in polar coordinates:
$$\int_0^{\frac{\pi}{4}} \int_0^{\sec \theta} \left\{\frac{\cos\left(\theta - \frac{\pi}{4}\right)}{\cos\left(\theta + \frac{\pi}{4}\right)}\right\}r\:dr\:d\theta = \frac{1}{2}\int_{\frac{\pi}{4}}^{\frac{\pi}{2}}\{\tan\theta\}\csc^2\theta\:d\theta$$
Now let $u = \cot\theta$ to get
$$\frac{1}{2}\int_0^1 \left\{\frac{1}{u}\right\}\:du = \frac{1-\gamma}{2}$$
which is a somewhat fairly known result, depending on who you ask. The top integral should follow similarly, I'll let you try it on your own (aka left to the reader as an exercise!)
$\textbf{EDIT}$: Given that the inside of the curly braces is negative on the upper piece as @EricTowers noted, we can use
$$\{-x\} = 1 - \{x\}$$
for noninteger $x$ to get that the integral of the upper piece equals
$$\iint_{\text{upper}} \left\{\frac{x+y}{x-y}\right\}\:dx\:dy = \iint_{\text{lower}}1 - \left\{\frac{x+y}{x-y}\right\}\:dx\:dy$$
which implies the sum of the two pieces is simply
$$\iint_{\text{triangle}} dA = \frac{1}{2}$$
Best Answer
Substitute $y=z^2/[4x(z-x)]$. Then $0<x<z/2$ maps to $1<y<\infty$, and $$f_Z(z)=\frac{8}{z^3}\int_1^\infty\left(1-\frac1y\right)^{-1/2}\exp\left(-\frac{4y}{z}\right)\,{\rm d}y=-\frac{4}{z^3}g'\left(\frac2z\right),$$ where $$g(a)=\int_1^\infty\frac{e^{-2ay}\,{\rm d}y}{\sqrt{y(y-1)}}\ \underset{y=(1+t)/2}{\phantom{\Big[}=\phantom{\Big]}}\ \int_1^\infty\frac{e^{-a(1+t)}}{\sqrt{t^2-1}}\,{\rm d}t=e^{-a}K_0(a)$$ uses the modified Bessel function $K_0$. Hence $f_Z(z)=(4/z^3)e^{-2/z}\big(K_0(2/z)+K_1(2/z)\big)$.