######################################################################
##SAP: Save this file as SAP. To use it, stay in the                 #
##same directory, get into Maple (by typing: maple <Enter> )         #
##and then type:  read SAP : <Enter>                                 #
##Then follow the instructions given there                           #
##                                                                   #
##Written by Doron Zeilberger, Temple University ,                   #
#zeilberg@math.temple.edu.                                           # 
#######################################################################
 
#Created: June 15, 1999
#This version:Dec. 16, 1999
#SAP: A Maple package to study SAPs
#It it is one of the packages that accompanies Doron Zeilberger's
#article: "Symbol-Crunching with the Transfer-Matrix Method with 
#Applications to the Enumeration of Skinny Physical Creatures
#
#Please report bugs to zeilberg@math.temple.edu
 
 
 
print(`Created: June 15, 1999.`):
print(`This version: June 17, 1999`):
print(` SAP: A Maple package to study SAPs `):
print(` It it is one of the packages that accompany Doron Zeilberger's `):
print(`article: "Symbol-Crunching with the Transfer-Matrix Method with `):
print(`Applications to the Enumeration of Skinny Physical Creatures`):
 
lprint(``):
print(`Written by Doron Zeilberger, zeilberg@math.temple.edu`):
lprint(``):
print(`Please report bugs to zeilberg@math.temple.edu`):
lprint(``):
 print(`The most current version of these  packages and paper`):
 print(` are  available from`):
 print(`http://www.math.temple.edu/~zeilberg/`):
 print(`For a list of the  MAIN procedures type: ezra();`):
 print(`For a list of ALL procedures type: ezra1();, for help with`):
 print(`a specific procedure, type ezra(procedure_name)`):
 print(``):
 
ezra1:=proc()
if args=NULL then
 print(`Contains the following procedures: Alphabet, DerivedPreLetters,`):
 print(` FillIns1, FindL, FindR , Followers, gf, Gf, GF, `):
 print(` gfBatchelorYung, gfSeries, GfSeries  `):
 print(` GFSeries, GfSeriesHeight, GFSeriesHeight, IncSeq, `):
 print(` IncSeqSand, IsLeftLetter`):
 print(` IsRightLetter, JoinLL, JoinRL, JoinRR, LebensRaum, LebensRaums`):
 print(` LeftLetters, MarCha, OneStepAll, `):
 print(` OneStepLL, OneStepRL `):
 print(` OneStepRR, Parenths, PlacesToInsert, PreFollowers, PrePreFollowers `):
 print(` RightLetters, SAP,SAPlhw, SAPseries, SolveMC3series, SolveMC3seriesS `):
 print(` SetToUnInt , SolveMC3, Stick1 `):
 print(` For help with a specific procedure type : ezra(proc_name);`):
else
ezra(args):
fi:
end:
 
ezra:=proc()
if args=NULL then
 print(`the MAIN procedures are : `):
 print(`  GF, gfBatchelorYung, GFSeries , SAPseries `):
 print(``):
 print(` For help with a specific procedure type : ezra(proc_name);`):
 print(` For example with help with procedure GF, type: ezra(GF); `):
fi:
 
 
if nops([args])=1 and op(1,[args])=gfBatchelorYung then 
print(`gfBatchelorYung(n,t): The generating function of `):
print(`self-avoiding polygons`):
print(`immersed the strip [0,n] and such that the left-most letter`):
print(`has a 0 in it. It gives the quantity described in`):
print(`M.T. Batchelor, and C.M. Yang, J. Physics A: Math. gen. 27 (1994)`):
print(`4055-4067 (the denominators on p. 4059 and 4066 and the `):
print(` numerators on p. 4064) `):
fi:
 
if nops([args])=1 and op(1,[args])=SAPseries then 
print(`SAPseries(L): The sequence whose i^th term is`):
print(`the number of self-avoiding polygons of length 2i`):
print(`for i=1...L`):
fi:
 
 
if nops([args])=1 and op(1,[args])=SAP then 
print(`SAP(l): The number of self-avoiding polygons of`):
print(`lentgh 2*l`):
fi:
 
if nops([args])=1 and op(1,[args])=SAPlhw then 
print(`SAPlhw(l,h,w): the number of self-avoiding polygons`):
print(`with length 2l, height h, and width w`):
fi:
 
if nops([args])=1 and op(1,[args])=GfSeriesHeight then 
print(`GfSeriesHeight(n,M,s): The sequence whose i^th term is the`):
print(`g.f. of self-avoiding polygons of length 2i of width=n`):
print(`for 1<=i<=M with weigth s^(height)`):
fi:
 
if nops([args])=1 and op(1,[args])=GFSeriesHeight then 
print(`GFSeriesHeight(n,M,s): The sequence whose i^th term is the`):
print(`g.f. of self-avoiding polygons of length 2i of width<=n`):
print(`for 1<=i<=M with weigth s^(height)`):
fi:
 
 
if nops([args])=1 and op(1,[args])=SolveMC3seriesS then 
print(`SolveMC3seriesS(MC,M,s): Given a type-III Markov Chain finds the`):
print(` sequence`):
print(`whose i^th member is the generating function  of paths of weight i`):
print(`with weight s^(length) for 0<=i<=M `):
 
fi:
 
if nops([args])=1 and op(1,[args])=GfSeries then 
print(`GfSeries(n,M): The sequence whose i^th term is the`):
print(`number of self-avoiding polygons of length 2i`):
print(` of width=n `):
fi:
 
if nops([args])=1 and op(1,[args])=GFSeries then 
print(`GFSeries(n,M): The sequence whose i^th term is the`):
print(`number of self-avoiding polygons of length 2i`):
print(` of width<=n `):
fi:
 
if nops([args])=1 and op(1,[args])=gfSeries then 
print(`gfSeries(n,M): The sequence whose i^th term is the`):
print(`number of self-avoiding polygons of length 2i immersed`):
print(`in the strip 0<=y<=n, for 1<=i<=M `):
fi:
 
if nops([args])=1 and op(1,[args])=SolveMC3series then 
print(`SolveMC3series(MC,M): Given a profile of a type-III `):
print(`Markov Chain, and a positive integer M, finds the sequence`):
print(`whose i^th member is the number of paths of weight i`):
print(`for 0<=i<=M`):
fi:
 
if nops([args])=1 and op(1,[args])=gf then 
print(`gf(n,t): The generating function of self-avoiding polygons`):
print(`immersed in a fixed strip`):
fi:
 
 
if nops([args])=1 and op(1,[args])=GF then 
print(`GF(n,t): The generating function of self-avoiding polygons`):
print(`in the variable t, of width<=n`):
fi:
 
 
if nops([args])=1 and op(1,[args])=Gf then 
print(`Gf(n,t): The generating function of self-avoiding polygons`):
print(`in the variable t, of width=n`):
fi:
 
if nops([args])=1 and op(1,[args])=SolveMC3 then 
print(`SolveMC3(MC,t): Given a type-III Markov Chain`):
print(` [LeftLetters,RightLetters,`):
print(`TransMatrx], finds the sum of the weights of all words that`):
print(`start with a letter of LeftLetters and ends with a letter`):
print(`of RightLetters, as a generating function in the variable t`):
fi:
 
 
if nops([args])=1 and op(1,[args])=MarCha then 
print(`MarCha(n): The Markov chain of self-avoiding polygons`):
print(`that live in the strip [0,n]`):
print(`It is asummed that the vertices (states) are labelled`):
print(`by positive integers. The Markov Chain is given as a`):
print(`list of lengh 3: [LeftLetters,RightLetters,TransitionMatrix]`):
print(`The TransitionMatrix is a list of sets, where the`):
print(` i^th entry is the set of pairs [j,wt(i,j)], where `);
print(` j is the a neighbor and wt(i,j) is the weight of `):
print(` the arc from i to j `):
fi:
 
if nops([args])=1 and op(1,[args])=RightLetters then 
print(`RightLetters(a,b,L,R): All the right`):
print(`letters for a self-avoiding polygon`):
print(`bounded in the strip [a,b], using the letters`):
print(`L and R followed by their finishing-length`):
print(`(i.e. what it takes to close the polygon and`):
print(` terminate the word) `):
fi:
if nops([args])=1 and op(1,[args])=Followers then 
print(`Followers(LET1,a,b,L,R): All the letters that`):
print(`can follow the letter LET1 in the ambient strip, [a,b]`):
print(`where L and R denote ( and ) resp.`):
fi:
 
 
if nops([args])=1 and op(1,[args])=DerivedPreLetters then 
print(` DerivedPreLetters(LET1,L,R): Given a pre-letter, finds `):
print(` all pre-letters obtained from it by inserting any number `):
print(` of eligible intervals `):
fi:
 
if nops([args])=1 and op(1,[args])=FillIns1 then 
print(`FillIns1(LET1,L,R): Given a pre-letter, LET1, where`):
print(`L and R denote ( and ), finds all the possible pre-letters`):
print(`obtained from it by inserting ONE new interval in`):
print(`a legal place`):
fi:
 
if nops([args])=1 and op(1,[args])=Stick1 then 
print(`Stick1(LET1,L,R,a,b): given a pre-letter, and an interval`):
print(`[a,b] finds the resulting pre-letter obtained by inserting`):
print(`[L,R] in the appropirate place in LET1[1][1] and`):
print(`[a,b] in the appropirate place in LET1[1][2]`):
print(`and replacing LET1[2] by LET1[2] minus {a..b}`):
fi:
 
if nops([args])=1 and op(1,[args])=PlacesToInsert then 
print(`PlacesToInsert(UnInts): Given a set of free-space, UnInts`):
print(`in terms of a list of intervals, finds all the`):
print(`intervals [a,b] of length at least 2 that may be`):
print(` inserted `):
fi:
 
if nops([args])=1 and op(1,[args])=SetToUnInt then 
print(`SetToUnInt(kv): Given a set of integers, kv, converts`):
print(`it list-of-intervals notation [[a1,b1],[a2,b2],[a3,b3],...]`):
print(`s.t. kv={a1..b1} union {a2..b2} ....`):
fi:
 
if nops([args])=1 and op(1,[args])=LeftLetters then 
print(`LeftLetters(a,b,L,R): All the left`):
print(`letters for a self-avoiding polygon`):
print(`bounded in the strip [a,b], using the letters`):
print(`L and R`):
fi:
 
if nops([args])=1 and op(1,[args])=PreFollowers then 
print(`PreFollowers(LET1,a,b,L,R): All the pre-letters that`):
print(`can follow the letter LET1 in the ambient strip, [a,b]`):
print(`where L and R denote ( and ) resp.`):
fi:
 
if nops([args])=1 and op(1,[args])=IncSeqSand then 
 
print(`IncSeqSand(BOUNDS,Centers): given a list of pairs`):
print(` [[a1,b1],...,[ak,bk]]`):
print(` and a set of centers, Centers `):
print(` finds the set of increasing sequences x1<=x2<=..<=xk `):
print(` such that ai<=xi<=bi, for i=1...k `):
print(` and the extra teritory taken `):
 
fi:
 
if nops([args])=1 and op(1,[args])=LebensRaums then 
print(`LebensRaums(LET1,FS1): Given a pre-letter LET1, and `):
print(`a free space set FS1, finds the list of pairs`):
print(`indicating the LebensRaum of each component`):
fi:
 
if nops([args])=1 and op(1,[args])=LebensRaum then 
print(` LebensRaum(LET1,FS1,i): Given a Pre-Letter LET1 with its `):
print(` free-space set FS1, and an integer i, finds the set range `):
print(` where the i^th component can roam `):
print(` the output is the lowest and the highest that it can `):
print(` venture to `):
fi:
 
if nops([args])=1 and op(1,[args])=PrePreFollowers then 
print(`PrePreFollowers(LET1,a,b,L,R): All the pre-pre-letters that`):
print(`can follow the letter LET1 in the ambient strip, [a,b]`):
print(`where L and R denote ( and ) resp.`):
fi:
 
if nops([args])=1 and op(1,[args])=OneStepAll then 
print(`OneStepAll(LET1,FS1,L,R): Set of all pairs [LET2,FS2]`):
print(` that can be obtained `):
print(` from the letter LET1, with its accompanying free-space set `):
print(` FS1, by performing either ONE RL-operation or ONE LL- or one RR-`):
print(` operation, L an R denote ( and ) resp. `):
fi:
 
if nops([args])=1 and op(1,[args])=OneStepRL then 
print(`OneStepRL(LET1,FS1,L,R): Set of all pairs [LET2,FS2]`):
print(` that can be obtained `):
print(` from the letter LET1, with its accompanying free-space set `):
print(` FS1, by performing ONE RL-operation `):
print(` L an R denote ( and ) resp. `):
fi:
 
if nops([args])=1 and op(1,[args])=OneStepRR then 
print(`OneStepRR(LET1,FS1,L,R): Set of all pairs [LET2,FS2]`):
print(` that can be obtained `):
print(` from the letter LET1, with its accompanying free-space set `):
print(` FS1, by performing ONE RR-operation `):
print(` L an R denote ( and ) resp. `):
fi:
 
if nops([args])=1 and op(1,[args])=OneStepLL then 
print(`OneStepLL(LET1,FS1,L,R): Set of all pairs [LET2,FS2]`):
print(` that can be obtained `):
print(` from the letter LET1, with its accompanying free-space set `):
print(` FS1, by performing ONE LL-operation `):
print(` L an R denote ( and ) resp. `):
fi:
 
if nops([args])=1 and op(1,[args])=JoinRL then 
print(`JoinRL(LET1,FreeSpace,j,L,R): Given a letter LET1,`):
print(` (phrased in terms of L and R), `):
print(` a set of integers, FreeSpace, and a location j`):
print(` such that LET1[1][j-1]=R, LET1[1][j]=L does the`):
print(` operation pRLq->pq and its corresponding`):
print(` effect on LET1[2] (the list of places with the`):
print(` opening, and  also outputs the new FreeSpace `):
print(` obtained by removing all the points between the two letters of `):
print(` this RL i.e. LET1[2][j-1]+1..LET1[2][j]-1 `):
fi:
 
if nops([args])=1 and op(1,[args])=JoinRR then 
print(` JoinRR(LET1,FreeSpace,j,L,R): Given a letter LET1, `):
print(` (phrased in terms of L and R),  `):
print(` a set of integers, FreeSpace, and a location j `):
print(` such that LET1[1][j]=LET1[1][j-1]=R does the `):
print(` operation LpLqRR->LpRq and its corresponding `):
print(` effect on LET1[2] (the list of places with the `):
print(` opening, and  also outputs the new FreeSpace `):
print(` obtained by removing all the points between these `):
print(` two adjacent Rs i.e. LET1[2][j-1]+1..LET1[2][j]-1 `):
fi:
 
if nops([args])=1 and op(1,[args])=JoinLL then 
print(` JoinLL(LET1,FreeSpace,i,L,R): Given a letter LET1, `):
print(` (phrased in terms of L and R), `):
print(` a set of integers, FreeSpace, and a location i `):
print(` such that LET1[1][i]=LET1[1][i+1]=L does the `):
print(` operation LLpRqR->pLqR and its corresponding `):
print(` effect on LET1[2] (the list of places with the `):
print(` opening, and  also outputs the new FreeSpace `):
print(` obtained by removing all the points between these `):
print(` two adjacent Ls i.e. LET1[2][i]+1..LET1[2][i+1]-1 `):
fi:
 
if nops([args])=1 and op(1,[args])=FindL then 
print(` FindL(w,j,L,R): Given a legal word w (in {L,R}), and `):
print(` a place j  that is an R, find the location of its L-mate`):
fi:
 
if nops([args])=1 and op(1,[args])=FindR then 
print(` FindR(w,i,L,R): Given a legal word w (in {L,R}), and`):
print(` a place i that is an L, find the location of its R-mate`):
fi:
 
if nops([args])=1 and op(1,[args])=IsRightLetter then 
print(`IsRightLetter(LET1,L,R): Decides whether the letter`):
print(`LET1 is a right letter`):
fi:
 
if nops([args])=1 and op(1,[args])=IsLeftLetter then 
print(`IsLeftLetter(LET1,L,R): Decides whether the letter`):
print(`LET1 is a left letter`):
fi:
 
if nops([args])=1 and op(1,[args])=Alphabet then 
print(`Alphabet(n,L,R): All the letters for a self-avoiding polygon`):
print(`bounded in the strip [0,n], using the letters`):
print(` L and R `):
fi:
 
if nops([args])=1 and op(1,[args])=IncSeq then 
print(`IncSeq(m,n,k): The set of inreasing sequences [i1, ..., ik]`):
fi:
 
if nops([args])=1 and op(1,[args])=Parenths then
print(` Parenths(n,L,R): All legal parentheses with n L and n R as  `):
print(` left and right parantheses respectively `):
print(` e.g. Perenths(2,L,R) should give {[L,R,L,R],[L,L,R,R]} `):
fi:
 
end:
 
 
#Parenths(n,L,R): All legal parentheses with n L and n R as 
#left and right parantheses respectively
#e.g. Perenths(2,L,R) should give {[L,R,L,R],[L,L,R,R]}
 
Parenths:=proc(n,L,R)
local gu,gu1,gu2,k,i1,i2:
option remember:
if n=0 then
RETURN({[]}):
fi:
 
gu:={}:
for k from 0 to n-1 do
gu1:=Parenths(k,L,R):
gu2:=Parenths(n-k-1,L,R):
 
for i1 from 1 to nops(gu1) do
for i2 from 1 to nops(gu2) do
gu:=gu union {[L,op(gu1[i1]),R,op(gu2[i2])]}:
od:
od:
od:
gu:
end:
 
 
#IncSeq(m,n,k): The set of inreasing sequences [i1, ..., ik]
#of k integers, such that m<=i1< ...<ik<=n
 
 
IncSeq:=proc(m,n,k)
local gu,mu,i:
option remember:
if k<0 or m>n then
RETURN({}):
fi:
 
if k=0 then
RETURN({[]}):
fi:
 
if m=n then
if k=1 then
 RETURN({[m]}):
else
 RETURN({}):
fi:
fi:
 
 
gu:=IncSeq(m,n-1,k):
mu:=IncSeq(m,n-1,k-1):
 
for i from 1 to nops(mu) do
gu:=gu union {[op(mu[i]),n]}:
od:
 
gu:
end:
 
 
#Alphabet(n,L,R): All the letters for a self-avoiding polygon
#bounded in the strip [0,n], using the letters
#L and R
Alphabet:=proc(n,L,R)
local mu,lu,gu,k,i1,i2:
gu:={}:
 
for k from 1 to (n+1)/2 do
mu:=Parenths(k,L,R):
lu:=IncSeq(0,n,2*k):
 
for i1 from 1 to nops(mu) do
for i2 from 1 to nops(lu) do
gu:=gu union {[mu[i1],lu[i2]]}:
od:
od:
 
od:
 
gu:
 
end:
 
 
 
#IsLeftLetter(LET1,L,R): Decides whether the letter
#LET1 is a left letter
IsLeftLetter:=proc(LET1,L,R)
local mu,k,i:
mu:=LET1[1]:
k:=nops(mu)/2:
 
if mu=[seq(op([L,R]),i=1..k)] then
RETURN(true):
else
RETURN(false):
fi:
 
end:
 
 
 
#IsRightLetter(LET1,L,R): Decides whether the letter
#LET1 is a right letter
IsRightLetter:=proc(LET1,L,R)
local mu,k,i,mu1:
mu:=LET1[1]:
k:=nops(mu)/2:
mu1:=[L,seq(op([L,R]),i=1..k-1),R]:
if mu=mu1 then
RETURN(true):
else
RETURN(false):
fi:
 
end:
 
 
# FindR(w,i,L,R): Given a legal word w (in {L,R}), and
# a place i that is an L, find the location of its R-mate
FindR:=proc(w,i,L,R)
local w1,j,gu:
 
if w[i]<>L then
ERROR(` The `, i, `place of`, w, `should have been an `, L):
fi:
 
gu:=1:
w1:=subs({L=1,R=-1},w):
 
for j from i+1 while gu>0 do
gu:=gu+w1[j]:
od:
j-1:
end:
 
 
 
# FindL(w,j,L,R): Given a legal word w (in {L,R}), and
# a place j  that is an R, find the location of its L-mate
FindL:=proc(w,j,L,R)
local w1,i,gu:
 
if w[j]<>R then
ERROR(` The `, i, `place of`, w, `should have been an `, L):
fi:
 
gu:=1:
w1:=subs({L=-1,R=1},w):
 
for i from j-1 by -1 while gu>0 do
gu:=gu+w1[i]:
od:
i+1:
end:
 
 
#JoinLL(LET1,FreeSpace,extw,i,L,R): Given a letter LET1,
# (phrased in terms of L and R), 
# a set of integers, FreeSpace, and a location i
# such that LET1[1][i]=LET1[1][i+1]=L does the
# operation LLpRqR->pLqR and its corresponding
# effect on LET1[2] (the list of places with the
# opening, and  also outputs the new FreeSpace
# obtained by removing all the points between these
# two adjacent Ls i.e. LET1[2][i]+1..LET1[2][i+1]-1
JoinLL:=proc(LET1,FreeSpace,extw,i,L,R)
local w,v,w1,v1,FreeSpace1,j,j1:
w:=LET1[1]:
v:=LET1[2]:
 
if not (w[i]=L and w[i+1]=L) then
ERROR(` The locations `, i , ` and `, i+1, `should be `, L):
fi:
 
j:=FindR(w,i,L,R):
j1:=FindR(w,i+1,L,R):
 
w1:=[op(1..i-1,w),op(i+2..j1-1,w),L,op(j1+1..j-1,w),op(j..nops(w),w)]:
v1:=[op(1..i-1,v),op(i+2..nops(v),v)]:
 
FreeSpace1:=FreeSpace minus {seq(j,j=v[i]+1..v[i+1]-1)}:
 
[w1,v1],FreeSpace1,extw+v[i+1]-v[i]:
end:
 
 
 
#JoinRR(LET1,FreeSpace,extw,j,L,R): Given a letter LET1,
# (phrased in terms of L and R), 
# a set of integers, FreeSpace, and a location j
# such that LET1[1][j]=LET1[1][j-1]=R does the
# operation LpLqRR->LpRq and its corresponding
# effect on LET1[2] (the list of places with the
# opening, and  also outputs the new FreeSpace
# obtained by removing all the points between these
# two adjacent Rs i.e. LET1[2][j-1]+1..LET1[2][j]-1
JoinRR:=proc(LET1,FreeSpace,extw,j,L,R)
local w,v,w1,v1,FreeSpace1,i,i1:
w:=LET1[1]:
v:=LET1[2]:
 
if not (w[j]=R and w[j-1]=R) then
ERROR(` The locations `, j , ` and `, j-1, `should be `, R):
fi:
 
i:=FindL(w,j,L,R):
i1:=FindL(w,j-1,L,R):
 
w1:=[op(1..i1-1,w),R,op(i1+1..j-2,w),op(j+1..nops(w),w)]:
v1:=[op(1..j-2,v),op(j+1..nops(v),v)]:
 
FreeSpace1:=FreeSpace minus {seq(i,i=v[j-1]+1..v[j]-1)}:
 
[w1,v1],FreeSpace1,extw+v[j]-v[j-1]:
end:
 
 
 
 
#JoinRL(LET1,FreeSpace,extw,j,L,R): Given a letter LET1,
# (phrased in terms of L and R), 
# a set of integers, FreeSpace, and a location j
# such that LET1[1][j-1]=R, LET1[1][j]=L does the
# operation pRLq->pq and its corresponding
# effect on LET1[2] (the list of places with the
# opening, and  also outputs the new FreeSpace
# obtained by removing all the points between these
# the R and the L of the deleted RL
# i.e. LET1[2][j-1]+1..LET1[2][j]-1
JoinRL:=proc(LET1,FreeSpace,extw,j,L,R)
local w,v,w1,v1,FreeSpace1,i:
w:=LET1[1]:
v:=LET1[2]:
 
if not (w[j-1]=R and w[j]=L) then
ERROR(` The locations `, j-1 , ` and `, j, `should be `, R, L):
fi:
w1:=[op(1..j-2,w),op(j+1..nops(w),w)]:
v1:=[op(1..j-2,v),op(j+1..nops(v),v)]:
 
FreeSpace1:=FreeSpace minus {seq(i,i=v[j-1]+1..v[j]-1)}:
 
[w1,v1],FreeSpace1,extw+v[j]-v[j-1]:
end:
 
 
 
#OneStepLL(LET1,FS1,extw,L,R): Set of all triples [LET2,FS2,extw]
# that can be obtained
#from the letter LET1, with its accompanying free-space set
#FS1, and extra-weight, extw, by performing ONE `LL-operation'
#L an R denote ( and ) resp.
#
OneStepLL:=proc(LET1,FS1,extw,L,R)
local i,i1,w,v,gu:
w:=LET1[1]:
v:=LET1[2]:
gu:={}:
 
for i from 1 to nops(w)-1 do
if  w[i]=L and w[i+1]=L and 
  {seq(i1,i1=v[i]+1..v[i+1]-1)} minus FS1={} then
 gu:=gu union {[JoinLL(LET1,FS1,extw,i,L,R)]}:
fi:
od:
gu:
end:
 
 
 
#OneStepRR(LET1,FS1,extw,L,R): Set of all triples [LET2,FS2,exwt]
# that can be obtained
#from the letter LET1, with its accompanying free-space set
#FS1, by performing ONE `RR-operation'
#L an R denote ( and ) resp.
#
OneStepRR:=proc(LET1,FS1,extw,L,R)
local j,i1,w,v,gu:
w:=LET1[1]:
v:=LET1[2]:
gu:={}:
 
for j from 2 to nops(w) do
if  w[j]=R and w[j-1]=R and 
  {seq(i1,i1=v[j-1]+1..v[j]-1)} minus FS1={} then
 gu:=gu union {[JoinRR(LET1,FS1,extw,j,L,R)]}:
fi:
od:
gu:
end:
 
 
#OneStepRL(LET1,FS1,extw,L,R): Set of all pairs [LET2,FS2,extw]
# that can be obtained
#from the letter LET1, with its accompanying free-space set
#FS1, by performing ONE `RL-operation'
#L an R denote ( and ) resp.
#
OneStepRL:=proc(LET1,FS1,extw,L,R)
local j,i1,w,v,gu:
w:=LET1[1]:
v:=LET1[2]:
gu:={}:
 
for j from 2 to nops(w) do
if  w[j-1]=R and w[j]=L and 
  {seq(i1,i1=v[j-1]+1..v[j]-1)} minus FS1={} then
 gu:=gu union {[JoinRL(LET1,FS1,extw,j,L,R)]}:
fi:
od:
gu:
end:
 
 
 
#OneStepAll(LET1,FS1,extw,L,R): Set of all pairs [LET2,FS2,extw]
# that can be obtained
#from the letter LET1, with its accompanying free-space set
#FS1, by performing ONE `RL, or RR- or LL-operation'
#L an R denote ( and ) resp.
#
OneStepAll:=proc(LET1,FS1,extw,L,R):
OneStepLL(LET1,FS1,extw,L,R) union
OneStepRR(LET1,FS1,extw,L,R) union
OneStepRL(LET1,FS1,extw,L,R):
end:
 
#PrePreFollowers(LET1,a,b,L,R): All the pre-pre-letters that
#can follow the letter LET1 in the ambient strip, [a,b]
#where L and R denote ( and ) resp.
PrePreFollowers:=proc(LET1,a,b,L,R)
local gu,mu,mu1,i:
 
 
gu:={[LET1,{seq(i,i=a..b)} minus convert(LET1[2],set),0]}:
mu:={[LET1,{seq(i,i=a..b)} minus convert(LET1[2],set),0]}:
 
while mu<>{} do
mu1:={}:
 
for i from 1 to nops(mu) do
mu1:=mu1 union OneStepAll(op(mu[i]),L,R):
od:
 
gu:=gu union mu1:
mu:=mu1:
od:
 
gu:
end:
 
 
#LebensRaum(LET1,FS1,i): Given a Pre-Letter LET1 with its
#free-space set FS1, and an integer i, finds the set range
#where the i^th component can roam
#the output is the lowest and the highest that it can
#venture to
LebensRaum:=proc(LET1,FS1,i)
local v,mina,maxa,i1:
v:=LET1[2]:
 
for i1 from v[i] while member(i1,FS1 union {v[i]}) do
od:
maxa:=i1-1:
 
for i1 from v[i] by -1 while member(i1,FS1 union {v[i]}) do
od:
mina:=i1+1:
 
mina,maxa:
end:
 
 
 
#LebensRaums(LET1,FS1): Given a pre-letter LET1, and 
#a free space set FS1, finds the list of pairs
#indicating the LebensRaum of each component
LebensRaums:=proc(LET1,FS1)
local i,gu:
gu:=[]:
 
for i from 1 to nops(LET1[1]) do
 gu:=[op(gu),[LebensRaum(LET1,FS1,i)]]:
od:
gu:
end:
 
 
#IncSeqSand(BOUNDS,Centers): given a list of pairs [[a1,b1],...,[ak,bk]]
#and a set of centers, Centers
#finds the set of increasing sequences x1<=x2<=..<=xk
#such that ai<=xi<=bi, for i=1...k
#followed the extra teritory taken
IncSeqSand:=proc(BOUNDS,Centers)
local k,mu,gu,a,b,i,j,j1,BOUNDS1,Centers1,vec,vk:
option remember:
k:=nops(BOUNDS):
 
if k=0 then
RETURN({[]}):
fi:
 
a:=BOUNDS[k][1]:
b:=BOUNDS[k][2]:
vk:=Centers[k]:
if k=1 then
if not (a<=vk and vk<=b) then
 ERROR(`Bad input`):
fi:
 
RETURN({seq([[i], {seq(j,j=i..vk)},vk-i],i=a..vk)}
union {seq([[i], {seq(j,j=vk..i)},i-vk],i=vk+1..b)}):
fi:
 
BOUNDS1:=[op(1..k-1,BOUNDS)]:
Centers1:=[op(1..k-1,Centers)]:
a:=BOUNDS[k][1]:
b:=BOUNDS[k][2]:
mu:=IncSeqSand(BOUNDS1,Centers1):
gu:={}:
for i from 1 to nops(mu) do
vec:=mu[i][1]:
 
for j from max(vec[k-1]+1,a) to vk do
gu:=gu union {[[op(vec),j],mu[i][2] union  {seq(j1,j1=j..vk)},mu[i][3]+vk-j]}
od:
 
for j from vk+1 to b do
gu:=gu union {[[op(vec),j],mu[i][2] union  {seq(j1,j1=vk..j)},mu[i][3]+j-vk]}
od:
 
od:
 
gu:
end:
 
 
 
 
#PreFollowers(LET1,a,b,L,R): All the pre-letters that
#can follow the letter LET1 in the ambient strip, [a,b]
#where L and R denote ( and ) resp.
PreFollowers:=proc(LET1,a,b,L,R)
local gu,mu,BOUNDS,lu,i,j:
 
mu:=PrePreFollowers(LET1,a,b,L,R):
 
gu:={}:
 
for i from 1 to nops(mu) do
BOUNDS:=LebensRaums(op(mu[i])):
lu:=IncSeqSand(BOUNDS,mu[i][1][2]):
 
for j from 1 to nops(lu) do
 gu:=gu union {[[mu[i][1][1],lu[j][1]],mu[i][2] minus lu[j][2],
   mu[i][3]+lu[j][3]+nops(mu[i][1][1])]}:
od:
 
od:
 
gu:
 
end:
 
 
 
