#######################################################################
## Bearoff. Save this file as Bearoff   To use it, stay in the        #
## same directory, get into Maple (by typing: maple <Enter> )         #
## and then type:  read Bearoff:   <Enter>                            #
## Then follow the instructions given there                           #
##                                                                    #
## Written by Doron Zeilberger, Rutgers University ,                  #
##  zeilberg@math.rutgers.edu.                                        # 
#######################################################################


print(`First Written: Sept. 21, 2005: tested for Maple 10 `):
print(`Version of Sept. 21, 2005:  `):
print():
print(`This is Bearoff a Maple package`):
print(`accompanying the article `):
print(`"How to  Play Backgammon (if you must) and how to Research it",`):
print(` (if you have time)" `):
print(`By Shalosh B. Ekhad and Doron Zeilberger.`):
print():
print(`The most current version of this packageand article`):
print(`are available on WWW at:`):
print(` http://www.math.rutgers.edu/~zeilberg .`):
print(`Please report all bugs to: zeilberg at math dot rutgers dot edu .`):
print():
print(`For general help, and a list of the available functions,`):
print(` type "ezra();". For specific help type "ezra(procedure_name);" `):
print():
 
 
 
ezra:=proc()
if args=NULL then
print(`Bearoff`):
print(`A Maple package for  Winning Backgammon Horse-Race`):
print():
print(`For help with a specific procedure, type "ezra(procedure_name);"`):
print(`Contains procedures: `):
print():
print(`Ek, Ek2, FindP, G, GP, GP2, Moves , Moves2, P, P2 `):
print(`POSki, POSkis, Tkis , Tkis2,TkisE, TkisE2, Box `):
print(`Test1, Exc1, Exc2, Sipur, Sipur2, SipurC`):

 elif nargs=1 and args[1]=Box  then
print(`Box(m): all the vertices with the given range`):
print(`e.g. Box([[0,4],[1,5]]);`):

 elif nargs=1 and args[1]=Ek  then
print(`Ek(Li): Given a list of length k, finds the expected`):
print(`number of moves in optimal play with a k-die in`):
print(`backgammon with ONE die`):

 elif nargs=1 and args[1]=Ek2  then
print(`Ek2(Li): Given a list of length k, finds the expected`):
print(`number of moves in optimal play with a k-die `):
print(`backgammon with two dice, played with the usual convention`):
print(` that a double means four identical moves`):
print(`For example, try Ek2([1,1]);`):

 elif nargs=1 and args[1]=Exc1  then
print(`Exc1(r,K,i): tests the hypothesis that in`):
print(`r-die backgammon end-game, if the`):
print(`roll is i, it is always`):
print(`best to remove a checker from the i^th place`):
print(`if it is occupied for all positions with  K checkers`):
print(`outputs all the exceptions to the rule`):

 elif nargs=1 and args[1]=Exc2  then
print(`Exc2(r,K,i,j): tests the hypothesis that in`):
print(`r-die backgammon end-game (with TWO DICE), if the`):
print(`roll is {i,j}, it is always`):
print(`best to remove a checker from the i^th place and a checker from the jth place`):
print(`if they are both occupied for all positions with  K checkers`):
print(`outputs all the exceptions to the rule`):

 elif nargs=1 and args[1]=FindP then
print(`FindP(S,a,r): Given a list of exceptions and one exceptional move of the `):
print(`highest number of checkes, finds its parents (if any) at level r`):
print(`For example, try: gu:=Sipur(4,12): FindP(gu,gu[nops(gu)][1],nops(gu)-1)`):


 elif nargs=1 and args[1]=G  then
print(`G(Li,i): the expected duration of r(:=nops(Li))-faced `):
print(` die backgammon`):
print(`solitaire with position Li and die-roll i`):
print(`with optimal play followed by the optimal move(s)`):
print(`For example, try: G([2,3],1);`):

 elif nargs=1 and args[1]=G2  then
print(`G2(Li,S): the expected duration of r(:=nops(Li))-faced `):
print(` die backgammon with two dice (the usual way)`):
print(`solitaire with position Li and die-roll S, (S is a set of one or two elementes`):
print(` representing the roll-die `):
print(`with optimal play followed by the optimal move(s)`):
print(`For example, try: G2([2,3],{1,2});`):

 elif nargs=1 and args[1]=GP  then
print(`GP(Li1,Li2,i): the probability of Player one winning`):
print(`with position Li1, who just rolled i in k-faced die`):
print(`backgammon end-game`):
print(`with optimal play followed by the optimal move(s)`):

 elif nargs=1 and args[1]=GP2  then
print(`GP2(Li1,Li2,S): the probability of Player one winning`):
print(`with position Li1, who just rolled S in k-faced die`):
print(`backgammon end-game with TWO dice `):
print(`with optimal play followed by the optimal move(s)`):

 elif nargs=1 and args[1]=Moves  then
print(`Moves(Li,i): all the legal moves from Li`):
print(`with die-roll i (k=nops(Li) (in a k-die backgammon with ONE die)`):
print(`For example, try: Noves([2,3],2);`):

 elif nargs=1 and args[1]=Moves2  then
print(`Moves2(Li,S), Given a position Li, and a set, S, `):
print(` of one or two rolls finds all the possible legal moves`):
print(`where if S has one element it means  a double-roll`):
print(`For example, try:Moves2([1,1,1,1,1,1],{1});`):

 elif nargs=1 and args[1]=P  then
print(`P(Li1,Li2): Given two lists of length k, say,`):
print(`representing the positions `):
print(`of the two players in a k-faced die backgammon end-game`):
print(` (with ONE die) finds the probability of the first `):
print(` (Li1) player winning `):

 elif nargs=1 and args[1]=P2  then
print(`P2(Li1,Li2): Given two lists of length k, say,`):
print(`representing the positions `):
print(`of the two players in a k-faced die backgammon end-game`):
print(` (with TWO dice) finds the probability of the first `):
print(` (Li1) player winning `):

 elif nargs=1 and args[1]=POSki  then
print(`POSki(k,i): all positions in one-die backgammon end-game solitaire`):
print(`with i counters `):

 elif nargs=1 and args[1]=POSkis  then
print(`POSkis(k,i,s): all positions in one-die backgammon end-game solitaire`):
print(`with i counters whose pip-count is s`):
print(`For example, try: POSkis(4,3,6);`):

 elif nargs=1 and args[1]=Sipur  then
print(`Sipur(r,K): all the excepion in an r-faced Backgammon end-game`):
print(`up to K checkers`):
print(`For example, try: Sipur(6,11);`):

 elif nargs=1 and args[1]=Sipur2  then
print(`Sipur2(r,K): all the excepion in an r-faced Backgammon end-game`):
print(`up to K checkers with TWO DICE`):
print(`For example, try: Sipur2(6,11);`):

 elif nargs=1 and args[1]=SipurC  then
print(`SipurC(r,K,x,s,t): all the excepion in an r-faced Backgammon end-game`):
print(`up to K checkers, compact and silent in terms of a generating`):
print(` function where a term of the form s^i*t^j*x[1]^a1*...*x[r]^ar`):
print(` indicates that the position is [a1, ..., ar] and the optimal move is`):
print(`to move from position i to position j`):
print(`For example, try: SipurC(4,12,x,s,t);`):

 elif nargs=1 and args[1]=Test1  then
print(` Test1(r,K,i): tests the hypothesis that in`):
print(`r-die backgammon end-game, if the`):
print(`roll is i, it is always`):
print(`best to remove a checker from the i^th place`):
print(`if it is occopied for all others with <=K checkers`):
print(`outputs all the exceptions to the rule`):

 elif nargs=1 and args[1]=Tkis  then
print(`Tkis(k,i,s): a sorted table of expected optimal-play life for `):
print(`all the positions in k-die backgammon`):
print(`with ONE die`):
print(`end game with number of counters i and pip-count  s`):
print(`For example, try: Tkis(6,6,10);`)

 elif nargs=1 and args[1]=TkisE  then
print(`TkisE(k,i,s): a sorted table of expected optimal-play life for `):
print(`all the positions in k-die backgammon`):
print(`with ONE die`):
print(`end game with number of counters i and pip-count  s`):
print(`and also listing the corresponding number of moves`):
print(`For example, try: TkisE(6,6,10);`)


 elif nargs=1 and args[1]=Tkis2  then
print(`Tkis2(k,i,s): a sorted table of expected optimal-play life for `):
print(`all the positions in k-die backgammon`):
print(`with TWO dice`):
print(`end game with number of counters i and pip-count  s`):
print(`For example, try: Tkis2(6,6,10);`)

 elif nargs=1 and args[1]=TkisE2  then
print(`TkisE(k,i,s): a sorted table of expected optimal-play life for `):
print(`all the positions in k-die backgammon`):
print(`with TWO dice`):
print(`end game with number of counters i and pip-count  s`):
print(`and also listing the corresponding number of moves`):
print(`For example, try: TkisE2(6,6,10);`)

fi:
 
end:





#Ek(Li): Given a list of length k, finds the expected
#number of moves in optimal play with a k-die in
#backgammon with ONE die
Ek:=proc(Li) local k,i:
option remember:
k:=nops(Li):

if Li=[0$k] then
 RETURN(0):
fi:

add(G(Li,i)[1],i=1..k)/k:
end:



#Moves(Li,i): all the legal moves from Li
#with die-roll i (k=nops(Li) (in a k-die backgammon with ONE die)
#For example, try: Noves([2,3],2);
Moves:=proc(Li,i) local k,i1,gu:
k:=nops(Li):

if i>k or i< 1 then
ERROR(`Second argument must be beteen 1 and `, k):
fi:

gu:={}:

if {op(i..k,Li)}={0} then

for i1 from i-1 by -1 to 1 do
 if Li[i1]>0 then
  RETURN({[op(1..i1-1,Li),Li[i1]-1,op(i1+1..k,Li)]}):
 fi

od:

 RETURN({}):

fi:

if Li[i]>0 then
 gu:=gu union {[op(1..i-1,Li),Li[i]-1,op(i+1..k,Li)]}:
fi:

for i1 from i+1 to k do

 if Li[i1]>0 then
  gu:= gu union 
{[op(1..i1-i-1,Li),Li[i1-i]+1,op(i1-i+1..i1-1,Li),Li[i1]-1,op(i1+1..k,Li)]}:
 fi:

od:

RETURN(gu):
end:

#G(Li,i): the expected duration of the game
#with optimal play followed by the optimal move(s)
G:=proc(Li,i)
local gu,alufim,shi,i1,muam,nase:
gu:=Moves(Li,i):
if gu={} then
 RETURN(0,{}):
fi:
alufim:={gu[1]}:
shi:=Ek(gu[1]):

for i1 from 2 to nops(gu) do
muam:=gu[i1]:
nase:=Ek(muam):
 if nase<shi then
  shi:=nase:
  alufim:={muam}:
 elif nase=shi then
  alufim:=alufim union {muam}:
 fi:
od:
1+shi,alufim:
end:

#POSki(k,i): all positions in one-die backgammon end-game solitaire
#with i counters 
POSki:=proc(k,i)
local gu,mu,j,i1:
option remember:

if k=1 then
 RETURN({[i]}):
fi:

gu:={}:

for j from 0 to i do
mu:=POSki(k-1,i-j):
gu:=gu union {seq([op(mu[i1]),j],i1=1..nops(mu))}:
od:
gu:
end:


#POSkis(k,i,s): all positions in one-die backgammon end-game solitaire
#with i counters whose pip-count is s
#For example, try: POSkis(4,3,6);
POSkis:=proc(k,i,s)
local gu,mu,j,i1:
option remember:

if k=1 then
 if s=i then
  RETURN({[i]}):
 else
  RETURN({}):
 fi:
fi:

gu:={}:

for j from 0 to i do
mu:=POSkis(k-1,i-j,s-k*j):
gu:=gu union {seq([op(mu[i1]),j],i1=1..nops(mu))}:
od:
gu:
end:

#Tkis(k,i,s): a sorted table of expected optimal-play life for 
#all the positions in k-die backgammon
#end game with number of counters i and pip-count s
#For example, try: Tkis(6,6,10);
Tkis:=proc(k,i,s) local lu,gu1,j,T,mu,mu1:
lu:=POSkis(k,i,s):
gu1:={seq([lu[j],Ek(lu[j])],j=1..nops(lu))}:

mu:=[seq(gu1[j][2],j=1..nops(gu1))]:


for j from 1 to nops(mu) do
T[mu[j]]:={}:
od:

for j from 1 to nops(gu1) do
T[gu1[j][2]]:=T[gu1[j][2]] union {gu1[j][1]}:
od:
mu1:=sort(mu):

[seq(op(T[mu1[j]]),j=1..nops(mu1))]:
end:


#TkisE(k,i,s): a sorted table of expected optimal-play life for 
#all the positions in k-die backgammon
#end game with number of counters i and pip-count s
#For example, try: TkisE(6,6,10), together with the expectations;
TkisE:=proc(k,i,s) local lu,gu1,j,T,mu,mu1:
lu:=POSkis(k,i,s):
gu1:={seq([lu[j],Ek(lu[j])],j=1..nops(lu))}:

mu:=[seq(gu1[j][2],j=1..nops(gu1))]:


for j from 1 to nops(mu) do
T[mu[j]]:={}:
od:

for j from 1 to nops(gu1) do
T[gu1[j][2]]:=T[gu1[j][2]] union {gu1[j][1]}:
od:
mu1:=sort(mu):

mu1:=[seq(op(T[mu1[j]]),j=1..nops(mu1))]:
[seq([mu1[j],evalf(Ek(mu1[j]))] ,j=1..nops(mu1))]:
end:

#Box(m): all the vertices with the given range
#e.g. Box([[0,4],[1,5]]);
Box:=proc(m) local lu,r,m1,i,j,gu:
r:=nops(m):

 if r=0 then
  RETURN({[]}):
 fi:

m1:=[op(1..r-1,m)]:

lu:=Box(m1):
gu:={}:

for i from 1 to nops(lu) do
 for j from m[r][1] to m[r][2] do
 gu:=gu union {[op(lu[i]),j]}:
od:
od:
gu:
end:

#Test1(r,K,i): tests the hypothesis that in
#r-die backgammon end-game, if the
#roll is i, it is always
#best to remove a checker from the i^th place
#if it is occopied for all others with <=K checkers
#outputs all the exceptions to the rule
Test1:=proc(r,K,i)
local gu,i1,lu:
if not 1<=i and i<=r then
 ERROR(`Out of range`):
fi:
gu:=Box([[0,K]$(i-1),[1,K],[0,K]$(r-i)]):

lu:={}:
for i1 from 1 to nops(gu) do

if gu[i1][i]-G(gu[i1],i)[2][1][i]<>1 then
lu:= lu union {gu[i1]}:
 print(`For position `, gu[i1], `and roll`, i , `the best moves are`):
 print(G(gu[i1],i)[2]):
fi:
od:
lu:
end:

#Exc1(r,K,i): tests the hypothesis that in
#r-die backgammon end-game, if the
#roll is i, it is always
#best to remove a checker from the i^th place
#if it is occupied for all positions with  K checkers
#outputs all the exceptions to the rule
Exc1:=proc(r,K,i)
local gu,i1,lu:
if not 1<=i and i<=r then
 ERROR(`Out of range`):
fi:
gu:=POSki(r,K):

lu:={}:
for i1 from 1 to nops(gu) do

if gu[i1][i]>0 then
if gu[i1][i]-G(gu[i1],i)[2][1][i]<>1 then
lu:= lu union {gu[i1]}:
 print(`For position `, gu[i1], `and roll`, i , `the best moves are`):
 print(G(gu[i1],i)[2]):
fi:
fi:
od:
lu:
end:

#Exc1Silent(r,K,i): tests the hypothesis that in
#r-die backgammon end-game, if the
#roll is i, it is always
#best to remove a checker from the i^th place
#if it is occupied for all positions with  K checkers
#outputs all the exceptions to the rule
Exc1Silent:=proc(r,K,i)
local gu,i1,lu:
if not 1<=i and i<=r then
 ERROR(`Out of range`):
fi:
gu:=POSki(r,K):

lu:={}:
for i1 from 1 to nops(gu) do

if gu[i1][i]>0 then
if gu[i1][i]-G(gu[i1],i)[2][1][i]<>1 then
lu:= lu union {gu[i1]}:
fi:
fi:
od:
lu:
end:


#Sipur(r,K): all the excepion in an r-faced Backgammon end-game
#up to K checkers
Sipur:=proc(r,K) local K1,i,lu,LU,i1,LU1:

LU:=[]:
for K1 from 1 to K do
LU1:=[]:
 for i from 1 to r do
  lu:=Exc1(r,K1,i):
   if lu<>{} then

     print(`When the roll was`,i, `and there were`, K1, `checkers , the exceptions to greedy play were`):
      lu:={seq([i,lu[i1],G(lu[i1],i)[2][1]],i1=1..nops(lu))}:
       print(lu):
     LU1:=[op(LU1),op(lu)]:
   fi:
 od:
if LU1<>[] then
 LU:=[op(LU),LU1]:
fi:
od:
LU:
end:






#P(Li1,Li2): Given two lists of length k, say,
#representing the positions 
#of the two players in a k-faced die backgammon end-game
#(with ONE die) finds the probability of the first
#(Li1) player winning
P:=proc(Li1,Li2) local k,i:
option remember:
k:=nops(Li1):
 if nops(Li2)<>k then
   ERROR(`Bad input`):
 fi:

if Li2=[0$k] then
 RETURN(0):
fi:

add(GP(Li1,Li2,i)[1],i=1..k)/k:
end:

#GP(Li1,Li2,i): the probability of Player one winning
#with position Li1, who just rolled i in k-faced die
#backgammon end-game
#with optimal play followed by the optimal move(s)
GP:=proc(Li1,Li2,i)
local gu,alufim,shi,i1,muam,nase:
gu:=Moves(Li1,i):
if gu={} then
 RETURN(0,{}):
fi:
alufim:={gu[1]}:
shi:=P(Li2,gu[1]):

for i1 from 2 to nops(gu) do
muam:=gu[i1]:
nase:=P(Li2,muam):
 if nase<shi then
  shi:=nase:
  alufim:={muam}:
 elif nase=shi then
  alufim:=alufim union {muam}:
 fi:
od:
1-shi,alufim:
end:





#FindP(S,a,r): Given a list of exceptions and one exceptional move of the 
#highest number of checkes, finds its parents (if any) at S[r]
FindP:=proc(S,a,r)
local lu,k,P1,M2,P2,i,j,x:
lu:={}:
k:=a[1]:
P1:=a[2]:
for i from 1 to nops(S[r]) do
 M2:=S[r][i]:

 if k=M2[1]  then
  P2:=M2[2]:
   if coeff(add(x[evalb(P1[j]-P2[j]=0)],j=1..nops(P1)),x[true],1)=nops(P1)-1 then
     lu:=lu union {M2}:
   fi:
 fi:
od:
lu:
end:


ez:=proc():print(`Cone(V,S,k), SipurC(r,K,x,t) `):end:

#Cone(V,S,k): all the vectors starting at V active along S with extra sum k
#For example try: Cone([2,3],{1},4);
Cone:=proc(V,S,k) local mu,i,j,x,t:
mu:=convert([seq(x[i]^V[i],i=1..nops(V))],`*`):
mu:=coeff(expand(mu*convert([seq(add(x[S[j]]^i*t^i,i=0..k),j=1..nops(S))],`*`)),t,k):

if not type(mu,`+`) then
 ERROR(`k too small`):
fi:

{seq([seq(degree(op(i,mu),x[j]),j=1..nops(V))],i=1..nops(mu))}:

end:


#SipurC(r,K,x,s,t): all the excepion in an r-faced Backgammon end-game
#up to K checkers, compact and silent in terms of a generating function
#where a term of the form s^i*t^j*x[1]^a1*...*x[r]^ar indicates
#that the position is [a1, ..., ar] and the optimal move is
#to move from position i to position j
SipurC:=proc(r,K,x,s,t) local K1,i,lu,LU,i1,LU1,j,ku:

LU:=[]:
for K1 from 1 to K do
LU1:=[]:
 for i from 1 to r do
  lu:=Exc1Silent(r,K1,i):
   if lu<>{} then

      lu:={seq([i,lu[i1],G(lu[i1],i)[2][1]],i1=1..nops(lu))}:
     LU1:=[op(LU1),op(lu)]:
   fi:
 od:
if LU1<>[] then
 LU:=[op(LU),LU1]:
fi:
od:


LU1:=[]:
for i from 1 to nops(LU) do

lu:=LU[i]:
  ku:=0:
 for j from 1 to nops(lu) do
  ku:= ku +GF1(lu[j],x,t,s):
 od:
 LU1:=[op(LU1),ku]:
od:
LU1:

end:



GF1:=proc(V,x,t,s) local k,vu,j,vu1,vu2:
vu:={}:
for j from 1 to nops(V[2]) do
if V[2][j]<>V[3][j] then
   vu:=vu union {j}:
fi:
od:

if nops(vu)<>2 then
 print(`vu is`,vu):
 ERROR(`Something is wrong`):
fi:

vu1:=min(op(vu)):vu2:=max(op(vu)):

 if vu2-vu1<>V[1] then
  ERROR(`Something is wrong`):
 fi:

s^vu2*t^V[1]*mul(x[k]^V[2][k],k=1..nops(V[2])):

end:


#IsGoodMono(LU,x,t,s,M1,i1,k):Given a list of sets of Exceptionals Polynomials
#and a monomial M1 in LU[i1], finds out whether all descendats of M1 are
#exceptionals (in k-faced one die Backgammon end-game), also returns the left-overs
IsGoodMono:=proc(LU,x,t,s,M1,i1,k) local z,vu,i,j,gu,LU1:
i:=degree(M1,t):j:=degree(M1,s):
vu:={j,j-i}:

gu:=M1:

for i from 1 to k do
 if not member(i,vu) then
   gu:=gu/(1-z*x[i]):
 fi:
od:
gu:=taylor(gu,z=0,nops(LU)-i1+1):

LU1:=LU:
for j from 1 to nops(LU)-i1 do
if {op(expand(coeff(gu,z,j)))} minus {op(LU[i1+j])}<>{} then
  RETURN(false):
else
 LU1:=[op(1..i1+j-1,LU1), {op(LU[i1+j])} minus {op(expand(coeff(gu,z,j)))},op(i1+j+1..nops(LU1),LU1 )]:
fi:
od:
true,LU1:

end:


###Starts Real-Life Backgammon

#Moves2 Given a set of one or two rolls finds all the possible legal moves
Moves2:=proc(Li,S) local i,j,gu,i1,i2,gu1,k:
k:=nops(Li):
if nops(S)=2 then
i:=max(op(S)):
j:=min(op(S)):
   gu1:=Moves(Li,i):
 if member([0$k],gu1) then
  RETURN({[0$k]}):
 fi:
  gu:= {seq(op(Moves(gu1[i1],j)),i1=1..nops(gu1))}:
 if member([0$k],gu) then
  RETURN({[0$k]}):
 fi:

  RETURN(gu):
fi:

if nops(S)=1 then
    i:=S[1]:
  gu:={Li}:
 for i1 from 1 to 4 do
  gu:= {seq(op(Moves(gu[i2],i)),i2=1..nops(gu))}:
if member([0$k],gu) then
  RETURN({[0$k]}):
 fi:

 od:
RETURN(gu):
fi:
print(`The second argument should have one or two elements`):
 FAIL:
end:


#G2(Li,S): the expected duration of the real game
#with optimal play followed by the optimal move(s)
G2:=proc(Li,S)
local gu,alufim,shi,i1,muam,nase:

if Li=[0$nops(Li)] then
 RETURN(0,{}):
fi:
gu:=Moves2(Li,S):
if gu={} then
 RETURN(1,{}):
fi:
alufim:={gu[1]}:
shi:=Ek2(gu[1]):

for i1 from 2 to nops(gu) do
muam:=gu[i1]:
nase:=Ek2(muam):
 if nase<shi then
  shi:=nase:
  alufim:={muam}:
 elif nase=shi then
  alufim:=alufim union {muam}:
 fi:
od:
1+shi,alufim:
end:

#Ek2(Li): Given a list of length k, finds the expected
#number of moves in optimal play with a k-die in
#backgammon with ONE die
Ek2:=proc(Li) local k,i,j:
option remember:
k:=nops(Li):

if Li=[0$k] then
 RETURN(0):
fi:

2*convert([seq(seq(G2(Li,{i,j})[1],i=j+1..k),j=1..k)],`+`)/k^2
+add(G2(Li,{i})[1],i=1..k)/k^2:
end:

#P2(Li1,Li2): Given two lists of length k, say,
#representing the positions 
#of the two players in a k-faced die backgammon end-game
#(with ONE die) finds the probability of the first
#(Li1) player winning
P2:=proc(Li1,Li2) local k,i,j:
option remember:
k:=nops(Li1):
 if nops(Li2)<>k then
   ERROR(`Bad input`):
 fi:

if Li2=[0$k] then
 RETURN(0):
fi:

2*add(add(GP2(Li1,Li2,{i,j})[1],j=i+1..k),i=1..k)/k^2+
add(GP2(Li1,Li2,{i})[1],i=1..k)/k^2:
end:

#GP2(Li1,Li2,S): the probability of Player one winning
#with position Li1, who just rolled S in k-faced die
#backgammon end-game with two dice
#with optimal play followed by the optimal move(s)
GP2:=proc(Li1,Li2,S)
local gu,alufim,shi,i1,muam,nase:
gu:=Moves2(Li1,S):
if gu={} then
 RETURN(0,{}):
