High frequency wave propagation is well described even at caustics by Gaussian beams and the complex eikonal equation. In contrast to the real eikonal equation, the complex eikonal equation is elliptic and not well posed as an initial value problem. We develop a new model that approximates the 2D complex eikonal equation but is well posed as an initial value problem. This model consists of a coupled system of partial and ordinary differential equations. We prove that there exists a local solution to this new system by a Picard iteration method and show uniqueness under certain constraints. Different numerical approximations are then developed based on direct finite difference approximations or the method of characteristics. Numerical simulations with a variety of velocity profiles are presented and compared with solutions to the corresponding Helmholtz equation. / text