The computation of Schrödinger equations in the semiclassical regime presents several enduring challenges due to the presence of the small semiclassical parameter. Standard approaches for solving these equations commence with spatial discretisation followed by exponentiation of the discretised Hamiltonian via exponential splittings.

In this thesis we follow an alternative strategy$-$we develop a new technique, called the symmetric Zassenhaus splitting procedure, which involves directly splitting the exponential of the undiscretised Hamiltonian. This technique allows us to design methods that are highly efficient in the semiclassical regime. Our analysis takes place in the Lie algebra generated by multiplicative operators and polynomials of the differential operator.

This Lie algebra is completely characterised by Jordan polynomials in the differential operator, which constitute naturally symmetrised differential operators. Combined with the $Z2$-graded structure of this Lie algebra, the symmetry results in skew-Hermiticity of the exponents for Zassenhaus-style splittings, resulting in unitary evolution and numerical stability.

The properties of commutator simplification and height reduction in these Lie algebras result in a highly effective form of $\textit{asymptotic splitting:} $exponential splittings where consecutive terms are scaled by increasing powers of the small semiclassical parameter. This leads to high accuracy methods whose costs grow quadratically with higher orders of accuracy.

Time-dependent potentials are tackled by developing commutator-free Magnus expansions in our Lie algebra, which are subsequently split using the Zassenhaus algorithm. We present two approaches for developing arbitrarily high-order Magnus--Zassenhaus schemes$-$one where the integrals are discretised using Gauss--Legendre quadrature at the outset and another where integrals are preserved throughout.

These schemes feature high accuracy, allow large time steps, and the quadratic growth of their costs is found to be superior to traditional approaches such as Magnus--Lanczos methods and Yoshida splittings based on traditional Magnus expansions that feature nested commutators of matrices.

An analysis of these operatorial splittings and expansions is carried out by characterising the highly oscillatory behaviour of the solution.