print( `This is FPL, Version of April 8, 2008`):


ezra1:=proc()
if args=NULL then
 print(` The supporting procedures are: AllInterviews`):
 print(` ASM , AtF, Conns, DrawFPL1t, DrawFPL1tGif,`):
 print(` DrawHefresh, ddGif, eipi, EligibleMatches, ShortestEligibleMathces, FPLcAll`):
 print(`Hefresh, Interview, IsSimplePath, Mate, Milon, PathFrom`) :
else
  ezra(args):
fi:
end:

ezra:=proc()
if args=NULL then
print(`This is FPL, to study the Razumov-Stroganov conj. `):
 print(` The main procdures are:  CheckRS`):
print(` Conn, DrawFPL , DrawFPLs, FPL, FPLc, RS, Nei, NeiL, NeiS, Yafe`):


elif nargs=1 and args[1]=AllInterviews then
print(`AllInterviews(e,pi): given a connectivity pi, does all the`):
print(`interviews with all the fpls of all its neighbors.e is just a symbol.`):
print(`For example, do: AllInterviews(e,{{1,2},{3,4}});`):

elif nargs=1 and args[1]=ASM then
 print(`ASM(n): the set of ASMs of size n`):

elif nargs=1 and args[1]=AtF then
print(`AtF(a): inputs an asm a, outputs an FPL as a set of edges`):
print(`in [0,k+1]^2`):

elif nargs=1 and args[1]=CheckRS then
print(`CheckRS(k): verifies the R-S conj. for fpl's of length k`):
print(`for example,try: CheckRS(3);`):

elif nargs=1 and args[1]=Conn then
print(`Conn(F,k): the connectivity of the FPL F of size k`):
print(`For example, try: Conn(FPL(3)[1],3);`):

elif nargs=1 and args[1]=Conns then
print(`Conns(k): all the connectivities of {1, ..., 2k}`):

elif nargs=1 and args[1]=DrawHefresh then
print(`DrawHefresh(k,gu,a,b): Given a pair [Black,Red] of`):
print(`edges, draws them (red also dashed).` ):
print(`for example, try: DrawHefresh(3,Hefresh(FPL(3)[1],FPL(3)[2]),0,0 );`):

elif nargs=1 and args[1]=ddGif then
print(`ddGif(k,F1,F2,a,b,FN): Given two FPLs`):
print(`draws their difference (red are non-edges).` ):
print(`and produces a gif file called FN.gif`):
print(`for example, try: `):
print(`ddGif(3,FPL(3)[1],FPL(3)[2],0,0,out1 );`):

elif nargs=1 and args[1]=DrawFPLno then
print(`DrawFPLno(gu,a,b): draws an fpl with origin at (a,b)`):
print(`For example DrawFPLno(AtF(ASM(2)[1]),0,0);`):

elif nargs=1 and args[1]=DrawFPL then
print(`DrawFPL(gu,a,b,k): draws an fpl with origin at (a,b)`):
print(`For example DrawFPL(FPL(3)[1],0,0,3);`):

elif nargs=1 and args[1]=DrawFPL1 then
print(`DrawFPL1(gu,a,b,k): draws an fpl with origin at (a,b)`):
print(`but does not display it`):
print(`For example DrawFPL(FPL(3)[1],0,0,3);`):

elif nargs=1 and args[1]=DrawFPL1t then
print(`DrawFPL1t(gu,a,b,k): draws an fpl with origin at (a,b), with thickness`):
print(`but does not display it`):
print(`For example DrawFPL1t(FPL(3)[1],0,0,3);`):

elif nargs=1 and args[1]=DrawFPL1tGif then
print(`DrawFPL1tGif(gu,a,b,k,FN): draws an fpl with origin at (a,b), with thickness`):
print(`into a file FN.gif`):
print(`For example DrawFPL1tGif(FPL(3)[1],0,0,3,out2);`):

elif nargs=1 and args[1]=DrawFPLs then
print(`DrawFPLs(S,k,w): draws a set of FPLs S of size k`):
print(`w is the spacing`):
print(`For example, try: DrawFPLs(FPL(3),3,6);`):


elif nargs=1 and args[1]=eipi then
print(`eipi(pi,i): applies e_i to connectivity pi`):
print(`For example, try:`):
print(`eipi({{1,2},{3,4}},2);`):


elif nargs=1 and args[1]=EligibleMatches then
print(`EligibleMatches(k,f,i): given an fpl f of size k, and i in [1,2k], decides`):
print(`who are the  eligible matches, followed by their lengths`):
print(`For example, try:`):
print(`EligibleMathces(3,FPL(3)[1],1);`):

elif nargs=1 and args[1]=ShortestEligibleMatches then
print(`ShortestEligibleMatches(k,f,i): given an fpl f of size k, and i in [1,2k], decides`):
print(`who are the  shortest eligible matches`):
print(`For example, try:`):
print(`ShortestEligibleMathces(3,FPL(3)[1],1);`):

elif nargs=1 and args[1]=FPL then
print(`FPL(n): The set of fpls of size n`):

elif nargs=1 and args[1]=FPLc then
print(`FPLc(k,C): The set of fpls of size k with connectivity C`):
print(`Try: FPLc(3,Conns(3)[1]);`):

elif nargs=1 and args[1]=FPLcAll then
print(`FPLcAll(n): The table of fpls of size n according to connectivity`):
print(`followed by the set of connectivities`):


elif nargs=1 and args[1]=Hefresh then
print(`Hefresh(F1,F2): the set of black edegs in F1 that are not in F2`):
print(`followed by the set of black edges in F2 not in F1`):
print(`for example, try: Hefresh(FPL(3)[1],FPL(3)[2]);`):

elif nargs=1 and args[1]=Interview then
print(`Interview(k,f,C): given an fpl f and a connectivity C, finds the`):
print(`all the the Hefresh (differences) between f and the members`):
print(`of FPLc(k,C). For example, try:`):
print(`Interview(3,FPL(3)[1],Conns(3)[1]);`):

elif nargs=1 and args[1]=IsSimplePath then
print(`IsSimplePath(F1,F2): Given two FLPs, decides whether their`):
print(`symmetric difference is a simple path, and if yes`):
print(`returns its length. For example, try`):
print(`IsSimplePath(FPL(4)[1],FPL(4)[2]);`):

elif nargs=1 and args[1]=Mate then
print(`Mate(F,k,i): Given an FPL of size k and one of the starting`):
print(`vertices i, (1<=i<=2k) finds its mate`):
print(`For example, try:`):
print(`Mate(FPL(3)[1],3,1);`):

