I am facing a simple (at first glance) problem. I need to implement a numerical scheme for the solution of the first order wave propagation equation with chromatic dispersion included. My original problem is (for a forward propagating wave):

\begin{equation} \frac{1}{c} \frac{\partial u(x,t)}{\partial t} = -\frac{ \partial u}{ \partial x} - \frac{i \beta_2}{2} \frac{ \partial^2 u}{ \partial t^2}, \end{equation} where $c$ is the velocity of light, $u$ is the (complex) envelope of the field, $\beta_2$ is the 2nd order dispersion coefficient. Assume also that the wave is propagating inside a ring cavity of length , say, L where I take periodic boundary conditions: $u(x+L,t) = u(x,t)$ and also that at $t=0$ we know $u(x,0) $ and $u_t(x,0)$.

I am trying to implement a time-stepping numerical scheme and in the process I tried the following:

1) MOL approach, where I do semidiscretization along $x$, reduce the set of equations to a system of first order ODEs (by setting $v = \dot{u}$) and I establish a system: \begin{equation} \begin{bmatrix} \dot{v} \\ \dot{u} \end{bmatrix} = A \begin{bmatrix} v \\ u \end{bmatrix} . \end{equation}

When I solve the corresponding ODEs via 4th Runge-Kutta, Crank-Nicholson , or simply precomputing the matrix exponential, unfortunatelly, all my solutions eventually blow up to Inf. I implemented the periodic boundary conditions by modifying the matrix $A$ as $A \leftarrow PA$ , where $P$ is the identity matrix with the first row identical copy of the last row.

I also tried a simple finite differences approach where the spatial derivative is approximated via an upwind FD (first order) but to no avail.

Lastly I tried a strang splitting approach based on the two equations: \begin{align} \frac{1}{c}\dot{u} &= -\frac{\delta u }{\delta x} \\ \frac{1}{c}\dot{u} &= -\frac{i\beta_2}{2}\frac{\delta^2 u }{\delta t^2} , \end{align}

where this time the solution does not blow up but it looks somehow unphysical.

Does someone here know a stable and possibly higher than first order time-stepping scheme for this equation? Please, note that a solution based on a Fourier transform in $x$ is also not a good option for me because I would like to have the flexibility to implement different non-periodic boundary conditions. I would also dislike substituting the second order derivative in time with one in space due to the fact that this complicates the implementation of boundary conditions.

Thanks.