I am looking for an explicit finite difference scheme for modelling a transverse wave propagation in discontinued 2D media. Let me pose the problem (little bit simplified one) $$ \begin{cases} u_{tt} = a^2(u_{xx} + u_{yy}), & (x, y) \in (-\infty, +\infty) \times (0, +\infty) \\ u_{tt} = b^2(u_{xx} + u_{yy}), & (x, y) \in (-\infty, +\infty) \times (-\infty, 0) \end{cases} $$
$$ u(x, -0, t) = u(x, +0, t) $$ $$ u_t(x, -0, t) = u_t(x, +0, t) $$ $$ u_x(x, -0, t) = u_x(x, +0, t) $$ $$ t \geq 0 $$
With boundary conditions $$ u(x, y, 0) = u_0(x, y) $$ $$ u_t(x, y, 0) = v_0(x, y) $$
The conservativity of scheme is an important factor, because excitation and absorption of energy occur in different media (real problem contains damping).
Ideally I need a single formula that can be applied to any point including joint ones and any time. This is important because I would like to perform calculations using CUDA in real time and there are much more layers in real problem. Independent calculation of them would be a performance hit. I need to calculate whole time slice from two previous ones in a single pass.