elif nargs=1 and args[1]=Milon then
print(`Milon(k): all the terminal vertices labelled 1,..., 2k`):
print(`in a size-k FPL, listed in order. Try: Milon(3);`):


elif nargs=1 and args[1]=PathFrom then
print(`PathFrom(k,P1,f): Given a point P1 in [0,k+1]^2`):
print(`with degree 1, and an fpl f,`):
print(`finds the path that starts with P1.`):
print(`For example, try: `):
print(`PathFrom(3,[1,0],FPL(3)[1]);`):

elif nargs=1 and args[1]=pASM then
 print(`pASM(n): prints the set of ASMs of size n`):

elif nargs=1 and args[1]=RS then
print(`RS(pi): verifies the Razumov-Strogranov Conj. for `):
print(`connectivity pi`):
print(`For example, try: RS(Conns(3)[1]);`):

elif nargs=1 and args[1]=Nei then
print(`Nei(pi): all the pairs [i,sig] such that e_i(sig)=pi`):
print(`For example, try: Nei({{1,2},{3,4}});`):

elif nargs=1 and args[1]=NeiL then
print(`NeiL(pi): all the pairs [i,j] such that e_i(ConnsL(k)[j])=pi`):
print(`For example, try: NeiL([[1,2],[3,4]]);`):

elif nargs=1 and args[1]=NeiS then
print(`NeiS(k): prints out all the neighbors of all connectivities`):
print(`of {1, ..., 2k}. For example, try: NeiS(2);`):

elif nargs=1 and args[1]=Yafe then
print(`Yafe(k,w): makes gif files for all FPL's of size k`):
print(`according to the connectivity outputs .gif files`):
print(`names oka.gif where a ranges from 1 to the  Catalan(k)`):
print(`w is the spacing between fpls`):
print(`For example, try: Yafe(3,w);`):


else
print(`No ezra for`, args):
fi:
end:


#Tkn gives the set T_k(n;a), where a=[a_1, ..., a_k] defined in the
#proof of 1.2.1.1
 
Tkn:=proc(k,n,a)
local nakh,i1,b,i,gu,mu:
 
mu:=DECSEQ(k,n-1):
gu:={}:
 
for i from 1 to nops(mu) do
  b:=op(i,mu):
 
nakh:=1:
 
 
 
  for i1 from 1 to 1 do
    if not ( k-i1+1<=op(i1,b) and op(i1,b)<=op(i1,a) )
        then
            nakh:=0:
    fi:
 
  od:
 
  for i1 from 2 to k do
    if not ( k-i1+1<=op(i1,b) and op(i1,b)<=min( op(i1,a),op(i1-1,a)-1 )   )
        then
            nakh:=0:
    fi:
 
  od:
 
 
if nakh=1 then
  gu:=gu union {b}:
fi:
 
od:
 
gu:
 
end:



#GOGa(k,n,a)  gives  the   set of  k by n Gog-trapezoids  such  that
#the rightmost border is the  vector a
 
GOGa:=proc(k,n,a)
local pip,kvu,firow,b,mu,gu,i,j,l,trap,trap1:
 
if  not k>=1 or not n>=k  or  not nops(a)=k  then
  ERROR(`Improper intput`):
fi:
 
 
if n=k and k=1 then
   if not op(1,a)=1 then
     RETURN({}):
   else
     RETURN({[[1]]}):
   fi:
fi:
 
if n=k then
   if not op(1,a)=k then
      ERROR(`Wrong input`):
   fi:
   mu:=GOGa(k-1,k,[op(2..k,a)]):
    gu:={}:
 
    for i from 1 to nops(mu) do
     trap1:=op(i,mu):
     firow:=op(1,trap1):
    
     gu:=gu union {[[op(firow),k],op(2..k,trap1)]}:
        od:
 RETURN(gu):
fi:
 
gu:={}:    
 
kvu:=Tkn(k,n,a):
 
 
 for pip from 1 to nops(kvu) do
   b:=op(pip,kvu):
 
  mu:=GOGa(k,n-1,b):
 
   for j from 1 to nops(mu) do
    trap:=op(j,mu):
    trap1:=[op(1..n-k,trap)]:
 
     for l from n-k+1 to n-1 do
       trap1:=[op(trap1),[op(op(l,trap)),op(l-(n-k),a)]]:
     od:
      
      trap1:=[op(trap1),[op(k,a)]]:
     gu:=gu union {trap1}:
 
   od:
 
 
od:
 
gu:
 
end:



DECSEQ:=proc(k,n)
local gu,gu1,i1:
option remember:
 
if n=1 and  k=1 then
 RETURN({[1]}):
fi:
 
if k=0  and  n>=1 then
    RETURN({[]}):
fi:
 
if n<1 or  k<1  then
   RETURN({}):
fi:
 
gu:=DECSEQ(k,n-1):
gu1:=DECSEQ(k-1,n):
 
for i1 from 1   to  nops(gu1)  do
  gu:=gu  union   {[n,op(op(i1,gu1))]}:
od:
 
gu:
 
end:
 
 
LOGOG:=proc(k,n)
local gu,gu1,i,i1,vec,nakh:
 
if not k>=1  or  not n>=k then
  RETURN({}):
fi:
 
gu1:=DECSEQ(k,n):
 
gu:={}:
 
 
for i from 1 to nops(gu1) do
 vec:=op(i,gu1):
  nakh:=1:
    for  i1  from 1 to k  do
       if  not op(i1,vec)>=k-i1+1 then
           nakh:=0:

           exit:

       fi:
    od:
 
      if nakh=1 then
          gu:=gu union {vec}:
        fi:
od:
 
 
gu:
 
end:

GOGset:=proc(k,n)
local i,gu,lu:
 
lu:=LOGOG(k,n):
gu:={}:
 
for i from 1 to nops(lu) do
 gu:=gu union GOGa(k,n,op(i,lu)):
od:
 
gu:
 
gu:
 
end:

GOGTOASM:=proc(mt)
local k,mat,mat1,ro,i,j:
 
k:=nops(mt):
 
mat:=array(1..k,1..k):
 
for i from 1 to k do
 for j from 1 to k do
   mat[i,j]:=0:
 od:
od:
 
 
  for i from 1 to k do
   ro:=op(i,mt):
    for j from 1 to nops(ro) do
      mat[i,op(j,ro)]:=1:
    od:
  od:
 
for i from 1 to k-1 do
 for j from 1 to k do
  mat1[i,j]:=mat[i,j]-mat[i+1,j]:
 od:
 
 for j from 1 to k do
  mat1[k,j]:=mat[k,j]:
 od:
 
