Arbitrary-Order Trigonometric Fourier Collocation Methods for Multi-Frequency Oscillatory Systems

Abstract

We rigorously study a novel type of trigonometric Fourier collocation methods for solving multi-frequency oscillatory second-order ordinary differential equations (ODEs) $q^{\prime \prime }(t)+Mq(t)=f(q(t))$ with a principal frequency matrix $M\in \mathbb {R}^{d\times d}$. If $M$ is symmetric and positive semi-definite and $f(q) = -\nabla U(q)$ for a smooth function $U(q)$, then this is a multi-frequency oscillatory Hamiltonian system with the Hamiltonian $H(q,p)=p^{T}p/2+q^{T}Mq/2+U(q),$ where $p = q'$. The solution of this system is a nonlinear multi-frequency oscillator. The new trigonometric Fourier collocation method takes advantage of the special structure brought by the linear term $Mq$, and its construction incorporates the idea of collocation methods, the variation-of-constants formula and the local Fourier expansion of the system. The properties of the new methods are analysed. The analysis in the paper demonstrates an important feature, namely that the trigonometric Fourier collocation methods can be of an arbitrary order and when $M\rightarrow 0$, each trigonometric Fourier collocation method creates a particular Runge–Kutta–Nyström-type Fourier collocation method, which is symplectic under some conditions. This allows us to obtain arbitrary high-order symplectic methods to deal with a special and important class of systems of second-order ODEs in an efficient way. The results of numerical experiments are quite promising and show that the trigonometric Fourier collocation methods are significantly more efficient in comparison with alternative approaches that have previously appeared in the literature.

Appendices

Appendix 1: Proof of Theorem 4.3

By virtue of Lemma 4.1, (27) and (28), one has

$$\begin{aligned}&H(\omega (h))-H(\mathbf {y}_0) =h \int _{0}^{1} \nabla H(\omega (\xi h))^T\omega '(\xi h)\mathrm{d}\xi \\&\quad =h \int _{0}^{1} \Big ( \big (Mv(\xi h)-\sum \limits _{j=0}^ {\infty }\widehat{P}_j(\xi )\gamma _j(v)\big )^T,\ u(\xi h)^T\Big )\\&\qquad \cdot \left( \begin{array}{c} u(\xi h) \\ -Mv(\xi h)+\sum \limits _{j=0}^ {r-1}\widehat{P}_j(\xi )\sum \limits _{l=1}^ {k}b_l\widehat{P}_j(c_l)f(v(c_l h)) \end{array} \right) \mathrm{d}\xi \\&\quad =h \int _{0}^{1} u(\xi h)^T \Big (\sum \limits _{j=0}^ {r-1}\widehat{P}_j(\xi )\sum \limits _{l=1}^ {k}b_l\widehat{P}_j(c_l)f(v(c_l h))-\sum \limits _{j=0}^ {\infty }\widehat{P}_j(\xi )\gamma _j(v) \Big )\mathrm{d}\xi \\&\quad =h \int _{0}^{1} u(\xi h)^T \Big (-\sum \limits _{j=0}^ {r-1}\widehat{P}_j(\xi )\varDelta _j(h)-\sum \limits _{j=r}^ {\infty }\widehat{P}_j(\xi )\gamma _j(v) \Big )\mathrm{d}\xi \\&\quad =-h \sum \limits _{j=0}^ {r-1}\int _{0}^{1} u(\xi h)^T \widehat{P}_j(\xi )\mathrm{d}\xi \varDelta _j(h)-h\sum \limits _{j=r}^ {\infty }\int _{0}^{1} u(\xi h)^T \widehat{P}_j(\xi )\mathrm{d}\xi \gamma _j(v)\\&\quad =h \sum \limits _{j=0}^ {r-1}\mathcal {O}(h^{j}\times h^{m-j})+h\sum \limits _{j=r}^{\infty } \mathcal {O}(h^{j}\times h^{j}) =\mathcal {O}(h^{m+1})+\mathcal {O}(h^{2r+1}). \end{aligned}$$

$\square$

Appendix 2: Proof of Theorem 4.4

From $Q(\mathbf {y})=q^TDp$ and $D^T=-D$, it follows that

$$\begin{aligned}&Q(\omega (h))-Q(\mathbf {y}_0) =h \int _{0}^{1} \nabla Q(\omega (\xi h))^T\omega '(\xi h)\mathrm{d}\xi \\&\quad =h \int _{0}^{1} \Big (- u(\xi h)^TD,\ v(\xi h)^TD\Big )\left( \begin{array}{c} u(\xi h) \\ -Mv(\xi h)\!+\!\sum \limits _{j=0}^ {r-1}\widehat{P}_j(\xi )\sum \limits _{l=1}^ {k}b_l\widehat{P}_j(c_l)f(v(c_l h)) \end{array} \right) \mathrm{d}\xi . \end{aligned}$$

Since $q^TD(f(q)-Mq)=0$ for any $q\in \mathbb {R}^{d}$, we obtain

$$\begin{aligned}&Q(\omega (h))-Q(\mathbf {y}_0) =h \int _{0}^{1} v(\xi h)^TD \Big (-Mv(\xi h)+\sum \limits _{j=0}^ {r-1}\widehat{P}_j(\xi )\sum \limits _{l=1}^ {k}b_l\widehat{P}_j(c_l)f(v(c_l h))\Big )\mathrm{d}\xi \\&\quad =h \int _{0}^{1} v(\xi h)^TD \Big (\sum \limits _{j=0}^ {r-1}\widehat{P}_j(\xi )\sum \limits _{l=1}^ {k}b_l\widehat{P}_j(c_l)f(v(c_l h))-\sum \limits _{j=0}^ {\infty }\widehat{P}_j(\xi )\gamma _j(v)\Big )\mathrm{d}\xi \\&\quad =h \int _{0}^{1} v(\xi h)^TD \Big ( -\sum \limits _{j=0}^ {r-1}\widehat{P}_j(\xi )\varDelta _j(h)-\sum \limits _{j=r}^ {\infty }\widehat{P}_j(\xi )\gamma _j(v)\Big )\mathrm{d}\xi \\&\quad =-h \sum \limits _{j=0}^ {r-1}\int _{0}^{1} v(\xi h)^TD \widehat{P}_j(\xi )\mathrm{d}\xi \varDelta _j(h)-h\sum \limits _{j=r}^ {\infty }\int _{0}^{1} v(\xi h)^TD \widehat{P}_j(\xi )\mathrm{d}\xi \gamma _j(v) \\&\quad =h \sum \limits _{j=0}^ {r-1}\mathcal {O}(h^{j}\times h^{m-j})+h\sum \limits _{j=r}^{\infty } \mathcal {O}(h^{j}\times h^{j}) =\mathcal {O}(h^{m+1})+\mathcal {O}(h^{2r+1}). \end{aligned}$$

$\square$