(*       Series for the energies of the linearly bound states in the                   *)
(*                       Shastry-Sutherland model                                      *)
(*                                                                                     *)
(*                           (MATHEMATICA version)                                     *)
(*                                                                                     *) 
(*       C. Knetter, A. Buehler, E. Mueller-Hartmann,  and  G. Uhrig                   *)
(*                                                                                     *)
(*                    Phys. Rev. Lett. 85, 3958 (2000)                                 *)
(*                                                                                     *)
(*    "Dispersion and Symmetry of Bound States in the Shastry-Sutherland Model"        *)
(*                                                                                     *)
(*                       University of Cologne                                         *)
(*                             May 2000                                                *) 
(*                                                                                     *)
(*                                                                                     *)
(*                                                                                     *)
(*      ......... Read this file in Mathematica with:  In[i] <<filename  .........     *)
(*                                                                                     *)
(*                                                                                     *)
(*                                                                                     *)
(*                                                                                     *)
(*                              SUMMARY                                                *)
(*                                                                                     *)
(*       The expansion parameter is x=J2/J1. There are 4 bound states                  *)
(*       with a linear term O(x) in the binding energy. For these                      *)
(*       4 states for S=0,1 the dispersion is given in 5th order and the               *)
(*       series up to x^14 at k=(0,0), (0,Pi).                                         *)
(*                                                                                     *)
(*               Auxiliary variables                                                   *)
(*                                                                                     *)

u = Cos[k1]*Cos[k2]
v = Cos[k1]+Cos[k2]
sq= Sqrt[1+(Cos[k1]-Cos[k2])^2]

(*                                                                                     *)
(*                       S = 0                                                         *)
(*                                                                                     *)
(*               Dispersion                                                            *)

w1s0= 2-x-2*x^2-x^3/2*(1+u)+x^4/8*(7-v^2+2*u^2-v*sq)+x^5/16*(10-2*v^2+23*u+2*u*(v^2+2*u)-4*u^3+v*(-6-2*v^2+9*u+2*u*(v^2-4*u))/sq)

w2s0= 2-x-2*x^2-x^3/2*(1+u)+x^4/8*(7-v^2+2*u^2+v*sq)+x^5/16*(10-2*v^2+23*u+2*u*(v^2+2*u)-4*u^3-v*(-6-2*v^2+9*u+2*u*(v^2-4*u))/sq)

w3s0= 2-x-x^2+x^3/2*(u-1)+x^4/8*(-5+6*u-2*u^2+2*Cos[k1]^2)+x^5*(-3/8+u/16+Cos[k2]^2/4*(3-u)-u^2/2+u^3/4)

w4s0= 2-x-x^2+x^3/2*(u-1)+x^4/8*(-5+6*u-2*u^2+2*Cos[k2]^2)+x^5*(-3/8+u/16+Cos[k2]^2/4*(3-u)-u^2/2+u^3/4)

(* at |k1|=|k2|=q there is a degeneracy between states 3 and 4 in low orders           *)
(* which is broken in x^5 by the additional term                                       *)

w3qqs0= w3s0 - Cos[q]*Sin[q]^2*x^5/8 /. {k1 -> q, k2 -> q} 

w4qqs0= w4s0 + Cos[q]*Sin[q]^2*x^5/8 /. {k1 -> q, k2 -> q}

(*               k=(0,0)                                                               *)

e100s0= 2-x-2*x^2-x^3+3/8*x^4+23/16*x^5-37/64*x^6-3593/768*x^7-12109/4608*x^8+62261/5184*x^9+1392949559/79626240*x^10-31969538107/1194393600*x^11-31775088956269/382205952000*x^12+637813109001179/22932357120000*x^13+2060594364294311131/6421059993600000*x^14

e200s0= 2-x-2*x^2-x^3+7/8*x^4+43/16*x^5-25/64*x^6-2435/256*x^7-12659/1536*x^8+98741/3456*x^9+949469191/15925248*x^10-39313719143/597196800*x^11-371668244803807/1146617856000*x^12-31302358410107/2548039680000*x^13+9400924284794617963/6421059993600000*x^14

e300s0= 2-x-x^2+1/8*x^4-1/16*x^5-53/192*x^6-13/768*x^7+377/6912*x^8-2341/13824*x^9-2366921/5308416*x^10-415612943/2388787200*x^11-33486987683/127401984000*x^12-4533115819027/5733089280000*x^13-92584505133810481/57789539942400000*x^14

e400s0= 2-x-x^2+1/8*x^4-1/16*x^5-53/192*x^6-13/768*x^7+377/6912*x^8-2341/13824*x^9-2366921/5308416*x^10-415612943/2388787200*x^11-33486987683/127401984000*x^12-4533115819027/5733089280000*x^13-92584505133810481/57789539942400000*x^14

(*               k=(0,Pi)                                                              *)

e10ps0= 2-x-2*x^2+9/8*x^4-5/16*x^5-449/192*x^6+1645/2304*x^7+70427/13824*x^8-452555/165888*x^9-364407181/26542080*x^10+20932260889/2388787200*x^11+43212785310413/1146617856000*x^12-1082384121586093/34398535680000*x^13-6423538821754896467/57789539942400000*x^14

e20ps0= e10ps0

e30ps0= 2-x-x^2-x^3-11/8*x^4-3/16*x^5+149/64*x^6+13129/2304*x^7+17267/6912*x^8-100291/6144*x^9-3919582531/79626240*x^10-103864493141/2388787200*x^11+14394828883589/127401984000*x^12+15990456557778707/34398535680000*x^13+33608754413732880413/57789539942400000*x^14

e40ps0= e30ps0

(*                                                                                     *)
(*                       S = 1                                                         *)
(*                                                                                     *)
(*               Dispersion                                                            *)

w1s1= 2-x/2-5*x^2/4+x^3*(u-2)/4+x^4/16*(-v^2+13*u+2*u^2-19-v*sq)-x^5/64*(118-119*u+20*v^2+4*u*(v^2-14*u)-8*u^3+v*(11+20*v^2-2*u*(39-2*v^2+8*u))/sq)

w2s1= 2-x/2-5*x^2/4+x^3*(u-2)/4+x^4/16*(-v^2+13*u+2*u^2-19+v*sq)-x^5/64*(118-119*u+20*v^2+4*u*(v^2-14*u)-8*u^3-v*(11+20*v^2-2*u*(39-2*v^2+8*u))/sq)

w3s1= 2-x/2-3*x^2/4-x^3*(u+1)/4+x^4*(-35+2*Cos[k1]^2-15*u-2*u^2)/16+x^5*(-175/64*u-137/32+Cos[k1]^2/8*(7+u)-3*u^2/4-u^3/8)

w4s1= 2-x/2-3*x^2/4-x^3*(u+1)/4+x^4*(-35+2*Cos[k2]^2-15*u-2*u^2)/16+x^5*(-175/64*u-137/32+Cos[k1]^2/8*(7+u)-3*u^2/4-u^3/8)

(* at |k1|=|k2|=q there is a degeneracy between states 3 and 4 in low orders          *)
(* which is broken in x^5:                                                            *)

w3qqs1= w3s1 - Cos[q]*Sin[q]^2*x^5/8 /. {k1 -> q, k2 -> q}

w4qqs1= w4s1 + Cos[q]*Sin[q]^2*x^5/4 /. {k1 -> q, k2 -> q}

(*               k=(0,0)                                                              *)

e100s1= 2-1/2*x-5/4*x^2-1/4*x^3-5/8*x^4-57/64*x^5-337/384*x^6-1087/1152*x^7-49745/27648*x^8-1230671/663552*x^9-84134083/31850496*x^10-71961335353/19110297600*x^11-6931489716539/1146617856000*x^12-41618542755193/5503765708800*x^13-134624403747285691/9631589990400000*x^14

e200s1=2-1/2*x-5/4*x^2-1/4*x^3-3/8*x^4-5/64*x^5-15/128*x^6-1/384*x^7-8677/27648*x^8-248747/663552*x^9-61964869/53084160*x^10-28161170261/19110297600*x^11-1286215153651/382205952000*x^12-722520360562511/137594142720000*x^13-694140469564491133/57789539942400000*x^14

e300s1= 2-1/2*x-3/4*x^2-1/2*x^3-25/8*x^4-441/64*x^5-3023/192*x^6-121721/4608*x^7-34775/6144*x^8+11768317/55296*x^9+38294096717/31850496*x^10+83630268980713/19110297600*x^11+650630140801309/57330892800*x^12+348272518978158487/22932357120000*x^13-4990505545003648955371/115579079884800000*x^14

e400s1= e300s1

(*               k=(0,Pi)                                                            *)

e10ps1= 2-1/2*x-5/4*x^2-3/4*x^3-15/8*x^4-189/64*x^5-71/24*x^6+16555/4608*x^7+521443/18432*x^8+57239671/663552*x^9+23845114583/159252480*x^10+62475412451/3822059520*x^11-1111196843333209/1146617856000*x^12-61012927483426099/15288238080000*x^13-18791746943890357267/2063912140800000*x^14

e20ps1= e10ps1

e30ps1= 2-1/2*x-3/4*x^2-5/4*x^4-91/64*x^5-243/128*x^6-575/288*x^7-6245/3072*x^8-388771/165888*x^9-29632859/10616832*x^10-11285534159/6370099200*x^11-2412392495713/1146617856000*x^12-251320198387873/45864714240000*x^13-1951253385625429681/57789539942400000*x^14

e40ps1= e30ps1