The numerical solution of the invariant imbedding equation, describing time domain, one dimensional direct scattering from a slab in which the material properties are spatially varying, is considered. It is proven that the equation discretised by the Trapezoidal rule has an asymptotic expansion for the global error involving only even powers of h. This expansion is utilised to generate a high order integration method by use of polynomial extrapolation. The method is suitable for adaptation to parallel computation, and by virtue of this together with its higher order integration, it constitutes a fast algorithm when compared with the current methods of solution of this equation.