Integral $\int_0^\infty\Big[\log\left(1+x^2\right)-\psi\left(1+x^2\right)\Big]dx$ I found this intriguing integral:
$$\int_0^\infty\Big[\log\left(1+x^2\right)-\psi\left(1+x^2\right)\Big]dx\approx0.84767315533332877726581...$$
where $\psi(z)=\partial_z\log\Gamma(z)$ is the digamma.
A lookup in ISC+ (in "advanced" mode) returned a possible closed form $\frac\pi2\Big[\zeta\!\left(\frac12\right)+2\Big]$ that seems to agree with the numeric value of the integral when calculated to a higher precision. Any ideas how to prove it?
 A: This can be done by a kind of residue calculation. First of all by parity the integral is $\frac12 \int_{-\infty}^{\infty}$. Considering this as a contour integral in the $x$-plane, we would like to pull the contour to $i\infty$. There will be two obstacles:


*

*an infinite number of poles of $\psi(x)$ given by $x_n=i\sqrt n$ with $n=1,2,3,\ldots$ Note that 
$$\operatorname{res}_{x=x_n}\psi\left(1+x^2\right)=\color{blue}{\frac{i}{2\sqrt n}}$$

*a logarithmic branch cut $(i,i\infty)$ of the function $\ln(1+x^2)$. Note that on this branch cut the logarithm jumps by $\color{red}{2\pi i}$.
So we deform the contour to the following: a horizontal ray $r_L=(-\infty+i\sqrt N,-0+i\sqrt N)$, a contour $C$ going counterclockwise around the logarithmic branch cut from $-0+i\sqrt N$ to $-0+i\sqrt N$, and another horizontal ray $r_R=(+0+i\sqrt N,\infty+i\sqrt N)$. Now we have 
$$\int_C= \color{red}{2\pi i} \cdot i \left(\sqrt N-1\right)-2\pi i\sum_{1\leq n<N}
\color{blue}{\frac i{2\sqrt n}}=\pi\left[2+\sum_{1\leq n<N}\frac1{\sqrt n}-2\sqrt{N}\right].$$
Using the known asymptotic behavior $\psi(z\to\infty)=\ln z+\frac1{2z}+O\left(z^{-2}\right)$, we can show that integrals over the rays $r_L$ and $r_R$ vanish in the limit $N\to\infty$.
Therefore the integral we are interested in is equal to
$$\mathcal I=\frac12\lim_{N\to\infty}\int_C=\pi+\frac{\pi}{2}\lim_{N\to\infty}\left[\sum_{1\leq n<N}\frac1{\sqrt n}-2\sqrt{N}\right].$$
Finally, that the remaining limit is equal to $\zeta\left(\frac12\right)$ can be shown relatively easily: see, for instance, this approach.
A: There is also a quick derivation available by utilizing Binet's first integral for the digamma function
$$ \psi(z) =\log z -\int_0^{\infty} \left(\frac{1}{e^t-1}+1-\frac{1}{t} \right) e^{-zt} \, dt, \hspace{0.5cm} \text{Re}(z)>0 $$
and the integral
$$ \zeta(s)-\frac{1}{s-1} = \frac{1}{\Gamma(s)} \int_0^{\infty} \left(\frac{1}{e^t-1}+1-\frac{1}{t} \right) t^{s-1} e^{-t} \, dt, \hspace{0.5cm} \text{Re}(s)>0 $$
where setting $s=\frac{1}{2}$ gives
$$ \sqrt{\pi} \left(\zeta \left(\tfrac{1}{2} \right)+2 \right) = \int_0^{\infty} \left(\frac{1}{e^t-1}+1-\frac{1}{t} \right) \frac{e^{-t}}{\sqrt{t}} \, dt $$
Letting $z=1+x^2$ and applying the Fubini theorem gives
$$ \int_0^{\infty} \left(\psi(1+x^2)-\log(1+x^2) \right) dx = -\int_0^{\infty} \int_0^{\infty} \left(\frac{1}{e^t-1}+1-\frac{1}{t} \right) e^{-(1+x^2)t} \, dt \, dx $$
$$ = - \int_0^{\infty} \left(\frac{1}{e^t-1}+1-\frac{1}{t} \right) e^{-t} \int_0^{\infty} e^{-tx^2} \, dx \, dt = - \frac{\sqrt{\pi}}{2} \int_0^{\infty} \left(\frac{1}{e^t-1}+1-\frac{1}{t} \right) \frac{e^{-t}}{\sqrt{t}} , dt = -\frac{\pi}{2} \left(\zeta \left(\tfrac{1}{2} \right)+2 \right) $$
This method admits considerable generalization. For example, we can use it to show that
$$ \int_0^{\infty} x^{r-1} \left(\psi(1+x^s)-\log(1+x^s) \right) \, dx = - \frac{\pi}{s \sin \left(\frac{\pi r}{s} \right)} \left(\zeta \left(1-\tfrac{r}{s} \right) + \frac{s}{r} \right), \hspace{0.5cm} 0<r<s. $$
A: To supplement Dave's awesome answer, I'm adding a proof of the integral $$\zeta(s)-\frac{1}{s-1} = \frac{1}{\Gamma(s)} \int_0^{\infty} \left(\frac{1}{e^x-1}+1-\frac{1}{x} \right) x^{s-1} e^{-x} \, \mathrm dx \ , \quad \Re (s)>0. $$
All we need are the well-known integrals $$\Gamma(s) = \int_{0}^{\infty} x^{s-1}e^{-x} \, \mathrm dx \ ,  \quad \Re(s)>0 \tag{1} $$ and $$\Gamma(s) \zeta(s) = \int_{0}^{\infty} \frac{x^{s-1}}{e^x-1} \, \mathrm dx \ , \quad \Re(s)>1. \tag{2}$$
Subtracting $(1)$ from $(2)$, we get $$ \Gamma(s) \zeta(s) - \Gamma(s)=\int_{0}^{\infty} x^{s-1} \left(\frac{1}{e^{x}-1} - e^{-x} \right) \, \mathrm dx = \int_{0}^{\infty} \frac{x^{s-1}e^{-x}}{e^{x}-1} \, \mathrm dx \ , \quad \Re(s)>1.$$
The integral on the right side of the equation can be rewritten as
$$ \begin{align} \ \int_{0}^{\infty} \frac{x^{s-1}e^{-x}}{e^{x}-1} \, \mathrm dx &= \int_{0}^{\infty} \left(\frac{1}{e^x-1}-\frac{1}{x}  \right) x^{s-1} e^{-x} \,  \mathrm dx+ \int_{0}^{\infty} x^{s-2} e^{-x} \, \mathrm dx  \\ &= \int_{0}^{\infty} \left(\frac{1}{e^x-1}-\frac{1}{x} \right) x^{s-1} e^{-x} \,  \mathrm dx + \Gamma(s-1)  , \end{align} $$  where $\frac{1}{x}$ is the first term of the Laurent series expansion of $\frac{1}{e^x-1}$ at the origin.
Dividing both sides of the equation by $\Gamma(s)$, we get $$\zeta(s) - 1  = \frac{1}{\Gamma(s)}  \int_0^{\infty} \left(\frac{1}{e^x-1}- \frac{1}{x}\right) x^{s-1} e^{-x} \, \mathrm dx + \frac{1}{s-1} , \quad \Re(s) >0. \tag{3}$$
To get the integral in the form given in Dave's answer, add $1$ to both sides of the equation and use the fact that $$1 = \frac{1}{\Gamma(s)} \int_{0}^{\infty} x^{s-1}e^{-x} \, \mathrm dx.$$
