We formulate a unitary perturbation theory for quantum mechanics inspired by the Lie-Deprit formulation of canonical transformations. The original Hamiltonian is converted into a solvable one by a transformation obtained through a Magnus expansion. This ensures unitarity at every order in a small parameter. A comparison with the standard perturbation theory is provided. We work out the scheme up to order ten with some simple examples.
Historically, classical mechanics played a crucial role in the early formulation of quantum mechanics. Afterwards, the new ideas led to a quantum formalism that loosened contact with the classical language. But nevertheless, every now and then, mutual relationships between them have been profitably considered.
Fixing our attention on the treatment of the effect of a perturbation on a solvable physical system, matrix mechanics adapted, with daring modifications and reinterpretations, the classical lore of canonical transformations, generating functions, action and angle variables and the like. This is apparent, e.g. in the famous paper by Born et al. (1926). Nowadays, both time-independent and time-dependent quantum perturbation theories are presented in the Hilbert space formalism and make no mention of their classical counterparts.
Suppose one is interested in a system with Hamiltonian H(t,ϵ)=H0+H′(t,ϵ), where H0 is time-independent and H′(t,ϵ) is a time-dependent perturbation, depending on a small parameter ϵ that controls the intensity of the perturbation in such a way that H′(t,ϵ=0)=0. The aim is then to determine the quantum time-evolution operator U(t,t0), which takes the state of the system from time t0 to time t. One of its most salient features is its unitary character. The evolution operator satisfies the Schrödinger equation which, in the conventional approach, is solved by series expansion in ϵ. After truncation at a given order, it allows approximation of the different observable magnitudes of physical interest. The main drawback of the method is that any finite order approximation of U(t,t0) ceases to be unitary.
By contrast, the procedure followed in the founding years of quantum mechanics and directly inherited from the classical perturbation theory was more geometric in character. The power series expansions were also the basic tool, but the goal was to construct a geometric (unitary) transformation that renders solvable the Hamiltonian H(t,ϵ) (e.g. Birtwistle 1928). The reader acquitted with the Hamilton–Jacobi equation will recognize the idea of finding a canonical transformation (namely, a change of phase space coordinates) that transforms a Hamiltonian H into a Hamiltonian K easily solvable in the new coordinates. The problem is then translated to the search of the corresponding generating function.
The original canonical transformations theory used in classical mechanics relies on the effective use of generating functions expressed in terms of mixed, old and new, variables (e.g. Arnold et al. 2006). It means that, once the equations of motion have been integrated, one has to solve implicit functional equations so as to express the results in terms of either the old or the new variables. In general, this is difficult to do.
This classical canonical theory, however, fitted in the late 1960s into a Lie-algebraic setting, mainly in connection with problems of celestial mechanics. In that respect, the paper (Deprit 1969) Canonical transformations depending on a small parameter played a seminal role by introducing what it is known nowadays as the Lie-Deprit perturbation method. Essentially, a canonical transformation is generated by a new ‘Hamiltonian’ W in such a way that the shift by ‘time’ ϵ along the trajectories of W produces the required transformation. In this way, no mixed variables appear in the formalism. The algorithm can be readily implemented in a computer, thus allowing carrying out of perturbative analysis up to high order.
Since then, techniques and procedures of classical perturbation theory have been used in a quantum mechanical context from time to time (e.g. Kummer 1971; Eckhardt 1986; Scherer 1994; Daems et al. 2004). In fact, this is quite natural, given the existing parallelism between perturbation theory in classical and quantum systems when both are formulated in terms of geometric transformations (canonical maps and unitary transformations, respectively). The common formal setting is indeed provided by the language of Lie algebras. From this perspective, the perturbation problem defined by H0+H′(t,ϵ) can be solved by changing through a (canonical or unitary) transformation to a Hamiltonian K whose dynamics is easier to solve, at least up to a certain order in the ϵ parameter. When the perturbation expansion converges, the original problem is completely solved. Even when the convergence problem is not settled, the procedure very often provides very accurate results.
The purpose of the present paper consists in adapting and developing an implementation of the Lie-Deprit method that gives results amenable to analyse time-dependent problems in quantum mechanics. For this reason, we have adapted the title of the 1969 Deprit paper. The procedure leads to operator linear differential equations, which we solve by use of the Magnus expansion. This will ensure unitarity at any order of truncation in ϵ. This is an important feature of the algorithm proposed here. Thus, it is guaranteed that the constructed approximation shares several important qualitative features with the exact solution. In addition, in this formalism, several approximation schemes are possible, better suited to the mathematical nature of the problem at hand. We will explain how to take advantage of this. The connection with the standard time-dependent perturbation theory will also be provided. Here, we will proceed at a purely formal level, keeping the technical requisites to a minimum.
The plan of the paper is the following: in §2, we briefly summarize the standard time-dependent perturbation theory in quantum mechanics in the interaction picture defined by H0. In §3, we develop the new unitary formalism, paying special attention to the construction of the transformation and the derivation of the basic equation. The procedure is highly flexible in the sense that several possible choices of the new Hamiltonian K are possible. Also the time-independent case fits naturally into the algorithm. These issues are analysed in §4, whereas §5 is devoted to a pair of examples to illustrate the main features of the scheme. Finally, §6 contains an additional discussion and further comparisons with the so-called quantum averaging approach.
2. Time-dependent perturbation theory in quantum mechanics
Let us consider the dynamics governed by the Hamiltonian 2.1where we have assumed that H′(t,ϵ) is analytic in ϵ. In most situations, the power series expansion of H′(t,ϵ) collapses to the first term: H′(t,ϵ)=ϵH1(t).
The corresponding time-evolution operator U(t,t0) is unitary, i.e. U(t,t0)U†(t,t0)=U†(t,t0)U(t,t0)=I, where the dagger means adjoint. This property guarantees, in physical terms, probability conservation. Although for simplicity we have not shown it explicitly, U(t,t0) depends also on ϵ. In turn, it obeys the Schrödinger equation 2.2
When the Hamiltonian is independent of time, i.e. H(t)=H0, the solution of (2.2) reads . Otherwise, the solution is by no means so simple, unless and H(t) do commute.
For the Hamiltonian (2.1), the problem can be approximately solved by the conventional time-dependent perturbation theory, yielding the solution of (2.2) in power series of ϵ. In this setting, it is very important to isolate the piece of the problem to be treated approximately and allow the solvable part of the Hamiltonian to be exactly integrated. A way to deal with this question consists in transforming the problem into the interaction picture, i.e. we factorize 2.3where the unknown operator UI obeys the Schrödinger equation 2.4with the new Hamiltonian given by 2.5In this way, the interaction picture defined in (2.3) provides a unitary transformation that allows the mathematical integration of the solvable piece H0 in H(t,ϵ) and focuses the approximations just on H′(t,ϵ). Explicitly, one deals with the expansion 2.6where Un(t,t0) stands for the contribution of order n.
Substitution of (2.6) in (2.4) and then equating terms of the same power ϵn yields the successive terms Un. In particular, the first two read 2.7Unfortunately, neither the approximant UI≃I+U1 nor any of higher finite order are unitary. This, of course, has undesirable consequences in practice. For instance, the computed transition probability between different quantum states may exceed unity. Only when the infinite series is summed up, the unitary character is restored.
3. Unitary time-dependent perturbation theory
(a) The formalism
We will start by rephrasing, in modern quantum mechanics language, the main idea of perturbation theory such as it was first formulated in matrix mechanics. It will help us to introduce the central issue of the present article.
Suppose we have a hard time-independent quantum problem defined by the Hamiltonian H(ϵ)=H0+ϵH1+ϵ2H2+⋯ , and assume the mathematical problem defined by H0 has been solved. This is equivalent to finding a basis of the Hilbert space such that H0 is diagonal. Then the idea is to look for a unitary transformation T(ϵ) (or ‘canonical transformation’ in the language of matrix mechanics) such that the new Hamiltonian K=T†HT is diagonal. The transformation T can be obtained recursively as a power series in the perturbation parameter ϵ. Notice that the starting point of the method is rather geometric in nature, in clear contrast with time-dependent perturbation theory that is algebraic from the very beginning. Of course, an appropriate choice of K is essential.
In the following, we consider the more general, time-dependent problem defined by (2.1). Assume that the dynamics corresponding to H0 has been solved, i.e. we have determined Let us now define 3.1where T(t,ϵ) is a unitary transformation we want to construct in order to solve our problem.
The Schrödinger equation for UK(t,t0) reads 3.2where the new Hamiltonian is 3.3
Our goal is to determine T(t,ϵ) in such a way that (3.2) is easier to solve than the original equation (2.2). Notice that the factorization (3.1) may be considered as a generalization of (2.3). Very often we will take T(t0,ϵ)=I in (3.1) and t0=0.
In perturbation theory, the perturbed and the non-perturbed Hamiltonians are always very close to each other in some mathematical sense. It is right that what makes meaningful perturbative computations. It is then natural to determine T(t,ϵ) as a unitary near-identity transformation Two different ways of guaranteeing these features are the following:
— Introduce a skew-Hermitian operator such that .
— Alternatively, introduce a skew-Hermitian operator L(t,ϵ) such that T(t,ϵ) is the solution of the operator differential equation 3.4
This second prescription is precisely the approach proposed by Deprit in classical mechanics, and the one we follow here. Somehow, the perturbed Hamiltonian may be regarded as the evolution, with respect to the parameter ϵ, of the non-perturbed Hamiltonian (which may be thought of as an initial condition). Deprit's idea consists in materializing such an abstract ϵ-evolution in a manner similar to the usual t-evolution, since this way the unitary character of T is ensured. To this end, an evolution law similar to the one in (2.2) is introduced. In (3.4), L(t,ϵ) is the (unknown) generator of the transformation T, namely, of the ϵ-evolution. Equivalently, one has 3.5
The assessment of the situation is the following. Initially, we are given a time-dependent Hamiltonian H(t,ϵ)=H0+H′(t,ϵ), and the problem is to solve the Schrödinger equation. However, we are only able to solve the piece H0. Then we introduce the idea of an ϵ-evolution in such a way that a (unknown) continuous transformation T(t,ϵ) takes H(t,ϵ) to a new Hamiltonian K(t,ϵ) to be determined in such a way that the corresponding Schrödinger equation (3.2) can be solved. The inverse transformation T†(t,ϵ) may be used then to transform the solution of K(t,ϵ) into the one for H(t,ϵ). However, by introducing the law of the ϵ-evolution, we have added a further unknown: the generator L(t,ϵ). Thus, the situation may now appear worse than at the beginning: neither T nor L are known!
This apparently adverse scheme will provide, however, an excellent starting point to deal with power series expansions and preserve the unitary character even if the series are truncated.
We first notice that L and T† are related by the linear differential equation (3.5). Thus, once L has been determined, we can obtain T and T† by formally applying the Magnus expansion (Magnus 1954) to (3.5). As is well known, his proposal is to solve a generic operator linear equation dY (t)/dt=A(t)Y (t) with Y (0)=I in the form The exponent Ω(t) is given by an infinite series of terms that are linear combinations of integrals and nested commutators involving m-fold products of the operator A at different times (Blanes et al. 2009).
Although explicit formulae for Ωm exist, it is sometimes more convenient to generate them by a recursive procedure. The one proposed in Klarsfeld & Oteo (1989) is particularly well suited from a computational point of view (Blanes et al. 2009): and 3.6where Bj stand for Bernoulli numbers and [⋅,⋅] represents the commutator, [A,B]≡AB−BA. Notice that this approach can be directly applied to obtain approximate formal solutions of equation (3.5), so that we can write 3.7where Ω is a skew-Hermitian operator. Thus, the two above enumerated alternatives to render T unitary are, in fact, equivalent.
Next, we derive equation (3.3) with respect to ϵ, 3.8Since this expression involves both the operators T and L, it is convenient to insert (3.7) here. In fact, it is straightforward to show that 3.9Here we have defined the operator adΩB≡[Ω,B], with . In consequence, we have 3.10an equation formulated only in terms of L,H and K that is easier to handle than (3.8).
In the present formalism, then, three different problems have to be studied:
Choose the new Hamiltonian K such that equation (3.2) is easy to solve.
Compute the skew-Hermitian generator L of the required transformation to the new picture.
Construct the unitary transformation T from the generator L, or equivalently, the operator Ω in .
The first two problems enumerated earlier can be solved perturbatively with equation (3.10), whereas the third one can be tackled independently. We treat these problems in reverse order.
(b) Determining T
We next proceed to compute Ω, once L has been fixed. To do that, we introduce, in addition to the expression (2.1) for the Hamiltonian H, the following series expansions: 3.11and 3.12It turns out that Ω(t,ϵ) in (3.7) can also be determined as a power series in ϵ: 3.13with vn(t) expressed in terms of Lj(t). This can be done algorithmically with the Magnus expansion up to any given order. In fact, the procedure already used to construct numerical integrators based on Magnus expansion (Blanes et al. 2009) is of use here also. Specifically, if the series (3.12) is inserted into the recurrence (3.6), we get up to order ϵ5 Reordering according with the power of ϵ, we have the following expressions for the first terms in the series (3.13): 3.14This procedure can be carried out for higher orders with a symbolic algebra package. On the other hand, if the series (3.13) is inserted into (3.9), then 3.15and, again, this series can be determined algorithmically. The first terms read 3.16It is worth stressing here that the series (3.13) and (3.15) have to be determined only once. In fact, computing these series up to n=10 takes only a few minutes on a personal computer. Terms up to n=10 and n=9, respectively, can be found at the website (http://www.gicas.uji.es/Research/Perturbation.html) as Mathematica files.
(c) Determining L
When the series expansions for K, L and (3.15) are substituted into (3.10), collecting terms of the same power in ϵ yields a system of recurrence equations. Specifically, we arrive at K0=H0 and 3.17with 3.18The next step in the procedure is to propose a suitable Kn and solve (3.17) to determine Ln and finally the transformation T, which is unitary at any order by construction.
Equation (3.17) has the following formal solution: 3.19where n>0, and we have taken t0=0 for simplicity. Very often we will fix Ln(0)=0, so that 3.20
4. On the choice of K
The choice of a particular K is a degree of freedom of the method but it may also constitute a nuisance. It might be the case that two particular H and K are associated in a natural way. Here we will comment about a number of possible choices of K in the general case.
Perhaps, this is the most immediate option. It means that systematically Equivalently, one tries to construct a unitary transformation in such a way that in the new picture there is no perturbation at all. In that case, and 4.1Notice that with this choice, the whole ϵ dependency is contained in T.
(b) K diagonal
Suppose that H0 has a pure non-degenerate point spectrum. In that case, we can choose Kn to be diagonal. In other words, we perturbatively construct a unitary transformation to a new picture where the Hamiltonian K is diagonal and therefore equation (3.2) is easy to solve. Specifically, and finally
(c) Time-independent case
The formalism developed earlier works also when H′ is time-independent. In that case, it makes sense to construct a time-independent unitary transformation T(ϵ) so that A natural choice is then to take the new Hamiltonian so that [H0,K]=0. This is analogous to the standard procedure of constructing the Birkhoff–Gustavson normal form, both in the classical and quantum setting (Ali 1985; Eckhardt 1986). From (3.17), the equations to solve read then 4.2It has been shown in Jauslin et al. (2000) that if these equations are defined in the space of bounded skew-Hermitian multiplication operators, then the formal solution can be written in terms of averages as and even if H0 has degenerate eigenvalues. As a matter of fact, these expressions reproduce the standard results (Birtwistle 1928; Wu 1986). If |ν,j〉 denotes an eigenbasis of H0, where ν labels the distinct eigenvalues and j distinguishes different basis vectors of the degeneracy subspace, then whereas The operator evolution in this setting reads being unitary even upon truncation of the series. This coincides with the formal solution proposed by Kummer (1971) to avoid secular terms in perturbative approximations. In this sense, the formalism developed here could be considered as a generalization to time-dependent problems of Kummer's approach.
In this section, we apply the previous formalism to two simple examples with the aim of illustrating its main features.
First, we consider the Hamiltonian 5.1where σj denote the Pauli matrices and ϵ, ω0, ω≠ω0 are real parameters. The exact solution of the Schrödinger equation (2.2) with U(0)=I is then Let 1 and 2 denote the spin up and down states, respectively. Then the exact transition probability from state 1 to state 2 is given by and may serve as a test to analyse the various results in the following.
Our perturbative scheme starts by writing H(t)=H0+ϵH1(t), with We have applied the procedure up to n=10 twice. First by choosing Kn=0 for n≥1 and second by taking Kn as a diagonal matrix. In either case, we get and the corresponding evolution operator UK(t), whereas leads to Ω(t,ϵ) up to this order. Finally,
In figure 1, we depict the difference between the exact and approximate transition probabilities as a function of time computed with this procedure up to n=4 and n=10. Solid lines correspond to choosing K=H0, whereas dot-dashed curves are generated by taking Kn as a diagonal matrix. Notice that increasing the number of terms in L(t,ϵ) in both cases provides a more accurate approximation.
The second example illustrates the time-independent case. We take the two-level spin Hamiltonian Then the corresponding equations (4.2) provide for the first orders the following results: Carrying out the computation up to order ϵ10, we get 5.2On the other hand, the exact eigenvalues of the Hamiltonian H are in this case and its expansion in power series of ϵ agrees with (5.2). Finally, Ω(ϵ)=iωaσ2, with and thus the approximate solution is given by .
(a) Connection with the standard time-dependent perturbation theory
The formalism we have developed here contains as a limiting case the usual time-dependent perturbation theory summarized in §2. This indeed corresponds to choosing K=H0 as a new Hamiltonian and expanding the operator T(t,ϵ)=e−Ω(t,ϵ) in power series of ϵ. Obviously, the resulting scheme is no longer unitary when the series is truncated.
Let us illustrate this point by computing explicitly the first order in ϵ. Since , it is clear from (3.20) that whereas But then with This, as we know, is precisely the result achieved by standard perturbation theory in the interaction picture defined by H0.
The previous treatment can be generalized to any order in ϵ. In fact, by following this approach, we have explicitly reproduced up to order ϵ6 the results achieved by the usual time-dependent perturbation theory for example 5.1.
(b) Quantum averaging
In a series of papers published in the mid-1990s, Scherer used the averaging method to construct quantum mechanical analogues of the classical perturbation expansions by Poincaré and Kolmogorov, both in the autonomous case and for time-dependent operators (Scherer 1994, 1997, 1998). The idea was further elaborated by Daems et al. (2003a, 2004) and Sugny et al. (2004), and successfully applied in the perturbative treatment of pulse-driven quantum problems. This quantum averaging algorithm can also be reproduced with our formalism as follows.
Scherer's basic proposal also consists of transforming the original Hamiltonian (2.1) with the help of a unitary transformation T(t,ϵ) such that the problem of finding the time evolution of the transformed Hamiltonian K(t,ϵ) is easier than the original one. In particular, he chooses the terms Kn in the series defining K so as to satisfy 6.1Then, it can be shown that the time evolution generated by K(t,ϵ) is given by (t0=0) In other words, the time evolution is reduced to a time-independent problem if one introduces an interaction picture defined by H0 for the new Hamiltonian K. Finally, 6.2whereas equation (6.1) can be solved by averaging, generalizing the autonomous case: 6.3The unitary transformation T(t,ϵ) is then computed order by order from Ln(t). In fact, a recursive procedure is proposed for determining the terms Tn(t) in the series 6.4If not all terms Ln can be computed, the corresponding truncated series (6.4) provides approximations that are no longer unitary. It is argued, however, that these violations of unitarity do not grow with time, since they do not contain secular terms.
Notice that Scherer's algorithm fits quite naturally in the present scheme. Moreover, with the choice (6.3), we can determine Ln(t) by applying (3.20) and then compute the transformation as with the series (3.13). Now, contrarily to the original formulation, the resulting scheme is unitary upon truncation.
We have applied this modified scheme again to the Hamiltonian (5.1) and compared with the previous choices of Kn=0 and Kn diagonal (n>0). In figure 2, we represent the corresponding errors in the transition probability obtained with the three procedures when n=10 terms are considered in the series. The results illustrate the fairly different character of each approximation. It should be stressed again that the approximation (6.2) with (6.3) considered here is a unitary modification of that originally proposed by Scherer (1997, 1998), i.e. with instead of (6.4).
(c) Open problems and further comments
There are several issues only touched upon in this work that deserve a much more detailed analysis. First, we have been only concerned here with the formal aspects of the transformations and expansions. It is clear, then, that a more rigorous treatment is necessary, with special attention to the convergence of the expansion. Second, the procedure developed here is fairly flexible since it enables different choices of the new Hamiltonian. We have seen on a simple example that different options lead indeed to different approximations; so the natural question is: for any particular problem, which is the optimal choice? In other words, which particular K allows one to get the best approximation once the expansion is truncated at a given order? In this sense, it is worth mentioning that there are several free parameters in the formalism that eventually may be used to tune the algorithm and improve the accuracy (e.g. Daems et al. 2004).
In dealing with quantum time-dependent problems, one also finds in the literature the methods that are iterative in character. So, for example, when the Hamiltonian is slowly varying with time, sequences of unitary transformations intended to diagonalize instantaneously the Hamiltonian lead to non-perturbative approaches (Berry 1987, 1990). In turn, the so-called Fer expansion (Blanes et al. 2009) is based on the iterative use of the transformation that defines the interaction picture. A perturbation theory can also be used in an iterative way: the Hamiltonian obtained in each step is used as unperturbed Hamiltonian for the next step. This quantum mechanical version of the KAM scheme has been shown to produce very accurate results in treatments, different to the one proposed in this paper, both for time-independent and time-dependent quantum perturbations (Scherer 1997; Daems et al. 2003a,b, 2004). The same idea could be implemented with the procedure developed here. In this work, by following closely the Lie-Deprit method in classical mechanics, we have developed a unitary perturbative algorithm well adapted to quantum problems. The technique provides a, in practice approximate, single unitary transformation designed to simplify the dynamics as much as possible.
This work has been partially supported by MICINN (Spain), under grant no. AYA2010-22111-C03-02 (J.A.O. and J.R.) and MTM2010-18246-C03-02 (F.C.), co-financed by FEDER Funds of the European Union.
- Received June 20, 2011.
- Accepted October 4, 2011.
- This journal is © 2011 The Royal Society