Operator splitting methods provide efficient approximations to evolution equations of the form $u'(t)=(A+B)u(t)$, particularly when the flows generated by the operators $A$ and $B$ can be computed separately. In applications to PDEs, however, these operators are typically unbounded: differential operators such as the Laplacian are defined only on suitable domains and generate strongly continuous semigroups rather than bounded operator exponentials. Consequently, arguments based on norm-convergent Baker–Campbell–Hausdorff series or direct Taylor expansions are generally no longer applicable, or require unnecessarily strong regularity of the initial data. In this talk, we analyze symmetric Zassenhaus splittings in both the bounded and unbounded settings. Introduced in 2014 for the time-dependent semiclassical Schrödinger equation, this family of methods has proved particularly effective in the semiclassical regime and has become an important tool for the numerical simulation of quantum dynamics. We present a rigorous error analysis of symmetric Zassenhaus splittings in the presence of unbounded operators $A$ and $B$. To overcome the limitations of formal infinite expansions, we employ an integral representation of the local error together with finite commutator expansions whose remainders are controlled explicitly. This approach identifies sufficient regularity conditions for recovering the local convergence order known from the bounded setting. It also reveals an “overshoot” phenomenon: controlling the remainder may require commutators of greater length — and hence stronger regularity — than the formal truncation order alone would suggest.