fi:
alufim:={gu[1]}:
shi:=P2(Li2,gu[1]):

for i1 from 2 to nops(gu) do
muam:=gu[i1]:
nase:=P2(Li2,muam):
 if nase<shi then
  shi:=nase:
  alufim:={muam}:
 elif nase=shi then
  alufim:=alufim union {muam}:
 fi:
od:
1-shi,alufim:
end:


#Tkis2(k,i,s): a sorted table of expected optimal-play life for 
#all the positions in k-die backgammon
#with TWO DICE
#end game with number of counters i and pip-count s
#For example, try: Tkis(6,6,10);
Tkis2:=proc(k,i,s) local lu,gu1,j,T,mu,mu1:
lu:=POSkis(k,i,s):
gu1:={seq([lu[j],Ek2(lu[j])],j=1..nops(lu))}:

mu:=[seq(gu1[j][2],j=1..nops(gu1))]:


for j from 1 to nops(mu) do
T[mu[j]]:={}:
od:

for j from 1 to nops(gu1) do
T[gu1[j][2]]:=T[gu1[j][2]] union {gu1[j][1]}:
od:
mu1:=sort(mu):

[seq(op(T[mu1[j]]),j=1..nops(mu1))]:
end:


#TkisE2(k,i,s): a sorted table of expected optimal-play life for 
#all the positions in k-die backgammon
#with TWO DICE
#end game with number of counters i and pip-count s
#For example, try: TkisE(6,6,10), together with the expectations;
TkisE2:=proc(k,i,s) local lu,gu1,j,T,mu,mu1:
lu:=POSkis(k,i,s):
gu1:={seq([lu[j],Ek2(lu[j])],j=1..nops(lu))}:

mu:=[seq(gu1[j][2],j=1..nops(gu1))]:


for j from 1 to nops(mu) do
T[mu[j]]:={}:
od:

for j from 1 to nops(gu1) do
T[gu1[j][2]]:=T[gu1[j][2]] union {gu1[j][1]}:
od:
mu1:=sort(mu):

mu1:=[seq(op(T[mu1[j]]),j=1..nops(mu1))]:
[seq([mu1[j],evalf(Ek2(mu1[j]))] ,j=1..nops(mu1))]:
end:




#Exc2(r,K,i,j): tests the hypothesis that in
#r-die backgammon end-game, with TWO dice, if the
#roll is {i,j}, it is always
#best to remove a checker from the i^th place
#if it is occupied for all positions with  K checkers
#outputs all the exceptions to the rule
Exc2:=proc(r,K,i,j)
local gu,i1,lu,ObvMove:
if i=j then
 print(`Two last argumentes must be different`):
  RETURN(FAIL):
fi:
if i>j then
print(` Third argument muts be less than Fourth argument`):
  RETURN(FAIL):
fi:

if not 1<=i and i<=r then
 ERROR(`Out of range`):
fi:
gu:=POSki(r,K):

lu:={}:
for i1 from 1 to nops(gu) do

if gu[i1][i]>0 and gu[i1][j]>0 then
ObvMove:=
[op(1..i-1,gu[i1]),gu[i1][i]-1,
op(i+1..j-1,gu[i1]),gu[i1][j]-1,
op(j+1..nops(gu[i1]),gu[i1]) ]:

if not member(ObvMove, G2(gu[i1],{i,j})[2]) then
lu:= lu union {gu[i1]}:
 print(`For position `, gu[i1], `and roll`, {i,j} , `the best moves are`):
 print(G2(gu[i1],{i,j})[2]):
fi:
fi:
od:
lu:
end:

#Sipur2(r,K): all the excepion in an r-faced Backgammon end-game
#up to K checkers with TWO dice
Sipur2:=proc(r,K) local K1,i,lu,LU,i1,LU1,j:

LU:=[]:
for K1 from 1 to K do
LU1:=[]:
 for i from 1 to r do
 for j from i+1 to r do
  lu:=Exc2(r,K1,i,j):
   if lu<>{} then

     print(`When the roll was`,{i,j}, `and there were`, K1, `checkers , the exceptions to greedy play were`):
      lu:={seq([{i,j},lu[i1],G2(lu[i1],{i,j})[2][1]],i1=1..nops(lu))}:
       print(lu):
     LU1:=[op(LU1),op(lu)]:
   fi:
 od:
od:
if LU1<>[] then
 LU:=[op(LU),LU1]:
fi:
od:
LU:
end: