---
abstract: 'Following our earlier research, we use the method introduced by the author in [@prevost1996] named Remainder Padé Approximant in [@rivoalprevost], to construct approximations of the Hurwitz zeta function. We prove that these approximations are convergent on the positive real line. Applications to new rational approximations of $\ze(2)$ and $\ze(3)$ are given.'
address: 'LMPA Joseph Liouville, Centre Universitaire de la Mi-Voix, Bât H. Poincaré, 50 rue F. Buisson, BP 699, 62228 Calais Cedex'
author:
- Marc Prévost
bibliography:
- 'bibfile3.bib'
title: Remainder Padé approximants for the Hurwitz zeta function
---

Introduction
============

In [@prevost1996], we gave a new proof of the irrationality of $\zeta(2)=\sum_{k=1}^{\infty}1/k^2$ (and also of $\zeta(3)=\sum_{k=1}^{\infty}1/k^3$) based on an explicit construction of some Padé approximants of the remainder term $R_2(1/n)=\sum_{k=n}^{\infty}1/k^2$. Pore precisely, we have that $\zeta(2)=\sum_{k=1}^{n-1}1/k^2+R_2(1/n)$ with $$R_2(z)=\sum_{k=0}^{\infty}\dfrac{z^2}{(z k+1)^2}.$$ The function $R_2$ is meromorphic in $\C\setminus \{0,-1,-\dfrac{1}{2},-\dfrac{1}{3},\cdots\}$ and so is not holomorphic at $0$.

However, it is $\cC^\infty$ at $z=0$ and admits a Taylor series $\hat R_2(z)=\sum_{k=0}^{\infty}B_k z^{k+1}$ which is an asymptotic development convergent only at $z=0$, (the radius of convergence is $0$), where $B_k$ is the $k$th Bernoulli number.

The idea was to compute an explicit representation of diagonal Padé approximant $P_n(z)/Q_n(z)=[n/n]_ {\hat R_2}(z)
\in \Q(z)$ of the series $\hat R_2(z)$ with a good estimation of the error term $E_n(z)=R_2(z)-[n/n]_ {\hat R_2}(z).$

At last, we find that $$Q_n(1/n)E_n(1/n)=Q_n(1/n)\ze(2)-Q_n(1/n)\sum_{k=1}^{n-1}1/k^2-P_n(1/n)=q_n\ze(2)-p_n$$ provides a sequence of rational approximation $P_n/q_n$ which proves tje irrationlaity of $\ze(2).$ The surprise was to find that it is exactly the sequences used by Apéry [@apery] for the same purpose.

The same method applied to the series $\zeta(3)=\sum_{k=1}^{\infty}1/k^3$ provides also the Apéry numbers given in [@apery].

The name Remainder Padé Approximant (RPA) has been used later in a paper written with T. Rivoal [@rivoalprevost] for the exponential function.

In this paper, we apply this method and prove that some sequence of RPA for Hurwitz zeta function is convergent on the real line.

For $a\in \C, \Re(a)>0$ , the Hurwitz zeta function is defined as

$$\begin{aligned}
\ze (s,a)& = &\dis{\sum ^\ii _{k=0}} \dis{\frac{1}{(k+a)^s}}
\label{1}\end{aligned}$$

where the Dirichlet series on the right hand side of (\[1\]) is convergent for $\Re(s) > 1$ and uniformly convergent in any finite region where $\Re(s) \geq 1 + \de $, with $\de  > 0$. It defines an analytic function for $\Re(s) > 1$. The Riemann $\ze$ function is $\ze(s,1)$.

An integral representation of $\ze(s,a)$ is $$\begin{aligned}
\ze (s,a) = \dis{\frac{1}{\Ga(s)}} \int ^\infty _0
\displaystyle{\frac{x^{s-1}e^{-a x}}{1-e^{-x}}}dx \quad \Re(s) > 1,\Re(a)>0,\end{aligned}$$

[where]{} $\Ga (s) =  \int ^\ii_0 y^{s-1} e^{-y}dy$ is the gamma function.

If $\Re(s) \leq 1$ the Hurwitz zeta function can be defined by the equation : $$\ze (s,a)=\Ga(1-s)\frac{1}{2 \pi i} \int_\cC \frac{z^{s-1}e^{a z}}{1-e^z}dz$$

where ${\cC}$ is some path in $\C$, which provides the analytic continuation of $\ze (s,a)$ over the whole $s$-plane.

If we write the formula (\[1\]) as :

$$\begin{aligned}
\ze (s,a) = \dis{\sum ^{n-1}_{k=0} \frac{1}{(k+a)^s}}+
\dis{\sum ^\ii _{k=0}} \dis{\frac{1}{(n+a+k)^s}}
 \end{aligned}$$

and set $$\Psi(s,x) : =  \dis t^{1-s}\sum ^\ii_{k=0}\({\frac{t}{1+kt}}\)^s$$

then $$\begin{aligned}
 \ze (s,a) = \dis{\sum ^{n-1} _{k=0}}
\dis{\frac{1}{(k+a)^s}} + (n+a)^{1-s}
\Psi\(s,\frac{1}{n+a}\).
\label{formzeta}
 \end{aligned}$$

Another expression of $\Psi(s,t)$ is [@watson-witt chapt. XIII]

$$\begin{aligned}
\Psi(s,t) = \dis{\frac{t^{1-s}}{\Ga (s)}}\int ^\ii _0 u^{s-1}\dis{
\frac{e^{-u/t}}{1-e^{-u}}}du\end{aligned}$$

whose Taylor series is

$$\begin{aligned}
\Psi(s,t) = \dis{\sum ^\ii_{k=0}}
\dis{\frac{B_k}{k!}} (s)_{ k-1} (-1)^k t^{k},\label{DAS}
 \end{aligned}$$

convergent only at $t=0$, where $B_k$ is the $k$th Bernoulli number and $(s)_k$ is the Pochammer symbol ($(s)_p:=s(s+1)(s+2)\cdots(s+p-1),(s)_{-1}=1/(s-1)$).

We replace in the remainder term in (\[formzeta\]), $\dis\Psi\(s,\frac{1}{n+a}\)$ by some Padé approximant $[m_1/m_2]_{\Psi(s,t)}=\widetilde
 Q_{m_1}(t)/\widetilde
 P_{m_2}(t)$ to the function $\Psi(s,t)
$ computed at $ t={1}/(n+a)$. We obtain the RPA approximant $$\begin{aligned}
  RPA{(n,m_1,m_2)}:={\sum ^{n-1} _{k=0}}
{\frac{1}{(k+a)^s}} + (n+a)^{1-s}\dfrac{\widetilde Q_{m_1}({1}/{(n+a)})}{\widetilde P_{m_2}({1}/{(n+a)})}.\end{aligned}$$

The crucial point is the convergence of this sequence $RPA{(n,m_1,m_2)}$ when one or more of the parameter $n,m_1,m_2$ tends to infinity.\
The convergence of Padé approximants is proved for many functions as meromorphic function, Stieltjes function, etc...([@bgm1996]). In this paper, we prove that the function $\Psi\(s,t\)$ is a Stieltjes function in the variable $t$ when $s$ is a real positive number. First we have to recall the definition and some properties of Padé approximants

Padé approximants
=================

Let $f$ be a function whose Taylor expansion about $t = 0$ is

$$f(t)=\sum ^\infty_{i=0} c_i  t^i.$$

The Padé approximant $[m_1/m_2]_f$ to $f$ is a rational fraction $\widetilde
 Q_{m_1}/\widetilde P_{m_2}$ whose denominator has degree $m_2$ and whose numerator has degree $m_1$ so that its expansion in ascending powers of $t$ coincides with the expansion of $f$ up to the degree $m_1+m_2$ , i.e.

$$\begin{aligned}
 \widetilde P_{m_2}(t) f(t) - \widetilde Q_{m_1}(t) = 0(t^{m_1+m_2+1})\;,\quad t
 \rightarrow 0 \end{aligned}$$

The theory of Padé approximation is linked with the theory of orthogonal polynomials as following [@PTA]: let us define the linear functional $c$ acting on ${\mathcal P}$, space of polynomials by

$$\begin{aligned}
 &c : &{\mathcal P}\rightarrow \R  (or \; \C )\\
 &&x^i \rightarrow <c, x^i> = c_i \quad i = 0, 1, 2, \ldots
 \end{aligned}$$

and

$$\begin{aligned}
 \mathrm{if} \; p \in {\Z}\;\;  <c^{(p)}, x^i> := <c,x^{i+p}> = c_{i+p} \quad i = 0, 1, 2,
 \ldots (c_i=0, i<0)
 \end{aligned}$$

then the denominators of the Padé approximants $[m_1/m_2]_f$ satisfy the following orthogonality property :

$$<c^{(m_1-m_2+1)}, x^i P_{m_2}(x)> = 0 \quad i = 0,1,2,\ldots,n-1$$ where $P_{m_2}(x) = x^{m_2} \widetilde P_{m_2}(x^{-1})$.

Let us define the associated polynomials :

$$\begin{aligned}
 R_{m_2-1}(t) := <c^{(m_1-m_2+1)},\frac{P_{m_2}(x)-P_{m_2}(t)}{x-t}>, R_{{m_2}-1}\in {\mathcal P}_{{m_2}-1}\label{associe}
 \end{aligned}$$

where $c^{(m_1-{m_2}+1)}$ acts on the variable $x$ . Then

$$\begin{aligned}
 \widetilde Q_{m_1}(t) =\( {\sum ^{{m_2}-m_1}_{i=0}} c_it^i\)
 \widetilde P_{m_2}(t) + t^{m_1-{m_2}+1}\widetilde R_{{m_2}-1}(t)
 \end{aligned}$$

where $\widetilde R_{{m_2}-1}(t) = t^{{m_2}-1} R_{{m_2}-1}(t^{-1})$ and $c_j = 0$ for $j < 0.$

If $c$ admits an integral representation with a function $\alpha $ non decreasing, with bounded variation,

$$\begin{aligned}
 c_i = \int _{L} x^i d\alpha (x),\label{integralrepr}
 \end{aligned}$$

then the function $f$ is a Stieltjes function and the theory of Gaussian quadrature shows that the polynomials $P_n$ orthogonal with respect to $c$, have all their roots in the support of the function $\alpha $. Moreover, the convergence is proved if the coefficients $c_i$ satisfy the Carleman condition [@bgm1996 p.240].

The aim of the paper is to find, in the case of Hurwitz zeta function, the weight function $d\alpha$ depending on $s$ and prove that it is a positive function.

The error is defined by $$\begin{aligned}
 \mbox{\underline{error}}\; : f(t) - [m_1/m_2]_f(t) =
 \displaystyle{\frac{t^{m_1+m_2+1}}{\widetilde P^2_{m_2}(t)}}c^{(m_1-{m_2}+1)} \left (
 \displaystyle{\frac{P^2_{m_2}(x)}{1-xt}}\right ).
 \end{aligned}$$

The above expression of the error is understood as a formal one if $c$ is only a formal linear functional [@PTA chapt. 3], but if $c$ admits the integral representation (\[integralrepr\]) then the error becomes :

$$\begin{aligned}
 f(t) - [m_1/m_2]_f(t) = f(t) - {\frac{\widetilde Q_{m_1}(t)}{\widetilde P_{m_2}(t)}} = {\frac{t^{m_1+{m_2}+1}}{\widetilde P^2_{m_2}(t)}} 
 \int _{L} 
 x^{m_1-m_2+1} \frac{P^2_{m_2}(x)}{1-x t}d\alpha (x)
 \end{aligned}$$

In the particular case $m_1 = m_2-1$,

$$\begin{aligned}
 f(t) - [m_2-1/m_2]_f(t) = \displaystyle{\frac{t^{2m_2}}{\widetilde
        P^2_{m_2}(t)}}\int _{L}
 {\frac{P^2_{m_2}(x)}{1-xt}} d\alpha (x).\label{errorpade}
 \end{aligned}$$

Note that

$$[{m_2}/{m_2}]_f(t) = c_0+t [{m_2}-1/{m_2}]_{f_1}$$

and

$$[m_2+p/m_2]_f(t) = c_0+c_1 t+\cdots +c_{p} t^{p}+t^{p+1} [m_2-1/m_2]_{f_p}(t)$$

where $$f_p(t)=\sum ^\infty_{i=0} c_{i+p+1 }  t^i.$$

The orthogonal polynomial can be expressed with determinant:

Let $(e_k)_{k\in \N}$ a basis of the space of polynomials. If $$\left\langle c,e_i\;e_j\right\rangle =\alpha_{i,j}$$ then $$\begin{aligned}
 P_n(x)&=&k_n\left|  \begin{array}{cccc}
\alpha_{0,0}&\alpha_{1,0}&\cdots&\alpha_{n,0}\\
\alpha_{0,1}&\alpha_{1,1}&\cdots&\alpha_{n,1}\\
\hdots&\hdots&\hdots&\hdots\\
\alpha_{0,n-1}&\alpha_{1,n-1}&\cdots&\alpha_{n,n-1}\\
e_0(x)&e_1(x)&\cdots&e_n(x)
 \end{array}\right|  \label{hadamard}
 \end{aligned}$$ where $k_n$ is some constant.

This property will be used in section \[particular case\].

Statements of the results
=========================

In this section, we will prove that it is possible to compute the weight function underlying the coefficients $\frac{B_k}{k!} (s)_{ k-1} (-1)^k$ of the function $\Psi(s,t) = {\sum ^\ii_{k=0}}
 {\frac{B_k}{k!}} (s)_{ k-1} (-1)^k t^{k}.$ The weight function $w_s$ depends on $s$. Thus, for particular values of $s$, $\Psi(s,t)$ is a Stieltjes function in the variable $t$ and the convergence of RPA will be proved for $t=1/(n+a)$.

For any complex numbers $ s,a $, $\Re(s)>-1,\Re(a)> 0$ and any positive integer $m$, $m>\Re(s)-1$, then

$$\begin{aligned}
\zeta(s,a)=\frac{1}{ a^{s-1}}\left(\frac{1}{s-1}+\frac{1}{2 a}+
\int_0^\infty
\frac{  1}{a^2+ x^2} w_s(x)dx\right)\label{zeta}\end{aligned}$$

where the weight $w_s$ is defined by:

$$w_s(x):=\frac{2 (-1)^{m}x^s}{\Gamma(s)\Gamma(m+1-s)} \int_x^\infty (t-x)^{m-s}\frac{d^m}{dt^m}\left(\frac{1}
{e^{2\pi t}-1}\right)dt.$$

\[theo1\]

We will prove in the section \[section3\] that this formula (\[zeta\]) is a consequence of Hermite’s formula for the function $\ze(s,a)$.

As found by Touchard [@touch], Bernoulli numbers satisfy $$B_k=-i \frac{\pi}{2}\int_L x^k \frac{dx}{\sin^2(\pi x)},$$ where $L$ is the line $L:=-\frac{1}{2}+i \R$ and thus can be viewed as moments of order $k\geq 0$ of the positive weight function $1/{\sin^2(\pi x)} $ on the line $L:=-\frac{1}{2}+i \R$: So for $s=2$, the coefficients ${\frac{B_k}{k!}} (s)_{ k-1} =B_k$ of (\[DAS\]) are moments of a positive weight. For $s=3$, these coefficients appear as their derivative $(k+1)B_k$. The derivative of $1/{\sin^2(\pi x)}$ equals is a weight function symmetric around $-1/2$ whose support is also the line $L$. But for integer value of $s$ greater than 4, it is no more true for all $k\geq 0$. A much more difficult case is the case $s$ is real.

From the previous Theorem, we obtain an integral representation of $\dis{\frac{B_k}{k!}} (s)_{ k-1} $ but for $k\geq 2$.

Integral representation of Bernoulli numbers.\[theo2\]

If $s$ a complex number such that $\Re(s)>-1$, then

$$\begin{aligned}
B_{k+2}\frac{(s)_{k+1}}{(k+2)!}&=&(-1)^{k/2}
\int_{0}^{+\infty}
x^{k} w_s(x)dx \hskip 1 cm (k \;even)  \label{bernoulli}\\
&=& 0  \hskip 5 cm\;(k\; odd)\end{aligned}$$

with

$$\dis w_s(x):=\frac{2 (-1)^{m}x^s}{\Gamma(s)\Gamma(m+1-s)} \int_x^\infty (t-x)^{m-s}\frac{d^m}{dt^m}\left(\frac{1}{e^{2\pi t}-1}\right)dt$$

[where $m$ is any integer satisfying]{} $ m > \Re(s)-1.$

\[teo2\]

**Proof** We use the expansion of $\zeta(s,a)$ in terms of Bernoulli numbers [@srivastava-choi p.160]. $$\zeta(s,a)=a^{-s}+\sum_{k=0}^n (s)_{k-1}\frac{B_k}{k!}a^{-k-s+1}+\frac{1}{\Gamma(s)}\int_0^\infty \left(\frac{1}{e^t-1}-\sum_{k=0}^n \frac{B_k}{k!}t^{k-1}\right)e^{-a\;t}\;t^{s-1}\;dt$$ valid for $\Re(s)>-2n+1,\Re(a)>0$ and $n$ a positive integer.

By identification with formula (\[zeta\]), we get the formula (\[bernoulli\]).

[**Remark**]{}: The positivity of the weight function is important because it will imply the convergence of Padé approximants. In the case where $s$ is a positive real number then $w_s$ is positive on its support. This gives the following main theorem.

For all positive real number $s$, for all complex number $a, \Re(a)>0$, for all integers $n\geq 0, p\geq 1$, the following sequence

$$\begin{aligned}
RPA{(n,m+p,m)}:={\sum ^{n-1} _{k=0}}
    {\frac{1}{(k+a)^s}}  +
    (n+a)^{1-s}\restriction{[m+p/m]_{\Psi(s,t)}}{t=1/(n+a)}  .\label{RPA approx}
    \end{aligned}$$

converges to $\ze(s,a)$ when $m$ tend to infinity.\[teo4\]

These Theorems will be proved in Section \[section3\].

For $n=0,p=-1$, this result is proved only for $s=2, 3$ in [@prevost1996].

If $n=0$, then formula (\[RPA approx\]) reduces to

$$\forall s >0 ,\forall a \in \C, \Re(a)>0, \forall p\geq 1,  \lim_{m\rightarrow \infty} a^{1-s}\restriction{[m+p/m]_{\Psi(s,t)}}{t=1/a}=\ze(s,a).$$

Moreover, for the particular case $n=0$, $a=1$, then $RPA{(0,m+p,m)}$ is a rational function depending on the variable $s$ and which converges to $\ze(s,1)=\ze(s)$ when $m$ tends to $\ii$ (by Theorem \[teo4\]).

For example, for $p=1$,\
$RPA(0,2,1)=\frac{(s+2) (s+3)}{12 (s-1)} ,\\RPA(0,3,2)=\frac{(s+2) \left(s^2+12 s+31\right)}{2 (s-1) \left(s^2+3 s+62\right)},\\RPA(0,4,3)=\frac{(s+2) (s+4) \left(3 s^2+88 s+405\right)}{120 (s-1) \left(s^2+7 s+54\right)},\\RPA(0,5,4)=\frac{(s+2) (s+4) \left(s^5+124 s^4+2644 s^3+23730 s^2+92939 s+122130\right)}{2 (s-1) \left(s^6+101 s^5+1755 s^4+21415 s^3+153884 s^2+657564 s+977040\right)} \\
RPA(0,6,5)=\frac{(s+2) (s+4) (s+6) \left(s^5+746 s^4+28162 s^3+498112 s^2+3612925 s+8457750\right)}{84 (s-1) \left(s^6+429 s^5+10387 s^4+134511 s^3+1044772 s^2+4891020 s+9666000\right)}$

If $n\geq 1$ and $a=1$, then $RPA(n,m+p,m)$ is no more a rational fraction since it contains powers of $(n+a)$ with $s$ as exponent:\
$RPA(1,2,1)=\frac{2^{1-s} \left(s^2+11 s+36\right)}{48 (s-1)}+1,\\
RPA(1,3,2)=\frac{2^{-s-1} \left(s^3+26 s^2+231 s+726\right)}{(s-1) \left(s^2+3 s+242\right)}+1,\\
RPA(1,4,3)=\frac{2^{-s-4} \left(3 s^4+166 s^3+2937 s^2+22334 s+64800\right)}{15 (s-1) \left(s^2+7 s+180\right)}+1,
\\
RPA(2,2,1)=\frac{3^{1-s} \left(s^2+17 s+90\right)}{108 (s-1)}+2^{-s}+1,\\
RPA(2,3,2)=\frac{3^{-s} \left(s^3+38 s^2+527 s+2710\right)}{2 (s-1) \left(s^2+3 s+542\right)}+2^{-s}+1.$

Nevertheless, if $s$ is an integer, if $a=1$ and $p\geq 1$, we obtain a sequence of rational approximations of $\ze(s)$ by computing $RPA(n,m+p,m)$, depending on the two integers parameters $n$ and $m$.

Proof of Theorems {#section3}
=================

Proof of Theorem \[theo1\]
--------------------------

We start from Hermite’s formula for $\zeta(s,a),\Re(a)>0,$ which is a consequence of Plana’s summation formula: $$\zeta(s,a)=\frac{1}{2}\;a^{-s}+\frac{a^{1-s}}{s-1}+2\int_0^\infty
(a^2+y^2)^{-s/2}\sin\left(s \arctan{\frac{y}{a}}\right)\frac{dy    }{e^{2\pi y}-1}.
\label{hermite}$$ Note that this integral converges for all complex number $s\neq1$.

Let us define $$I_s:= \int _0^t\frac{x^s(t-x)^{-s}}{1+x^2}dx \;\;\;(t\geq 0,-1<\Re(s)<1).$$

Using the identity $$\int_0^\infty\frac{v^{s-1}}{1+\alpha\; v}={\alpha^{-s}}{\Gamma(s)\Gamma(1-s)}$$ we get, with the change of variable $x=\dis\frac{t\;v}{1+v}$,

$$\begin{aligned}
I_s&=&t\int_0^\infty \frac{v^s}{(1+v)^2+t^2v^2}dv\\&=&\frac{1}{2i}\int_0^\infty v^{s-1}\left( \frac{1    }{1+v(1-i t)}-
\frac{1    }{1+v(1+i t)}\right)dv\\&=&
\Gamma(s)\Gamma(1-s)
\frac{1}{2i}\left((1-i t)^{-s}-(1+i t)^{-s}\right).\end{aligned}$$

Setting $t=\tan \theta $, we obtain $$I_s=\Gamma(s)\Gamma(1-s)\cos^s( \theta) \sin(s \theta)=\Gamma(s)\Gamma(1-s)(1+t^2)^{-s/2}\sin (s \arctan t).$$

Thus, the following identity $$\frac{1}{\Gamma(s)\Gamma(1-s)}\int _0^t\frac{x^s(t-x)^{-s}}{1+x^2}dx=(1+t^2)^{-s/2}\sin (s \arctan t)$$

holds for all real $t\geq 0$ and complex $s$, $-1<\Re(s)<1$.

Now, for $s$ complex satisfying $\Re(s)>-1$ and $m$ an integer such that $m>\Re(s)-1$, by recurrence on $m$, it is easy to prove the following formula,

$$\frac{\sin(s \arctan t)}{(1+t^2)^{s/2}}=\frac{1}{\Gamma(s)\Gamma(m+1-s)}\frac{d^m}{dt^m}\left(\int_0^t
\frac{x^s(t-x)^{m-s}    }{1+x^2}\;dx\right).
\label{sin}$$

To prove Theorem \[theo1\], we replace $\displaystyle \frac{\sin(s \arctan (t/a))}{(a^2+t^2)^{s/2}}$ in Hermite’s formula (\[hermite\]), then apply $m$ integration by parts and permutation of the integrals.

$$\begin{aligned}
\zeta(s,a)&=&\!\!\!\frac{a^{-s}}{2}\!+\!\frac{a^{1-s}}{s\!-\!1}\!+\!\frac{2\;a^{1-s}}{\Gamma(s)\Gamma(m\!+\!1\!-\!s)}\int_0^\infty
\frac{d^m}{dt^m}\left(\int_0^t
\frac{x^s(t\!-\!x)^{m-s}}{a^2+x^2}dx\right)\frac{dt    }{e^{2\pi t}\!-\!1}\\
&=&\!\!\!\frac{a^{-s}}{2}\!+\!\frac{a^{1-s}}{s\!-\!1}\!+\!\frac{2(-1)^ma^{1-s}}{\Gamma(s)\Gamma(m\!+\!1\!-s\!)}\int_0^\infty\!\!\!
\left(\int_0^t
\!\!\frac{x^s(t\!-\!x)^{m-s}}{a^2\!+\!x^2}dx\right)\!\frac{d^m}{dt^m}\!\!\left(\!\!\frac{1 }{e^{2\pi t}\!-\!1}\!\!\right)dt\\
&=&\!\!\!
\frac{a^{-s}}{2}\!+\!\frac{a^{1-s}}{s-1}\!+\!\frac{2(-1)^ma^{1-s}}{\Gamma(s)\Gamma(m\!+\!1\!-\!s)}\int_0^\infty\!\!
\!\!\!\frac{x^sdx}{a^2\!+\!x^2}
\int_x^\infty\!\!(t\!-\!x)^{m-s}
\frac{d^m}{dt^m}\!\!\!\left(\!\!\frac{1    }{e^{2\pi t}\!-\!1}\!\!\right)dt.\label{zeta(s,a)}\end{aligned}$$

The permutation in the last integral is valid since the function $
\frac{x^s}{a^2+x^2}
 (t-x)^{m-s}
\frac{d^m}{dt^m}\left(\frac{1    }{e^{2\pi t}-1}\right)$ is less than $ \frac{t^m}{a^2}\frac{d^m}{dt^m}\left(\frac{1    }{e^{2\pi t}-1}\right)$ which integrable on $[0,\infty[$.

**Remark**

If we consider $a=1$, we get an integral representation of the Riemann zeta function:

$$\begin{aligned}
\zeta(s)&=&\frac{1}{2}+\frac{1}{s-1}+\int_0^\infty\frac{ 1}{1+ x^2} w_s(x)dx\end{aligned}$$

with

$$\begin{aligned}
\dis w_s(x):=\frac{2 (-1)^{m}x^s}{\Gamma(s)\Gamma(m+1-s)} \int_x^\infty (t-x)^{m -s}\frac{d^m}{dt^m}\left(\frac{1}{e^{2\pi t}-1}\right)dt\end{aligned}$$

where $m$ is an integer satisfying $ m > \Re(s)-1>0.$

So, if $s$ is an integer, we can take $m=s$, and

$$\begin{aligned}
\zeta(s)&=&\frac{1}{2}+\frac{1}{(s-1)}+ \frac{2 (-1)^{s-1}}{\Gamma(s)}
\int_0^\infty
\frac{  x^s}{1+ x^2} \frac{d^{s-1}}{dx^{s-1}}\left(\frac{1}{e^{2\pi x}-1}\right)dx\\
&=&\frac{1}{2}+\frac{1}{(s-1)}+ \frac{2}{\Gamma(s)}
\int_0^\infty
\frac{ 1}{e^{2\pi x}-1}\; \frac{d^{s-1}}{dx^{s-1}}\left(\frac{x^s}{1+x^2}\right)dx
\label{eq:}\end{aligned}$$

$$\begin{aligned}
\zeta(2)&=&\frac{3}{2} +\pi
\int_0^\infty
\frac{  x^2}{1+ x^2}  \frac{1}{\sinh ^2 \pi x}dx= \frac{3}{2}+4\int_0^\infty
\frac{1}{e^{2\pi x}-1}\;\frac{x}{(1+x^2)^2}dx\\
\zeta(3)&=& 1+ \pi^2
\int_0^\infty
\frac{  x^3}{1+ x^2} \frac{\cosh \pi x}{\sinh ^3\pi x}dx=
1-2\int_0^\infty \frac{1}{e^{2\pi y}-1}\;\frac{x(x^2-3)}{(1+x^2)^3}dx\\
\zeta(4)&=&\frac{5}{6} + \frac{\pi^3}{3}
\int_0^\infty
\frac{  x^4}{1+ x^2}\frac{2+\cosh(2\pi x)}{\sinh^4\pi x}dx=\frac{5}{6}-8 \int_0^\infty \frac{1}{e^{2\pi x}-1}\;\frac{x(x^2-1)}{(1+x^2)^4}dx.\end{aligned}$$

Proof of Theorem \[theo3\]
--------------------------

We consider the function $\Psi(s,t)$ given in (\[DAS\]) written as $$\begin{aligned}
\Psi(s,t) =\frac{B_0}{s-1}-B_1\, t+t^{2}\Psi_2(s,t)\end{aligned}$$ where $$\Psi_2(s,t):={\sum ^\ii_{k=0}}
{\frac{B_{k+2}}{(k+
    2)!}} (s)_{ k+1} (-1)^k t^k.$$

Thus the Padé approximant $[m+p/m]_{\Psi(s,t)}$ satisfies $$[m+p/m]_{\Psi(s,t)}=\frac{B_0}{s-1}-B_1 \,t+t^{2}[m+p-2,m]_{\Psi_2(s,t)}.$$

For $s$ positive real number, the weight $w_s$ is positive (Theorem \[theo2\]) and so $\Psi_2(s,t)$ is a Stieltjes function. Its coefficients $c_k:=\dis{\frac{B_{k+2}}{(k+1)!}} (s)_{ k+1} (-1)^k$ are positive and satisfy the Carleman condition, i.e. the series $$\sum_{k=0}^{\ii}c_k^{-1/2k} \mathrm{diverges}.$$ Actually, the Bernoulli numbers of even index satisfy (Bernoulli of odd index greater than 2 are zero) $$\left| B_{2k}\right| \sim 4 \sqrt{\pi k}\(\frac{k}{\pi e}\)^{2k}$$ and so $$c_{2k}^{-1/4k}\sim \(\frac{k}{\pi e}\)^{-1/2}$$ which is the general term of a divergent series.

We now apply the Theorem of [@bgm p.240] to conclude that for all integer $p, p\geq1,$ $$\lim_{m\rightarrow\ii}t^{s-1}[m+p/m]_{\Psi(s,t)}=t^{s-1}\Psi(s,t)=\ze(s,1/t)\mbox{ for all complex number  } t, \Re(t)>0$$

Particular case: $s$ is an integer {#particular case}
==================================

If $s$ is an integer, we can improve the approximation of $\ze(s,a)$. In (\[formzeta\]), we replace $\Psi\(s,\frac{1}{n+a}\)$ by some Padé approximant $[m+p/m]_{\Psi(s,t)}$ to obtain an approximation of $\ze (s,a)$:

$$\ze (s,a) = \dis{\sum ^{n-1} _{k=0}}
 \dis{\frac{1}{(k+a)^s}} + (n+a)^{1-s}\lim_{m \rightarrow \infty}
 \restriction{[m+p/m]_{\Psi(s,t)}}{t=1/(n+a)}  .$$

For $s=2$ and $s=3$, we will find in this section the formal expression of these approximants and the expression of the error.

The case for $s\geq 4$ remains an open problem.

Case $s=2$
----------

For $s=2$, the weight in the expression (\[zeta\]) is $$\begin{aligned}
w_2(x):=2 x^2  \int_x^\infty  \frac{d^2}{dt^2}\left(\frac{1}
{e^{2\pi t}-1}\right)dt=\frac{\pi x^2  }{\sinh^2 (\pi x)}.\end{aligned}$$ Another expression of $w_2(x)$ is $$w_2(x)=\frac{1}{\pi}\left| \Gamma(1+i x)\Gamma(1-i x)\right|  ^2.$$ Generalizing a result by Carlitz [@carlitz], Askey and Wilson [@askeywilson] gave an explicit expression: the orthogonal polynomial $P_n$ with respect to the weight function $w_2$ satisfies: $$\begin{aligned}
P_m(x)&=&(m+1)(m+2)_3F_2\(\begin{array}{ccccc}
-m&m+3&1-x&\\
&&&;1\\
&2&2&
\end{array}\)\\&=&
\sum_{k=0}^{m}\(\begin{array}{c}
m+1\\
k+1
\end{array}\)
\(\begin{array}{c}
m+k+2\\
k+1
\end{array}\)
\(\begin{array}{c}
x-1\\
k
\end{array}\),\label{polorths=2}\end{aligned}$$

and

$$\begin{aligned}
\int_{i \R}P_n(x)P_m(x)\frac{\pi x^2}{\sin^2\pi x}dx&=&0 \;\;(n\neq m),\\
&=&\frac{(-1)^n 2(n+1)(n+2)}{2n+3}\;\;(n=m).\end{aligned}$$

The roots of the polynomials $P_n$ are located on the imaginary axis since the weight $\dis\frac{\pi x^2}{\sin^2\pi x}$ is positive on this line (see \[integralrepr\] ).

The three terms recurrence relation is $$\begin{aligned}
P_{m+1}(x)&=&\frac{2(2m+3)}{(m+1)(m+2)}\,x\,P_m(x)+P_{m-1}(x)\end{aligned}$$ with initial conditions: $P_{-1}=0,P_0(x)=2.$

[**Associated polynomials**]{} (see \[associe\])

Before computing the associated polynomials, we need the modified moments, i.e., the moment of the binomial $ 
{\left(\!\!\!\begin{array}{c}x-1\\#2\end{array}\!\!\!\right)}.$

Let us define the linear functional $c^{(s)}$ acting on the space of polynomials as$$\left\langle c^{(s)},x^j\right\rangle :=\frac{B_j}{j!}(s)_{j-1}(-1) ^j,\;\;j\in \N,$$ and $x^2c^{(s)}$ by $$\left\langle x^2c^{(s)},x^j\right\rangle :=\left\langle c^{(s)},x^{j+2}\right\rangle.$$ By recurrence, it is not difficult to prove that $$\begin{aligned}
\left\langle x^2c^{(2)},
\(\begin{array}{c}
x-1\\
k
\end{array}\)\right\rangle &=&(-1)^k\frac{(k+1)(k+1)!}{(k+3)!}\end{aligned}$$ and $$\begin{aligned}
\left\langle x^2c^{(2)},
\(\begin{array}{c}
x-1\\
k
\end{array}\)\(\begin{array}{c}
x-1\\
j
\end{array}\)\right\rangle &=&(-1)^{k+j}\frac{(k+1)(k+1)!(j+1)(j+1)!}{(k+j+3)!}.\end{aligned}$$

If we define the following polynomial basis $$e_k(x):=\frac{1}{(k+1)(k+1)!}\(\begin{array}{c}
x-1\\
k
\end{array}\),k\in \N,$$

then $$\left\langle x^2c^{(2)},
e_k\;e_j\right\rangle =\frac{(-1)^{k+j}}{(k+j+3)!}.$$ These moments are related with the coefficients of the exponential function $$g(x):=\frac{e^{-x}-(1-x+x^2/2)}{x^3}.$$ So we can recover the expression (\[polorths=2\]) of the orthogonal polynomials for the functional $x^2c^{(2)}$ by substituting in the orthogonal polynomials for the function $g$ which are $$\sum_{k=0}^{m}\(\begin{array}{c}
m+1\\
k+1
\end{array}\)
\(\begin{array}{c}
m+k+2\\
k+1
\end{array}\)
x^k$$ the monomials $x^k$ by the $e_k's$ (see \[hadamard\]).

The associated polynomials $R_{m-1}$ (\[associe\]) are defined as $$R_{m-1}(t) := 
<x^ 2c ^ {(2)},\frac{P_m(x)-P_m(t)}{x-t}>$$ where the variable is $x$. From the expression (\[polorths=2\]) for $P_m$, we get the following formula for $R_{m-1}$ $$\begin{aligned}
R_{m-1}=\sum_{k=0}^{m}\(\begin{array}{c}
m+1\\
k+1
\end{array}\)
\(\begin{array}{c}
m+k+2\\
k+1
\end{array}\)\left\langle x^2c^{(2)},\frac{\(\begin{array}{c}
    x-1\\
    k
    \end{array}\)-\(\begin{array}{c}
    t-1\\
    k
    \end{array}\)}{x-t}\right\rangle.\end{aligned}$$

Using the expression of the polynomial $\dis\frac{\(\begin{array}{c}
    x-1\\
    k
    \end{array}\)-\(\begin{array}{c}
    t-1\\
    k
    \end{array}\)}{x-t}$ on the Newton basis on $0,1,2,\cdots,k-1$

$$\frac{{\left(\!\!\!\begin{array}{c}x-1\\#2\end{array}\!\!\!\right)}-{\left(\!\!\!\begin{array}{c}t-1\\#2\end{array}\!\!\!\right)}}{x-t}={\left(\!\!\!\begin{array}{c}t-1\\#2\end{array}\!\!\!\right)}\sum_{i=1} ^k
\frac{{\left(\!\!\!\begin{array}{c}x-1\\#2\end{array}\!\!\!\right)}}{i{\left(\!\!\!\begin{array}{c}t-1\\#2\end{array}\!\!\!\right)}},$$ we can write a compact formula for $R_{n-1}$: $$R_{m-1}(t)=\sum_{k=0}^m
{\left(\!\!\!\begin{array}{c}m+1\\#2\end{array}\!\!\!\right)}{\left(\!\!\!\begin{array}{c}m+k+2\\#2\end{array}\!\!\!\right)}{\left(\!\!\!\begin{array}{c}t-1\\#2\end{array}\!\!\!\right)}\sum_{i=1}^k
\frac{(-1)^{i-1}}{{\left(\!\!\!\begin{array}{c}t-1\\#2\end{array}\!\!\!\right)}(i+1)(i+2)}
.$$

So, we get the $[m+1/m]$ Padé approximant to the function $\Psi(2,t)$: $$[m+1/m]_{\psi(2,t)}=B_0-B_1t+t^2\frac{\widetilde{R}_{m-1}(t)}{\widetilde{P}_m(t)}=
B_0-B_1t+t\frac{{R}_{m-1}(1/t)}{{P}_m(1/t)}$$ and an approximation of $\ze(2,a)$

$$\begin{aligned}
\ze (2,a) &\approx &\dis{\sum ^{n-1} _{k=0}}
\frac{1}{(k+a)^2}+ \frac{1}{n+a}[m+1/m]_{\Psi(2,t)}\vert_{t=1/(n+a)}  \\&=&
\dis{\sum ^{n-1} _{k=0}}
\frac{1}{(k+a)^2} + \frac{1}{n+a}+\frac{1}{2(n+a)^2}+\frac{1}{(n+a)^2}\varepsilon_m(a)=:\dfrac{v_{n,m}(a)}{u_{n,m}(a)}\end{aligned}$$

where $$\varepsilon_m(a)=\frac{\sum_{k=0}^m
    {\left(\!\!\!\begin{array}{c}m+1\\#2\end{array}\!\!\!\right)}{\left(\!\!\!\begin{array}{c}m+k+2\\#2\end{array}\!\!\!\right)}\sum_{i=1}^k
    \frac{{\left(\!\!\!\begin{array}{c}n+a-i-1\\#2\end{array}\!\!\!\right)}(-1)^{i-1}}{{\left(\!\!\!\begin{array}{c}k\\#2\end{array}\!\!\!\right)}(i+1)(i+2)}}{\sum_{k=0}^{m}{\left(\!\!\!\begin{array}{c}m+1\\#2\end{array}\!\!\!\right)}
    {\left(\!\!\!\begin{array}{c}m+k+2\\#2\end{array}\!\!\!\right)}
    {\left(\!\!\!\begin{array}{c}n+a-1\\#2\end{array}\!\!\!\right)}}.$$

[**Irrationality of $\ze(2)$**]{}

A consequence of the previous formula is another proof of the well known irrationality of $\ze(2)$ since it equals to $\pi^2/6.$

Actually, if $d_k:=LCM[1,2,\cdots,k]$, then ([@prevost1996]) $$\frac{d_k }{i{\left(\!\!\!\begin{array}{c}k\\#2\end{array}\!\!\!\right)}}\in \N.$$ Then, the numerator of $\varepsilon_m(1)$ multiplied by $d_{m+2}^2$ is an integer and for all integers $n,m$, $d_{m+2}^2v_{n,m}(1)\in \N,d_{m+2}^2u_{n,m}(1) \in \N.$

The error (formula \[errorpade\]) applied to the function $\Psi(2,t)$ becomes $$\Psi(2,t) - [m+1/m]_{\Psi(2,t)}(t) = \displaystyle{\frac{t^{2m}}{\widetilde
        P^2_m(t)}}\int _{i\R}
{\frac{P^2_m(x)}{1-xt}} \frac{x^2}{\sin^2(\pi x)}dx$$ and the error term satisfies $$\begin{aligned}
\left| \ze(2)-\frac{v_{n,m}(1)}{u_{n,m}(1)}\right| &\leq& \frac{1}{P^2_m(n+1)}\int _{i\R}
\left| {\frac{P^2_m(x)}{1-x/(n+1)}} \right|  \frac{x^2}{\sin^2(\pi x)} dx\\
&\leq&\frac{1}{P^2_m(n+1)}\int _{i\R}
\left| {P^2_m(x)} \right|  \frac{x^2}{\sin^2(\pi x)} dx\\
&\leq&\frac{1}{P^2_m(n+1)}\frac{2(m+1)(m+2)}{\pi(2m+3)}\end{aligned}$$

Now, we consider $r \in \Q$ such that $m=r n\in\N$. Applying the Stirling formula to the expression (\[polorths=2\]) for orthogonal polynomials $P_n$, we get $$\begin{aligned}
  \limsup_n(P_{rn}(n+1)^{1/n}&=&\max_{t\in [0,1]}\frac{(r+t)^{r+t}}{t^{3t}(r-t)^{r- t}(1-t)^{1-t}}\\
  &=&\frac{(r+\sigma(r))^r}{(1-\sigma(r))(r-\sigma(r))^r}=\rho(r)
  \end{aligned}$$ where $\sigma(r)=\frac{-r^2+\sqrt{r^4+4r^2}}{2}$ is a zero of $t^2+r^2t-r^2=0$.

So

$$\begin{aligned}
   \limsup_n
   \left| 
   \ze(2)
   d_{r n+2}^2u_{n,rn}(1)
   -d_{rn+2}^2v_{n,rn}(1)
   \right|
   ^{1/n} 
   &\le & 
   \limsup_n (d_{rn+2}^2u_{n,rn}(1))^{1/n}
   \limsup_n\frac{1}{(u_{n,rn}(1))^{2/n}}\\
   &=&
  e^{2\max(r,1)}/\rho(r)
    \end{aligned}$$

The inequation $e^{2\max(r,1)}/\rho(r)<1$ is satisfied for $r\in]0.74,1.53[.$

The irrationality of $\ze(2)$ follows from the following limit $$\lim_n\( \ze(2)d_{rn+2}^2u_{n,rn}(1)-d_{rn+2}^2v_{n,rn}(1)\)=0.$$

Case $s=3$
----------

For $s=3$, the weight in the expression (\[zeta\]) is $$\begin{aligned}
w_3(x):=- \frac{x^3}{3}  \int_x^\infty  \frac{d^3}{dt^3}\left(\frac{1}
{e^{2\pi t}-1}\right)dt=\frac{\pi^2 x^3 \cosh \pi x  }{3\sinh^3 (\pi x)}.\end{aligned}$$ Another expression of $x^2w_3(x)$ is $$x^2w_3(x)=\frac{1}{12\pi}\frac{\left| \Gamma(1+i x)\right|^8}{\left| \Gamma(2i x) \right|  ^2}.$$ This weight has been investigated by Wilson [@wilson]. The orthogonal polynomial $P_n$ satisfies: $$\begin{aligned}
P_{m}(x^2)&=&(m+1)(m+2)_4F_3\(\begin{array}{cccccc}
-m&m+3&1-x&1+x\\
&&&&;1\\
&2&2&2
\end{array}\)\\&=&
\sum_{k=0}^{m}{\left(\!\!\!\begin{array}{c}m+1\\#2\end{array}\!\!\!\right)}
{\left(\!\!\!\begin{array}{c}m+k+2\\#2\end{array}\!\!\!\right)}
{\left(\!\!\!\begin{array}{c}x-1\\#2\end{array}\!\!\!\right)}{\left(\!\!\!\begin{array}{c}x+k\\#2\end{array}\!\!\!\right)}/(k+1) \label{polorths=3}\end{aligned}$$

and $$\begin{aligned}
\int_{i \R}P_{n}(x^2)P_{m}(x^2)\frac{\pi^2}{3}\frac{x^5\cos \pi x}{\sin^3\pi x}dx&=&0 \;\;(n\neq m)\\
&=&\frac{1}{3}\frac{(-1)^n (n+1)(n+2)}{2n+3}\;\;(n=m)\end{aligned}$$

The roots of the polynomials $P_n$ are located on the imaginary axis because the weight $\dis\frac{x^5\cos \pi x}{\sin^3\pi x} $is positive on this line.

We set $\Pi_m(x):=P_n(x^2).$

The three terms recurrence relation is $$\begin{aligned}
\Pi_{m+1}(x)&=&\(\frac{2(2m+3)}{(m+1)(m+2)^2}x^2+\frac{2m+3}{m+2}\)\,\Pi_m(x)-\frac{m+1}{m+2}\Pi_{m-1}(x)\end{aligned}$$ with initial conditions: $\Pi_{-1}=0,\Pi_0=2.$

[**Associated polynomials**]{}

Before computing the associated polynomials, we need the modified moments, i.e., the moment of the product of binomials $ 
{\left(\!\!\!\begin{array}{c}x-1\\#2\end{array}\!\!\!\right)}{\left(\!\!\!\begin{array}{c}x+k\\#2\end{array}\!\!\!\right)}.$

In the previous subsection, we have define the linear functional $c^{(s)}$ acting on the space of polynomials as $$\left\langle c^{(s)},x^j\right\rangle :=\frac{B_j}{j!}(s)_{j-1}(-1) ^j,\;\;j\in \N,$$ and $x^4c^{(s)}$ by $$\left\langle x^4c^{(s)},x^j\right\rangle :=\left\langle c^{(s)},x^{j+4}\right\rangle.$$ By recurrence, we can prove $$\begin{aligned}
\left\langle x^4c^{(3)},
{\left(\!\!\!\begin{array}{c}x-1\\#2\end{array}\!\!\!\right)}{\left(\!\!\!\begin{array}{c}x+k\\#2\end{array}\!\!\!\right)}\right\rangle &=&(-1)^{k+1}\frac{(k+1)^2}{2(k+3)(k+2)}\end{aligned}$$ and $$\begin{aligned}
\left\langle x^4c^{(3)},
{\left(\!\!\!\begin{array}{c}x-1\\#2\end{array}\!\!\!\right)}{\left(\!\!\!\begin{array}{c}x+k\\#2\end{array}\!\!\!\right)}{\left(\!\!\!\begin{array}{c}x-1\\#2\end{array}\!\!\!\right)}{\left(\!\!\!\begin{array}{c}x+j\\#2\end{array}\!\!\!\right)}\right\rangle &=&(-1)^{k+j+1}\frac{(k+1)^2(j+1)^2}{2(k+j+3)(k+j+2)}.\end{aligned}$$

If we define the following polynomial basis $$\theta_k(x):=\frac{(-1)^k}{(k+1)^2}{\left(\!\!\!\begin{array}{c}x-1\\#2\end{array}\!\!\!\right)}{\left(\!\!\!\begin{array}{c}x+k\\#2\end{array}\!\!\!\right)}$$

then $$\left\langle x^4c^{(3)},
            \theta_k\;\theta_j\right\rangle =\frac{(-1)}{2(k+j+3)(k+j+2)}=\frac{1}{2(k+j+3)}-\frac{1}{2(k+j+2)}.$$ These moments are those of the weight function $x(1-x)$ on the interval $[0,1]$.

So we can recover the expression of orthogonal polynomials for the functional $x^4c^{(3)}$ by substituting in the Jacobi orthogonal polynomials with parameters $\alpha=1,\beta=1$ .

$$P_m^{(1,1)}(2x-1)=\sum_{k=0}^m{\left(\!\!\!\begin{array}{c}m+1\\#2\end{array}\!\!\!\right)}{\left(\!\!\!\begin{array}{c}m+k+2\\#2\end{array}\!\!\!\right)} (-1)^kx^k$$ the monomials $x^k$ by the $\theta_k's$(see \[hadamard\]).

The associated polynomials $\Theta_{m-1}$ (\[associe\]) are defined as $$\Theta_{m-1}(t) := 
<x^ 4c ^ {(3)},\frac{\Pi_m(x)-\Pi_m(t)}{x-t}>$$ where the variable is $x$. From the expression (\[polorths=3\]) for $\Pi_m$, we get the following formula for $\Theta_{m-1}\in \mathcal{P}_{2m-2}
$ $$\begin{aligned}
\Theta_{m-1}=\sum_{k=0}^{m}\(\begin{array}{c}
m+1\\
k+1
\end{array}\)
{\left(\!\!\!\begin{array}{c}m+k+2\\#2\end{array}\!\!\!\right)}
\frac{1}{k+1}
\left\langle x^4c^{(3)},
\frac{
    {\left(\!\!\!\begin{array}{c}x-1\\#2\end{array}\!\!\!\right)}{\left(\!\!\!\begin{array}{c}x+k\\#2\end{array}\!\!\!\right)}-{\left(\!\!\!\begin{array}{c}t-1\\#2\end{array}\!\!\!\right)}{\left(\!\!\!\begin{array}{c}t+k\\#2\end{array}\!\!\!\right)}
    }
    {x-t}
    \right\rangle.\end{aligned}$$

Using the expression of the polynomial $\dis
\frac{
    {\left(\!\!\!\begin{array}{c}x-1\\#2\end{array}\!\!\!\right)}{\left(\!\!\!\begin{array}{c}x+k\\#2\end{array}\!\!\!\right)}-{\left(\!\!\!\begin{array}{c}t-1\\#2\end{array}\!\!\!\right)}{\left(\!\!\!\begin{array}{c}t+k\\#2\end{array}\!\!\!\right)}
    }
    {x-t}$ on the Newton basis on $0,1,-1,2,-2\cdots,k,-k$

$$\dis\frac{
    {\left(\!\!\!\begin{array}{c}x-1\\#2\end{array}\!\!\!\right)}{\left(\!\!\!\begin{array}{c}x+k\\#2\end{array}\!\!\!\right)}-{\left(\!\!\!\begin{array}{c}t-1\\#2\end{array}\!\!\!\right)}{\left(\!\!\!\begin{array}{c}t+k\\#2\end{array}\!\!\!\right)}}{x-t}
=
\sum_{j=0}^{2k-1}
\frac
{\nu_{2k}(t)}
{\nu_{j+1}(t)}
\frac
{\nu_j(x)}
{[(j+2)/2]}$$

where $$\begin{aligned}
\nu_{2j+1}(x)&=&{\left(\!\!\!\begin{array}{c}x-1\\#2\end{array}\!\!\!\right)}{\left(\!\!\!\begin{array}{c}x+j\\#2\end{array}\!\!\!\right)}\\
\nu_{2j}(x)&=&{\left(\!\!\!\begin{array}{c}x-1\\#2\end{array}\!\!\!\right)}{\left(\!\!\!\begin{array}{c}x+j\\#2\end{array}\!\!\!\right)}.\end{aligned}$$

Using $$\begin{aligned}
    \left\langle x^4c^{(3)},
\nu_{2j+1}(x)\right\rangle &=-\left\langle x^4c^{(3)},
\nu_{2j}(x)\right\rangle=&(-1)^{j}\frac{(j+1)^2}{2(j+3)(j+2)},\end{aligned}$$

we can write a compact formula for $\Theta_{m-1}$: $$\begin{aligned}
\Theta_{m-1}(t)&\!\!\!\!\!=&\!\!\!\!\!\sum_{k=0}^m
{\left(\!\!\!\begin{array}{c}m+1\\#2\end{array}\!\!\!\right)}{\left(\!\!\!\begin{array}{c}m+k+2\\#2\end{array}\!\!\!\right)}\frac{1}{k+1}\sum_{j=0}^{2k-1}
\frac{\nu_{2k}(t)}{\nu_{j+1}(t)}\left\langle x^4c^{(3)},\nu_j(x)\right\rangle\\
=&&\sum_{k=0}^m
{\left(\!\!\!\begin{array}{c}\!m\!+\!1\!\\#2\end{array}\!\!\!\right)}{\left(\!\!\!\begin{array}{c}\!m\!+\!k\!+\!2\!\\#2\end{array}\!\!\!\right)}\frac{\nu_{\!2k\!}(t)}{\!k\!+\!1\!}\sum_{p=1}^{k}
\frac{t}{{\left(\!\!\!\begin{array}{c}\!t\!-\!1\!\\#2\end{array}\!\!\!\right)}{\left(\!\!\!\begin{array}{c}\!t\!+\!p\!\\#2\end{array}\!\!\!\right)}}\frac{(-1)^p}{2(\!p\!+\!1\!)
    (\!p\!+\!2\!)}\end{aligned}$$

So, for $m\geq 2$, we get the $[2m+2/2m]$ Padé approximant to the function $\Psi(3,t)$, $$[2m+2/2m]_{\Psi(3,t)}=\frac{B_0}{2}-B_1t+\frac{3}{2}t^2B_2+t^5\frac{\widetilde{\Theta}_{m-1}(t)}{\widetilde{\Pi}_m(t)}=
\frac{B_0}{2}-B_1t+\frac{3}{2}t^2B_2+t^4\frac{     {\Theta}_{m-1}(1/t)}{ {\Pi}_m(1/t)}$$ and an approximation of $\ze(3,a)$ $$\ze (3,a) \approx \dis{\sum ^{n-1} _{k=0}}
\frac{1}{(k+a)^3} + \frac{1}{2(n+a)^2}+\frac{1}{2(n+a)^3}+\frac{1}{4(n+a)^4}+\frac{1}{(n+a)^5}\epsilon_m(a)=:\dfrac{f_{n,m}(a)}{g_{n,m}(a)}$$

where $$\begin{aligned}
\epsilon_m(a)&=&
    \frac{\Theta_{m-1}(n+a)}
{\Pi_m(n+a)}\\
&=&
\frac{
\sum_{k=0}^m{\left(\!\!\!\begin{array}{c}m+1\\#2\end{array}\!\!\!\right)}{\left(\!\!\!\begin{array}{c}m+k+2\\#2\end{array}\!\!\!\right)}\frac{n+a}{k+1}
    \sum_{p=1}^k
\frac{  {\left(\!\!\!\begin{array}{c}n+a-p-1\\#2\end{array}\!\!\!\right)}{\left(\!\!\!\begin{array}{c}n+a+k\\#2\end{array}\!\!\!\right)}(-1)^p}
{{\left(\!\!\!\begin{array}{c}k\\#2\end{array}\!\!\!\right)}^2 2(p+1)(p+2)}
}
{
    \sum_{k=0}^{m}{\left(\!\!\!\begin{array}{c}m+1\\#2\end{array}\!\!\!\right)}{\left(\!\!\!\begin{array}{c}m+k+2\\#2\end{array}\!\!\!\right)}{\left(\!\!\!\begin{array}{c}n+a-1\\#2\end{array}\!\!\!\right)}{\left(\!\!\!\begin{array}{c}n+a+k\\#2\end{array}\!\!\!\right)}/(k+1)  }.\end{aligned}$$

[**Irrationality of $\ze(3)$**]{}

The irrationality of $\ze(3)$ has been proved by Apery in a celebrated paper [@apery]. A little later, a particular straightforward and elegant proof of this irrationality was given by Beukers [@beukers]. The author gave another proof in [@prevost1996] using the RPA. Actually, $RPA(m,2m-1,2m)$ provides exactly the Apery’s numbers. In the following, we show that the previous approximation also proves irrationality of $\ze(3)$ with RPA’s of various degree.

Actually, $$\forall n,a \in \N,d_{m+1}^3\Theta_{m-1}(n+a)\in \Z$$ and $d_{m+1}^3 f_{n,m}(n+a)\in \Z,d_{m+1}^3g_{n,m}(n+a)\in \Z.$

The proof is similar than in previous subsection.

We take $a=1$.

The error (formula \[errorpade\]) applied to the function $\Psi(3,t)$ becomes $$\Psi(3,t) - [2m+2/2m]_{\Psi(3,t)} = \displaystyle{\frac{t^{2m}}{\widetilde
        \Pi^2_m(t)}}\int _{i\R}
{\frac{\Pi^2_m(x)}{1-xt}} \frac{\pi ^2 x^5\cos\pi x}{3 \sin^3(\pi x)}dx$$ and the error term $$\ze(3)-\frac{f_{n,m}(1)}{g_{n,m}(1)}$$ satisfies

$$\begin{aligned}
\left| \ze(3)-\frac{f_{n,m}(1)}{g_{n,m}(1)}\right| &\leq&
 \dfrac{1}{\Pi^2_m(n+1)} \int _{i\R}
 \frac{\Pi^2_m(x)}{\left| 1-x/(n+1)\right| } \frac{\pi ^2 x^5\cos\pi x}{3 \sin^3(\pi x)}dx  \\
 &\leq&\dfrac{1}{\Pi^2_m(n+1)} \int _{i\R}
 \Pi^2_m(x) \frac{\pi ^2 x^5\cos\pi x}{3 \sin^3(\pi x)}dx  \\
 &\leq &
\dfrac{1}{\Pi^2_m(n+1)}\frac{1}{3}\frac{ (m+1)(m+2)}{2m+3}\end{aligned}$$

Now, we consider $r \in \Q$ such that $m=r n\in\N$. Using the Stirling formula in the expression of the orthogonal polynomials (\[polorths=3\]), $$\begin{aligned}
    \limsup_n(\Pi_{rn}(n+1)^{1/n}&=&\max_{t\in [0,1]}\frac{(r+t)^{r+t}(1+t)^{1+t}}{t^{3t}(r-t)^{r- t}(1-t)^{1-t}}\\
    &=&\frac{(1+\mu(r))(r+\mu(r))^r}{(1-\mu(r))(r-\mu(r))^r}=\eta(r)\end{aligned}$$ where $\dis\mu(r)=\frac{r}{\sqrt{1+r^2}}$ is a zero of $(1-t^2)(r^2-t^2)=t^4$.

So $$\begin{aligned}
    \limsup_n
    \left| 
    \ze(3)
    d_{r n+1}^3 g_{n,rn}(1)
    -d_{rn+1}^3 f_{n,rn}(1)
    \right|
    ^{1/n} 
    &\le & 
    \limsup_n (d_{rn+1}^3g_{n,rn}(1))^{1/n}
    \limsup_n\frac{1}{(g_{n,rn}(1))^{2/n}}\\
    &=&
    e^{3\max(r,1)}/\eta(r)\end{aligned}$$

The inequation $e^{3\max(r,1)}/\eta(r)<1$ is satisfied for $r\in]0.74,1.36[.$

The irrationality of $\ze(3)$ then follows from the following limit $$\lim_n\( \ze(3)d_{rn+1}^3g_{n,rn}(1)-d_{rn+1}^3f_{n,rn}(1)\)=0.$$
