{VERSION 5 0 "IBM INTEL NT" "5.0" }
{USTYLETAB {CSTYLE "Maple Input" -1 0 "Courier" 0 1 255 0 0 1 0 1 0 0 
1 0 0 0 0 1 }{CSTYLE "2D Math" -1 2 "Times" 0 1 0 0 0 0 0 0 2 0 0 0 0 
0 0 1 }{CSTYLE "2D Input" 2 19 "" 0 1 255 0 0 1 0 0 0 0 0 0 0 0 0 1 }
{CSTYLE "" -1 256 "" 1 16 0 0 0 0 0 1 0 0 0 0 0 0 0 0 }{CSTYLE "" -1 
257 "" 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 0 }{PSTYLE "Normal" -1 0 1 
{CSTYLE "" -1 -1 "Times" 1 12 0 0 0 1 2 2 2 2 2 2 1 1 1 1 }1 1 0 0 0 
0 1 0 1 0 2 2 0 1 }}
{SECT 0 {EXCHG {PARA 0 "" 0 "" {TEXT 256 39 "Taylor expansion of an op
erator product" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" 
{TEXT -1 208 "We check the order of the split-operator time propagatio
n technique. First, we analyze a popular version of the time splitting
 technique (used in Ehrenfest.mws). Then, we show why the simple opera
tor product " }{TEXT 19 44 "exp(-I*(T+V)*dt) = exp(-I*T*dt) exp(-I*V*d
t)" }{TEXT -1 16 "  results in an " }{TEXT 19 7 "O(dt^2)" }{TEXT -1 
10 " error if " }{TEXT 19 1 "T" }{TEXT -1 5 " and " }{TEXT 19 1 "V" }
{TEXT -1 37 " do not commute. Then, we analyze an " }{TEXT 19 7 "O(dt^
4)" }{TEXT -1 73 " accurate scheme which involves a large amount of al
gebraic computations." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 
"" {TEXT -1 103 "This is, of course, related to the Campbell-Baker-Hau
sdorff formula, which reads ([A,B] : = A B - B A):" }}{PARA 0 "" 0 "" 
{TEXT -1 88 "exp(A) exp(B) = exp(O), where O = A + B + 1/2 [A,B] + 1/1
2 ([[A,B],B] + [A,[A,B]]) + ..." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}
{PARA 0 "" 0 "" {TEXT -1 110 "This worksheet can be used in maple6 and
 maple7, but it gobbles lots of memory there. It works fine in maple8.
" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 8 "restart;" }}}{EXCHG 
{PARA 0 "> " 0 "" {MPLTEXT 1 0 66 "unassign(`&*`); \ndefine(`&*`,multi
linear,zero=0,identity=1,flat);\n" }}}{EXCHG {PARA 0 "> " 0 "" 
{MPLTEXT 1 0 5 "A&*B;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 5 "a&*
A;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 11 "gamma &* A;" }}}
{EXCHG {PARA 0 "" 0 "" {TEXT -1 73 "If we need a variable to be treate
d as a constant, we add it to the list:" }}}{EXCHG {PARA 0 "> " 0 "" 
{MPLTEXT 1 0 24 "constants:=constants,dt;" }}}{EXCHG {PARA 0 "> " 0 "
" {MPLTEXT 1 0 11 "&*(A,dt,B);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 
1 0 12 "simplify(%);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 31 "tay
lor(exp(-I*A&*B*dt),dt=0,3);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 8 "Se
tting " }{TEXT 19 2 "dt" }{TEXT -1 112 " to a constant is problematic,
 can't Taylor expand about it. So the trick is to Taylor expand about \+
a variable (" }{TEXT 19 1 "h" }{TEXT -1 228 " in our case), and later \+
to make substitutions, so that one may make use of the constant proper
ty. The above Taylor expansion attempt failed, because we can't expand
 about a variable that has just been declared to be a constant." }}}
{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 29 "taylor(exp(-I*A&*B*h),h=0,3)
;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 13 "subs(h=dt,%);" }}}
{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 32 "No:=5; #order to which we ex
pand" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 52 "Ut:=convert(taylor(
exp(-I*(A+B)*h),h=0,No),polynom);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 
66 "We need to do something here: if we leave it as is, the powers of \+
" }{TEXT 19 3 "A,B" }{TEXT -1 46 " get evaluated assuming a commutativ
e product." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 97 "Ut:=subs((A+B
)^2=(A+B)&*(A+B),(A+B)^3=(A+B)&*(A+B)&*(A+B),(A+B)^4=(A+B)&*(A+B)&*(A+
B)&*(A+B),Ut);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 53 "UtA:=conv
ert(taylor(exp(-I*(A)*h/2),h=0,No),polynom);" }}}{EXCHG {PARA 0 "> " 
0 "" {MPLTEXT 1 0 51 "UtB:=convert(taylor(exp(-I*(B)*h),h=0,No),polyno
m);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 45 "Uspo:=simplify(subs(
h=dt,UtA &* UtB &* UtA));" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 
35 "Udif:=simplify(subs(h=dt,Ut)-Uspo);" }}}{EXCHG {PARA 0 "> " 0 "" 
{MPLTEXT 1 0 34 "coeff(Udif,dt,0),coeff(Udif,dt,1);" }}}{EXCHG {PARA 
0 "> " 0 "" {MPLTEXT 1 0 21 "C2:=coeff(Udif,dt,2);" }}}{EXCHG {PARA 0 
"> " 0 "" {MPLTEXT 1 0 21 "C3:=coeff(Udif,dt,3);" }}}{EXCHG {PARA 0 ">
 " 0 "" {MPLTEXT 1 0 21 "C4:=coeff(Udif,dt,4);" }}}{EXCHG {PARA 0 "> \+
" 0 "" {MPLTEXT 1 0 13 "simplify(C2);" }}}{EXCHG {PARA 0 "> " 0 "" 
{MPLTEXT 1 0 112 "definemore(`&*`,\n&*(1,A::symbol,B::symbol)=A&*B,&*(
A::symbol,1,B::symbol)=A&*B,&*(A::symbol,B::symbol,1)=A&*B);\n" }}}
{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 13 "simplify(C2);" }}}{EXCHG 
{PARA 0 "> " 0 "" {MPLTEXT 1 0 47 "definemore(`&*`,&*(B::symbol,B::sym
bol) = B^2);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 77 "This one was tric
ky, some other obvious versions didn't work (putting in the " }{TEXT 
19 8 "::symbol" }{TEXT -1 16 " was essential)!" }}}{EXCHG {PARA 0 "> \+
" 0 "" {MPLTEXT 1 0 17 "simplify(A &* A);" }}}{EXCHG {PARA 0 "> " 0 "
" {MPLTEXT 1 0 13 "simplify(C2);" }}}{EXCHG {PARA 0 "> " 0 "" 
{MPLTEXT 1 0 13 "simplify(C3);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 77 
"Will this vanish?  Make sure the original expansions are to sufficien
t order!" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 
64 "Note that obvious simplifications were not carried out! Such as " 
}{TEXT 19 23 "&*(A,B^2,1) = A &* B^2 " }{TEXT -1 39 ", despite the for
mal prescription. Why?" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 14 "w
hattype(B^2);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 130 "definemor
e(`&*`,\n&*(1,A::algebraic,B::algebraic)=A&*B,&*(A::algebraic,1,B::alg
ebraic)=A&*B,&*(A::algebraic,B::algebraic,1)=A&*B);\n" }}}{EXCHG 
{PARA 0 "> " 0 "" {MPLTEXT 1 0 13 "simplify(C3);" }}}{EXCHG {PARA 0 "
" 0 "" {TEXT -1 118 "OK, so we managed to get the triple product reduc
ed when a constant (1) is inside, while a power of operators appears.
" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 124 "definemore(`&*`,\n&*(A
::algebraic,A::algebraic,B::algebraic)=&*(A^2,B),&*(A::algebraic,B::al
gebraic,B::algebraic)=&*(A,B^2));" }}}{EXCHG {PARA 0 "> " 0 "" 
{MPLTEXT 1 0 13 "simplify(C3);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 
1 0 141 "definemore(`&*`,\n&*(A::algebraic,A::algebraic,A::algebraic)=
A^3,&*(A::algebraic,(A::algebraic)^2)=A^3,&*((A::algebraic)^2,A::algeb
raic)=A^3);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 13 "simplify(C3)
;" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 68 "The third-order error term c
an be expressed in terms of commutators:" }}}{EXCHG {PARA 0 "> " 0 "" 
{MPLTEXT 1 0 25 "CM:=(A,B)->A &* B-B &* A;" }}}{EXCHG {PARA 0 "> " 0 "
" {MPLTEXT 1 0 30 "C3f:=I*(CM(A+2*B,CM(A,B)))/24;" }}}{EXCHG {PARA 0 "
> " 0 "" {MPLTEXT 1 0 17 "simplify(C3-C3f);" }}}{EXCHG {PARA 0 "> " 0 
"" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 54 "Now show t
hat this is better than the naive splitting:" }}}{EXCHG {PARA 0 "> " 
0 "" {MPLTEXT 1 0 51 "UtA:=convert(taylor(exp(-I*(A)*h),h=0,No),polyno
m);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 51 "UtB:=convert(taylor(
exp(-I*(B)*h),h=0,No),polynom);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 
1 0 38 "Uspo:=simplify(subs(h=dt,UtA &* UtB));" }}}{EXCHG {PARA 0 "> \+
" 0 "" {MPLTEXT 1 0 35 "Udif:=simplify(subs(h=dt,Ut)-Uspo);" }}}
{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 34 "coeff(Udif,dt,0),coeff(Udif,
dt,1);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 21 "C2:=coeff(Udif,dt
,2);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 59 "The simple product is onl
y second-order accurate, i.e., an " }{TEXT 19 7 "O(dt^2)" }{TEXT -1 
17 " term survive if " }{TEXT 19 1 "A" }{TEXT -1 5 " and " }{TEXT 19 
1 "B" }{TEXT -1 15 " don't commute." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }
}{PARA 0 "" 0 "" {TEXT 257 9 "Exercise:" }}{PARA 0 "" 0 "" {TEXT -1 
105 "Look at the Campbell-Baker-Hausdorff expansion with a smallness p
arameter dt to generate a Taylor series." }}{PARA 0 "" 0 "" {TEXT 19 
31 "exp(A*dt) exp(B*dt) = exp(O*dt)" }{TEXT -1 8 ", where " }{TEXT 19 
58 "O = A + B + 1/2 [A,B] + 1/12 ([[A,B],B] + [A,[A,B]]) + ..." }}
{PARA 0 "" 0 "" {TEXT -1 94 "Use the Taylor expansion to verify the qu
oted result, and calculate the next order, i.e., the " }{TEXT 19 7 "O(
dt^3)" }{TEXT -1 12 " expression." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}
{EXCHG {PARA 0 "" 0 "" {TEXT -1 139 "The next issue: can we generate a
 more complicated split-operator technique which will agree with the e
xact answer to a higher order in dt?" }}{PARA 0 "" 0 "" {TEXT -1 87 "T
he answer is yes. We first define the basic split-operator technique a
s a single step:" }}{PARA 0 "" 0 "" {TEXT 19 47 "U2 = exp(-I*A*dt/2) e
xp(-I*B*dt) exp(-I*A*dt/2)" }{TEXT -1 18 " as a function of " }{TEXT 
19 8 "A, B, dt" }}{PARA 0 "" 0 "" {TEXT 19 48 "U3 = U2(A,B,g*dt) U2(A,
B,(1-2*g)*dt) U2(A,B,g*t)" }{TEXT -1 38 " will work for a particular c
hoice of " }{TEXT 19 1 "g" }{TEXT -1 1 "." }}}{EXCHG {PARA 0 "> " 0 "
" {MPLTEXT 1 0 59 "U2:=(A,B,h)->exp(-I*h/2*A) &* exp(-I*h*B) &* exp(-I
*h/2*A);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 68 "g:='g': # a par
ameter to subdivide the time step dt to be optimized." }}}{EXCHG 
{PARA 0 "> " 0 "" {MPLTEXT 1 0 52 "U3:=U2(A,B,g*h) &* U2(A,B,(1-2*g)*h
) &* U2(A,B,g*h);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 35 "conver
t(taylor(U3,h=0,No),polynom);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 155 
"The direct approach won't work. We need to do what was done before, n
amely generate the Taylor expansions for the individual pieces, and pu
t them together." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 53 "UtA:=co
nvert(taylor(exp(-I*(A)*h/2),h=0,No),polynom);" }}{PARA 0 "> " 0 "" 
{MPLTEXT 1 0 51 "UtB:=convert(taylor(exp(-I*(B)*h),h=0,No),polynom);" 
}{TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 38 "U2:=unappl
y(UtA &* UtB &* UtA, A,B,h):" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 40 "T
he naive implementation is given below." }}}{EXCHG {PARA 0 "> " 0 "" 
{MPLTEXT 1 0 54 "# U3:=U2(A,B,g*h) &* U2(A,B,(1-2*g)*h) &* U2(A,B,g*h)
:" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 23 "constants:=constants,g
;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 76 "#U3:=simplify( subs(h=
dt,U2(A,B,g*h) &* U2(A,B,(1-2*g)*h) &* U2(A,B,g*h)) ):" }}}{EXCHG 
{PARA 0 "" 0 "" {TEXT -1 56 "The above line would at least make use of
 the fact that " }{TEXT 19 2 "dt" }{TEXT -1 242 " is a constant, and s
implify operator products. It is, however, too cumbersome for Maple to
 evaluate within a reasonable amount of memory and time. The trick is \+
to build up the expression gradually, and to keep only the lowest four
 orders in " }{TEXT 19 2 "dt" }{TEXT -1 86 " for each factor as one bu
ilds the expression. Note that back-and-forth substitutions " }{TEXT 
19 2 "dt" }{TEXT -1 4 " vs " }{TEXT 19 1 "h" }{TEXT -1 45 " are needed
. Taylor expansion will work with " }{TEXT 19 1 "h" }{TEXT -1 15 ", an
d not with " }{TEXT 19 2 "dt" }{TEXT -1 50 ". Non-commutative product \+
simplification requires " }{TEXT 19 2 "dt" }{TEXT -1 41 ", so that it \+
is recognized as a constant." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 
0 81 "U3a:=convert(taylor(subs(dt=h,simplify(subs(h=dt,U2(A,B,g*h)))),
h=0,No),polynom);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 87 "U3b:=c
onvert(taylor(subs(dt=h,simplify(subs(h=dt,U2(A,B,(1-2*g)*h)))),h=0,No
),polynom);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 48 "U3p:= expand
(subs(h=dt,U3a) &* subs(h=dt,U3b) ):" }}}{EXCHG {PARA 0 "> " 0 "" 
{MPLTEXT 1 0 52 "U3p:=convert(taylor(subs(dt=h,U3p),h=0,No),polynom):
" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 48 "U3p:= expand(subs(h=dt,
U3p) &* subs(h=dt, U3a)):" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 
51 "U3:=convert(taylor(subs(dt=h,U3p),h=0,No),polynom):" }}}{EXCHG 
{PARA 0 "" 0 "" {TEXT -1 132 "Note that this relatively fast computati
on relies on the Taylor truncations. Try it without them and see memor
y and time gobbled up." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 14 "c
oeff(U3,h,1);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 3 "Ut;" }}}
{EXCHG {PARA 0 "" 0 "" {TEXT -1 88 "Now we verify that the time evolut
ion operator expressions agree in their orders 0 to 2:" }}}{EXCHG 
{PARA 0 "> " 0 "" {MPLTEXT 1 0 21 "DIF:=simplify(Ut-U3):" }}}{EXCHG 
{PARA 0 "> " 0 "" {MPLTEXT 1 0 15 "coeff(DIF,h,0);" }}}{EXCHG {PARA 0 
"> " 0 "" {MPLTEXT 1 0 15 "coeff(DIF,h,1);" }}}{EXCHG {PARA 0 "> " 0 "
" {MPLTEXT 1 0 15 "coeff(DIF,h,2);" }}}{EXCHG {PARA 0 "> " 0 "" 
{MPLTEXT 1 0 21 "DIF3:=coeff(DIF,h,3);" }}}{EXCHG {PARA 0 "" 0 "" 
{TEXT -1 37 "Now show that a particular choice of " }{TEXT 19 1 "g" }
{TEXT -1 22 " turns this into zero:" }}}{EXCHG {PARA 0 "> " 0 "" 
{MPLTEXT 1 0 37 "simplify(subs(g=1/(2-2^(1/3)),DIF3));" }}}{EXCHG 
{PARA 0 "" 0 "" {TEXT -1 91 "This completes the proof that a fourth-or
der accurate time evolution operator can be found." }}{PARA 0 "" 0 "" 
{TEXT -1 192 "In the numerical implementation one needs to ask the que
stion whether this algorithm is economical. The improved accuracy come
s at the price of a substantial number of additional calculations." }}
}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "> " 0 "
" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}}
{MARK "0 0 0" 0 }{VIEWOPTS 1 1 0 1 1 1803 1 1 1 1 }{PAGENUMBERS 0 1 2 
33 1 1 }