######################################################################
##NumberCrunch.txt: Save this file as  NumberCrunch.txt              #
## To use it, stay in the                                            #
##same directory, get into Maple (by typing: maple <Enter> )         #
##and then type:  read  NumberCrunch.txt<Enter>                      #
##Then follow the instructions given there                           #
##                                                                   #
##Written by Doron Zeilberger, Rutgers University ,                  #
#DoronZeil at gmail dot com                                          #
######################################################################
 
#Created: July 7, 2019
 
print(`Created:  July 7, 2019`):
print(` This is NumberCrunch.txt`):
print(`To create and solve Number Crunch puzzles that appear in the daily Australian newpaper "The Australian"`):
print(`It is the package that created the webook `):
print(`100 Number Crunch Puzzles by`):
print(`by Shalosh B. Ekhad `):

print(`and also available from Zeilberger's website`):
print(``):
print(`Please report bugs to DoronZeil at gmail dot com `):
print(``):
 print(`The most current version of this  package and paper`):
 print(` are  available from`):
 print(`http://sites.math.rutgers.edu/~zeilberg/  .`):

print(`---------------------------------------`):
 print(`For a list of the Supporting procedures type ezra1();, for help with`):
 print(`a specific procedure, type ezra(procedure_name);   .`):
 print(``):
print(`---------------------------------------`):

print(`---------------------------------------`):
 print(`For a list of the MAIN procedures type ezra();, for help with`):
 print(`a specific procedure, type ezra(procedure_name);   .`):
 print(``):
print(`---------------------------------------`):

with(combinat):

ezra1:=proc()

if args=NULL then
 print(` The supporting procedures are: `):
 print(``):

else
ezra(args):
fi:

end:

ezra:=proc()

if args=NULL then
 print(`The main procedures are: DrawPu, DrawSol, MakePuz, RP `):
 print(` `):



elif nops([args])=1 and op(1,[args])=DrawPu then
print(`DrawPu(P): Draws a Number Crunch puzzle, P . Try:`):
print(`Try: `):
print(`gu:=RP(1000,X); `):
print(`DrawPu(gu[1])`):
print(``):

elif nops([args])=1 and op(1,[args])=DrawSol then
print(`DrawSol(P): Draws a Number Crunch puzzle, P . Try:`):
print(`Try: `):
print(`gu:=RP(1000,X); `):
print(`DrawSol(gu[2])`):
print(``):

elif nops([args])=1 and op(1,[args])=MakePuz then
print(`MakePuz(L,K,X): makes L puzzles RP(K,X). Try:`):
print(`MakePuz(10,100,X);`):

elif nops([args])=1 and op(1,[args])=RP then
print(`RP(K,X): Makes a random Number Cruncher problem with a unique solution, the unknowns indicated by X, trying K times (or returns FAIL).`):
print(`It also outputs the solution`):
print(`Try: `):
print(`RP(100,x);`):
print(``):


else
print(`There is no ezra for`,args):
fi:
 
end:



ez:=proc(): print(`RPS(K), RPS1() , EvalTri(L) , Trios(S) , SolTri(L,S)`): 
print(` SolTri([X,`+`,X,`+`,X,10],{seq(i,i=1..9)}), RP(K,X)  `):
print("SolP([[X,`/`,X,`*`,X,10],[`*`,0,`+`,0,`+`,0],[X,`-`,X,`/`,X,4],[`+`,0,`-`,0,`-`,0],[X,`*`,X,`/`,X,3],[20,0,6,0,0,0]]);"):
print("SolP1([[X,`/`,X,`*`,X,10],[`*`,0,`+`,0,`+`,0],[X,`-`,X,`/`,X,4],[`+`,0,`-`,0,`-`,0],[X,`*`,X,`/`,X,3],[20,0,6,0,0,0]]);");
print("Cont1s1(op(SolP1([[X,`/`,X,`*`,X,10],[`*`,0,`+`,0,`+`,0],[X,`-`,X,`/`,X,4],[`+`,0,`-`,0,`-`,0],[X,`*`,X,`/`,X,3],[20,0,6,0,0,0]])[1]));"):
print("SolP2([[X,`/`,X,`*`,X,10],[`*`,0,`+`,0,`+`,0],[X,`-`,X,`/`,X,4],[`+`,0,`-`,0,`-`,0],[X,`*`,X,`/`,X,3],[20,0,6,0,0,0]]);");
print("Cont2s1(op(SolP2([[X,`/`,X,`*`,X,10],[`*`,0,`+`,0,`+`,0],[X,`-`,X,`/`,X,4],[`+`,0,`-`,0,`-`,0],[X,`*`,X,`/`,X,3],[20,0,6,0,0,0]])[1]));"):
print("SolP3([[X,`/`,X,`*`,X,10],[`*`,0,`+`,0,`+`,0],[X,`-`,X,`/`,X,4],[`+`,0,`-`,0,`-`,0],[X,`*`,X,`/`,X,3],[20,0,6,0,0,0]]);");
print(`IsGood(P); `):
print("SolP([[X,`/`,X,`*`,X,10],[`*`,0,`+`,0,`+`,0],[X,`-`,X,`/`,X,4],[`+`,0,`-`,0,`-`,0],[X,`*`,X,`/`,X,3],[20,0,6,0,0,0]]);");
#print("SolP([[X,"/",X,"*",X,10],["*",0,`+`,0,`+`,0],[X,`-`,X,`/`,X,4],[`+`,0,`-`,0,`-`,0],[X,`*`,X,`/`,X,3],[20,0,6,0,0,0]]);");
end:

