On the method of Neumann series for highly oscillatory equations (Q1770921)

The purpose of the author is to describe and analyse techniques for the numerical solution of highly oscillatory ordinary differential equations by exploiting a Neumann expansion. In the first part he presents the theoretical notions which are used in the sequel. Thus, given the formal equation \[ Y'= A(t)Y,\quad t\geq 0,\quad Y(0)=Y_0,\tag{1} \] where \(A\) might be a matrix, a differential operator or an integral operator, its solution is given by the series \[ Y(y)= (I- K)^{-1}Y_0= \sum^\infty_{m=0} K^m Y_0,\tag{2} \] where \[ K^0X= X(t),\quad K^{m+1}X= \int^t_0 A(x)K^m(x)\,dx,\quad m\geq 0.\tag{3} \] The Neumann expansions is known in physics as Dyson series and in matrix analysis as the Peano series. The second part computes the size of \(\| K\|\), where \(\|\cdot\|\) is the Euclidian norm and explores how it is attenuated by rapid oscillation. The third part describes the practical implementation of the Neumann series for the solution of a linear ordinary differential equation \[ y'= A(t)y,\quad t\geq 0,\quad y(0)= y_0\in\mathbb{R}^N\tag{4} \] having a fast oscillating solution. The author preconditions the solution of (4) locally by the known solution of a linear system with constant coefficients and solves \[ x'= B(t)x,\quad t\geq 0,\quad x(0)= y_n\in \mathbb{R}^N,\quad y_n= y(t_n),\quad n\geq 0, \] where \(B(t)= e^{-t\overline A}[A(t_n+ t)-\overline A] e^{t\overline A}\) with \(\overline A= A(t_n+{1\over 2} h)\), \(h= t_{n+1}- t_n\). Thus truncating (2) he gets \[ \begin{multlined} x_{n+1}= y_n+ \int^h_0 B(x)\,dxy_n+ \int^h_0\int^x_0 B(x_1)B(x_2)\,dx_2 dx_1y_n+\cdots\\ +\int^h_0 \int^{x_1}_0\cdots \int^{x_{r-1}}_{01} B(x_1) B(x_2)\cdots B(x_r)\,dx_r\cdots dx_1 y_n\end{multlined} \] which needs the replacement of multivariate integrals by quadrature. He describes how to accomplish this to high accuracy with a modest computational effort. The fourth section contains the application to the wave equation with a forcing term.