A second order adjoint method is described for calculating directional derivatives of stiff {ODE} embedded functionals. The derivation of the general directional second order adjoint equations for point- and integral-form functionals is presented. A numerical procedure for calculating these directional derivatives that is relatively insensitive to the number of parameters is described and showcased. By combining automatic differentiation ({AD}) to obtain the adjoint and sensitivity equations with the staggered corrector method to solve the sensitivity systems, we achieve computational costs noticeably lower than directional finite differences based on a first order adjoint code.