Prove that $ \int_0^1 x^2 \psi(x) \, dx = \ln\left(\frac{A^2}{\sqrt{2\pi}} \right) $ Basically what the title says:

Prove that  $ \displaystyle \int_0^1 x^2 \psi(x) \, dx = \ln\left(\dfrac{A^2}{\sqrt{2\pi}} \right). $

where $A\approx 1.2824$ denotes the Glaisher–Kinkelin constant and $\psi(x) $ denote the digamma function.
I've tried reflection formula and some Feynmann method on the integrals in the wikipedia but nothing works.
(I've copied this from AoPS because there's no response there)
 A: Let $(A_n)$ and $(B_n)$ by
$$ A_n = \frac{1^1 2^2 \cdots n^n}{n^{n^2/2+n/2+1/12} e^{-n^2/4}}
\quad \text{and} \quad
B_n = \frac{n!}{e^{-n}n^{n+1/2}}. $$
In view of the Wikipedia page, we know that $A_n \to A$ and $B_n \to \sqrt{2\pi}$. We can also equivalently state these definitions as follows:
\begin{align*}
\sum_{k=1}^{n} k \log k
&= \log A_n + \left( \frac{n^2}{2} + \frac{n}{2} + \frac{1}{12} \right) \log n - \frac{n^2}{4}, \\
\sum_{k=1}^{n} \log k
&= \log B_n + \left(n + \frac{1}{2}\right) \log n - n.
\end{align*}
Now using the partial fractional decomposition of $\psi(x)$, we find that
\begin{align*}
\int_{0}^{1} x^2 \psi(x) \, dx
&= \lim_{n\to\infty} \int_{0}^{1} x^2 \left( -\gamma + H_{n+1} - \frac{1}{x} - \sum_{k=1}^{n} \frac{1}{k+x} \right) \, dx \\
&= \lim_{n\to\infty} \left( \frac{H_{n+1}-\gamma}{3} + \frac{n^2}{2} - \frac{1}{2} + \sum_{k=1}^{n} k^2 (\log k - \log(k+1)) \right) \\
&= \lim_{n\to\infty} \left( \frac{1}{3}\log n + \frac{n^2}{2} - \frac{1}{2} + \sum_{k=1}^{n} (2k-1)\log k - n^2 \log(n+1) \right) \\
&= \lim_{n\to\infty} \left(2\log A_n - \log B_n - n^2 \log\left(1+\frac{1}{n}\right) + n - \frac{1}{2} \right) \\
&= 2\log A - \log\sqrt{2\pi}.
\end{align*}
A: Integrating by parts twice leads to $$\int x^2 \psi(x) \, dx =x^2 \,\text{log$\Gamma $}(x)-2 \,x\, \psi ^{(-2)}(x)+2\,\psi ^{(-3)}(x)$$ For $x=0$ the expression is $0$. For $x=1$, the first term is $0$ which makes the result to be $$\int_0^1 x^2 \psi(x) \, dx =2\left(\psi ^{(-3)}(1)-\psi ^{(-2)}(1)\right)$$ Now, I am stuck but I found this which reports the value of $\psi ^{(-n)}(1)$ for $n=1,\cdots,10$. Hence the result.
Based on the results given by Wolfram page, we can play with $$I_n=\int_0^1 x^n \psi(x) \, dx$$ and get $$I_1=-\frac{1}{2} \log (2 \pi )$$ $$I_2=2 \log (A)-\frac{1}{2} \log (2 \pi )$$ $$I_3=3 \log (A)-\frac{3 \zeta (3)}{4 \pi ^2}-\frac{1}{2} \log (2 \pi )$$ $$I_4=4 \log (A)-4 \zeta '(-3)-\frac{3 \zeta (3)}{2 \pi ^2}-\frac{11}{180}-\frac{1}{2} \log
   (2 \pi )$$ $$I_5=5 \log (A)-10 \zeta '(-3)+\frac{5 \left(3 \zeta (5)-2 \pi ^2 \zeta (3)\right)}{4 \pi
   ^4}-\frac{11}{72}-\frac{1}{2} \log (2 \pi )$$ $$I_6=6 \log (A)-6 \zeta '(-5)-20 \zeta '(-3)-\frac{15 \left(\pi ^2 \zeta (3)-3 \zeta
   (5)\right)}{4 \pi ^4}-\frac{211}{840}-\frac{1}{2} \log (2 \pi )$$
A: A little late, but here is yet another solution not shown previously.
I will first show that
$$\int_0^1 x \ln\left(\Gamma(x)\right)\,dx=\frac{1}{4} \ln \frac{2 \pi}{A^4} \tag{1}$$
And then from $(1)$ derive the desired result.

Recall Kummer´s fourier expansion for LogGamma $0<x<1$
$$\ln\left(\Gamma(x)\right)=\frac{\ln\left(2 \pi\right)}{2}+ \sum_{k=1}^{\infty}\frac{1}{2k}\cos\left(2 \pi k x\right) + \sum_{k=1}^{\infty}\frac{\gamma+\ln\left(2 \pi k\right)}{ \pi k}\sin\left(2 \pi k x \right) \tag{2}$$
Then, if we plug $(2)$ in $(1)$
$$
\begin{aligned}
\int_0^1 x \ln\left(\Gamma(x)\right)\,dx&=\frac{\ln 2 \pi}{2}\int_0^1 x\,dx+ \sum_{k=1}^{\infty}\frac{1}{2k}\int_0^1 x \cos\left(2 \pi k x\right)\,dx + \sum_{k=1}^{\infty}\frac{\gamma+\ln\left(2 \pi k\right)}{ \pi k}\int_0^1 x \sin\left(2 \pi k x \right)\,dx\\
&=\frac{\ln 2 \pi}{4}+ \sum_{k=1}^{\infty}\frac{\gamma+\ln\left(2 \pi k\right)}{ \pi k}\left(-\frac{1}{2 \pi k}  \right)\\
&=\frac{\ln 2 \pi}{4}-\frac{\gamma}{2 \pi^2}\sum_{k=1}^{\infty}\frac{1}{k^2}-\frac{1}{2\pi^2}\sum_{k=1}^{\infty}\frac{\ln\left(2 \pi k\right)}{k^2}\\
&=\frac{\ln 2 \pi}{4}-\frac{\gamma}{12}-\frac{\ln 2 \pi}{2\pi^2}\sum_{k=1}^{\infty}\frac{1}{k^2}-\frac{1}{2\pi^2}\sum_{k=1}^{\infty}\frac{\ln  k}{k^2}\\
&=\frac{\ln 2 \pi}{4}-\frac{\gamma}{12}-\frac{\ln 2 \pi}{12}-\frac{1}{2\pi^2}\sum_{k=1}^{\infty}\frac{\ln k}{k^2}\\
&=\frac{\ln 2 \pi}{6}-\frac{\gamma}{12}+\frac{1}{2\pi^2}\zeta^\prime(2)\\
&=\frac{\ln 2 \pi}{6}-\frac{\gamma}{12}+\frac{\zeta(2)}{2\pi^2} \left(-12 \ln A + \gamma + \ln2 \pi \right)\\
&=\frac{\ln 2 \pi}{6}- \ln A +\frac{\ln2 \pi}{12} \\
&=\frac{\ln 2 \pi}{4}- \ln A  \\
&=\frac{1}{4} \ln \frac{2 \pi}{A^4}  \qquad \blacksquare\\
\end{aligned}
$$
We used that $\zeta^\prime(2) =\zeta(2)\left(-12 \ln A + \gamma + \ln2 \pi \right)$

On the other hand, if we integrate by parts the L.H.S. of $(1)$ we have that ($du=\psi(x)$ and $v=\frac{x^2}{2}$):
$$\int_0^1 x^2 \psi(x)\,dx =-2\int_0^1 x \ln\left(\Gamma(x)\right)\,dx\tag{3}$$
Which imiediatly gives the desired result
$$\boxed{\int_0^1 x^2 \psi(x) \, dx = \ln\left(\dfrac{A^2}{\sqrt{2\pi}} \right)}$$

Appendix
$$
\begin{aligned}
\int_0^1 x \cos(2 \pi kx)\,dx&=\frac{1}{(2 \pi k)^2}\int_0^{2 \pi k} x \cos(x)\,dx \qquad (2 \pi kx \to x)\\
&=\frac{1}{(2 \pi k)^2}\left(x\sin(x)\Big|_0^{2 \pi k} -\int_0^{2 \pi k}\sin(x)\,dx\right)\\
&=\frac{1}{(2 \pi k)^2}\left(-\cos(x)\Big|_0^{2 \pi k}\right)\\
&=\frac{1}{(2 \pi k)^2}\left(-1+1\right)\\
&=0 \qquad \blacksquare
\end{aligned}
$$

$$
\begin{aligned}
\int_0^1 x \sin(2 \pi kx)\,dx&=\frac{1}{(2 \pi k)^2}\int_0^{2 \pi k} x \sin(x)\,dx \qquad (2 \pi kx \to x)\\
&=\frac{1}{(2 \pi k)^2}\left(-x\cos(x)\Big|_0^{2 \pi k} +\int_0^{2 \pi k}\cos(x)\,dx\right)\\
&=\frac{1}{(2 \pi k)^2}\left(-2 \pi k\right)\\
&=-\frac{1}{2 \pi k}\qquad \blacksquare
\end{aligned}
$$