od:
 
[seq([seq(mat1[i,j],j=1..k)],i=1..k)]:
 
end:

pASM:=proc(k)
local i,gu,asm:
 
gu:=GOGset(k,k):
 
print(`There are`, nops(gu),`Alternating Sign Matrices of size`,k):
print(`Here they all are:`):
 
for i from 1 to nops(gu) do
  asm:=GOGTOASM(op(i,gu)):
 print(op(asm)):
od:
end:

ASM:=proc(k)
local i,gu,asm,lu:
 
gu:=GOGset(k,k):
 
lu:={}:
for i from 1 to nops(gu) do
  asm:=GOGTOASM(op(i,gu)):
 lu:=lu union {asm}:
od:
lu:
end:


with(plots):
#DrawFPL(gu,a,b): draws an fpl with origin at (a,b)
#For example DrawFPL(AtF(ASM(2)[1]),0,0);
DrawFPLno:=proc(gu,a,b) local pic,g,gu1:
g:=gu[1]:
pic:=plot([[a+g[1][2],b-g[1][1]],[a+g[2][2],b-g[2][1]]],axes=none,
scaling=constrained):
gu1:=gu minus {g}:
for g in gu1 do
pic:=pic,plot([[a+g[1][2],b-g[1][1]],[a+g[2][2],b-g[2][1]]],axes=none,
scaling=constrained):
od:
display(pic):
end:




#AtF(a): inputs an asm a, outputs an FPL as a set of edges
#in [0,k+1]^2
AtF:=proc(a) local i,j,gu,k,Deg,NonEdges:

k:=nops(a):

for i from 0 to k+1 do
for j from 0 to k+1 do
 Deg[i,j]:=0:
od:
od:

gu:={}:

NonEdges:={}:

if k mod 2=1 then

for i from 1 by 2 to k do
gu:=gu union {{[i,1],[i,0]}}:
Deg[i,1]:=Deg[i,1]+1:
Deg[i,0]:=Deg[i,0]+1:
NonEdges:=NonEdges union {{[i+1,1],[i+1,0]}}:
gu:=gu union {{[i,k],[i,k+1]}}:
Deg[i,k]:=Deg[i,k]+1:
Deg[i,k+1]:=Deg[i,k+1]+1:
NonEdges:=NonEdges union {{[i+1,k],[i+1,k+1]}}:
od:

for j from 2 by 2 to k do
gu:=gu union {{ [1,j],[0,j]}}:
Deg[1,j]:=Deg[1,j]+1:
Deg[0,j]:=Deg[0,j]+1:

NonEdges:=NonEdges union {{[1,j-1],[0,j-1]}}:

gu:=gu union {{ [k,j],[k+1,j]}}:

Deg[k,j]:=Deg[k,j]+1:
Deg[k+1,j]:=Deg[k+1,j]+1:


NonEdges:=NonEdges union {{[k,j-1],[k+1,j-1]}}:
od:

NonEdges:=NonEdges union {{[k,k],[k+1,k]}}:
NonEdges:=NonEdges union {{[0,k],[1,k]}}:

else

for i from 1 by 2 to k-1 do
gu:=gu union {{[i,1],[i,0]}}:
Deg[i,1]:=Deg[i,1]+1:
Deg[i,0]:=Deg[i,0]+1:

NonEdges:=NonEdges union {{[i+1,1],[i+1,0]}}:

gu:=gu union {{[i+1,k],[i+1,k+1]}}:
Deg[i+1,k]:=Deg[i+1,k]+1:
Deg[i+1,k+1]:=Deg[i+1,k+1]+1:

NonEdges:=NonEdges union {{[i,k],[i,k+1]}}:
od:

for j from 2 by 2 to k do
gu:=gu union {{ [1,j],[0,j]}}:
Deg[1,j]:=Deg[1,j]+1:
Deg[0,j]:=Deg[0,j]+1:

NonEdges:=NonEdges union {{[1,j-1],[0,j-1]}}:

gu:=gu union {{ [k,j-1],[k+1,j-1]}}:
Deg[k,j-1]:=Deg[k,j-1]+1:
Deg[k+1,j-1]:=Deg[k+1,j-1]+1:

NonEdges:=NonEdges union {{[k,j],[k+1,j]}}:
od:

fi:




for i from 1 to k do
 for j from 1 to k do
  if (i+j mod 2=1 and a[i][j]=1) or (i+j mod 2=0 and a[i][j]=-1) then
  Deg[i,j]:=2:

    if not member({[i,j],[i-1,j]},gu) then
      Deg[i-1,j]:=Deg[i-1,j]+1:
     fi:

    if not member({[i,j],[i+1,j]},gu) then
      Deg[i+1,j]:=Deg[i+1,j]+1:
     fi:

   gu:=gu union {{[i,j],[i-1,j]},{[i,j],[i+1,j]}}:

 elif
 (i+j mod 2=0 and a[i][j]=1) or (i+j mod 2=1 and a[i][j]=-1) then



  Deg[i,j]:=2:

    if not member({[i,j-1],[i,j]},gu) then
      Deg[i,j-1]:=Deg[i,j-1]+1:
     fi:

    if not member({[i,j],[i,j+1]},gu) then
      Deg[i,j+1]:=Deg[i,j+1]+1:
     fi:

   gu:=gu union {{[i,j-1],[i,j]},{[i,j],[i,j+1]}}:

 fi:
od:
od:



for i from 1 to k do
for j from 1 to k do
 if Deg[i,j]=0 then
   gu:=gu union {{[i,j],[i+1,j]}, {[i,j],[i,j+1]}}:

    Deg[i,j]:=2:
    Deg[i+1,j]:=Deg[i+1,j]+1:
    Deg[i,j+1]:=Deg[i,j+1]+1:
 
 elif Deg[i,j]=1 then

 if member({[i,j],[i,j-1]},gu) or member({[i,j],[i,j+1]},gu)  then
   gu:=gu union {{[i,j],[i+1,j]}}:
   Deg[i,j]:=2:
   Deg[i+1,j]:=Deg[i+1,j]+1:
 
elif member({[i,j],[i-1,j]},gu) or member({[i,j],[i+1,j]},gu)  then

   gu:=gu union {{[i,j],[i,j+1]}}:

   Deg[i,j]:=2:
   Deg[i,j+1]:=Deg[i,j+1]+1:
else
  RETURN(FAIL):
 fi:

fi:
od:
od:

 
gu:

end:



#DrawFPL(gu,a,b): draws an fpl with origin at (a,b)
#For example DrawFPL(AtF(ASM(2)[1]),0,0);
DrawFPL:=proc(gu,a,b,k) local pic,g,gu1,i,kulam,mu,vu,v:
kulam:={seq(seq({[i,j],[i,j+1]},j=1..k-1),i=1..k),
seq(seq({[i,j],[i+1,j]},i=1..k-1),j=1..k)}:

vu:={}:
for g in gu do
 if {op(g[1]),op(g[2])} intersect {0,k+1}<>{} then  
   vu:=vu union {g}:
 fi:
od:

mu:=kulam minus gu:

g:=gu[1]:
if not member(g,vu) then
pic:=plot([[a+g[1][2],b-g[1][1]],[a+g[2][2],b-g[2][1]]],axes=none):
else
pic:=plot([[a+g[1][2],b-g[1][1]],[a+g[2][2],b-g[2][1]]],axes=none,color=green,linestyle=4):
fi:


gu1:=gu minus {g}:
for g in gu1 do
if not member(g,vu) then
pic:=pic,plot([[a+g[1][2],b-g[1][1]],[a+g[2][2],b-g[2][1]]],axes=none):
else
pic:=pic,plot([[a+g[1][2],b-g[1][1]],[a+g[2][2],b-g[2][1]]],axes=none,color=green, linestyle=4):
fi:
od:


for g in mu do
pic:=pic,plot([[a+g[1][2],b-g[1][1]],[a+g[2][2],b-g[2][1]]],color=blue,
linestyle=2,axes=none):
od:

pic:=pic,textplot([a-1/2,b-1,1]):

if k mod 2=1 then

for i from 2 to (k+1)/2 do
 pic:=pic,textplot([a+2*(i-1),b+1/2,i]):
od:

for i from (k+1)/2+1 to k+1 do
 pic:=pic,textplot([a+k+3/2,b-2*(i-(k+1)/2)+1,i]):
od:




for i from k+2  to k+(k+1)/2 do
 pic:=pic,textplot([a+k+1-2*(i-k-1),b-k-3/2,i]):
od:


for i from k+(k+1)/2+1 to 2*k do
 pic:=pic,textplot([a-1/2,b-k-2+2*(i-k-(k+1)/2),i]):
od:



else



for i from 2 to k/2+1 do
 pic:=pic,textplot([a+2*(i-1),b+1/2,i]):
od:

for i from k/2+2 to k+1 do
 pic:=pic,textplot([a+k+3/2,b-2*(i-(k+1)/2)+1,i]):
od:


for i from k+2  to k+k/2+1 do
 pic:=pic,textplot([a+k+1-2*(i-k-1),b-k-3/2,i]):
od:


for i from k+k/2+2 to 2*k do
 pic:=pic,textplot([a-1/2,b-k-2+2*(i-k-(k+1)/2),i]):
od:







fi:

display(pic):
end:


#DrawFPL1(gu,a,b): draws an fpl with origin at (a,b)
#For example DrawFPL(AtF(ASM(2)[1]),0,0);
DrawFPL1:=proc(gu,a,b,k) local pic,g,gu1,i,kulam,mu,vu,v:
kulam:={seq(seq({[i,j],[i,j+1]},j=1..k-1),i=1..k),
seq(seq({[i,j],[i+1,j]},i=1..k-1),j=1..k)}:

vu:={}:
for g in gu do
 if {op(g[1]),op(g[2])} intersect {0,k+1}<>{} then  
   vu:=vu union {g}:
 fi:
od:

mu:=kulam minus gu:

g:=gu[1]:
if not member(g,vu) then
pic:=plot([[a+g[1][2],b-g[1][1]],[a+g[2][2],b-g[2][1]]],axes=none,scaling=constrained):
else
pic:=plot([[a+g[1][2],b-g[1][1]],[a+g[2][2],b-g[2][1]]],axes=none,color=green,scaling=constrained,linestyle=4):
fi:


gu1:=gu minus {g}:
for g in gu1 do
if not member(g,vu) then
pic:=pic,plot([[a+g[1][2],b-g[1][1]],[a+g[2][2],b-g[2][1]]],axes=none):
else
pic:=pic,plot([[a+g[1][2],b-g[1][1]],[a+g[2][2],b-g[2][1]]],axes=none,color=green, linestyle=4):
fi:
od:


for g in mu do
pic:=pic,plot([[a+g[1][2],b-g[1][1]],[a+g[2][2],b-g[2][1]]],color=blue,
linestyle=2,axes=none):
od:

pic:=pic,textplot([a-1/2,b-1,1]):

if k mod 2=1 then

for i from 2 to (k+1)/2 do
 pic:=pic,textplot([a+2*(i-1),b+1/2,i]):
od:

for i from (k+1)/2+1 to k+1 do
 pic:=pic,textplot([a+k+3/2,b-2*(i-(k+1)/2)+1,i]):
od:




for i from k+2  to k+(k+1)/2 do
 pic:=pic,textplot([a+k+1-2*(i-k-1),b-k-3/2,i]):
od:


for i from k+(k+1)/2+1 to 2*k do
 pic:=pic,textplot([a-1/2,b-k-2+2*(i-k-(k+1)/2),i]):
od:



else



for i from 2 to k/2+1 do
 pic:=pic,textplot([a+2*(i-1),b+1/2,i]):
od:

for i from k/2+2 to k+1 do
 pic:=pic,textplot([a+k+3/2,b-2*(i-(k+1)/2)+1,i]):
od:


for i from k+2  to k+k/2+1 do
 pic:=pic,textplot([a+k+1-2*(i-k-1),b-k-3/2,i]):
od:


for i from k+k/2+2 to 2*k do
 pic:=pic,textplot([a-1/2,b-k-2+2*(i-k-(k+1)/2),i]):
od:







fi:

pic:
end:







#FPL(n): The set of fpls of size n
FPL:=proc(n) local gu,g:
gu:=ASM(n):
{seq(AtF(g),g in gu)}:
end:




#DrawFPLs(S,k,w): draws a set of FPLs S of size k
#For example, try: DrawFPLs(FPL(3),3); (w is the spacing)
DrawFPLs:=proc(S,k,w)
local pic,godel,i,j,co:

godel:=trunc(sqrt(nops(S))):


pic:=DrawFPL1(S[1],0,0,k) :

co:=2:

 for i from 1 to godel-1 while co<=nops(S) do
 pic:=pic,DrawFPL1(S[co],i*(k+w),0,k) :
  co:=co+1:
 od:

for j from 1 to godel+1 while co<=nops(S) do
 for i from 0 to godel-1 while co<=nops(S) do
 pic:=pic,DrawFPL1(S[co],i*(k+w),-j*(k+w),k) :
  co:=co+1:
 od:
od:

display(pic):

end:





#Milon(k): all the terminal vertices labelled 1,..., 2k
#in a size-k FPL, listed in order. Try: Milon(3);
Milon:=proc(k) local i,gu:

gu:=[[1,0]]:
if k mod 2 =1 then
gu:=[op(gu),seq([0,2*i],i=1..(k-1)/2)]:
gu:=[op(gu),seq([1+2*i,k+1],i=0..(k-1)/2)]:
gu:=[op(gu),seq([k+1,k-1-2*i],i=0..(k-1)/2-1)]:
gu:=[op(gu),seq([k-2*i,0],i=0..(k-1)/2-1)]:
RETURN(gu):
else

gu:=[op(gu),seq([0,2*i],i=1..k/2)]:
gu:=[op(gu),seq([2*i,k+1],i=1..k/2)]:
gu:=[op(gu),seq([k+1,k-1-2*i],i=0..k/2-1)]:
gu:=[op(gu),seq([k-1-2*i,0],i=0..k/2-2)]:
RETURN(gu):
fi:



end:



#Mate(F,k,i): Given an FPL of size k and one of the starting
#vertices i, (1<=i<=2k) finds its mate
#For example, try:
#Mate(FPL(3)[1],3,1);
Mate:=proc(F,k,i) local M,AM,V,Shakhen,sta,khaveryashan,lu,i1,
khaver,f,v:
M:=Milon(k):

for i1 from 1 to nops(M) do
 AM[M[i1]]:=i1:
od:

V:={seq(op(f),f in F)}:

for v in V do
 Shakhen[v]:={}:
od:

for f in F do
 Shakhen[f[1]]:= Shakhen[f[1]] union {f[2]}:
 Shakhen[f[2]]:= Shakhen[f[2]] union {f[1]}:
od:

sta:=M[i]:


khaver:=Shakhen[sta][1]:
khaveryashan:=sta:
while nops(Shakhen[khaver])=2 do
 lu:=Shakhen[khaver] minus {khaveryashan}:
  khaveryashan:=khaver:
  khaver:=lu[1]:
od:

AM[khaver]:
 end:


#Conn(F,k): the connectivity of the FPL F of size k
#For example, try: Conn(FPL(3)[1],3);
Conn:=proc(F,k) local i,gu,S,s,t:
option remember:
gu:={}:
S:={seq(i,i=1..2*k)}:

while S<>{} do
 s:=S[1]:
 t:=Mate(F,k,s):
 gu:=gu union {{s,t}}:
 S:=S minus {s,t}:
od:
gu:
end:


#Conns(k): all the connectivities of {1, ..., 2k}
Conns:=proc(k) local gu1,gu2,i,lu1,j,g1,g2,gu:
if k=0 then
 RETURN({{}}):
fi:

gu:={}:
for i from 2 by 2 to 2*k do
lu1:={1,i}:
gu1:=Conns((i-2)/2):
gu2:=Conns(k-i/2):

gu1:=subs({seq(j=j+1,j=1..i-2)},gu1):
gu2:=subs({seq(j=j+i,j=1..2*k-i)},gu2):

gu:= gu union { seq(seq({lu1,op(g1),op(g2)},g1 in gu1),g2 in gu2)}:

od:
gu:
end:


ConnsD:=proc(k) local gu,g:
gu:=FPL(k):
{seq(Conn(g,k),g in gu)}:
end:



#FPLcAll(k): the table of all FPLs of size k according to
#connectivity, followed by the list of connectivities
FPLcAll:=proc(k) local gu,mu,T,g,m:
option remember:
gu:=FPL(k):
mu:=Conns(k):

for m in mu do
T[m]:={}:
od:

for g in gu do
 T[Conn(g,k)]:=T[Conn(g,k)] union {g}:
od:

{seq([m,T[m]],m in mu)}:
end:



#FPLc(k,C): Given a k and a connectivity of {1, ..., 2k}
#outputs the set of all FPLs with that connectivity
#for example, try:
#FPLc(3,Conns(3)[1]);
FPLc:=proc(k,C) local gu,g:
gu:=FPLcAll(k):

for g in gu do
 if g[1]=C then
  RETURN(g[2]):
 fi:
od:
FAIL:
end:

#eipi(pi,i): applies e_i to connectivity pi
#For example, try:
#eipi({{1,2},{3,4}},2);
eipi:=proc(pi,i) local ba,k,zug,L,zug1,zug2:
k:=nops(pi):

if i<2*k then
 ba:=i+1:
elif i=2*k then
 ba:=1:
fi:


if member(ba,pi) then
 RETURN(pi):
fi:


for zug in pi while not member(i,zug) do od:
zug1:=zug:
k:=(zug1 minus {i})[1]:

for zug in pi while not member(ba,zug) do od:
zug2:=zug:
L:=(zug2 minus {ba})[1]:

(pi minus {zug1, zug2}) union {{i,ba},{k,L}}:

end:


#RSold(pi): verifies the Razumov-Strogranov Conj. for 
#connectivity pi
RSold:=proc(pi) local k,sig,gu,C,vu,v,i:
k:=nops(pi):
C:=Conns(k):
gu:={}:

for i from 1 to 2*k do
for sig in C do
 if eipi(sig,i)=pi then
   vu:=FPLc(k,sig):
   gu:=gu union {seq([i,sig,v], v in vu)}:
 fi:
od:
od:
gu:
end:

#CheckRSold(k): verifies the R-S conj. for fpl's of length k
#for example,try: CheckRSold(3);
CheckRSold:=proc(k) local gu,pi:
gu:=Conns(k):

evalb(
{true}={seq(evalb(nops(RSold(pi))=2*k*nops(FPLc(k,pi))), pi in gu) }):
end:

 









#Nei(pi): all the pairs [i,sig] such that e_i(sig)=pi
#For example, try: Nei({{1,2},{3,4}});
Nei:=proc(pi) local k,gu,cu,sig,i:
k:=nops(pi):
cu:=Conns(k):
gu:={}:
for sig in cu do
 for i from 1 to 2*k do
  if eipi(sig,i)=pi then
      gu:=gu union {[i,sig]}:
  fi:
 od:
od:
gu:

end:


#NeiNT(pi): all the non-trivial pairs [i,sig] such that e_i(sig)=pi
#with sig not pi
#For example, try: Nei({{1,2},{3,4}});
NeiNT:=proc(pi) local k,gu,cu,sig,i:
k:=nops(pi):
cu:=Conns(k):
gu:={}:
for sig in cu do
 for i from 1 to 2*k do
  if eipi(sig,i)=pi and sig<>pi then
      gu:=gu union {[i,sig]}:
  fi:
 od:
od:
gu:

end:



#NeiL(pi): all the pairs [i,sig] such that e_i(sig)=pi
#For example, try: NeiL([[1,2],[3,4]]);
NeiL:=proc(pi) local k,gu,cu,sig,i,cuL,T,pi1,sig1:
k:=nops(pi):
cu:=Conns(k):
cuL:=ConnsL(k):

for i from 1 to nops(cuL) do
T[cuL[i]]:=i:
od:

pi1:={seq(convert(pi[i1],set),i1=1..nops(pi))}:

gu:={}:
for sig in cuL do
 for i from 1 to 2*k do
  sig1:={seq(convert(sig[i1],set),i1=1..nops(sig))}:
  if eipi(sig1,i)=pi1 then
      gu:=gu union {[i,T[sig]]}:
  fi:
 od:
od:
gu:

end:


#NeiS(k): prints out all the neighbors of all connectivities
#of {1, ..., 2k}. For example, try: NeiS(2);
NeiS:=proc(k) local pi,gu:
gu:=Conns(k):

for pi in gu do
print(`The neighbors of`, pi):
print():
print(Nei(pi)):
od:

end:



RS:=proc(pi) local k,gu:
k:=nops(pi):
gu:=Nei(pi):
evalb(2*k*nops(FPLc(k,pi))=add(nops(FPLc(k,g[2])), g in gu)):
end:

CheckRS:=proc(k) local gu,pi:
gu:=Conns(k):
evalb({seq(RS(pi),pi in gu)}={true}):
end:






#ConnsL(k): all the connectivities of {1, ..., 2k}
ConnsL:=proc(k) local gu1,gu2,i,lu1,j,g1,g2,gu:
if k=0 then
 RETURN([[]]):
fi:

gu:=[]:
for i from 2 by 2 to 2*k do
lu1:=[1,i]:
gu1:=ConnsL((i-2)/2):
gu2:=ConnsL(k-i/2):

gu1:=subs([seq(j=j+1,j=1..i-2)],gu1):
gu2:=subs([seq(j=j+i,j=1..2*k-i)],gu2):

gu:= [op(gu),  seq(seq([lu1,op(g1),op(g2)],g1 in gu1),g2 in gu2)] :

od:
gu:
end:



#Sader(pi): converts pi (connectivity) to canonical form
Sader:=proc(pi) local T,p,katan,mu,i:

mu:=[]:
for p in pi do

katan:=min(op(p)):
mu:=[op(mu),katan]:
T[katan]:=(p minus {katan})[1]:
od:

mu:=sort(mu):

[seq([mu[i],T[mu[i]] ],i=1..nops(mu))]:

end:




#DrawFPLs1(S,k,koteret,koteret1,koteret2, koteret3,w): the plot-structre for FPLs S of size k
#For example, try: DrawFPLs1(FPL(3),3,Hi,Bye,4);
DrawFPLs1:=proc(S,k,koteret,koteret1,koteret2,koteret3,w)
local pic,i,j,co:



pic:=DrawFPL1(S[1],0,0,k) :

pic:=pic,textplot([w,w,koteret]):
pic:=pic,textplot([w,w+6,koteret1]):
pic:=pic,textplot([w,w+2,koteret2]):
pic:=pic,textplot([w,w+4,koteret3]):


co:=2:

 for i from 1 to 2 while co<=nops(S) do
 pic:=pic,DrawFPL1(S[co],i*(k+w),0,k) :
  co:=co+1:
 od:

for j from 1 to trunc(nops(S)/3)+1  while co<=nops(S) do
 for i from 0 to 2 while co<=nops(S) do
 pic:=pic,DrawFPL1(S[co],i*(k+w),-j*(k+w),k) :
  co:=co+1:
 od:
od:

pic:

end:

#Yafe(k,w): makes gif files for all FPL's of size k
#according to the connectivity outputs .gif files
#names oka.gif where a ranges from 1 to the  Catalan(k)
#w is the spacing between fpls
#For example, try: Yafe(3,4);

Yafe:=proc(k,w) local cu, i,gu,lu,lu1,K,j,K1,i1,mu,ku,ku1,Ku1,T,kvu,i2,
Ku1a,Ku2a:

cu:=ConnsL(k):

for i from 1 to nops(cu) do

lu:=cu[i]:

ku:=NeiL(lu):

Ku1:={}:

for ku1 in ku do
 if ku1[2]<>i then
  Ku1:=Ku1 union {ku1}:
 fi:
od:

kvu:={}:


for ku1 in Ku1 do
 kvu:=kvu union {ku1[2]}:
od:


for ku1 in kvu do
T[ku1]:={}:
od:


for ku1 in Ku1 do
 T[ku1[2]]:=T[ku1[2]] union {ku1[1]}:
od:
kvu:=sort(convert(kvu,list)):

Ku1:=[seq([T[kvu[i2]],kvu[i2]],i2=1..nops(kvu))]:

if nops(Ku1)<=4 then
 Ku1a:=Ku1:
 Ku2a:=Hi:
else
Ku1a:=[op(1..4,Ku1)]:
Ku2a:=[op(5..nops(Ku1),Ku1)]:
fi:


lu1:={seq(convert(lu[i1],set),i1=1..nops(lu))}:

gu:=FPLc(k,lu1):


K:=trunc(nops(gu)/9):

if nops(gu)=9*K then
 K1:=K:
else
 K1:=K+1:
fi:


for j from 1 to K do
plotsetup(gif,plotoutput= cat(o,k,i,`.`,j, `.gif`),plotoptions=`portrait,noborder`):

mu:=[i,[j,K1]]:
print(display(DrawFPLs1({op((j-1)*9+1..j*9,gu)},k,lu,mu,Ku1a,Ku2a,w)));
od:


if 9*K<nops(gu) then
plotsetup(gif,plotoutput= cat(o,k,i,`.`,K+1, `.gif`),plotoptions=`portrait,noborder`):
lu:=cu[i]:

mu:=[i,[K1,K1]]:
print(display(DrawFPLs1({op(K*9+1..nops(gu),gu)},k,lu,mu,Ku1a,Ku2a,w)));
fi:


od:


end:




#Hefresh(F1,F2): the set of black edegs in F1 that are not in F2
#followed by the set of black edges in F2 not in F1
#for example, try: Hefresh(FPL(3)[1],FPL(3)[2]);
Hefresh:=proc(F1,F2) :
[F1 minus F2, F2 minus F1]:
end:

#Interview(k,f,C): given an fpl f and a connectivity C, finds the
#all the the Hefresh (differences) between f and the members
#of FPLc(k,C). For example, try:
#Interview(3,FPL(3)[1],Conns(3)[1]);

Interview:=proc(k,f,C) local fu,f1,gu:
fu:=FPLc(k,C):

gu:={}:

for f1 in fu do
gu:=gu union {[f1,Hefresh(f,f1)]}:
od:
gu:
end:



#AllInterviews(e,pi): given a connectivity pi, does all the
#interviews with all the fpls of all its neighbors
#For example, do: AllInterviews(e,{{1,2},{3,4}});
AllInterviews:=proc(e,pi) local k, gu,g,Gug,i:
k:=nops(pi):
gu:=NeiNT(pi):

print(`The connectivity`, pi, `has the following non-trivial neigbors`):
print(gu):
print(`let's look, for each of them, the differences it makes`):
print(`with all the members of the FPLs with that connectivity `):

for g in gu do
 print(`for the neighbor`, g):
 print(`which means that `, e[g[1]], `has to be applied to all the fpls of`):
 print(`connectivity `, g[2]):
 print(`to get an fpl whose connectivity is`, pi):
  Gug:=FPLc(k,g[2]):
 print(`there are`, nops(Gug), `fpls `):
 
 for i from 1 to nops(Gug) do
 print(`The differences of the `, i, `member  are as follows`):

  print(Interview(k,Gug[i],pi)):
 od:

od:

end:

 
  
#DrawHefresh(k,gu,a,b): Given a pair [Black,Red] of
#edges, draws them (red also dashed).

DrawHefresh:=proc(k,gu,a,b) local pic,g,kulam,i,j:

kulam:={seq(seq({[i,j],[i,j+1]},j=1..k-1),i=1..k),
seq(seq({[i,j],[i+1,j]},i=1..k-1),j=1..k)}:


pic:=NULL:

for g in kulam do
pic:=pic,plot([[a+g[1][2],b-g[1][1]],[a+g[2][2],b-g[2][1]]],color=green,
axes=none,linestyle=2):
od:

for g in gu[1] do
pic:=pic,plot([[a+g[1][2],b-g[1][1]],[a+g[2][2],b-g[2][1]]],color=black,
axes=none):
od:


for g in gu[2] do
pic:=pic,plot([[a+g[1][2],b-g[1][1]],[a+g[2][2],b-g[2][1]]],color=red,
linestyle=2,axes=none):
od:



display(pic):
end:

#ddGif(k,F1,F2,a,b,FN): Given a pair [Black,Red] of
#edges, draws them (red also dashed) and produces a 
#gif file called FN

ddGif:=proc(k,F1,F2,a,b,FN) local pic,g,kulam,i,j,gu:

gu:=Hefresh(F1,F2) :
plotsetup(gif,plotoutput= cat(FN, `.gif`),plotoptions=`portrait,noborder`):

kulam:={seq(seq({[i,j],[i,j+1]},j=1..k-1),i=1..k),
seq(seq({[i,j],[i+1,j]},i=1..k-1),j=1..k)}:


kulam:=kulam minus (gu[1] union gu[2]):


pic:=NULL:

for g in kulam do
pic:=pic,plot([[a+g[1][2],b-g[1][1]],[a+g[2][2],b-g[2][1]]],color=green,
axes=none, linestyle=4, thickness=0, scaling=constrained):
od:

for g in gu[1] do
pic:=pic,plot([[a+g[1][2],b-g[1][1]],[a+g[2][2],b-g[2][1]]],color=red,
axes=none, scaling=constrained, thickness=6):
od:


for g in gu[2] do
pic:=pic,plot([[a+g[1][2],b-g[1][1]],[a+g[2][2],b-g[2][1]]],color=blue,
linestyle=2,axes=none,scaling=constrained, thickness=6):
od:



display(pic):
end:







#IsSimplePath(F1,F2): Given two FLPs, decides whether their
#symmetric difference is a simple path, and if yes
#returns its length. For example, try
#IsSimplePath(FPL(4)[1],FPL(4)[2]);

IsSimplePath:=proc(F1,F2) local gu,V,g,T,v,Sa,Sak,Sak1:

gu:=(F1 minus F2) union (F2 minus F1):

V:={seq(op(g),g in gu)}:


for v in V do
T[v]:=0:
Sa[v]:={}:
od:

for g in gu do
T[g[1]]:=T[g[1]]+1:
T[g[2]]:=T[g[2]]+1:

Sa[g[1]]:=Sa[g[1]] union {g[2]}:
Sa[g[2]]:=Sa[g[2]] union {g[1]}:

od:

if {seq(T[v],v in V)}<>{2} then
 RETURN(false,0):
fi:


v:=V[1]:

Sak:={v}:
Sak1:=Sak union Sa[v]:

while Sak<>Sak1 do
Sak:=Sak1:
Sak1:=Sak1 union {seq(op(Sa[v]),v in Sak1)}:
od:



if Sak1<>V then
 RETURN(false,0):
fi:

true,nops(gu)/2:

end:




#EligibleMatches(k,f,i): given an fpl f of size k, and i in [1,2k], decides
#who are the  eligible matches, followed by their lengths
#For example, try:
#EligibleMathces(3,FPL(3)[1],1);
EligibleMatches:=proc(k,f,i) local pi,sig,cands,f1,gu,vu:
sig:=Conn(f,k):
pi:=eipi(sig,i):
if pi=sig then
 RETURN({[f,0]}):
fi:

cands:=FPLc(k,pi):

gu:={}:

for f1 in cands do

vu:=IsSimplePath(f,f1):

if vu[1] then
 gu:=gu union {[f1,vu[2]]}:
fi:

od:

gu:

end:








#ShortestEligibleMatches(k,f,i): given an fpl f of size k, and i in [1,2k], finds the
#shortest  eligible match to be e_i(f)
#For example, try:
#ShortestEligibleMatches(3,FPL(3)[1],1);
ShortestEligibleMatches:=proc(k,f,i) local gu,shi,mu,g:
gu:=EligibleMatches(k,f,i):

if gu={} then 
RETURN(FAIL):
fi:

shi:=min(seq(g[2], g in gu)):

mu:={}:

for g in gu do
 if g[2]=shi then
  mu:=mu union {g[1]}:
 fi:
od:
mu:


end:




#DrawFPL1t(gu,a,b,k): draws an fpl with origin at (a,b) using thickness
#For example DrawFPL1t(FPL(3)[1],0,0,3);
DrawFPL1t:=proc(gu,a,b,k) local pic,g,gu1,i,kulam,mu,vu,v,NBE:
kulam:={seq(seq({[i,j],[i,j+1]},j=1..k-1),i=1..k),
seq(seq({[i,j],[i+1,j]},i=1..k-1),j=1..k)}:

#top border
NBE:={seq({[1,2*i-1],[0,2*i-1]},i=1..trunc((k+1)/2) )}:

#left border
NBE:=NBE union {seq({[2*i,1],[2*i,0]},i=1..trunc(k/2) )}:

#right border
if k mod 2= 0 then
NBE:=NBE union {seq({[2*i-1,k],[2*i-1,k+1]},i=1..k/2 )}:
 else
NBE:=NBE union {seq({[2*i,k],[2*i,k+1]},i=1..(k-1)/2 )}:
fi:


#bottom border
if k mod 2= 0 then
NBE:=NBE union {seq({[k,2*i],[k+1,2*i]},i=1..k/2 )}:
 else
NBE:=NBE union {seq({[k,2*i-1],[k+1,2*i-1]},i=1..(k+1)/2 )}:
fi:


vu:={}:
for g in gu do
 if {op(g[1]),op(g[2])} intersect {0,k+1}<>{} then  
   vu:=vu union {g}:
 fi:
od:

mu:=kulam minus gu :
mu:=mu union NBE:

pic:=NULL:

for g in gu do
 pic:=pic,plot([[a+g[1][2],b-g[1][1]],[a+g[2][2],b-g[2][1]]],axes=none,scaling=constrained, color=red,
thickness=4):
od:


for g in mu do
 pic:=pic,plot([[a+g[1][2],b-g[1][1]],[a+g[2][2],b-g[2][1]]],axes=none,scaling=constrained, color=blue,
thickness=4):

od:

pic:=pic,textplot([a-1/2,b-1,1]):

if k mod 2=1 then

for i from 2 to (k+1)/2 do
 pic:=pic,textplot([a+2*(i-1),b+1/2,i]):
od:

for i from (k+1)/2+1 to k+1 do
 pic:=pic,textplot([a+k+3/2,b-2*(i-(k+1)/2)+1,i]):
od:




for i from k+2  to k+(k+1)/2 do
 pic:=pic,textplot([a+k+1-2*(i-k-1),b-k-3/2,i]):
od:


for i from k+(k+1)/2+1 to 2*k do
 pic:=pic,textplot([a-1/2,b-k-2+2*(i-k-(k+1)/2),i]):
od:



else



for i from 2 to k/2+1 do
 pic:=pic,textplot([a+2*(i-1),b+1/2,i]):
od:

for i from k/2+2 to k+1 do
 pic:=pic,textplot([a+k+3/2,b-2*(i-(k+1)/2)+1,i]):
od:


for i from k+2  to k+k/2+1 do
 pic:=pic,textplot([a+k+1-2*(i-k-1),b-k-3/2,i]):
od:


for i from k+k/2+2 to 2*k do
 pic:=pic,textplot([a-1/2,b-k-2+2*(i-k-(k+1)/2),i]):
od:

fi:

pic:
end:





#DrawFPL1tGif(gu,a,b,k,FN): draws an fpl with origin at (a,b) using thickness
#into a file FN.gif
#For example DrawFPL1tGif(FPL(3)[1],0,0,3,out2);
DrawFPL1tGif:=proc(gu,a,b,k,FN) local pic:
plotsetup(gif,plotoutput= cat(FN, `.gif`),plotoptions=`portrait,noborder`):

pic:=DrawFPL1t(gu,a,b,k):

display(pic):

end:




#ClosestEligibleMatches(k,f,i): 
#given an fpl f of size k, and i in [1,2k], finds the
#shortest  eligible match to be e_i(f)
#For example, try:
#ClosestEligibleMatches(3,FPL(3)[1],1);
ClosestEligibleMatches:=proc(k,f,i) local gu,shi,mu,g:
gu:=EligibleMatches(k,f,i):

if gu={} then 
RETURN(FAIL):
fi:

shi:=min(seq(g[2], g in gu)):

mu:={}:

for g in gu do
 if g[2]=shi then
  mu:=mu union {g[1]}:
 fi:
od:
mu:


end:



#PathFrom(k,P1,f): Given a point P1 in [0,k+1]^2
#with degree 1, and an fpl f,
#finds the path that starts with P1.
#For example, try: 
#PathFrom(3,[1,0],FPL(3)[1]);
PathFrom:=proc(k,P1,f) local T,edg, i,j,lu,gu,luP,luN:
option remember:


for i from 0 to k+1 do
 for j from 0 to k+1 do
  T[[i,j]]:={}:
 od:
od:

for edg in f do
 T[edg[1]]:=T[edg[1]] union {edg[2]}:
 T[edg[2]]:=T[edg[2]] union {edg[1]}:
od:

if nops(T[P1])<>1 then
 RETURN(FAIL):
fi:


gu:=[P1]:

luP:=P1:
lu:=T[P1][1]:

while nops(T[lu])=2 do
 gu:=[op(gu),lu]:
 luN:=(T[lu] minus {luP})[1]:
 luP:=lu:
 lu:=luN:
od:
[op(gu),lu]:
end:



#WhereIsi(i,k): where is point labelled i in [0,k+1]^2
WhereIsi:=proc(i,k)

if k mod 2=0 then

if i=1 then
 RETURN([1,0]):
fi:

if 1<i and i<=k/2+1 then
 RETURN([0,2*(i-1)]):
elif i>k/2+1 and i<=k+1 then
  RETURN([2*(i-k/2-1),k+1]):
elif  i>=k+2 and i<=3*k/2+1 then
   RETURN([k+1,k-1-2*(i-k-2)]):
else
  RETURN([k-1-2*(i-(3*k/2+2)),0]):
fi:

else


if i=1 then
 RETURN([1,0]):
fi:

if 1<i and i<=(k+1)/2 then
 RETURN([0,2*(i-1)]):
elif i>(k+1)/2 and i<=k+1 then
  RETURN([1+2*(i-(k+1)/2),k+1]):
elif  i>=k+2 and i<=(3*k+1)/2 then
   RETURN([k+1,k-1-2*(i-k-2)]):
else
  RETURN([k-2*(i-(3*k+1)/2-1),0]):
fi:

fi:

end: