(* pu2.m:  S. A. Fulling, Nov. 1995 (first version summer 1994) *)

Off[General::spell1]
Share[]
orderRules = {
   x1^(n1_)*x2^(n2_) :> x1^n2*x2^n1 /; n2 > n1
   }

qmax = 9
max = 2(qmax+1)

link = Sum[t^(2j) c[j] x1^j x2^j, {j, 0, qmax}] +O[t]^max
point1 = Sum[t^(2j) c[j] x1^(2j), {j,0, qmax}] + O[t]^max
point2 = point1 /. x1 -> x2
pointpair = Sum[t^(4j) c[j]^2 x1^(2j) x2^(2j), 
                                {j,0, Floor[qmax/2]}] + O[t]^max

identity = link * point1 * point2
interchange = link * pointpair

pu2series = (1/2)(identity + interchange)
pu2temp = (Expand[Normal[pu2series] * x1^2 x2^2] //.
                                   orderRules)/ (x1^2 x2^2);
pu2 = Collect[pu2temp, Prepend[Table[c[qmax-i], {i,0,qmax}], t]]

pu2point = pu2 //. c[n_] -> 1
pu2line = pu2 //. {x1->1, x2->1}
pu2total = pu2line //. c[n_] -> 1

Save["pu2.out", pu2, pu2point, pu2line, pu2total]
Quit[]