useful for solving the TDME. In practice an efficient implementation of the first-order scheme is all that is needed to construct the higher-order algorithms Eqs.(\[secordapp\]) and (\[fouordapp\]).

To summarize: Suzuki’s product-formula approach provides the formal machinery to define algorithms that are unconditionally stable by construction. The accuracy of these algorithms can be improved systematically, to any desired order [@Suzuki8591]. The only assumption made so far is that the real, skew-symmetric matrix $H$ representing the TDME can be written as a sum of $p$ real, skew-symmetric matrices $H_i$. The next step is to choose the $H_i$’s such that the matrix exponentials $\exp(\tau H_1)$, ..., $\exp(\tau H_p)$ can be calculated efficiently. This will turn the formal expressions for $U_2(\tau)$ and $U_4(\tau)$ into efficient algorithms to solve the TDME.

Implementation {#sec4}
==============

In this section we present the details of our implementation of unconditionally stable algorithms to solve the TDME based on the Suzuki product-formula approach. For pedagogical reasons we start by considering the simplest case: A one-dimensional (1D) system. Then we show that the strategy adopted for 1D readily extends to higher spatial dimensions. The implementation we describe below is by no means unique, leaving a lot of room for further improvements. 