  Discrete updates of numerical partial differential equations (PDEs) rely on
two branches of temporal integration. The first branch is the widely-adopted,
traditionally popular approach of the method-of-lines (MOL) formulation, in
which multi-stage Runge-Kutta (RK) methods have shown great success in solving
ordinary differential equations (ODEs) at high-order accuracy. The clear
separation between the temporal and the spatial discretizations of the
governing PDEs makes the RK methods highly adaptable. In contrast, the second
branch of formulation using the so-called Lax-Wendroff procedure escalates the
use of tight couplings between the spatial and temporal derivatives to
construct high-order approximations of temporal advancements in the Taylor
series expansions. In the last two decades, modern numerical methods have
explored the second route extensively and have proposed a set of
computationally efficient single-stage, single-step high-order accurate
algorithms. In this paper, we present an algorithmic extension of the method
called the Picard integration formulation (PIF) that belongs to the second
branch of the temporal updates. The extension presented in this paper furnishes
ease of calculating the Jacobian and Hessian terms necessary for third-order
accuracy in time.
