I am trying to find the closed form for the following integral: $$\int_0^1 \left\{ \frac{m}{x}\right\}^n\,dx, \quad \forall m,n \in \mathbb{N}$$ where $\{x\}$ denotes the fractional part of $x$. This integral is a generalization of some of the fractional part integrals found in Ovidiu Furdui's book Limits, Series, and Fractional Part Integrals: Problems in Mathematical Analysis.
My attempt: $$I(m,n):= \int_0^1 \left\{ \frac{m}{x}\right\}^n\,dx, \quad \forall m,n \in \mathbb{N}.$$
\begin{align*} I(m,n) &\stackrel{x \mapsto \frac{m}{x}}{=} m\int_m^\infty \frac{\{x\}^n}{x^2}\,dx\\ &= m\sum_{k = m}^\infty \int_k^{k+1} \frac{(x-k)^n}{x^2}\,dx\\ &= m\sum_{k = m}^\infty \sum_{\ell=0}^n (-1)^\ell \binom{n}{\ell} k^n\int_k^{k+1} \frac{x^{n - \ell}}{x^2}\,dx\\ &= m\sum_{k = m}^\infty\left[ \sum_{\ell=0}^{n-2}\left[ (-1)^\ell \binom{n}{\ell} k^n\int_k^{k+1} \frac{x^{n - \ell}}{x^2}\,dx\right] + (-1)^{n-1}n k^{n-1} \int_0^1 \frac{1}{x}\,dx + (-1)^n k^n\int_k^{k+1} \frac{1}{x^2}\,dx\right]\\ &= m\sum_{k = m}^\infty\left[ \sum_{\ell=0}^{n-2}\left[ (-1)^\ell \binom{n}{\ell} k^n\frac{(k+1)^{n-\ell-1}-k^{n-\ell-1}}{n-\ell-1}\right] + (-1)^{n-1}n k^{n-1} \log\left(\frac{k+1}{k} \right) + (-1)^n \frac{k^{n-1}}{k+1}\right]\\ &= \lim_{N \to \infty} m\sum_{k = m}^N \left[ \sum_{\ell=0}^{n-2}\left[ (-1)^\ell \binom{n}{\ell} k^n\frac{(k+1)^{n-\ell-1}-k^{n-\ell-1}}{n-\ell-1} \right] + (-1)^{n-1}n k^{n-1} \log\left(\frac{k+1}{k} \right) + (-1)^n \frac{k^{n-1}}{k+1}\right]\\ &= \lim_{N \to \infty} m\left[ \sum_{\ell=0}^{n-2}\sum_{k = m}^N\left[ (-1)^\ell \binom{n}{\ell} k^n\frac{(k+1)^{n-\ell-1}-k^{n-\ell-1}}{n-\ell-1} \right] + (-1)^{n-1}n \sum_{k = m}^N \log\left[\left(\frac{k+1}{k} \right)^{k^{n-1}}\right] + (-1)^n \sum_{k = m}^N \frac{k^{n-1}}{k+1}\right]. \end{align*}
From here I began to struggle. I believe I I was able to make some good progress with the sum with the logarithm. To do so I defined something I will call the generalized hyperfactorials (I do not know if such a function exists in literature): $$\mathrm{H}(k,n) := \prod_{j = 1}^n j^{j^k}.$$ We see that if $k = 0$, we get the regular factorial. If we let $k = 1$, we get the hyperfactorial. With this function, I believe I the partial sum with the logarithms has a closed form: \begin{align*} \sum_{k = m}^N \log\left[\left(\frac{k+1}{k} \right)^{k^{n-1}}\right] &= \log\left[\prod_{k = m}^N\left(\frac{k+1}{k} \right)^{k^{n-1}}\right]\\ &= \log\left[\frac{(m+1)^{m^{n-1}} \times (m+2)^{(m+1)^{n-1}} \times \ldots \times N^{(N-1)^{n-1}} \times (N+1)^{N^{n-1}}}{m^{m^{n-1}} \times (m+1)^{(m+1)^{n-1}} \times (m+2)^{(m+2)^{n-1}} \times \ldots \times N^{N^{n-1}}}\right]\\ &= \log\left[\frac{(N+1)^{N^{n-1}}}{m^{m^{n-1}}} \left((m+1)^{m^{n-1} - (m+1)^{n-1}} \times (m+2)^{(m+1)^{n-1} - (m+2)^{n-1}} \times N^{(N - 1)^{n-1} - N^{n-1}}\right)\right]\\ &= \log\left[\frac{(N+1)^{N^{n-1}}}{m^{m^{n-1}}} \left((m+1)^{((m+1) - 1)^{n-1} - (m+1)^{n-1}} \times (m+2)^{((m+2) - 1)^{n-1} - (m+2)^{n-1}} \times N^{(N - 1)^{n-1} - N^{n-1}}\right)\right]\\ &= \log\left[\frac{(N+1)^{N^{n-1}}}{m^{m^{n-1}}} \left((m+1)^{\sum_{p=1}^{n-1} (-1)^p \binom{n-1}{p}(m+1)^{n - p -1}} \times (m+2)^{\sum_{p=1}^{n-1} (-1)^p \binom{n-1}{p}(m+2)^{n - p -1}} \times \ldots \times N^{\sum_{p=1}^{n-1} (-1)^p \binom{n-1}{p} N^{n - p -1}}\right)\right]\\ &= \log\left[\frac{(N+1)^{N^{n-1}}}{m^{m^{n-1}}} \prod_{p=1}^{n-1}\left((m+1)^{(m+1)^{n - p -1}} \times (m+2)^{(m+2)^{n - p -1}} \times \ldots \times N^{ N^{n - p -1}}\right)^{(-1)^p \binom{n-1}{p}}\right]\\ &= \log\left[\frac{(N+1)^{N^{n-1}}}{m^{m^{n-1}}} \prod_{p=1}^{n-1}\left(\frac{\mathrm{H}(n -p-1,N)}{\mathrm{H}(n -p-1,m)}\right)^{(-1)^p \binom{n-1}{p}}\right]. \end{align*}
Hopefully I did not make any typos. I believe the integrals has a strong connection to the Bendersky-Adamchik constants. Here are some papers on the constants
Some new quicker approximations of Glaisher–Kinkelin's and Bendersky–Adamchik's constants
Closed-form calculation of infinite products of Glaisher-type related to Dirichlet series
From these papers, we can use the limit definition of the Bendersky-Adamchik constants to derive Stirling-like asymptotic formulas for the generalized hyperfactorials, which I conjecture will be handy when evaluating the limit. I also believe that the other sums left to be simplified could be simplified via the Euler-Maclaurin summation formula or even Faulhaber's formula somewhere down the line. Possibly the Bernoulli numbers from the definition of the limit definition of the Bendersky-Adamchik constants interact with Bernoulli numbers from the Euler-Maclaurin summation formula. Is it possible that I am approaching this integral the wrong way? If so, what would be the best way to tackle this integral?
We have \begin{align*} I(m, n) & = m\sum\limits_{k = m}^\infty {\int_k^{k + 1} {\frac{{(x - k)^n }}{{x^2 }}{\rm d}x} } = m\sum\limits_{k = m}^\infty {\int_0^1 {\frac{{s^n }}{{(s + k)^2 }}{\rm d}s} } \\ & = m\int_0^1 {s^n \sum\limits_{k = m}^\infty {\frac{1}{{(s + k)^2 }}} {\rm d}s} = m\int_0^1 {s^n \psi '(m + s){\rm d}s}, \end{align*} where $\psi'$ is the derivative of the digamma function.
Integrating by parts $n+1$ times gives $$ I(m, n) = ( - 1)^n n!m\left[ {\psi ^{( - n)} (m + 1) - \psi ^{( - n)} (m)} \right] + m\sum\limits_{k = 0}^{n - 1} {( - 1)^k \frac{{n!}}{{(n - k)!}}\psi ^{( - k)} (m + 1)} $$ where $\psi^{(-k)}$ is a polygamma function of negative order. Their values may be expressed in terms of Bernoulli polynomials, harmonic numbers, derivatives of the Riemann zeta function and the Hurwitz zeta function, see Proposition $2$ in the above linked paper.
Alternatively, integration by parts once yields $$ I(m, n) = m\psi (m + 1) - nm\int_0^1 {s^{n - 1} \psi (m + s){\rm d}s} $$ where $\psi$ is the digamma function. Then $$ I(m, 1) = m\psi (m + 1) - nm\log m, $$ and, using the binomial series, $$ I(m, n) = m\psi (m + 1) + n\sum\limits_{k = 0}^{n - 1} {\binom{n-1}{k}( - m)^{n - k} \left[ {\int_0^{m + 1} {s^k \psi (s){\rm d}s} - \int_0^m {s^k \psi (s){\rm d}s} } \right]} $$ for $n\ge 2$. The integrals may be expressed again in terms of Bernoulli polynomials, harmonic numbers, derivatives of the Riemann zeta function and the Hurwitz zeta function, see Proposition $3$ in the above linked paper. The result should be identical to that arising from the first approach.