{'fraction_non_alphanumeric': 0.04416956060040754, 'fraction_numerical': 0.015211089745429498, 'mean_word_length': 4.110589615723086, 'pattern_counts': {'":': 0, '<': 2, '<?xml version=': 0, '>': 6, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 4, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Algebraic theory for higher-order methods in computational quantum mechanics

We present the algebraic foundations of the symmetric Zassenhaus algorithm and some of its variants. These algorithms have proven effective in devising higher-order methods for solving the time-dependent Schr\"{o}dinger equation in the semiclassical regime. We find that the favourable properties of these methods derive directly from the structural properties of a Z2-graded Lie algebra. Commutators in this Lie algebra can be simplified explicitly, leading to commutator-free methods. Their other structural properties are crucial in proving unitarity, stability, convergence, error bounds and quadratic costs of Zassenhaus based methods. These algebraic structures have also found applications in Magnus expansion based methods for time-varying potentials where they allow significantly milder constraints for convergence and lead to highly effective schemes. The algebraic foundations laid out in this work pave the way for extending higher-order Zassenhaus and Magnus schemes to other equations of quantum mechanics.


Introduction
Recently devised symmetric Zassenhaus splittings (Bader, Iserles, Kropielnicka & Singh 2014) and Magnus-Zassenhaus methods , Bader, Iserles, Kropielnicka & Singh 2015 for computationally solving semiclassical time-dependent Schrödinger equations (TDSEs) have proven highly effective in achieving arbitrarily high orders of accuracy with costs growing quadratically in the desired order. This is in contrast to Yoshida based methods (Yoshida 1990), whose costs grow exponentially with order. These new schemes are commutator-free and allow much larger time steps than otherwise thought possible for convergence of Magnus and Zassenhaus expansions. Underlying the effectiveness of these new methods are some algebraic structures which form the crux of this work. In a break from tradition in numerical solutions of partial differential equations, the methods of (Bader et al. 2014) are devised by working directly in the free Lie algebra of undiscretised operators ∂ 2 x and V (x). The quadratic costs of symmetric Zassenhaus methods are achievable due to the property of height reduction in these Lie algebras. This property entirely fails to materialise when working with matrices corresponding to the discretisation of ∂ 2 x and V (x). Hence the structure of these free Lie algebras is crucial to the effectiveness of Zassenhaus schemes.
The characterisation of these algebras forms a major motivation of this work. It is this characterisation that will allow us to formally prove the quadratic growth in costs. We will also find that these algebras possess a Z 2 -graded structure that, along with the symmetry of its constituents, proves crucial for the numerical stability as well as the unitary evolution of the wave-packet. Instead of characterising the free Lie algebra, however, we will attempt the characterisation of a bigger free Lie algebra, FLA(S ∪ P(∂ x )), generated by all functions in the function space S and all polynomials of the differential operator ∂ x .
These algebraic characterisations will be introduced in a more abstract form-in the context of an associative algebra A with a commutative subalgebra C and its Lie idealiser I. The algebraic structures of the Zassenhaus schemes for the TDSE are seen to be special cases of the Lie algebras introduced here. To be more concrete, the associative algebra A in the case of the TDSE will consist of linear operators, the commutative subalgebra will be formed by multiplicative operators (whose action is multiplication by V , for instance), and ∂ x will be seen to be an element of the Lie idealiser I. The properties of these Zassenhaus schemes are traced to the structural properties of these algebras.
This algebraic formulation will prove of crucial importance in cutting the Zassenhaus algorithms free of their semiclassical TDSE origins, allowing their use in devising higher-order computational methods for other equations of quantum mechanics and possibly beyond.
In Section 2 we introduce the abstract context in which our algebra F is defined. In Section 3, we show that F is an associative algebra and study its structure. While it follows immediately that F is also a Lie algebra, in Section 4 we find out that it possesses a very interesting structure which has significant ramifications for numerical methods discussed in Section 5.
The tables of coefficients in Appendix A should aid direct applications of the results of this paper. The lengthier proofs of Lemma 2 and Lemma 3 from Section 3 have been confined to Appendices B and C, respectively, for ease of narration.

Some notations and definitions
Consider a commutative algebra C which is a subalgebra of the unital associative algebra (A, ·, +) over the field of rational numbers Q. The commutator on the associative product, acts as the canonical Lie product while the anticommutator, acts as a Jordan product. A along with the Lie (Jordan) product forms a Lie (Jordan) algebra which we can identify with A again. The Lie idealiser of C in A, is the largest subalgebra of A in which C is a Lie ideal.
In this work we typically use a, b, c for elements of the associative algebra A, x, y, z for elements of the commutative subalgebra C and reserve d for elements of the Lie idealiser I, except where it is clear from context. For elements x, y of the commutative algebra C, we write xy to denote x · y, dropping the explicit use of the multiplication operator. The letters k, l, m, n, p, r, s, i are used for non-negative integers.
Consider the adjoint map, where the action of ad a for a ∈ A is described by ad a is always a derivation on A and is called an inner derivation. However, it need not be a derivation on C in general. For d ∈ I, ad d is also a derivation on C, Definition 1 For x ∈ C, d ∈ I and a non-negative integer k we introduce the notation Definition 2 For any d ∈ I and non-negative integer k, we define the linear space and denote the direct sum of these linear spaces by In this work we will only use a single non-trivial element of I which is not in C, d ∈ I \ C, and can drop the superscript d, writing x k , F k and F in place of x d k , F d k and F d , respectively. To simplify notation further, we denote ad i d (x) by D i x. Since D is a derivation on A (and C), D(a·b) = Da·b+a·Db, and it distributes binomially on A (and C), Starting from d · x = Dx + x · d, a simple inductive procedure leads to a similar binomial identity, 3 F as an associative algebra Lemma 1 The linear space F is an associative algebra with

1)
where the terms are in C and π k,l n,i ∈ Q.
Proof A proof for Lemma 1 follows by expanding (3.1) using the binomial identity (2.2) and comparing powers of d. At first, we expand the left hand side, The right hand side of (3.1) is expanded in similar fashion, A sequence of terms z n for which all terms accompanying d k+l−a , a ∈ {0, . . . , k + l}, match on both sides will certainly satisfy (3.1) although it need not be the unique solution. This leads to the relation which can easily be fashioned into a recursive procedure for solving z a , starting from z 0 = xy (found by substituting a = 0). The second part of the lemma states that z n of the form (3.2) satisfy (3.1) for some π k,l n,i ∈ Q. After substituting this form in (3.3), the left side is comprised of where the inner term D a−n · D i x D n−i y has been expanded using (2.1). We now equate terms accompanying D p x D a−p y in (3.3), noting that this would give a solution which need not be unique. We arrive at equations of the form and δ i,j is the Kronecker delta function. The fact that a recursive procedure for finding πs can be designed, starting from π k,l 0,0 = 1, is the proof that such π k,l n,i ∈ Q exist, whereby some z n ∈ C of the form (3.2) satisfy (3.1). Hence F is an associative algebra with the prescribed properties (3.1) and (3.2). 2

Explicit form of the coefficients
The recursive procedure (3.8) suffices for the proof of Lemma 1 and for expanding products in the associative algebra F using a symbolic computation algorithm. Coefficients derived using this procedure for a few combinations of k and l are listed in Table 1 in Appendix A. Nevertheless, an explicit form for these coefficients remains highly desirable for a variety of reasons. Firstly, this could give us an explicit form for some analytic functions on F (or some subspace of F), such as the exponential map. Secondly, an analysis of growth of these coefficients could allow us to deduce estimates concerning the convergence of some of these functions. Lastly, the observations which lead us to the structural properties of the Lie algebra in Section 4-the main result of this paper-are not immediately evident in the recursive form.
We state here some results concerning the form of the coefficients, confining the proofs to appendix B and appendix C for ease of narration.

Lemma 2
The explicit form of the coefficients is given by (3.10) and P r are defined in terms of the Bernoulli numbers B r , Some properties of the coefficients π k,l n,i become evident once we study their generating functions. In this context we make an exception in notation, using u, w, x and y as variables in which the formal series of the generating function is specified.

F as a Lie algebra
Since F is an associative algebra, it is immediately obvious that it is also a Lie algebra with the commutator as a canonical Lie product. From (3.1), we know that the commutator can be expanded to where µ k,l n,i = π k,l n,i − π l,k n,n−i . However, the structure of F turns out to be more interesting than this.
Theorem 5 F is a Lie algebra where commutators can be solved explicitly using the rule where λ k,l n,i = 2π k,l 2n+1,i ∈ Q.
Proof The even indexed coefficients µ k,l 2n,i in (4.1) vanish due to Lemma 4, µ k,l 2n,i = π k,l 2n,i − (−1) 2n π k,l 2n,i = 0, while µ k,l 2n+1,i = π k,l 2n+1,i − (−1) 2n+1 π k,l 2n+1,i = 2π k,l 2n+1,i . We conveniently rename µ k,l 2n+1,i as λ k,l n,i . 2 We now look at some interesting properties of this Lie algebra. Consider the linear spaces, We note that, as a consequence of Theorem 5, we have a natural filteration, This property of height reduction is evident in (4.2) where the largest index k + l − 1 occurs for s = 0. We define the height of a term in G n as, ht n k=0 x k k = n, if x n is non-zero. We define the height of 0 as −1, this being the only term with negative height. Thus ht ( x k ) = k and, as a consequence of (4.2), Corollary 6 (Height Reduction) For any commutator C featuring the letters a i ∈ G k i , i = 1, . . . , n, ht (C (a i , . . . , a n )) ≤ n i=1 k i − n + 1. (4.4) Another useful property of this Lie algebra is that it doesn't mix terms of certain forms. We Thus, if each letter a i is either in e or in o, the commutator C (a i , . . . , a n ) is either in e or o and does not mix terms from the two. Moreover, commutators with the same letters will fall in the same space, e or o. This has implications for certain numerical methods where this property translates into a matrix being either symmetric or skew-symmetric but not a mix of the two which would have resulted in unfavourable structures. The property (4.5), in fact, says that F is a Z 2 -graded Lie algebra with the degree 0 component o and degree 1 component e. The linear space o is a Lie algebra in its own right, while structures of the form e are also called Lie triple systems which are closed under double commutation: [e, [e, e]] ⊆ e.
These notions are closely related to Lie groups and symmetric spaces which have found applications in linear algebra methods such as the generalised polar decomposition (Munthe-Kaas, Quispel & Zanna 2001). As anticipated, there is a close relation between F and the free Lie algebra generated by C and polynomials in d.

Lemma 7
The free Lie algebra generated by C and P(d) (the ring of polynomials in d with constant coefficients) is contained in F, The two are identical if an inverse of the mapping D = ad d : C → C exists.
Proof The containment is not difficult to prove. Since d k = 1 k , it is contained in F k and therefore P(d) ⊆ F. The algebra C is also contained in F since every x ∈ C can be written in the form x = x 0 . The free Lie algebra g generated by C ∪ P(d) is the intersection of all Lie algebras containing C and P(d) and is therefore contained in the Lie algebra F.
The two algebras g and F are identical if inverse of the map D exists. Any term in F 0 is of the form y 0 for some y ∈ C and therefore trivially resides in g. Take any x ∈ C and note that d 2 ∈ P(d). These are both contained in g and by definition of being a Lie algebra Consequently, any y 1 ∈ F 1 can be expressed as 1 2 [d 2 , D −1 (y)] so long as the inverse of D exists. Thus we have F 0 , F 1 ⊆ g. These two cases form the base case of our induction argument.
Assume that for all k < n, F 2k and F 2k+1 are contained in g. For completing the inductive proof, consider The first term on the right hand side is in g since d 2n+1 ∈ P(d), while the terms in the summation fall in g due to the induction hypothesis. Thus, we find that y 2n ∈ g for any y ∈ C. A similar proof shows that y 2n+1 also resides in g, proving that F k ⊆ g for every k. Since g is also a linear space, the direct sum of these spaces, F = k∈Z + 0 F k , is also contained in it. This completes our proof, g = F. 2 We note that, in general, the free Lie algebra g is not a subalgebra of the Lie idealiser I since d 2 need not be in I for every d ∈ I.
Of less immediate and practical interest to us is the fact that F is also a Jordan algebra. This follows directly and trivially from the fact that it is an associate algebra (Lemma 1). More interestingly, it is also a Z 2 -graded Jordan algebra since, along similar lines to (4.1), where γ k,l n,i = (π k,l n,i + π l,k n,n−i )/2 vanishes for odd values of n due to Lemma 4, while γ k,l 2n,i = π k,l 2n,i survives. Corresponding observations about the free Jordan algebra generated by C and P(d) can also be made along similar lines, however a property corresponding to height reduction of Corollary 6 does not follow.

Some applications
The applications of the algebras of Section 4 in solving partial differential equations such as the time-dependent Schrödinger equation will arise from treating function spaces as commutative algebras and linear differential operators as elements of the associative algebra of endomorphisms. With this motivation, we will restrict our attention to cases where the commutative algebra C is isomorphic to a function space. Let S be a commutative algebra isomorphic to C and Θ be the isomorphism between them. Anticipating the case where S is a function space, we will use f, g, h to denote its elements. For f ∈ S, we use the notation Θ f and Θ(f ) interchangeably for the corresponding element in C.
For an element d ∈ I in the Lie idealiser of C, D = ad d is a derivation on C. For any f ∈ S, Since D is a derivation on C and S is isomorphic to C,d is a derivation on S. We say thatd is a derivation induced by d. The element of the Lie idealiser that we will require in the following sections will be the differential operator d = ∂ x , which will end up coinciding with the induced derivationd.
Our first example of Θ will be M, the left multiplication map, which we encounter while solving the TDSE. In the case of the Wigner equation, the isomorphism Θ is more complicated and maps functions to pseudo-differential operators.

Semiclassical time-dependent Schrödinger equation (TDSE)
The linear time-dependent Schrödinger equation (TDSE), is a fundamental equation of quantum mechanics. The reduced Planck's constant, ≈ 1.054× 10 −34 J · s, is a truly minute number which would typically result in considerable difficulties as far as numerical solutions are concerned. However, when working in atomic units, we can re-scale = 1. In these units the mass of an electron is 1. The typical length scales in these units are 10 −11 m, while typical time and mass scales are 10 −17 s and 10 −30 kg, respectively. Consequently, working in the atomic units restricts us to extremely small spatio-temporal windows.
When computation for any larger spatio-temporal windows is required, we arrive at the TDSE under the semi-classical scaling, The semiclassical parameter ε plays a role similar to the reduced Planck's constant. Although it is very small, 0 < ε 1, the semiclassical parameter is considerably larger than the Planck's constant and the range 10 −8 ≤ ε ≤ 10 −2 isn't unrealistic.
The equation (5.1) also arises out of the Born-Oppenheimer approximation of the molecular Schrödinger equation. In this case it describes the evolution of the nuclei in the electric potential, V (x). Working in the atomic units, ε is the square root of the ratio of the mass of an electron and the heaviest nucleus. When larger spatio-temporal windows are required, ε can decrease further in size as before.
Typical computational methods for solving this equation commence with imposition of periodic boundaries at ±1 followed by a discretisation of space. The infinite dimensional operators ∂ x and V are replaced by matrices K and D V , respectively. The TDSE is now replaced by a system of ODEs, where u is the vector of values at grid points. The exact solution of this equation can be written as Methods for evaluating the exact matrix exponential for this matrix are prohibitively expensive and one resorts to splitting methods such as the Trotter splitting, or the Strang splitting, Methods with higher-order accuracy than the Strang splitting can be created using Yoshida composition (Yoshida 1990). However the number of exponentials grows exponentially with the order desired. Other methods are based on the symmetric Baker-Campbell-Hausdorff (sBCH) formula (Dynkin 1947, Casas & Murua 2009, These require the evaluation of nested commutators of the matrices K 2 and D V such as , which is prohibitive and therefore typically avoided. In (Bader et al. 2014) we develop arbitrarily high order methods by working directly within the free Lie algebra generated by the undiscretised operators ∂ 2 x and V . The cost of these methods grows quadratically in contrast to Yoshida methods. Additionally, these methods are commutator-free and preserve unitary evolution of the wave function-a fundamental principle of quantum mechanics. In this section we discover how the remarkable properties of symmetric Zassenhaus methods of (Bader et al. 2014) arise from the structural properties of (a special case of) the Lie algebra F.

Notation
Let S ⊆ T be real and complex valued function spaces, respectively. Endomorphisms on T form an associative algebra A = (End(T ), •, +) where • is operator composition. For convenience, we consider S = C ∞ p ([−1, 1]; R) and T = C ∞ p ([−1, 1]; C), the space of smooth periodic functions over [−1, 1] with values in R and C, respectively. Let M : S → End(T ) be the left multiplication map, The image of S under M , is a commutative algebra of multiplication operators with M acting as an isomorphism between S and C. We also write M f to denote the operator M (f ) whose action is that of multiplying by f . The partial differentiation operator ∂ x is an element of the Lie idealiser of C in the associative operator algebra (End(T ), •, +), where the commutator is the canonical Lie product as usual. Following our terminology, the operator ∂ x ∈ End(T ) induces the derivation∂ x = ∂ x on S. Here the induced derivationd overlaps with the element d, but this need not be the case in general. To be very precise,∂ x is an operator on S while ∂ x is an operator on T , but the two coincide on S. Lemma 7 directly allows us to conclude that Here D −1 is isomorphic tod −1 which is the inverse of differentiation on the space S. Assuming S is closed under integration, therefore, an appropriate inverse of D exists and the free Lie algebra of multiplicative operators M(S) and polynomials (with constant coefficients) in the differential operators, P(∂ x ), is characterised completely by F ∂x . For convenience, we abuse notation and write In this notation where 1 S is the constant function over S with value 1.

Consequences for exponential splitting schemes
The semiclassical TDSE (5.1), where the unit in 1 2 is understood to be 1 S . We can find the solution of this equation by formally exponentiating the Hamiltonian without discretisation, We may now perform splittings directly on this undiscretised Hamiltonian. Methods which featured nested matrix commutators such as [[D V , K 2 ], K 2 ] now feature corresponding commutators of operators 1 2 and V 0 , such as [[ V 0 , 1 2 ], 1 2 ]. Since the two operators reside in our Lie algebra, so do all of their commutators. Moreover, as a consequence of Theorem 5, commutators in this Lie algebra can be solved explicitly, where λs are the same as before. This allows us to design commutator-free methods by explicitly working out nested commutators such as [ appearing in the sBCH of itε 1 2 and −itε −1 V 0 using the rules, which can be read off Table 2 in Appendix A. Consequently, The sBCH of itε 1 2 and −itε −1 V 0 up to O t 5 , for instance, is Another favourable consequence of Theorem 5 is that all terms in the sBCH of itε 1 2 and whose elements are skew-Hermitian operators. Upon discretisation using pseudospectral methods, these are replaced by where D f is a diagonal matrix and K a skew-symmetric circulant differentiation matrix. Consequently, the discretised forms of elements of H are skew-Hermitian matrices, the exponentials of which are unitary matrices. Zassenhaus splittings (Bader et al. 2014) feature exponentials of terms in H, thereby guaranteeing unitary evolution of the wavefunction and resulting in unconditional stability of these numerical methods.
It is well known that solutions of the semiclassical TDSE (5.1)-irrespective of smooth initial conditions-develop oscillations of frequency O ε −1 in both space and time (Jin, Markowich & Sparber 2011). These oscillations make computational solutions very costly since one typically needs M = O ε −1 grid points to resolve spatial oscillations and time steps of size O (ε) to resolve temporal oscillations.
The spectral radius of the differentiation matrix K scales as O (M ), i.e. O ε −1 . Keeping eventual discretisation in mind and noting that we use the shorthand f k = O ε −k . It is now easy to see how the height of a term a ∈ G n acts as a proxy for the spectral radius of the eventual discretisation, a ∈ G n ⇐⇒ ht(a) ≤ n =⇒ a = O ε −n .

Corollary 8 Commutators of
Proof Consider a grade n commutator C of 1 2 and V 0 , featuring k occurrences of the letter 1 2 and n − k occurrences of V 0 . Since 1 2 ∈ G 2 and V 0 ∈ G 0 , using Corollary 6, ht(C) ≤ 2k − n + 1.
Thus C is O ε −2k+n−1 . However, the corresponding commutator of A and B is scaled by k occurrences of ε and n − k occurrences of ε −1 , bringing its size to O ε −1 . 2 When we scale the time step as t = O (ε σ ) for some σ > 0, a grade n commutator of A = itε 1 2 and B = itε −1 V 0 , scales as O ε nσ−1 in contrast to the O ε n(σ−1) scaling seen when working with commutators of matrices. Whereas the latter requires σ > 1 for convergence of sBCH (which corresponds to very small time steps), our approach does not have such restrictions and allows convergence for extremely large time steps. With time steps of size O ( √ ε), which corresponds to σ = 1/2, for instance, grade n commutators are O ε n/2−1 and the sBCH series easily converges.
The property of height reduction underlies the asymptotic splitting of the symmetric Zassenhaus kind where it results in inexpensive exponentiations via Lanczos iterations (Bader et al. 2014). The advantages of height reduction also extend to error analysis for existing methods such as Yoshida splittings since commutators discarded in these methods can be assumed to be undiscretised and analysed appropriately. An order six Yoshida splitting discards all commutators of grade seven and higher, for instance, thereby committing an error of O ε 7σ−1 .
Despite this, the quadratic cost of the Zassenhaus splitting, which results from a systematic exploitation of the structures of F and the property of height reduction, trumps the exponential cost of Yoshida.
Theorem 9 Cost of Zassenhaus splittings grows quadratically in the order desired.
The overall cost for an O ε (2n+3)σ−1 splitting is O C σ (n)ε −1 log ε −1 per time step. We remark that careful choices can allow us to reduce the exact number of FFTs required. However, we are largely concerned with asymptotic growths here. Since this proof does not assume a specific form of the exponents, W [k] s, and works solely on the basis of the structure of G k and e, the odd graded nature of the sBCH formula and the height reduction of Corollary 6, the cost estimates derived here also apply directly to Magnus-Zassenhaus approaches of , Bader et al. 2015 discussed in Section 5.2. 2

Some related algebras
The gain of powers of ε has also been observed by Gaim & Lasser (2014) using Moyal brackets in the phase space. However, the analysis is not as easily generalised as the algebraic approach here. These ideas are closely related to the Weyl algebra (Dixmier 1968, Coutinho 1997) which is the universal enveloping algebra of the Lie algebra of Heisenberg groups. The univariate Weyl algebra can be written in the form These are special cases of a more general form under the choice of polynomials in x, P(x), as the function space S. We note that S is contained within F as an associative algebra since M(f k ), ∂ k x ∈ F, while, in the other direction, F ⊆ S is evident using the Leibniz rule. Thus, the two are identical as associative algebras. It is not difficult to prove that a height reduction phenomenon also holds for S.
However, there is a significant advantage to working in the symmetrised form of F since a S lacks a Z 2 -graded structure. The Z 2 grading on F and the symmetric nature of its elements ensures that elements of H are skew-Hermitian after discretisation. This proves crucial for devising stable numerical schemes such as the symmetric Zassenhaus schemes of (Bader et al. 2014) once we start utilising nested commutators. In contrast, lacking a Z 2 grading, a clean separation of terms does not occur in the form S.
One could argue that any nested commutator of skew-Hermitian operators in S should be skew-Hermitian. This is indeed true prior to discretisation. However, unless some highly specialised differentiation matrices can be constructed, the discretised version of the commutator simplified in S possesses no such structure. What is more, even prior to discretisation we prefer to discard terms smaller than a certain size (while analysing in powers of ε). Working in S instead of F, it becomes difficult to discern which components of a term such as n k=0 f k (x)∂ k x can be discarded and which need to be kept despite their small size in order to preserve the skew-Hermiticity of the undiscretised operator.

Time-varying potentials in the TDSE
Consider the TDSE with time-varying electric potential, V (x, t), written in the form The solution of such Lie group equations can be expressed using the Magnus expansion (Magnus 1954), where the infinite series Ω(t) = ∞ k=1 Ω k (t), is an element of the underlying Lie algebra H. The exponent Ω(t) satisfies the dexpinv equation (Iserles & Nørsett 1999), , where B k are the Bernoulli numbers. The solution of (5.5), originally proposed as an infinite series using Picard iterations (Magnus 1954), has been widely analysed in (Iserles & Nørsett 1999, Iserles, Munthe-Kaas, Nørsett & Zanna 2000, Blanes, Casas, Oteo & Ros 2009). The graded version of the Magnus expansion (Iserles & Nørsett 1999) reads Once again, since A(t) ∈ H, all commutators can be simplified to commutator-free terms in H. This is in contrast to (Hochbruck & Lubich 2003) where the Magnus expansion features nested commutators of matrices. Bounds for convergence of the Magnus expansion analysed in (Iserles & Nørsett 1999) were improved upon in (Moan & Niesen 2008). Applied directly to the semiclassical TDSE, these would impose a highly stringent restriction of t = O ε 2 on the time step. Hochbruck & Lubich (2003) analysed the convergence of the Magnus expansion in the context of the TDSE, making significant improvements on the earlier and more general bounds. Nevertheless, under semiclassical scaling their analysis suggests a need for very small time steps of size O (ε) for convergence.
However, due to Corollary 8, a grade n term in the Magnus expansion of A(t) is O ε nσ−1 . Thus, asymptotically speaking in terms of ε, the terms in the expansion are decreasing in size with increasing n for any σ > 0, so that convergence of the Magnus expansion also occurs for much larger time steps such as t = O ( √ ε). This is a considerable improvement over earlier analysis.
Different versions of the Magnus expansion based on these observations have been successfully combined with Zassenhaus splittings to yield efficient methods for solving the TDSE with time-varying potentials , Bader et al. 2015. The skew-Hermitian nature of elements of H once again leads to unitary evolution of our methods. As noted in Theorem 9, the cost of these methods also grows quadratically in the order desired. This should be contrasted against Yoshida based methods-instead of requiring 3 exponentials, the Strang splitting features 5 or more terms for the first non-trivial Magnus expansion and consequently the cost of Yoshida methods (which are derived by composing Strang splittings) grows as O (5 n ) for an order-2n accuracy. This is considerably more expensive than the O (3 n ) cost for time-independent potentials.

A finite dimensional example
While the examples of the previous sections are based in infinite dimensional algebras and function spaces, it is also possible to construct finite dimensional examples of such structures. Consider the commutative subalgebra of 2n × 2n matrices, A = M 2n (R), where I n is the n × n identity matrix and O n is the n × n zero matrix. Consider any where E ij ∈ M n (R). Then for every x ∈ C, Thus, d is in the Lie idealiser of C, D = ad d is a derivation on C, and, consequently, is a Z 2 -graded matrix Lie algebra. Similar structures should be found in other commutative subalgebras of M n (R). Applications of these structural observations to linear algebra algorithms are yet to be explored.

Future work
1. Significantly relaxed bounds, improving on the bounds implied by (Suzuki 1977, Hochbruck & Lubich 2003, Moan & Niesen 2008, Casas, Murua & Nadinic 2012 for the convergence of the sBCH, Magnus and Zassenhaus expansions in the case of semiclassical TDSEs have been obtained due to the property of height reduction. Tighter bounds will require a careful analysis of the growth of the coefficients π k,l n,i and derivatives of functions in (5.4).
2. The TDSEs analysed here are one-dimensional. Extending these ideas to ∂ t u = iε∆u − iε −1 V (x)u will require an analysis of the free Lie algebra generated by polynomials in the commuting elements of the Lie idealiser, ∂ x , ∂ y and ∂ z .

The Wigner equation
is a pseudo-differential operator, features operators which are considerably more complicated (Jin et al. 2011). Yet, the Lie algebraic structure is similar to F with the role of the isomorphism played by Θ. In addition to ∂ x , multiplication by the frequency variable ξ is also an element of the Lie idealiser. While ∂ x induces the derivation ∂ x as usual, ξ induces the derivation i∂ y .
4. Matrix-valued potentials appear in TDSEs when we need to consider multiple energy levels at once. Magnus expansion based methods such as (Hochbruck & Lubich 2003) are used once these potentials start featuring time-varying components (Kormann, Holmgren & Karlsson 2008).
Unfortunately, matrix-valued potentials do not directly fall into the framework of F as a Lie algebra. However suitable extensions for these contexts are being actively sought. Our initial findings suggest that a modified version of the height reduction rule holds and it might be possible to extend the significantly milder time step restrictions and favourable computational complexity of (Bader et al. 2014, Bader et al. 2015 to these cases.

A Tables of coefficients
The coefficients π k,l n,i appearing in Lemma 1 for k + l ranging between 1 and 4 are presented in Table 1. The case x 0 · y 0 = xy 0 is trivial and not listed. From Lemma 4 we know that π k,l n,i = (−1) n π l,k n,n−i , so that specifying the rows (1, 3) and (0, 4), for instance, would be redundant.
Similarly, in the context of the function space S and differential operator d = ∂ x , (1, 0) 0 1 0 (2, 0) 0 2 0 (1, 1) 0 1 −1 (3, 0) 0  is found by substituting in (5.4) the coefficients from row (k, l) = (2, 1). Note that these brackets are linear, αx + βy k = α x k + β y k . To prove Lemma 2, we need to prove that A in (3.9) is such that R k,l a,p = L k,l a,p is satisfied. This certainly holds if q n=0 q i=0 S q (a,p),(n,i) A q (n,i),(s,j) = δ a,s δ p,j , (B.2) holds for any q and a, s, p, j ∈ {0, . . . , q}. To prove this we note that S and A, and therefore their product SA, are lower triangular. Thus we may concern ourselves solely with the case 0 ≤ s ≤ a ≤ q and 0 ≤ j ≤ p ≤ q. Denoting (SA) q (a,p),(s,j) as T for brevity, can be separated into four parts, δ a,n δ p,i δ n,s δ i,j = δ a,s δ p,j , In the case of T 4 we note that the binomial coefficients vanish except for s ≤ n ≤ a and j ≤ i ≤ p, when expanding them leads to the expression In order to reduce this expression, will need the following identities, where δ b>0 is 1 if b > 0 and 0 otherwise. Using these identities, which we prove at the end of this appendix, This allows us to complete the proof of Lemma 2 by proving (B.2), in proving which we have used P 1 = 1/2 and the fact that SA is lower triangular (thus the only case of a ≤ s that we need to consider is s = a). This completes the proof of Lemma 2. Proofs for the identities (B.3), (B.4) and (B.5) are given below.
which is well known, where we use the fact that P 0 = 0, and Two standard tricks for exchanging summations that we exploit are With these tools, we proceed to seek expressions for the generating functions h j , writing h j as shorthand for h j (u, w, y, x). where we have changed summation limits on n using (C.6), and for i since r ≤ n and r i = 0 for i > r ≥ 0. Using (C.4) and changing limits again using (C.5), y n n r=0 P r+1 (r + 1)!(n − r)! w n−r ((w + u)(1 + x)) r = − 1 1 − (w + u) (y(w + u)) n n r=0 P r+1 (r + 1)!(n − r)!
x j (1 + x) r w j (w + u) n−j , since r i−j vanishes unless j ≤ i ≤ r + j and r + j ≤ n. Once again, exchanging limits twice using (C.5), once for r and once for j, e z e z + 1 .