Quantum Statistical Calculations and Symplectic Corrector Algorithms

0.16 MB·English

by  Siu A. Chin

Quantum Statistical Calculations and
Symplectic Corrector Algorithms
Siu A. Chin
Department of Physics, Texas A&M University, College Station, TX 77843, USA
4
0 The quantumpartition function at finitetemperaturerequirescomputingthetraceof theimagi-
0 nary time propagator. For numerical and Monte Carlo calculations, the propagator is usually split
2 intoitskineticandpotentialparts. Ahigherordersplittingwillresult inahigherorderconvergent
algorithm. Atimaginarytime,thekineticenergypropagatorisusuallythediffusionGreensfunction.
n
Sincediffusion cannot besimulated backward in time, thesplittingmust maintain thepositivity of
a
all intermediate time steps. However, since the trace is invariant under similarity transformations
J
of thepropagator, one can usethis freedom to “correct” thesplit propagator tohigher order. This
1
useofsimilaritytransformsclassically giverisestosymplecticcorrectoralgorithms. Thesplitprop-
2
agatoristhesymplectickernelandthesimilarity transformation isthecorrector. Thiswork proves

a generalization of the Sheng-Suzukitheorem: no positive time step propagators with only kinetic
]
h and potential operators can be corrected beyond second order. Second order forward propagators
c canhavefourthordertracesonlywiththeinclusion ofanadditionalcommutator. Wegivedetailed
e derivations of four forward correctable second order propagators and their minimal correctors.
m
-
at I. INTRODUCTION
t
s
. The quantum partition function requires computing the trace
t
a
m Z =Tr(ρ)=Tr(e−βH), (1.1)
-
d whereρisimaginarytimepropagator,β =1/(k T)istheinversetemperatureandH =T+V istheusualHamiltonian
B
n
operator. Although specific forms of the kinetic and potential energy operators will not be used in the following, it
o is useful to keep in mind the many-body case where T = ( ¯h2/2m) 2 and V = v(r ). In numerical or
c − i=1∇i i<j ij
[ Monte Carlo calculations, the imaginary time propagator is first discretized as
P P

2 n
e−β(T+V) = eε(T+V) , (1.2)
v
1 h i
2 whereε= ∆β = β/n,andtheshort-timepropagatoreε(T+V)isthenapproximatedinvariousways. Onesystematic
0 − −
method is to decompose, or split the short-time propagatorinto the product form
2
1
N
3 eε(T+V) etiεTeviεV, (1.3)
0 ≈
/ iY=1
t
a
with coefficients t ,v determined by the required order of accuracy. For quantum statistical calculations, since
m r′ etiεT r e−({r′−ir)2i}/(2ti∆β) isthediffusionkernel,the coefficientt mustbepositiveinorderforitto besimulated
i
- h | | i∝
d or integrated. If ti were negative, the kernel is unbounded and unnormalizable, and no probabilistic based (Monte
n Carlo) simulation is possible. However, as first proved by Sheng1, and later by Suzuki2, beyond second order, any
o factorizationoftheform(1.3)mustcontainsomenegativecoefficientsintheset t ,v . GoldmanandKaper3 further
i i
c provedthatanyfactorizationoftheform(1.3)mustcontainatleastonenegativ{ecoeffi}cientforbothoperators. Thus,
:
v despitemyriadoffactorizationschemesoftheform(1.3)proposedintheclassicalsymplecticintegratorliterature4,5,6,7,
i nonecanbeusedfordoingquantumstatisticalcalculationsbeyondsecondorder. Itisonlyrecentlythatfourthorder,
X
all positive-coefficient factorization schemes have been found8,9 and applied to time-irreversible problems containing
r the diffusion kernel10,11,12,13,14. In order to bypass the Sheng-Suzuki’s theorem, one must include other operators in
a
the factorization (1.3), such as the double commutator [V,[T,V]], where [A,B] AB BA.
Incomputingthe quantumpartitionfunctionZ,onlythe traceofρ=e−βH is≡requir−ed. Since the traceisinvariant
under the similarity transformation
ρ˜=SρS−1, (1.4)
one is free to use any such ρ˜ to compute Z. This is immaterial if ρ is known exactly. However, if the short-time
propagatorisonlyknownapproximately,thenonemayuseacleverchoiceofS tofurtherimprovetheapproximation.
This is a well known idea in many areas of physics. For example, to calculate the exact quantum many-body ground
2
state using the Diffusion Monte Carlo algorithm, one can choose S = φ , where φ is a known trial function close
0 0
to the exact ground state. This is the idea of “importance sampling” as introduced by Kalos et al.15. Its operator
formulation as described above has been implemented by Chin16 some time ago. Similar ideas have been used to
improve path-integrals, as detailed by Kleinert17. If the short time propagator is approximated by the product form
(1.3), the error terms can be calculated explicitly and eliminated by S. When implemented classically, these are
knownassymplectic“corrector”,or“process”algorithms18,19,20,21,22,23. Inthiscontextthepropagatorρisthekernel
algorithm and S is the corrector. Since S disappears in the calculation of Z, there is no restriction on the form of S.
If S were also expanded in the product form (1.3), there is no restrictionon the sign of its coefficients. This suggests
that there may exist a product form (1.3) of ρ with only positive coefficients such that its trace is correct to higher
order. This would not be precluded by the existing Sheng-Suzuki theorem.
In this work, we show that this is not possible. If ρ is approximated by the product form (1.3) with positive
coefficients t , then ρ˜ cannot be corrected by S to higher than second order. The proof of this generalizes the
i
{ }
Sheng-Suzukitheorem. Thecorrectedpropagatorρ˜canbefourthorderonlyifadditionaloperators,suchas[V,[T,V]],
are used in the splitting of ρ. By understanding the “correctability” requirement, we can systematically deduce the
four fundamental correctable second order propagatorsand their correctors.
In the following Section, we recall some basic results of similarity transforms. Beyond second order, only a special
class of approximate ρ satisfying the “correctability” condition can be corrected to higher order. In Section III, we
computetheexplicitformoftheerrorcoefficientsrequiredbythecorrectabilitycriterion. InSectionIV,weshowthat
this requirement cannot be satisfied for propagatorsof the product form (1.3) with only positive t coefficients. In
i
{ }
Section V, based on our understanding of the correctability restriction, we deduce all four second order correctable
propagators and their minimal correctors. Some conclusions are given in Section VI.
II. SIMILARITY TRANSFORMS AND THE CORRECTABILITY CRITERION
Similarity transforms on approximate propagators of the product form (1.3) have been studied extensively in the
context of symplectic correctors18,19,20,21,22. However, not all use the language of operators and some are specific to
celestialmechanics. Here,we recallsome elementaryresultsandestablishthe fundamentalcorrectabilityrequirement
in the context of quantum statistical physics.
Since
n
SρS−1 = Seε(T+V)S−1 , (2.1)
h i
it is sufficient to study the similarity transformation of the approximate short-time propagator ρ . Let ρ approxi-
A A
mates eε(T+V) in the product form such that
N
ρ = etiεTeviεV =eεHA, (2.2)
A
i=1
Y
where H is the approximate Hamiltonian
A
H =T +V +ε(e [T,V])+ε2(e [T,[T,V]]+e [V,[T,V]])+O(ε3) (2.3)
A TV TTV VTV
witherrorcoefficients e , e ,e determinedby factorizationcoefficients t ,v . The transformedpropagator
TV TTV VTV i i
{ }
is
ρ˜ =Sρ S−1 =SeεHAS−1 =eε(SHAS−1) =eεHA, (2.4)
A A
where the last equality defines the transformed approximate Hamiltonian HA. Ief now we take
S =exp[εC] (2.5)
e
where C is the to-be-determined corrector,then we have the fundamental result
1 1
H =eεCH e−εC =H +ε[C,H ]+ ε2[C,[C,H ]]+ ε3[C,[C,[C,H ]]]+ . (2.6)
A A A A A A
2 3! ···
Let’s first consieder the case where the product form (2.2) for HA is left-right symmetric, i.e., either t1 = 0 and
v = v , t = t , or v = 0 and v = v , t = t . In this cases, the propagator is reversible,
i N−i+1 i+1 N−i+1 N i N−i i N−i+1
ρ (ε)ρ ( ε)=1, and H (ε) is an even function of ε with e =0. In this case
A A A TV
−
H = H +ε[C,H ]+ ,
A A A
···
= T +V +ε2(e [T,[T,V]]+e [V,[T,V]])+ε[C,T +V]+ , (2.7)
TTV VTV
e ···
3
and one immediately sees that the choice C = εC with C c [T,V] would eliminate either second order error
1 1 TV
≡
term with c =e or c =e . So, if H is constructed such that
TV TTV TV VTV A
e =e (2.8)
TTV VTV
thenbothcanbe simultaneouslyeliminatedbythecorrector. Thisisthe fundamental“correctability”requirementfor
correcting a second order ρ to fourth order. This observation can be generalized to higher order. At higher orders,
A
H will have error terms of the form [T,Q ] and [V,Q ] where Q are some higher order commutator generatedby T
A i i i
and V. If H is of order 2n in ε, then H can be of order 2n+2 only if H ’s errorcoefficients for [T,Q ] and [V,Q ]
A A A i i
are equal for all Q ’s. This fundamental correctorinsight is often obscured by the more general case where odd order
i
errors are allowed. e
Sheng1 and Suzuki2 independently proved that no ρ of the form (2.2) can have positive coefficients t beyond
A i
second order. More precisely, if ρ is of the product form (2.2) with positive t ’s such that e =0, then e and
A i TV TTV
e cannot both be zero. We will prove a more general theorem that the product form (2.2) with positive t ’s such
VTV i
that e = 0 cannot be corrected beyond second order, i.e., e can never equal to e . From this perspective,
TV TTV VTV
the Sheng-Suzuki theorem is a special case where the common value for both coefficients is zero.
In the general case where e =0, we have
TV
6
H =T +V +ε(e [T,V])+ε2(e [T,[T,V]]+e [V,[T,V]])
A TV TTV VTV
1
+ε[C,T +V]+ε2e [C,[T,V]]+ ε2[C,[C,T +V]]+O(ε3). (2.9)
e TV 2
Since [c T +c V,T +V]=(c c )[T,V], the liner term in ε can be eliminated if we choose C =C c T +c V
T V T V 0 T V
− ≡
such that
(c c )= e . (2.10)
T V TV
− −
This is the first order correctability condition. This means that with a suitable choice of c and c , a first order
T V
propagator can always be corrected to second order. Hence, the trace of any first order propagator is always second
order. For example, the trace Tr(eεTeεV) is second order despite its appearance.
With the first order correctability condition satisfied, the remaining commutators in (2.9) are either [T,[T,V]] or
[V,[T,V]], and can again be corrected by adding to C the term εC =εc [T,V]. Thus with
1 TV
C =C +εC =c T +c V +εc [T,V] (2.11)
0 1 T V TV
such that (c c )= e , we have
T V TV
− −
H =T +V +ε2(e [T,[T,V]]+e [V,[T,V]])
A TTV VTV
1
+ε2[C ,T +V]+ε2e [C ,[T,V]]+ ε2[C ,[C ,T +V]]+O(ε3),
e 1 TV 0 2 0 0
1
=T +V +ε2(e c + c e )[T,[T,V]]
TTV TV T TV
− 2
1
+ε2(e c + c e )[V,[T,V]]+O(ε3) (2.12)
VTV TV V TV
− 2
If we now choose c = e + 1c e to eliminate the error term [T,[T,V]], then the error term [V,[T,V]] can
TV TTV 2 T TV
vanish only if
1
e =e + (e )2. (2.13)
TTV VTV TV
2
This is the general second order correctability requirement for correcting any first order propagator beyond second
order. The major result of this work is to show that this condition cannot be satisfied for product decomposition of
the form (2.2) with only positive t coefficients.
i
III. DETERMINING THE ERROR COEFFICIENTS
To check whether the correctability requirement (2.13) can ever be satisfied by an approximate propagator of the
product form (2.2), we need to determine e , e and e in terms of t ,v . From the assumed equality
TV TTV VTV i i
{ }
N
etiεTeviεV =eεHA, (3.1)
i=1
Y
4
with H given by (2.2), we can expand both sides and compare terms order by order in powers of ε. The left hand
A
side of (3.1) can be expanded as
N N
eεt1Teεv1Veεt2Teεv2V eεtNTeεvNV =1+ε t T +ε v V + , (3.2)
i i
··· ! ! ···
i=1 i=1
X X
and the right hand side as
eεHA =1+ε(T +V)+ε2e [T,V]+ε3e [T,[T,V]]+ε3e [V,[T,V]]
TV TTV VTV
1 1
+ ε2(T +V)2+ ε3e (T +V)[T,V]+[T,V](T +V)
TV
2 2 { }
1
+ ε3(T +V)3+ (3.3)
3! ···
Matching the first order terms in ε gives the primary constraints
N N
t =1 and v =1. (3.4)
i i
i=1 i=1
X X
To determine the error coefficients, we “tag” a particular operator in (3.3) whose coefficient contains e , e or
TV TTV
e and match the same operator’s coefficients in the expansionof (3.2). For example, in the ε2 terms of (3.3), the
VTV
coefficient of the operator TV is (1 +e ) Equating this to the coefficients of TV from (3.2) gives
2 TV
N
1
+e = s v . (3.5)
TV i i
2
i=1
X
where we have introduced the variable
i
s = t . (3.6)
i j
j=1
X
ThiswayofcomputingTV from(3.2)correspondstofirstpickingoutaV operatorfromamongallthev terms,then
i
combineallthe t termsto its leftinthe exponentialto generateaT operator. Alternatively,the samecoefficientcan
i
also be expressed as
N
1
+e = t u . (3.7)
VT i i
2
i=1
X
where
N
u = v . (3.8)
i j
j=i
X
This way of computing TV corresponds to first picking out a T operator from among all the t terms, then combine
i
all the v terms to its right in the exponential to generate a V operator. To demonstrate how these variables are to
i
be used, we can directly prove the equality of (3.5) and (3.7). First, note that s = 1 and u = 1. Second, since
N 1
t = s s , at i = 1 we must consistently set s = 0. Similarly, since v = u u , we must set u = 0.
i i i−1 0 i i i+1 N+1
− −
Therefore we have
N N N N
s v = s (u u )= (s s )u = t u (3.9)
i i i i i+1 i i−1 i i i
− −
i=1 i=1 i=1 i=1
X X X X
The determination of error coefficients is simplified if we pick operators whose expansion coefficients are easy to
calculate. Matching the coefficients of operators TTV and TVV (note, not the operator VTV) yields
N N
1 1 1 1
+ e +e = s2v = (s2 s2 )u , (3.10)
6 2 TV TTV 2 i i 2 i − i−1 i
i=1 i=1
X X
N
1 1 1
+ e e = t u2. (3.11)
6 2 TV − TVT 2 i i
i=1
X
5
IV. PROVING THE MAIN RESULT
Using the expression for e from (3.11), the correctability requirement (2.13) reads
TVT
N
1
t u2 =a, (4.1)
2 i i
i=1
X
with
1 1 1
a= ( +e )2+ e (4.2)
TV TTV
2 2 24 −
and e , e given by (3.7), (3.10) respectively. In Suzuki’s proof2, he recognizes that in terms of the variable
TV TTV
√t u , (4.1) is a hypersphere and (3.7), (3.10) are hyperplanes. His proofis basedon a geometric demonstrationthat
i i
his hyperplane cannot intersect his hypersphere. While this geometric language is very appealing, it is cumbersome
when dealing with more than one hyperplane. We will use a different strategy.
If t are all positive, then the LHS of (4.1) is a positive-definite quadratic form in u . There would be no real
i i
solutions for u if the minimum of the quadratic form is greater than a. Our strategy is therefore to minimize the
i
quadratic form subject to constraints (3.7) and (3.10)
N
t u = b, (4.3)
i i
i=1
X
N
t (s +s )u = c, (4.4)
i i i−1 i
i=1
X
withb= 1+e ,c= 1+e +2e ,andshowthattheresultingminimumisalwaysgreaterthana. (Theprimary
2 VT 3 TV TTV
constraints (3.4) are just s =1 and u =1.)
N 1
For constrained minimization, one can use the method of Lagrange multiplier. Minimizing
N N N
1
F = t u2 λ t u b λ t (s +s )u c (4.5)
2 i i − 1 i i− !− 2 i i i−1 i− !
i=1 i=1 i=1
X X X
gives
u =λ +λ (s +s ). (4.6)
i 1 2 i i−1
Substituting this back to satisfy constraints (4.3) and (4.4) determines λ and λ :
1 2
λ +λ =b, (4.7)
1 2
λ +λ +gλ =c. (4.8)
1 2 2
The only non-trivial evaluation is N t (s +s )2 =1+g, where
i=1 i i i−1
P N
g = (s2s s s2 ). (4.9)
i i−1− i i−1
i=1
X
The minimum of the quadratic form is therefore
N
1
F = t [λ +λ (s +s )]2
i 1 2 i i−1
2
i=1
X
1
= [(λ +λ )2+gλ2]
2 1 2 2
1 1
= [b2+ (c b)2]. (4.10)
2 g −
6
To minimize F, one must maximize g. Solving ∂g/∂s = 0 gives s = (s +s )/2, which means that s is linear
i i i+1 i−1 i
in i. The normalization s =1 fixes s =i/N, giving
N i
1 1
g = (1 ). (4.11)
max 3 − N2
This is indeed a maximum since one can directly verify that ∂2g/∂s2 = 2(s s )<0. Hence, at any finite N,
i − i+1− i−1
1 1 1 3 1
F > [b2+3(c b)2]= ( +e )2+ (2e )2 =a+6e2 . (4.12)
2 − 2 2 TV 2 TTV − 6 TTV
Thus the minimum of the quadratic form is always higher than the value required by the correctability condition.
Hence, no real solutions for u are possible if t are all positive.
i i
Wenotethattheaboveproofisindependentofe . Fore =0,thecorrectabilityconditionisjuste =e .
TV TV TTV VTV
Hence forsymmetricdecompositionswithpositive t ’s,where e =0is automatic,we haveasacorollarythate
i TV TTV
can never equal to e .
VTV
V. CORRECTABLE FORWARD PROPAGATORS AND THEIR CORRECTORS
The last section is the main result of this work. Here, we show how the correctability criterion can be applied
systematically to deduce forward correctable second order propagatorsand their minimal correctors.
The proof of non-correctability is limited to the conventional product form (2.2), which factorizes the propagator
only in terms of operators T and V. As shown in the last section, symmetrically decomposed positive-time-step
propagators cannot be corrected beyond second order because e cannot be made equal to e . For example,
TTV VTV
the second order propagator
1 1
exp( εT)exp(εV)exp( εT) (5.1)
2 2
has t = t = 1/2, v = u = 1, s = 1/2 and e = 0. From (3.10) and (3.11), we can determine indeed that the
1 2 1 1 1 TV
two error coefficients are not equal:
2
1 1 1 1
e = 1 = ,
TTV
2 2 − 6 −24
(cid:18) (cid:19)
1 1 1 1
e = 1= . (5.2)
TVT
6 − 2 2 −12
(cid:18) (cid:19)
Asimple wayto forcethemequalisto directlyincorporateeither operator[T,[T,V]]or[V,[T,V]]inthe factorization
process. Since [V,[T,V]] = (h¯2/m) v(r )2 is just another potential function, Suzuki24 suggested that
i|∇i j6=i ij |
one should keep the operator [V,[T,V]]. If now we add 1 ε3[V,[T,V]] to εV in (5.1), we can change the coefficient
P P 24
e from 1/12 to 1/24,matching that of e . The result is still only a second order propagator
VTV TTV
− −
1 1 1
ρ =exp εT exp εV + ε3[V,[T,V]] exp εT , (5.3)
TI
2 24 2
(cid:18) (cid:19) (cid:18) (cid:19) (cid:18) (cid:19)
but now has a fourth order trace. This propagator was first obtained by Takahashi and Imada25,26 by directly
computingthetrace. Itisaremarkablefindgivenhowlittlethey hadtoworkwith. Thisderivationexplains,without
doing any trace calculation, why the propagatorworked.
The alternative of keeping [T,[T,V]] would require adding 1 ε3[T,[T,V]] to make e equal to e ’s value
−24 TTV VTV
of 1/12. This operator is too complicated for practical use, but in the case of the harmonic oscillator, it can be
−
combined with the kinetic energy operator:
1 1 1 1
ρ′ =exp εT ε3[T,[T,V]] exp(εV)exp εT ε3[T,[T,V]] . (5.4)
2B 2 − 48 2 − 48
(cid:18) (cid:19) (cid:18) (cid:19)
This can also be written in the form of
1 1 1
ρ =exp εV exp εT ε3[T,[T,V]] exp εV . (5.5)
2B
2 − 24 2
(cid:18) (cid:19) (cid:18) (cid:19) (cid:18) (cid:19)
7
In this case exp(1εV)exp(εT)exp(1εV) has e = 1/12 and e = 1/24 and propagator ρ corresponds to
2 2 TTV VTV 2B
changing e ’s value to match that of e . The Takahashi-Imada propagator (5.3) can also be written as
TTV VTV
1 1 1 1
ρ′ =exp εV + ε3[V,[T,V]] exp(εT)exp εV + ε3[V,[T,V]] , (5.6)
TI 2 48 2 48
(cid:18) (cid:19) (cid:18) (cid:19)
corresponding to changing e ’s value to match that of e . These are the four fundamental correctable second
VTV TTV
order propagators with a fourth order trace.
Forthecomputationofthetrace,itisunnecessarytoknowthecorrectorexplicitly. Inothercases,suchassymplectic
correctoralgorithms,onemaywishtoapplythecorrectoroccasionallytoseetheworkingofthecorrectedfourthorder
propagator ρ. We will give a detailed derivation of correctors for propagators (5.3)-(5.6), cumulating in a set of four
minimal correctors. These minimal correctors with analytical coefficients have not been previously described in the
literature18,19,20,21,22,23.
e
For the Takahashi-Imada propagator, we have e = e = e with e = 1/24. From (2.7), we see that a
TTV VTV 2 2
−
possible corrector is C =e ε[T,V]. This can be constructed in a straightforwardmanner as suggested by Wisdom et
2
al.18. Since
B(v ,t ) exp(εv V)exp(εt T)exp( εv V)exp( εt T)
1 1 1 1 1 1
≡ − −
1 1
= exp v t ε2[T,V] t2v ε3[T,[T,V]] t v2ε3[V,[T,V]]+O(ε4) , (5.7)
− 1 1 − 2 1 1 − 2 1 1
(cid:18) (cid:19)
by setting v t =(1/48), the following product is a workable corrector
1 1
1
B(v ,t )B( v , t )=exp ε2[T,V]+O(ε4) . (5.8)
1 1 1 1
− − −24
(cid:18) (cid:19)
Note that it is important to have the operator V before T to generate a negative e coefficient. However, without
2
fully determining both v and t , this corrector clearly under-utilizes B(v ,t ). It requires eight operators, which is
1 1 1 1
far from optimal. We will show below that four is sufficient.
LetH =T+V andG=[T,V]. SinceH =H+e ε2[H,G],wecanseefrom(2.7)thataddingatermc H toC,will
A 2 0
not affect the corrector term ε[C,T +V], but such a term will generate unwanted third order terms c e ε3[H,[H,G]]
0 2
from ε[C,H ] and 1c e ε3[H,[G,H]] from 1ε[C,[C,H ]]. To cancel them, we must add another term c ε2[H,G] to
A 2 0 2 2 A 2
the corrector such that c = 1c e . Thus the corrector can have the more general form
2 2 0 2
1
exp(εC) = exp c εH +e ε2G+ c e ε3[H,G] +O(ε4), (5.9)
0 2 0 2
2
(cid:18) (cid:19)
= exp(c εH)exp(e ε2G)+O(ε4), (5.10)
0 2
where the second line follows from the fundamental Baker-Campbell-Hausdorff formula, exp(A)exp(B) = exp(A+
B+(1/2)[A,B] ). To exploit the use of the free parameter c , we can approximate exp(c εH) by
0 0
···
c c
0 0
exp(ε V)exp(εc T)exp(ε V)
0
2 2
1 1
=exp c εH + c3ε3[T,[T,V]]+ c3ε3[V,[T,V]] +O(ε5), (5.11)
0 12 0 24 0
(cid:18) (cid:19)
and the term exp(e ε2G) by B(v ,t ). We can now choose c ,v ,t such that v t =1/24 and the third order terms
2 1 1 0 1 1 1 1
in (5.11) exactly cancel the third order terms in (5.7): 1t2v = 1 c3, 1t v2 = 1 c3. This gives c = 1/(2 31/6),
2 1 1 12 0 2 1 1 24 0 0 ·
v =1/(4√3) and t =1/(2√3). The result is a corrector with six operators:
1 1
c c
0 0
S =exp(ε V)exp(εc T)exp(ε( +v )V)exp(εt T)exp( εv V)exp( εt T). (5.12)
0 1 1 1 1
2 2 − −
Sincethiscorrectorhasmadegooduseofalltheparameters,itissurprisingthatonecanfindaevenshortercorrector.
Instead of B(v ,t ), consider just
1 1
exp(εd V)exp(εd T)
0 0
1 1 1
=exp d εH d2[T,V]+ d3ε3[T,[T,V]] d3ε3[V,[T,V]] +O(ε4). (5.13)
0 − 2 0 12 0 − 12 0
(cid:18) (cid:19)
8
The corrector
c c
0 0
S = exp(ε V)exp(εc T)exp(ε( +d )V)exp(εd T) (5.14)
TI 0 0 0
2 2
1 1 1
= exp (c +d )εH +( d2)ε2G+ ( d2)(c +d )[H,G]
0 0 −2 0 2 −2 0 0 0
(cid:18)
1 1
+ (c3+4d3)ε3[T,[T,V]]+ (c3+4d3)ε3[V,[T,V]] +O(ε4)
12 0 0 24 0 0
(cid:19)
will have the correct value for e if we take d2/2 = 1/24, fixing d = 1/(2√3). The corrector will also be of the
2 0 0
form (5.9) after both commutators have been eliminated by setting c3 = 4d3, giving c = 1/(21/3√3). This is the
0 − 0 0 −
minimal corrector for the Takahashi-Imada propagator.
The corrector of the form (5.9) is completely determined by a single number e . Its sign dictates the order of the
2
T and V operators and its value fixes their coefficients. For the alternative propagator ρ′ (5.4) with e = 1/12,
2B 2 −
its corrector is of the same form as (5.14), but now with d =1/√6 and c = 21/6/√3.
0 0
−
For positive values of e , the corrector is of the form
2
c c
0 0
S = exp(ε T)exp(εc V)exp(ε( +d )T)exp(εd V) (5.15)
0 0 0
2 2
1 1 1
= exp (c +d )εH +( d2)ε2G+ ( d2)(c +d )[H,G]
0 0 2 0 2 2 0 0 0
(cid:18)
1 1
(c3+4d3)ε3[T,[T,V]] (c3+4d3)ε3[V,[T,V]] +O(ε4).
−24 0 0 − 12 0 0
(cid:19)
Propagator ρ is dual to the TI propagator with e = 1/24. Its corrector is of the form (5.15) but with same
2B 2
coefficients d = 1/(2√3) and c = 1/(21/3√3). The ρ′ propagator (5.6) with e = 1/12 is dual to ρ′ . Its
0 0 − TI 2 2B
correctoris of the form (5.15) with d =1/√6 and c = 21/6/√3. These compact correctorsare fitting companions
0 0
−
to their equally compact propagators.
VI. CONCLUSIONS
In this work,we proveda fundamentalresult onthe correctabilityofforwardtime step propagators. We show that
ifρ=eε(T+V) weretobeapproximatedbytheproductform(2.2),thennoproductformwithpositivecoefficients t
i
{ }
iscorrectablebeyondsecondorder. Whereasaconventionalhigherorderpropagatorrequiresitserrortermstovanish,
a correctable propagator only require its error terms to satisfy the correctability condition. The latter requirement
seemedfarlessstringent. Asurprisingelementofthis workisthat,this isnotthe case. Forsymmetricdecomposition
with positive t , the two second order error coefficients cannot both vanish because, they can never be equal! The
i
{ }
correctability requirementitself is stringent enough. This proof of non-correctabilitygeneralizesthe previous work of
Sheng1 and Suzuki2.
From knowing correctability requirement, we derived systematically the four forward correctable second order
propagators and their minimal correctors. These minimal correctors follow from a more general form (5.10) of the
corrector with free parameters. Much of the existing literature on symplectic corrector is rather opaque, concerned
only with how to satisfy “order conditions” numerically22,23. This work suggests that a more analytical approach is
possible.
TheTakahashi-Imadatypeofpropagatorsconsideredhereareuniqueinthattheyaretheonlyknownsecond-order,
forward-time-step propagators with a fourth order trace. If one is willing to evaluate the potential at least twice,
then with the inclusion of [V,[T,V]], one can make both error coefficients e and e vanish8,9. The result is
TTV VTV
a whole family of positive time step fourth order propagators27,28,29,30 with a fourth order trace. While this class
of forward decomposition algorithms is indispensable for solving time-irreversible problems10,11,12,13,14, they are less
interesting from the point of view of calculating the trace. For correctable propagators, their key attraction is that
one can obtain a higher order trace without using a higher order propagator. Methods and results of this work can
be used to study ways of correcting these fourth order propagatorsto higher orders.
Acknowledgments
I wish to thank J. Boronat and J. Casulleras, for their invitation to lecture in Barcelona during the summer of
2003 which initiated this work, E. Krotscheck, for his interest and hospitality at Linz, H. Forbert, for discussing the
9
correctability requirement, and G. Chen, on the use of constrained minimization. This work was supported, in part,
by a National Science Foundation grant, No. DMS-0310580.
REFERENCES
1 Q. Sheng,IMA J. Num.Anaysis, 9, 199 (1989).
2 M. Suzuki,J. Math. Phys. 32, 400 (1991).
3 D. Goldman and T. J. Kaper, SIAMJ. Numer.Anal., 33, 349 (1996).
4 R. I.McLachlan and P. Atela, Nonlinearity, 5, 542 (1991).
5 R. I.McLachlan, SIAMJ. Sci. Comput. 16, 151 (1995).
6 P. V. Koseleff, in Integration algorithms and classical mechanics, Fields Inst.Commun., 10, Amer. Math. Soc., Providence,
RI,P.103, (1996).
7 R. I.McLachlan and G. R. W. Quispel, Acta Numerica, 11, 241 (2002).
8 M. Suzuki, Computer Simulation Studies in Condensed Matter Physics VIII, eds, D. Landau, K. Mon and H. Shuttler
(Springler, Berlin, 1996).
9 S.A. Chin, Physics Letters A 226, 344 (1997).
10 H. A.Forbert and S. A.Chin, Phys. Rev. E 63, 016703 (2001).
11 H. A.Forbert and S. A.Chin, Phys. Rev. B 63, 144518 (2001).
12 J. Auer,E. Krotscheck, and S. A.Chin, J. Chem. Phys. 115, 6841 (2001).
13 O. Ciftja and S.A. Chin, Phys. Rev. B 68, 134510 (2003).
14 S. Jang, S.Jang and G. A.Voth,J. Chem. Phys. 115 7832, (2001).
15 M. Kalos, D.Levesque,and L. Verlet,Phys. Rev. A 9, 2178 (1974).
16 S. A.Chin, Phys. Rev. A 42, 6991 (1990).
17 H. Kleinert, Path Integrals in Quantum Mechanics, Statistics and Polymer Physics, World Scientific, Singapore, 1990.
18 J. Wisdom, M. Holman, AND J. Touma, “Symplecticcorrectors”, in Integration Algorithms and Classical Mechanics, J. E.
Marsden, G. W. Patrick, and W. F. Shadwick, eds., American Mathematical Society,Providence, RI, 1996.
19 R.I.McLachan,“Moreonsymplecticcorrectors”, inIntegration Algorithms and Classical Mechanics, J.E.Marsden,G.W.
Patrick, and W. F. Shadwick,eds., American Mathematical Society, Providence, RI,1996.
20 M. A. Lopez-Marcos, J. M. Sanz-Serna, and R. D. Skeel, in Numerical Analysis 1995, D. F. Griffiths and G. A. Watson,
eds., Longman, Harlow, UK,1996, pp.107-122.
21 M. A.Lopez-Marcos, J. M. Sanz-Serna, and R.D. Skeel, SIAMJ. Sci. Comput., 18 223, (1997).
22 S. Blanes, F. Casas, and J. Ros, Siam J. Sci. Comput., 21, 711 (1999).
23 S. Blanes, App.Numer.Math 37, 289 (2001).
24 M. Suzuki,Phys. Lett. A 201, 425 (1995).
25 M. Takahashi and M. Imada, J. Phys.Soc. Jpn 53, 3765 (1984).
26 X.-P. Li, J. Q. Broughton, J. Chem. Phys., 86, 5094 (1987)
27 S. A.Chin and C. R.Chen, J. Chem. Phys. 117, 1409 (2002).
28 I. P. Omelyan,I. M. Mryglod and R. Folk,Phys. Rev. E66, 026701 (2002).
29 I. P. Omelyan,I. M. Mryglod and R. Folk,Comput. Phys. Commun. 151 272 (2003)
30 S. A. Chin, and C. R. Chen, “Forward Symplectic Integrators for Solving Gravitational Few-Body Problems”, arXiv,
astro-ph/0304223.