ezP:=proc(): 
print(`DrawPu(P,eps); `):
print(`DrawSol(P,eps); `):
print(`DrawPlus(pt,eps); `): 
print(`DrawMinus(pt,eps); `): 
print(`DrawTimes(pt,eps); `): 
print(`DrawDivide(pt,eps); `): 
print(`DrawSign(S,pt,eps);`):
end:

with(combinat):

#EvalTri(T): given a triple [a,op1,b,op2,c] where a, b,c are integers and op1 and op2 are in {`+`,`-`,`*`,`/`} putputs
#its outcome, or returns FAIL if it is not a positive integer. Try:
#EvalTri([4,`+`,5,`*`, 6];
EvalTri:=proc(T) local L,gu,t:

gu:=[T[2],T[4]];
L:=[T[1],T[3],T[5]]:


if gu=[`+`,`+`] then
 t:= L[1]+L[2]+L[3]:

elif gu=[`+`,`-`] then
 t:= L[1]+L[2]-L[3]:

elif gu=[`+`,`*`] then
 t:=L[1]+L[2]*L[3]:

elif gu=[`+`,`/`] then
 t:=L[1]+(L[2]/L[3]):

elif gu=[`-`,`+`] then
 t:= L[1]-L[2]+L[3]:

elif gu=[`-`,`-`] then
 t:= (L[1]-L[2])-L[3]:

elif gu=[`-`,`*`] then
 t:=L[1]-L[2]*L[3]:

elif gu=[`-`,`/`] then
 t:=L[1]-(L[2]/L[3]):

elif gu=[`*`,`+`] then
 t:= L[1]*L[2]+L[3]:

elif gu=[`*`,`-`] then
 t:= L[1]*L[2]-L[3]:

elif gu=[`*`,`*`] then
 t:= L[1]*L[2]*L[3]:

elif gu=[`*`,`/`] then
 t:= L[1]*L[2]/L[3]:

elif gu=[`/`,`+`] then
 t:= (L[1]/L[2])+L[3]:

elif gu=[`/`,`-`] then
 t:= (L[1]/L[2])-L[3]:

elif gu=[`/`,`*`] then
 t:= (L[1]/L[2])*L[3]:

elif gu=[`/`,`/`] then
 t:= (L[1]/L[2])/L[3]:

fi:

if not (type(t,integer) and t>=0) then
 RETURN(FAIL):
else
 RETURN(t):
fi:

end:


#RPS1(): A random puzzle-with-solution Number Cruncher problem. Try:
#RPS1();
RPS1:=proc() local ra,i,j,pi,co,L,T,A,t:
 pi:=randperm(9):
ra:=rand(1..4):

L:=[`+`,`-`,`*`,`/`]:

for i from 1 to 6 do
 for j from 1 to 6 do
  T[i,j]:=0:
 od:
od:


co:=1:
 for i from 1 to 3 do
  for j from 1 to 3 do
   T[2*i-1,2*j-1]:= pi[co]:
   co:=co+1:
  od:
 od:


 for i from 1 to 3 do
  for j from 1 to 2 do
   T[2*i-1,2*j]:= L[ra()]:
  od:
 od:


 for i from 1 to 2 do
  for j from 1 to 3 do
   T[2*i,2*j-1]:= L[ra()]:
  od:
 od:


for i from 1 to 3 do
 A:=[seq(T[2*i-1,j],j=1..5)]:
 t:=EvalTri(A):
 if t=FAIL then
   RETURN(FAIL):
  else
 T[2*i-1,6]:=t:
 fi:
od:



for j from 1 to 3 do
 A:=[seq(T[i,2*j-1],i=1..5)]:
 t:=EvalTri(A):
 if t=FAIL then
   RETURN(FAIL):
  else
 T[6,2*j-1]:=t:
 fi:
od:



[seq([seq(T[i,j],j=1..6)],i=1..6)]:

end:






#RPS(K): A good random puzzle-with-solution Number Cruncher problem using K trials
RPS:=proc(K) local i,gu:
 
for i from 1 to K do
 gu:=RPS1():
if gu<>FAIL then
  RETURN(gu):
fi:
od:
FAIL:
end:


#Trios(S): the set of distict triples from the set S. Try Trios({1,5,7,8}).
Trios:=proc(S) local a,b,c,gu:

gu:={}:

for a in S do
 for b in S minus {a} do
   for c in S minus {a,b} do
    gu:=gu union {[a,b,c]}:
   od:
 od:
od:
gu:
end:



#SolTri(L,S): Given a list [X,op1,X,op2,X,outcome] and a set S, finds the set of all Trios(S) such that if you plug-in
#[a,b,c] for the X's you would get the outcome. Try:
#SolTri([X,`+`,X,`+`,X,10],{seq(i,i=1..9)});
SolTri:=proc(L,S) local mu,mu1,gu,L1:
mu:=Trios(S):
gu:={}:

for mu1 in mu do
 L1:=[mu1[1],L[2],mu1[2],L[4],mu1[3]]:
 if EvalTri(L1)=L[6] then
  gu:=gu union {mu1}:
 fi:
od:
gu:
end:

#RP(K,X): Makes a random Number Cruncher problem with a unique solution, the unknowns indicated by X, trying K times (or returns FAIL).
#It also outputs the solution
#Try:
#RP(100,x);
RP:=proc(K,X) local mu,gu,i,j,gu1,lu:

mu:=RPS(K):

if mu=FAIL then
 RETURN(FAIL):
fi:

gu:=[]:

for i from 1 to 6 do

 if i mod 2=0 then
  gu:=[op(gu),mu[i]]:
 else
  
  gu1:=[]: 
 
   for j from 1 to 6 do
    
      if j mod 2=0 then
       gu1:=[op(gu1),mu[i][j]]:
      else
       gu1:=[op(gu1),X]:
      fi:
  od:

 gu:=[op(gu),gu1]:
  
 fi:



od:

lu:=SolP(gu): 

if nops(lu)>1 then
# print(gu, `has more than one solution`):
 RETURN(FAIL):
fi:

[gu,lu[1]]:

end:



#SolP1(P): Given a Number-Cruncher problem P , finds its set of partial solutions only satisfying its first row. For example, try:
#SolP1([[X,`/`,X,`*`,X,10],[`*`,0,`+`,0,`+`,0],[X,`-`,X,`/`,X,4],[`+`,0,`-`,0,`-`,0],[X,`*`,X,`/`,X,3],[20,0,6,0,0,0]]);
SolP1:=proc(P) local gu,gu1,i:
gu:=SolTri(P[1],{seq(i,i=1..9)}):
{seq([ [ [gu1[1],P[1][2],gu1[2],P[1][4],gu1[3],P[1][6]],op(2..6,P)],{seq(i,i=1..9)} minus {op(gu1)}],gu1 in gu)}:
end:

#Cont1s1(P,S): inputs a partial solution P, with the first row filled  and a set of available integers,
#outputs all the possible continuations to the  second row, followed by the still available integers for the
#third row. Try:
#
#Cont1s1(op(SolP1([[X,`/`,X,`*`,X,10],[`*`,0,`+`,0,`+`,0],[X,`-`,X,`/`,X,4],[`+`,0,`-`,0,`-`,0],[X,`*`,X,`/`,X,3],[20,0,6,0,0,0]])[1]));
Cont1s1:=proc(P,S) local gu,gu1:
gu:=SolTri(P[3],S):
{seq([ [op(1..2,P),    [gu1[1],P[3][2],gu1[2],P[3][4],gu1[3],P[3][6]], op(4..6,P)],S minus {op(gu1)}],gu1 in gu)}:
end:


#SolP2(P): Given a Number-Cruncher problem P , finds its set of partial solutions only satisfying its first and second row. For example, try:
#SolP2([[X,`/`,X,`*`,X,10],[`*`,0,`+`,0,`+`,0],[X,`-`,X,`/`,X,4],[`+`,0,`-`,0,`-`,0],[X,`*`,X,`/`,X,3],[20,0,6,0,0,0]]);
SolP2:=proc(P) local gu,gu1:
gu:=SolP1(P):
{seq(op(Cont1s1(op(gu1))),gu1 in gu)}:
end:


#Cont2s1(P,S): inputs a partial solution P, with the first and second row filled  and a set of available integers,
#outputs all the possible continuations to the  third row
#Cont2s1(op(SolP2([[X,`/`,X,`*`,X,10],[`*`,0,`+`,0,`+`,0],[X,`-`,X,`/`,X,4],[`+`,0,`-`,0,`-`,0],[X,`*`,X,`/`,X,3],[20,0,6,0,0,0]])[1]));
Cont2s1:=proc(P,S) local gu,gu1:
gu:=SolTri(P[5],S):
{seq([op(1..4,P),    [gu1[1],P[5][2],gu1[2],P[5][4],gu1[3],P[5][6]], P[6]], gu1 in gu)}:
end:

#SolP3(P): Given a Number-Cruncher problem P , finds its set of partial solutions only satisfying its rows, disregarding the columns. Try:
#SolP3([[X,`/`,X,`*`,X,10],[`*`,0,`+`,0,`+`,0],[X,`-`,X,`/`,X,4],[`+`,0,`-`,0,`-`,0],[X,`*`,X,`/`,X,3],[20,0,6,0,0,0]]);
SolP3:=proc(P) local gu,gu1:
gu:=SolP2(P):
{seq(op(Cont2s1(op(gu1))),gu1 in gu)}:
end:

#IsGood(P): inputs a tentative solution obeying the rows restrictions, checks whether the column restrictions also hold
IsGood:=proc(P) local i,j,T:

for j from 1 by 2 to 5 do
 T:=[seq(P[i][j],i=1..5)]:
 if EvalTri(T)<>P[6,j] then
  RETURN(false):
 fi:
od:
true:
end:


#SolP(P): Given a Number-Cruncher problem P , finds its set of solution. For example, try:
#SolP([[X,`/`,X,`*`,X,10],[`*`,0,`+`,0,`+`,0],[X,`-`,X,`/`,X,4],[`+`,0,`-`,0,`-`,0],[X,`*`,X,`/`,X,3],[20,0,6,0,0,0]]);
SolP:=proc(P) local mu,mu1,gu:
mu:=SolP3(P):
gu:={}:

for mu1 in mu do
 if IsGood(mu1) then
   gu:=gu union {mu1}:
 fi:
od:
gu:

end:

with(plots):





DrawPlus:=proc(pt,eps) local d:
d:=plot([
[pt[1]-eps,pt[2]],[pt[1]+eps,pt[2]]
        ],axes=none):
d:=d,plot([[pt[1],pt[2]-eps],[pt[1],pt[2]+eps]]):
d:
end:

DrawMinus:=proc(pt,eps) local d:
d:=plot([
[pt[1]-eps,pt[2]],[pt[1]+eps,pt[2]]
        ],axes=none):
d:
end:

DrawTimes:=proc(pt,eps) local d:
d:=plot([
[pt[1]-eps,pt[2]+eps],[pt[1]+eps,pt[2]-eps]
        ],axes=none):
d:=d,plot([[pt[1]-eps,pt[2]-eps],[pt[1]+eps,pt[2]+eps]]):
d:
end:

DrawDivide:=proc(pt,eps) local d:
d:=plot([
[pt[1]-eps,pt[2]],[pt[1]+eps,pt[2]]
        ],axes=none):
d:=d,plot([[pt[1],pt[2]-eps]], style=point,symbol=circle):
d:=d,plot([[pt[1],pt[2]+eps]], style=point,symbol=circle):
d:
end:

#DrawSign(S,pt,eps): draws the sign S with width eps
DrawSign:=proc(S,pt,eps)

if S=`+` then
 DrawPlus(pt,eps):
elif S=`-` then
 DrawMinus(pt,eps):
elif S=`*` then
 DrawTimes(pt,eps):
elif S=`/` then
 DrawDivide(pt,eps):
else
 FAIL:
fi:
end:

#DrawPu(P):
DrawPu	:=proc(P) local d,i,j:

d:=plot([[1,1],[6,1]],axes=none,scaling=constrained):

for i from 1 to 6 do
d:=d,plot([[1,i],[6,i]],axes=none):
od:

for j from 1 to 6 do
d:=d,plot([[j,1],[j,6]],axes=none):
od:

for i from 0 to 2 do
 d:=d,plot([[1+2*i,1],[1+2*i,0],[2+2*i,0],[2+2*i,1]]):
od:

for j from 0 to 2 do

 d:=d,plot([[6,1+2*j],[7,1+2*j],[7,2+2*j],[6,2+2*j]]):

od:

#textplot([C[1],C[2],NC]):

d:=d,textplot([1.5,0.5,P[6][1]]), textplot([3.5,0.5,P[6][3]]), textplot([5.5,0.5,P[6][5]]):
d:=d,textplot([6.5,1.5,P[5][6]]),textplot([6.5,3.5,P[3][6]]),textplot([6.5,5.5,P[1][6]]):

for i from 1 to 5 do
for j from 1 to 5 do
 if i+j mod 2=1 then
  d:=d,textplot([i+0.5,j+0.5,P[6-j][i]])
 fi:
od:
od:

display(d):
end:




#DrawSol(P):
DrawSol:=proc(P) local d,i,j,eps:
eps:=0.1:
d:=plot([[1,1],[6,1]],axes=none,scaling=constrained):

for i from 1 to 6 do
d:=d,plot([[1,i],[6,i]],axes=none):
od:

for j from 1 to 6 do
d:=d,plot([[j,1],[j,6]],axes=none):
od:

for i from 0 to 2 do
 d:=d,plot([[1+2*i,1],[1+2*i,0],[2+2*i,0],[2+2*i,1]]):
od:

for j from 0 to 2 do

 d:=d,plot([[6,1+2*j],[7,1+2*j],[7,2+2*j],[6,2+2*j]]):

od:

#textplot([C[1],C[2],NC]):

d:=d,textplot([1.5,0.5,P[6][1]]), textplot([3.5,0.5,P[6][3]]), textplot([5.5,0.5,P[6][5]]):
d:=d,textplot([6.5,1.5,P[5][6]]),textplot([6.5,3.5,P[3][6]]),textplot([6.5,5.5,P[1][6]]):

for i from 1 to 5 do
for j from 1 to 5 do

 if not (i mod 2=0 and j mod 2=0) then

  if i+j mod 2=1 and j<=5 then
    d:=d,DrawSign(P[6-j][i], [i+0.5,j+0.5],eps):
else
 
 d:=d,textplot([i+0.5,j+0.5,P[6-j][i]])
 fi:

 fi:
od:
od:

display(d):
end:

#DrawPu(P):
DrawPu:=proc(P) local d,i,j,eps:
eps:=0.1:

d:=plot([[1,1],[6,1]],axes=none,scaling=constrained):

for i from 1 to 6 do
d:=d,plot([[1,i],[6,i]],axes=none):
od:

for j from 1 to 6 do
d:=d,plot([[j,1],[j,6]],axes=none):
od:

for i from 0 to 2 do
 d:=d,plot([[1+2*i,1],[1+2*i,0],[2+2*i,0],[2+2*i,1]]):
od:

for j from 0 to 2 do

 d:=d,plot([[6,1+2*j],[7,1+2*j],[7,2+2*j],[6,2+2*j]]):

od:

#textplot([C[1],C[2],NC]):

d:=d,textplot([1.5,0.5,P[6][1]]), textplot([3.5,0.5,P[6][3]]), textplot([5.5,0.5,P[6][5]]):
d:=d,textplot([6.5,1.5,P[5][6]]),textplot([6.5,3.5,P[3][6]]),textplot([6.5,5.5,P[1][6]]):

for i from 1 to 5 do
for j from 1 to 5 do
 if not (i mod 2=0 and j mod 2=0) then

   if i+j mod 2=1 and j<=5 then
    d:=d,DrawSign(P[6-j][i], [i+0.5,j+0.5],eps):
 fi:
fi:
od:
od:

display(d):
end:

#MakePuz(L,K,X): makes L puzzles RP(K,X)
MakePuz:=proc(K,L,X) local gu,mu:

gu:=[]:

while nops(gu)<L do
 mu:=RP(K,X):
 if mu<>FAIL then
  gu:=[op(gu),mu]:
 fi:
od:

gu:

end: