A review on numerical schemes for solving a linear stochastic oscillator

Abstract

In recent years several numerical methods to solve a linear stochastic oscillator with one additive noise have been proposed. The usual aim of these approaches was to preserve different long time properties of the oscillator solution. In this work we collect these properties, namely, symplecticity, linear growth of its second moment and asymptotic oscillation around zero. We show that these features can be studied in terms of the coefficients of the matrices that appear in the linear recurrence obtained when the schemes are applied to the oscillator. We use this study to compare the numerical schemes as well as to propose new schemes improving some properties of classical methods.

Appendix

Consider the recurrence (2.1) with initial conditions $x_0,\,y_0$ obtained when the PEM method is applied to the linear stochastic oscillator (1.1). The scheme may be written

$$\begin{aligned} \begin{pmatrix} x_{n}\\ y_{n} \end{pmatrix}=A^n \begin{pmatrix} x_0\\ y_0 \end{pmatrix} +\alpha \Delta W_0 A^{n-1} \begin{pmatrix} 0\\ 1 \end{pmatrix}+\alpha \,\Delta W_1A^{n-2} \begin{pmatrix} 0\\ 1 \end{pmatrix}+\cdots +\alpha \,\Delta W_{n-1} \begin{pmatrix} 0\\ 1 \end{pmatrix} \end{aligned}$$

where $A$ is given in (3.1). Suppose $\Delta <2$ and denote by

$$\begin{aligned} \lambda =1-\frac{\Delta ^2}{2}+ i\,\Delta \,\frac{\sqrt{4-\Delta ^2}}{2} \end{aligned}$$

and $\bar{\lambda }$ the eigenvalues of $A$. Since $|\lambda |=1$ we write $\lambda =e^{i\theta }$, $0<\theta <\pi$. It can be seen that

$$\begin{aligned} A^n=\begin{pmatrix} \cos (n\theta )+\frac{\Delta }{\sqrt{4-\Delta ^2}}\sin (n\theta ) &{} \frac{2}{\sqrt{4-\Delta ^2}}\sin (n\theta )\\ -\frac{2}{\sqrt{4-\Delta ^2}}\sin (n\theta ) &{} \cos (n\theta )-\frac{\Delta }{\sqrt{4-\Delta ^2}}\sin (n\theta ) \end{pmatrix}. \end{aligned}$$

If $\Vert \cdot \Vert$ denotes the Euclidean vector norm, we have

$$\begin{aligned} \text {E}[x_n^2+y_n^2]=\text {E}\left\| \begin{pmatrix} x_n\\ y_n \end{pmatrix}\right\| ^2=\left\| A^n \begin{pmatrix} x_0\\ y_0 \end{pmatrix}\right\| ^2+\alpha ^2\Delta \sum _{j=0}^{n-1}\left\| A^j \begin{pmatrix} 0\\ 1 \end{pmatrix}\right\| ^2. \end{aligned}$$

(7.1)

Since

$$\begin{aligned} \begin{aligned} \left\| A^n\begin{pmatrix} x_0\\ y_0\end{pmatrix}\right\| ^2=&\ (x_0^2+y_0^2)\left( 1+\frac{2\Delta ^2}{4-\Delta ^2}\sin ^2(n\theta )\right) + (x_0^2-y_0^2)\frac{\Delta }{\sqrt{4-\Delta ^2}}\sin (2n\theta )\\&+x_0y_0\frac{8\Delta }{4-\Delta ^2}\,\sin ^2(n\theta ) \end{aligned} \end{aligned}$$

and

$$\begin{aligned} \left\| A^j\begin{pmatrix} 0\\ 1\end{pmatrix}\right\| ^2=1+\frac{2\Delta ^2}{4-\Delta ^2}\sin ^2(j\theta )-\frac{\Delta }{\sqrt{4-\Delta ^2}}\sin (2j\theta ), \end{aligned}$$

from (7.1) we get

$$\begin{aligned} \text {E}[x_n^2+y_n^2]&= x_0^2+y_0^2+\sin ^2(n\theta )\frac{2\Delta }{4-\Delta ^2}\left( (x_0^2+y_0^2)\Delta +4 x_0y_0\right) \\&+\frac{\Delta }{\sqrt{4-\Delta ^2}}(x_0^2-y_0^2)\sin (2n\theta )\\&+ \alpha ^2\Delta \sum _{j=0}^{n-1}\left( 1+\frac{2\Delta ^2}{4-\Delta ^2}\sin ^2(j\theta )-\frac{\Delta }{\sqrt{4-\Delta ^2}}\sin (2j\theta )\right) . \end{aligned}$$

Finally, using that

$$\begin{aligned}&\displaystyle \sum _{j=0}^{n-1}\sin ^2(j\theta )=\frac{n}{2}-\frac{\cos ((n-1)\theta )\sin (n\theta )}{2\sin (\theta )}, \\&\displaystyle \sum _{j=0}^{n-1}\sin (2j\theta )=\frac{\sin (n\theta )\sin \left( (n-1)\theta \right) }{\sin (\theta )}, \end{aligned}$$

we obtain

$$\begin{aligned} \text {E}[x_n^2+y_n^2]&= x_0^2+y_0^2+\frac{4}{4-\Delta ^2}\alpha ^2 n\Delta +\frac{2\Delta }{4-\Delta ^2}\sin ^2(n\theta )\left( \Delta (x_0^2+y_0^2)+4x_0y_0\right) \\&+\frac{\Delta }{\sqrt{4-\Delta ^2}}(x_0^2-y_0^2)\sin (2n\theta )\\&-\,\alpha ^2\left( \frac{\Delta ^2}{\sqrt{4-\Delta ^2}}\,\frac{\sin ((n-1)\theta )\sin (n\theta )}{\sin (\theta )}+\frac{\Delta ^3}{4-\Delta ^2}\frac{\cos ((n-1)\theta )\sin (n\theta )}{\sin (\theta )}\right) . \end{aligned}$$