Error Analysis due to finite time step in the GLE integration in a simple case

When we discretize the equation 5.13 we introduce an error due to the finite integration time step $\tau$ . Following the idea of Ref. (111) we can evaluate this error analytically in the case of a simple harmonic oscillator. Consider the equation:

$$x^{n+1}-(1+e^{-\lambda \tau })x^n+e^{-\lambda \tau} x^{n-1} = \frac{\tau}{\lambda}(1- e^{-\lambda \tau}) ( \omega^2 x^n +r^n).$$ (E.1)

We are interested in a statistically stationary process and so we proceed to evaluate mean average energies and correlation functions as function of $\tau$ . To do so we multiply the equation E.1 for $x^n$ and $x^{n-1}$ , respectively, and then take the average. We obtain the following pair of equations:

$$\langle x^{n+1} x^{n} \rangle + \left [ \frac{\tau \omega^2 }{\la... ... ] \langle x^{n} x^{n} \rangle +e^{-\lambda \tau} \langle x^{n-1} x^{n} \rangle =$$ 0

$$\langle x^{n+1} x^{n-1} \rangle -(1+e^{-\lambda \tau }) \langle x... ...\lambda \tau} - 1) + e^{-\lambda \tau} \right ] \langle x^{n-1} x^{n-1} \rangle =$$ 0

Becuase we are interested in the equilibrium distribution, we can assume that $\langle x^{n} x^{n-1} \rangle =\langle x^{n+1} x^{n} \rangle$ and $\langle x^{n+1} x^{n+1} \rangle =\langle x^{n} x^{n} \rangle = \langle x^{n-1} x^{n-1} \rangle$ . Thus we have three unknown quantities $\langle x^{n+1} x^{n} \rangle$ , $\langle x^{n} x^{n} \rangle$ and $\langle x^{n+1} x^{n-1} \rangle$ . To get a third relation among this quantities we square the Eq. E.1 to obtain the relation:

$$2 \left [ \frac{\tau \omega^2 }{\lambda} (e^{-\lambda \tau} -1 ) - ... ...tau}) \right ] \left(1+ e^{-\lambda \tau} \right) \langle x^{n+1} x^{n} \rangle$$

$$+ \left [ 1 + e^{-2\lambda \tau } + \left ( \frac{\tau \omega^2 }{\... ...da \tau} )+(1+e^{-\lambda \tau}) \right )^2 \right] \langle x^{n} x^{n} \rangle$$

$$+ 2 e^{-\lambda \tau} \langle x^{n+1} x^{n-1} \rangle = \frac{\tau^2}{\lambda^2}(1- e^{-\lambda \tau})^2\langle r^n r^n\rangle$$

After solving this equation system we can evaluate the potential and the kinetic energy:

$$\langle E_{pot} \rangle = \frac{1}{2} \omega^2 \langle x^{n} x^{n} \rangle$$ (E.2)

$$\langle E_{kin} \rangle = \frac{1}{2} \langle V^{n} V^{n} \rangle = \frac{1}{2} \frac{ \langle (x^{n}-x^{n-1})^2 \rangle }{\tau^2}$$ (E.3)

$$= \frac{1}{\tau^2} \left [ \langle x^{n} x^{n} \rangle - \langle x^{n} x^{n-1} \rangle \right]$$ (E.4)

It is easy to show that in the limit of small $\tau$ the potential and the kinetic energies converge to:

$$\langle E_{kin} \rangle = \frac{1}{2} k_b T \left ( 1 + \frac{1}{4} \omega^2 \tau^2 + O(\tau^4) \right)$$ (E.5)

$$\langle E_{pot} \rangle = \frac{1}{2} k_b T \left ( 1 + O( \tau^2) + O(\tau^3) \right)$$ (E.6)

This show that at least in this simple model the impulse integrator leads to a quadratic error in $\tau$ in both kinetic and poterntial energy.