#LeftLetters(a,b,L,R): All the left
#letters for a self-avoiding polygon
#bounded in the strip [a,b], using the letters
#L and R followed by their length
LeftLetters:=proc(a,b,L,R)
local mu,lu,gu,k,i1,i2,ru,lu2,i:
gu:={}:
 
for k from 1 to (b-a+1)/2 do
mu:=[seq(op([L,R]),i1=1..k)]:
lu:=IncSeq(a,b,2*k):
 
 
for i2 from 1 to nops(lu) do
lu2:=lu[i2]:
ru:=0:
for i from 1 to k do
ru:=ru+lu2[2*i]-lu2[2*i-1]+2:
od:
 
gu:=gu union {[[mu,lu2],ru] }:
od:
 
od:
 
gu:
 
end:
 
 
#SetToUnInt(kv): Given a set of integers, kv, converts
#it list-of-intervals notation [[a1,b1],[a2,b2],[a3,b3],...]
#s.t. kv={a1..b1} union {a2..b2} ....
SetToUnInt:=proc(kv)
local a1,b1,i,gu:
if kv={} then
RETURN([]):
fi:
 
a1:=min(op(kv)):
 
for i from a1 while member(i,kv) do
od:
b1:=i-1:
gu:=SetToUnInt(kv minus {seq(i,i=a1..b1)}):
[[a1,b1],op(gu)]:
end:
 
 
 
#PlacesToInsert1(a,b): Given two integers
#finds all the intervals [i,j] such that a<=i<j<=b
PlacesToInsert1:=proc(a,b)
local gu,i,j:
gu:={}:
for i from a to b-1 do
for j from i+1 to b do
gu:=gu union {[i,j]}:
od:
od:
gu:
end:
 
#PlacesToInsert(UnInts): Given a set of free-space, UnInts
#in terms of a list of intervals, finds all the
#intervals [a,b] of length at least 2 that may be
#inserted
PlacesToInsert:=proc(UnInts)
local k,ak,bk:
k:=nops(UnInts):
 
if k=0 then
 RETURN({}):
fi:
 
ak:=UnInts[k][1]:
bk:=UnInts[k][2]:
 
PlacesToInsert([op(1..k-1,UnInts)]) union PlacesToInsert1(ak,bk):
 
 
end:
 
 
#Stick1(LET1,L,R,a,b): given a pre-letter, and an interval
#[a,b] finds the resulting pre-letter obtained by inserting
#[L,R] in the appropirate place in LET1[1][1] and
#[a,b] in the appropirate place in LET1[1][2]
#and replacing LET1[2] by LET1[2] minus {a..b}
#and replacing LET1[3] by LET1[3]+b-a+2
Stick1:=proc(LET1,L,R,a,b)
local w,v,FS,i,ku,extw:
 
 
w:=LET1[1][1]:
v:=LET1[1][2]:
FS:=LET1[2]:
extw:=LET1[3]:
ku:={seq(i,i=a..b)}:
 
if ku minus FS<>{} then
ERROR(`something is wrong`):
fi:
 
 
 
if a>v[nops(v)] then
w:=[op(w),L,R]:
v:=[op(v),a,b]:
RETURN([[w,v],FS minus ku,extw+b-a+2]):
fi:
 
 
if b<v[1] then
w:=[L,R,op(w)]:
v:=[a,b,op(v)]:
RETURN([[w,v],FS minus ku,extw+b-a+2]):
fi:
 
 
 
for i from 1 while v[i]<a do
od:
i:=i-1:
 
w:=[op(1..i,w),L,R,op(i+1..nops(w),w)]:
v:=[op(1..i,v),a,b,op(i+1..nops(v),v)]:
[[w,v],FS minus ku,extw+b-a+2]:
end:
 
 
#FillIns1(LET1,L,R): Given a pre-letter, LET1, where
#L and R denote ( and ), finds all the possible pre-letters
#obtained from it by inserting ONE new interval in
#a legal place
FillIns1:=proc(LET1,L,R)
local gu,mu,i:
gu:={}:
mu:=PlacesToInsert(SetToUnInt(LET1[2])):
 
for i from 1 to nops(mu) do
gu:=gu union {Stick1(LET1,L,R,mu[i][1],mu[i][2])}:
 
od:
gu:
end:
 
#FillInsOneStep(SetLET1,L,R): given a set of pre-letters
#finds all the pre-letters obtained from their members
#by inserting ONE interval
FillInsOneStep:=proc(SetLET1,L,R)
local gu,i:
gu:={}:
 
for i from 1 to nops(SetLET1) do
gu:=gu union FillIns1(SetLET1[i],L,R):
od:
gu:
end:
 
#DerivedPreLetters(LET1,L,R): Given a pre-letter, finds
#all pre-letters obtained from it by inserting any number
#of eligible intervals
DerivedPreLetters:=proc(LET1,L,R)
local mu,gu:
mu:={LET1}:
gu:={}:
 
 
while mu<>{} do
gu:=gu union mu:
mu:=FillInsOneStep(mu,L,R):
od:
 
gu:
 
end:
 
 
 
 
#Followers(LET1,a,b,L,R): All the letters that
#can follow the letter LET1 in the ambient strip, [a,b]
#where L and R denote ( and ) resp.
#followed by the extra-weight
Followers:=proc(LET1,a,b,L,R)
local mu,gu1,gu,i:
mu:=PreFollowers(LET1,a,b,L,R):
 
gu1:={}:
 
for i from 1 to nops(mu) do
gu1:=gu1 union DerivedPreLetters(mu[i],L,R):
od:
 
gu:={}:
 
for i from 1 to nops(gu1) do
gu:=gu union {[gu1[i][1],gu1[i][3]]}:
od:
 
gu:
end:
 
 
 
#RightLetters(a,b,L,R): All the right
#letters for a self-avoiding polygon
#bounded in the strip [a,b], using the letters
#L and R followed by their finishing-length
#(i.e. what it takes to close the polygon and
#terminate the word)
#
RightLetters:=proc(a,b,L,R)
local mu,lu,gu,k,i1,i2,ru,lu2,i:
gu:={}:
 
for k from 0 to (b-a-1)/2 do
mu:=[L,seq(op([L,R]),i1=1..k),R]:
lu:=IncSeq(a,b,2*k+2):
 
 
for i2 from 1 to nops(lu) do
lu2:=lu[i2]:
ru:=0:
for i from 1 to k+1 do
ru:=ru+lu2[2*i]-lu2[2*i-1]:
od:
 
gu:=gu union {[[mu,lu2],ru] }:
od:
 
od:
 
gu:
 
end:
 
 
#MarCha(n): The Markov chain of self-avoiding polygons
#that live in the strip [0,n]
#It is asummed that the vertices (states) are labelled
#by positive integers. The Markov Chain is given as a
#list of lengh 3: [LeftLetters,RightLetters,TransitionMatrix]
#The TransitionMatrix is a list of lists, where the
#i^th entry is the set of pairs [j,wt(i,j)], where
#j is the a neighbor and wt(i,j) is the weight of
#the arc from i to j
MarCha:=proc(n)
local Vertices,L,R,T,i,lu,gu,S,mu1,mu1a,j,mu2:
Vertices:=Alphabet(n,L,R):
 
T[L]:=LeftLetters(0,n,L,R):
 
 
for i from 1 to nops(Vertices) do
 T[Vertices[i]]:=Followers(Vertices[i],0,n,L,R):
od:
 
lu:=RightLetters(0,n,L,R):
 
for i from 1 to nops(lu) do
 T[lu[i][1]]:=T[lu[i][1]] union {[R,lu[i][2]]}:
od:
 
T[R]:={}:
 
Vertices:=Vertices union {L,R}:
Vertices:=convert(Vertices,list):
 
for j from 1 to nops(Vertices) do
S[Vertices[j]]:=j:
od:
 
gu:=[]:
 
for i from 1 to nops(Vertices) do
 
mu1:=T[Vertices[i]]:
mu2:={}:
 
for j from 1 to nops(mu1) do
 mu1a:=op(j,mu1):
 mu2:=mu2 union {[S[mu1a[1]],mu1a[2]]}:
od:
 
gu:=[op(gu),mu2]:
od:
 
[{S[L]},{S[R]},gu]:
end:
 
 
 
#SolveMC3(MC,t): Given a type-III Markov Chain [LeftLetters,RightLetters,
#TransMatrx], finds the sum of the weights of all words that
#start with a letter of LeftLetters and ends with a letter
#of RightLetters, as a generating function in the variable t
#
SolveMC3:=proc(MC,t)
local eq,var,a,gu,LL1,RL1,TM,eq1,lu,i,j:
 
LL1:=MC[1]:
RL1:=MC[2]:
TM:=MC[3]:
eq:={}:
var:={}:
 
for i from 1 to nops(TM) do
var:=var union {a[i]}:
 
if member(i,RL1) then
eq1:=a[i]-1:
else
eq1:=a[i]:
fi:
 
lu:=TM[i]:
 
for j from 1 to nops(lu) do
 eq1:=eq1-a[lu[j][1]]*t^lu[j][2]:
od:
 
eq:=eq union {eq1}:
od:
 
var:=solve(eq,var):
 
if var=NULL then
RETURN(0):
fi:
 
gu:=0:
 
for i from 1 to nops(LL1) do
gu:=normal(gu+subs(var,a[LL1[i]])):
od:
 
gu:
 
end:
 
 
#gf(n,t): The generating function of self-avoiding polygons
#immersed in a fixed strip
gf:=proc(n,t) option remember: SolveMC3(MarCha(n),t): end:
 
#GF(n,t): The generating function for animals of width<=n
GF:=proc(n,t):sort(factor(normal(gf(n,t)-gf(n-1,t)))):end:
 
#Gf(n,t):The generating function for animals of width=n
Gf:=proc(n,t):normal(gf(n,t)-2*gf(n-1,t)+gf(n-2,t)):end:
 
 
 
 
#SolveMC3series(MC,M): Given a Markov Chain finds the sequence
#whose i^th member is the number of paths of weight i
#for 0<=i<=M
SolveMC3series:=proc(MC,M)
local a,i,SID,gu,gadol,LL1,RL1,TM,mispar,j,lu,i1,m,mu,sakhen,kha:
LL1:=MC[1]:
RL1:=MC[2]:
TM:=MC[3]:
mispar:=nops(TM):
 
gu:={}:
for i from 1 to mispar do
lu:=TM[i]:
 
for j from 1 to nops(lu) do
gu:=gu union {lu[j][2]}:
od:
 
od:
 
gadol:=max(op(gu)):
 
for i from 1 to mispar do
if member(i,RL1) then
a[i]:=[seq(0,i1=1..gadol-1),1]:
else
a[i]:=[seq(0,i1=1..gadol)]:
fi:
od:
 
 
SID:=[]:
for m from 1 to M do
 
for i from 1 to mispar do
mu:=TM[i]:
lu:=0:
 
for j from 1 to nops(mu) do
sakhen:=a[mu[j][1]]:
lu:=lu+sakhen[nops(sakhen)-mu[j][2]+1]:
od:
 
kha[i]:=lu:
od:
 
gu:=0:
 
for i from 1 to nops(LL1) do
gu:=gu+kha[LL1[i]]:
od:
 
SID:=[op(SID),gu]:
 
for i from 1 to mispar do
a[i]:=[op(2..nops(a[i]),a[i]),kha[i]]:
od:
 
od:
SID:
end:
 
 
 
 
#gfSeries(n,M): The sequence whose i^th term is the
#number of self-avoiding polygons of length 2i immersed
#in the strip 0<=y<=n, for 1<=i<=M 
gfSeries:=proc(n,M)
local i,mu:
option remember:
mu:=SolveMC3series(MarCha(n),2*M):
[seq(mu[2*i],i=1..M)]:
end:
 
 
#GFSeries(n,M): The sequence whose i^th term is the
#number of self-avoiding polygons of length 2i of width<=n
#for 1<=i<=M 
GFSeries:=proc(n,M)
local i,mu1,mu2:
mu1:=gfSeries(n,M):
 
if n=1 then
RETURN(mu1):
fi:
 
mu2:=gfSeries(n-1,M):
[seq(mu1[i]-mu2[i],i=1..M)]:
end:
 
 
#GfSeries(n,M): The sequence whose i^th term is the
#number of self-avoiding polygons of length 2i of width=n
#for 1<=i<=M 
GfSeries:=proc(n,M)
local i,mu1,mu2:
 
mu1:=GFSeries(n,M):
 
if n=1 then
RETURN(mu1):
fi:
 
 
mu2:=GFSeries(n-1,M):
[seq(mu1[i]-mu2[i],i=1..M)]:
end:
 
 
 
###########
 
#SolveMC3seriesS(MC,M,s): Given a Markov Chain finds the sequence
#whose i^th member is the generating function  of paths of weight i
#with weight s^(length)
#for 0<=i<=M
SolveMC3seriesS:=proc(MC,M,s)
local a,i,SID,gu,gadol,LL1,RL1,TM,mispar,j,lu,i1,m,mu,sakhen,kha:
LL1:=MC[1]:
RL1:=MC[2]:
TM:=MC[3]:
mispar:=nops(TM):
 
gu:={}:
for i from 1 to mispar do
lu:=TM[i]:
 
for j from 1 to nops(lu) do
gu:=gu union {lu[j][2]}:
od:
 
od:
 
gadol:=max(op(gu)):
 
for i from 1 to mispar do
if member(i,RL1) then
a[i]:=[seq(0,i1=1..gadol-1),1]:
else
a[i]:=[seq(0,i1=1..gadol)]:
fi:
od:
 
 
SID:=[]:
for m from 1 to M do
 
for i from 1 to mispar do
mu:=TM[i]:
lu:=0:
 
for j from 1 to nops(mu) do
sakhen:=a[mu[j][1]]:
lu:=expand(lu+s*sakhen[nops(sakhen)-mu[j][2]+1]):
od:
 
kha[i]:=lu:
od:
 
gu:=0:
 
for i from 1 to nops(LL1) do
gu:=gu+kha[LL1[i]]:
od:
 
SID:=[op(SID),gu]:
 
for i from 1 to mispar do
a[i]:=[op(2..nops(a[i]),a[i]),kha[i]]:
od:
 
od:
SID:
end:
 
 
 
 
 
#gfSeriesHeight(n,M,s): The sequence whose i^th term is the
#g.f. of self-avoiding polygons of weight 2i 
#with weight=s^(height) immersed
#in the strip 0<=y<=n, for 1<=i<=M and 
gfSeriesHeight:=proc(n,M,s)
local i,mu:
option remember:
mu:=SolveMC3seriesS(MarCha(n),2*M,s):
[seq(normal(mu[2*i]/s),i=1..M)]:
end:
 
 
#GFSeriesHeight(n,M): The sequence whose i^th term is the
#g.f. of self-avoiding polygons of length 2i of width<=n
#for 1<=i<=M with weigth s^(height)
GFSeriesHeight:=proc(n,M,s)
local i,mu1,mu2:
option remember:
mu1:=gfSeriesHeight(n,M,s):
 
if n=1 then
RETURN(mu1):
fi:
 
mu2:=gfSeriesHeight(n-1,M,s):
[seq(expand(mu1[i]-mu2[i]),i=1..M)]:
end:
 
 
#GfSeriesHeight(n,M,s): The sequence whose i^th term is the
#g.f. of self-avoiding polygons of length 2i of width=n
#with weight s^(height)
#for 1<=i<=M 
GfSeriesHeight:=proc(n,M,s)
local i,mu1,mu2:
option remember:
mu1:=GFSeriesHeight(n,M,s):
 
if n=1 then
RETURN(mu1):
fi:
 
mu2:=GFSeriesHeight(n-1,M,s):
[seq(expand(mu1[i]-mu2[i]),i=1..M)]:
end:
 
 
#SAPlhw(l,h,w): the number of self-avoiding polygons
#with length 2l, height h, and width w
SAPlhw:=proc(l,h,w)
local gu,M,s:
option remember:
M:=l+h:
gu:=GfSeriesHeight(w,M,s):
coeff(gu[l],s,h):
end:
 
#SAP(l): The number of self-avoiding polygons of
#lentgh 2*l
SAP:=proc(l)
local h,w,gu:
option remember:
gu:=0:
 
 
for w from 1 to l/2 do
for h from w+1 to l-w do
gu:=gu+2*SAPlhw(l,h,w):
od:
od:
 
 
for w from 1 to l/2 do
gu:=gu+SAPlhw(l,w,w):
od:
 
gu:
 
end:
 
#SAPseries(L): The sequence whose i^th term is
#the number of self-avoiding polygons of length 2i
#for i=1...L
SAPseries:=proc(L)
local gu,i:
gu:=[]:
for i from 1 to L do
gu:=[op(gu),SAP(i)]:
print(gu):
od:
gu:
end:
 
########Batchelor-Yung stuff
 
 
#LeftLettersBY(a,b,L,R): 
#Like LeftLetters(a,b,L,R): but must have
#0 in them
LeftLettersBY:=proc(a,b,L,R)
local i,mu,gu:
mu:=LeftLetters(a,b,L,R):
gu:={}:
 
for i from 1 to nops(mu) do
if mu[i][1][2][1]=0 then
gu:=gu union {mu[i]}:
fi:
od:
gu:
end:
 
 
#MarChaBY(n): Like MarCha(n) BUT
#The Left letters all start with [L,...],[0,...]
MarChaBY:=proc(n)
local Vertices,L,R,T,i,lu,gu,S,mu1,mu1a,j,mu2:
Vertices:=Alphabet(n,L,R):
 
T[L]:=LeftLettersBY(0,n,L,R):
 
 
for i from 1 to nops(Vertices) do
 T[Vertices[i]]:=Followers(Vertices[i],0,n,L,R):
od:
 
lu:=RightLetters(0,n,L,R):
 
for i from 1 to nops(lu) do
 T[lu[i][1]]:=T[lu[i][1]] union {[R,lu[i][2]]}:
od:
 
T[R]:={}:
 
Vertices:=Vertices union {L,R}:
Vertices:=convert(Vertices,list):
 
for j from 1 to nops(Vertices) do
S[Vertices[j]]:=j:
od:
 
gu:=[]:
 
for i from 1 to nops(Vertices) do
 
mu1:=T[Vertices[i]]:
mu2:={}:
 
for j from 1 to nops(mu1) do
 mu1a:=op(j,mu1):
 mu2:=mu2 union {[S[mu1a[1]],mu1a[2]]}:
od:
 
gu:=[op(gu),mu2]:
od:
 
[{S[L]},{S[R]},gu]:
end:
 
 
#gfBatchelorYungBY(n,t): The generating function of 
#self-avoiding polygons
#immersed the strip [0,n] and such that the left-most letter
#has a 0 in it. It gives the quantity described in
#M.T. Batchelor, and C.M. Yang, J. Physics A: Math. gen. 27 (1994)
#4055-4067 (the denominators on p. 4059 and 4066 and the 
#numerators on p. 4064)
gfBatchelorYung:=proc(n,t) option remember: 
sort(factor(SolveMC3(MarChaBY(n),t))): end: