or equivalently with θ n instead of Θ n , where θn and φn are expressed in Eq. 113a-113b and Θn and φn in Eq. 104a-104b. In theory it is possible to have a singular value for θn (Eq. 113a) or φn (Eq. 104b) due to a possible zero value of the denominator. In practice, for the calculations performed in this study, this situation rarely occurs for the sampled initial condition for s = W approach, but may happen for the first time-step of the focused initial conditions with the choice of s = Q. Nevertheless, for the time-step where this situation occurs, one can switch back to Eq. 98c to avoid these rare numerical singularities in θn . In this case, one can switch to use a Verlet algorithm with Ωi and Ωi = 1 N 2 -1 j,k=1 f ijk H j Ωk , or to the MMST variables by propagating Eq. 117.