interface(echo = 0, quiet = true, screenwidth = infinity):

# The purpose of this script is to take each triple of weights from a
# list of such triples and to produce the LattE input file for the type
# B_r Berenstein-Zelevinsky polytope corresponding to that triple of
# weights.  All the weights in the list must have the same rank.  Recall
# that, given the tiple of weights (la, mu, nu), the corresponding
# BZ-polytope contains a number integer lattice points equal to the
# Clebsch--Gordan coefficient C_{la, nu}^mu (note the position of each
# weight).
#
# To use this program, create a file called "weights" containing triples
# of weights, each preceded by a NULL line.  The file "weights" should
# be in the same directory as this Maple script.  This program assumes
# that your weights are written in the basis of fundamental weights for
# B_r and that they all have the same rank.
#
# Here is an example of the contents of an input file "weights"
# containing two triples of weights, each with rank 4:
#
# NULL
# [17, 21, 20, 18]
# [34, 18, 33, 10]
# [20, 34, 18, 12]
# NULL
# [36, 9, 13, 37]
# [26, 32, 37, 32]
# [28, 21, 19, 0]
#
# Make sure (1) to place a NULL line before every triple -- even the
# first one, and (2) *not* to put a new line (LF) at the end of the last
# triple in your file.
#
# After you have created the file "weights", create a folder called
# "LattEPolytopes" in the same directory as this Maple file.
#
# You are now ready to run the maple script by typing
# 
#    maple<Br_BZpolytopesForLattE
#    
# in a unix terminal.
#
# This will produce all of the LattE input files for the triples of
# weights in the file "weights".  The LattE files will appear in the
# directory "LattEPolytopes".  Each file's name will express the triple
# of weights which produce the polytope described by that file.
#
# The code for the script follows.


with(linalg):
with(LinearAlgebra):
with(StringTools):

# Define r to be the rank of the weights from the file "weights".
#
readline(weights):

r := nops(parse(readline(weights))):

close(weights):

# Generate the simple roots alpha_0, ... , alpha_{r-1} for B_r (with
# respect to the basis of fundamental weights).
#
al[0] := [2, -1, seq(0, n=2..r)]:

al[1] := [-2, 2, -1, seq(0, n=4..r)]:

for i from 2 to r-2 do
	al[i] := [seq(0, n=1..i-1), -1, 2, -1, seq(0, n=i+3..r)];
od:

al[r-1] := [seq(0, n=1..r-2), -1, 2]:

# Generate the coroots alpha_0^\vee, ... , alpha_{r-1}^\vee corresponding
# to the simple roots generated above.
#
for i from 0 to r-1 do
	alcheck[i] := [seq(0, n=1..i), 1, seq(0, n=i+2..r)];
od:

# So long as there are more weights to be read from the file weights, keep
# generating files.
#
while readline(weights) <> 0 do
	
	# Set the interface to the ideal settings for writing to
	# BZfromMaple
	#
	interface(prettyprint = false):
	
	# Define la, mu, and nu to be the current weights we will use to
	# generate the BZ polytope.
	#
	la := parse(readline(weights)):
	mu := parse(readline(weights)):
	nu := parse(readline(weights)):
 	
	# Unassign the values of the t[i,j] assigned the last time through
	# the while loop.
	#
	unassign('t');
	
	#Set t[i,j] = 0 unless 0 <= |i| <= j < r.
	#
	for j from -2*r to 2*r do
		for i from -2*r to 2*r do
			if not (0 <= abs(i) and abs(i) <= j and j < r) then
				t[i,j] := 0;
			fi;
		od;
	od:
    
   	# We now generate the equalities and inequalities.  These
   	# equalities and inequalities are divided into four types.  For
   	# Types (1), (3), and (4), the equal signs in the equations will
   	# all be interpreted as inequalites <= by LattE, but they are
   	# treated as equalities for now, since this makes them easier to
   	# translate into matrices in Maple.
	
	# Type (1) inequalities.
	#
	type1 := [ -t[0,0] = 0 ]:
	 
	for j from 1 to r-1 do
		for i from -j to -2 do
			type1 := [ op(type1), -2*t[i,j] + 2*t[i+1,j] = 0 ]:
		od:
		
		type1 := [ op(type1), -2*t[-1,j] + t[0,j] = 0 ]:
		
		type1 := [ op(type1), -t[0,j] + 2*t[1,j] = 0 ]:
		
		for i from 1 to j-1 do
			type1 := [ op(type1), -2*t[i,j] + 2*t[i+1,j] = 0 ]:
		od:
		
		type1 := [ op(type1), -2*t[j,j] = 0 ]:
	od:
	
	# Type (2) equalities.  (These will be interpreted as actual
	# equalities by LattE.)
	#
	type2 := []:
	
	for k from 1 to r do
		kthlinesum := 0:
		
		for j from 0 to r-1 do
			for i from -j to j do
				kthlinesum := kthlinesum + t[i,j] * al[abs(i)][k]:
			od:
		od:
	
		type2 := [ op(type2), kthlinesum = la[k] + nu[k] - mu[k] ]:
	od:
	
	# Type (3) inequalities.
	#
	# But first, convert the weights and coroots to capital-V Vectors
	# so that we can do easy dot products with the LinearAlgebra
	# package.
	#
	la := convert(la, Vector):
	mu := convert(mu, Vector):
	nu := convert(nu, Vector):
	
	for i from 0 to r-1 do
		alcheck[i] := convert(alcheck[i], Vector):
	od:
	
	type3 := [ t[0,0] = la.alcheck[0] ]:
	
	for j from 1 to r-1 do
		type3 := [ op(type3), t[j,j] = la.alcheck[j] ]:  
		type3 := [ op(type3), t[0,j] - t[1,j-1] - t[-1,j] = la.alcheck[j] ]:
		type3 := [ op(type3), t[1,j-1] + t[-1,j] - t[0,j-1] = la.alcheck[j] ]:
		
		for i from 1 to j-1 by 1 do
			type3 := [ op(type3), t[i,j] + t[-i,j] - t[i+1,j-1] - t[-i-1,j] = la.alcheck[j] ]:
			type3 := [ op(type3), t[i+1,j-1] + t[-i-1,j] - t[-i,j-1] - t[i,j-1] = la.alcheck[j] ]:
			type3 := [ op(type3), t[-i,j] - t[-i,j-1] = la.alcheck[j] ]:
			type3 := [ op(type3), t[i,j] - t[i,j-1] = la.alcheck[j] ]:
		od:
	od:
	
	# Type (4) inequalities.
	#
	type4 := []:
	
	for j from 0 to r-1 do
		sum1 := t[0,j]:
		
		for k from j+1 to r do
			sum1 := sum1 + 2*(t[0,k] - t[-1,k] - t[1,k-1]):
		od:
		
		type4 := [ op(type4), sum1 = nu.alcheck[0] ]:
	od:
	
	for j from 1 to r-1 do
		sum2 := t[-1,j]:
		sum3 := t[1,j]:
		
		for k from j+1 to r do
			sum2 := sum2 + 2*t[-1,k] + 2*t[1,k-1] - t[0,k-1] - t[-2,k] - t[2,k-1]:
			sum3 := sum3 + 2*t[-1,k] + 2*t[1,k] - t[0,k] - t[-2,k] - t[2,k-1]:
		od;
		
		type4 := [ op(type4), sum2 = nu.alcheck[1], sum3 = nu.alcheck[1] ]:
	od:
	
	for i from 2 to r-1 do
		for j from i to r-1 do
			sum4 := t[-i,j]:
			sum5 := t[i,j]:
			
			for k from j+1 to r do
				sum4 := sum4 + 2*t[-i,k] + 2*t[i,k-1] - t[-i+1,k-1] - t[i-1,k-1] - t[-i-1,k] - t[i+1,k-1]:
				sum5 := sum5 + 2*t[-i,k] + 2*t[i,k] - t[-i+1,k] - t[i-1,k] - t[-i-1,k] - t[i+1,k-1]:
			od;
			
			type4 := [ op(type4), sum4 = nu.alcheck[i], sum5 = nu.alcheck[i] ]:
		od;
	od:
	
	# Produce the matrix A and b such that Ax = b corresponds to the
	# linear equalities produced above
	#
	AllIneq := [ op(type2), op(type1), op(type3), op(type4) ]:

	variables := []:
	
	for j from 0 to r-1 do
		for i from -j to j do
			variables := [ op(variables), t[i,j] ];
		od:
	od:
	
	A := genmatrix(AllIneq, variables, b):
	
	# Produce the array (b, -A), denoted bnegA.
	#
	bnegA := augment(b,-A):
	
	bnegA := convert(bnegA, Matrix):
	
	# Create a Latte-readable file (modulo square brackets, commas, and
	# other such junk) for this system of inequalities
	#
	writeto("BZfromMaple"):
	print(Dimension(bnegA));
	print(convert(bnegA, listlist));
	print([linearity, nops(type2), seq(i, i = 1..nops(type2))]);
	close(BZfromMaple):
	
	# Now make that file really LattE-readable, and print the results
	# to a file whose file-name expresses la, mu, and nu.
	#
	interface(indentamount = 0 , prettyprint = 1):
	
	s := readline(BZfromMaple):
	news := SubstituteAll(s, ",", " "):
	
	t := readline(BZfromMaple):
	newt := SubstituteAll(t, "\[\[", ""):
	newt := SubstituteAll(newt, ",", " "):
	newt := SubstituteAll(newt, "\]  \[", "\n"):
	newt := SubstituteAll(newt, "\]\]", ""):
	
	u := readline(BZfromMaple):
	newu := SubstituteAll(u, "\[", ""):
	newu := SubstituteAll(newu, ",", " "):
	newu := SubstituteAll(newu, "\]", ""):
	
	new := cat(news, "\n", newt, "\n", newu):
		
	# Generate the filename for the LattE file that will contain this
	# BZ-polytope.
	#
	filename := "LattEPolytopes/Br-":
	
	for i from 1 to r-1 do
		filename := cat(filename, la[i], ".");
	od:
	
	filename := cat(filename, la[r], "_"):
	
	for i from 1 to r-1 do
		filename := cat(filename, mu[i], ".");
	od:
	
	filename := cat(filename, mu[r], "_"):
	
	for i from 1 to r-1 do
		filename := cat(filename, nu[i], ".");
	od:
	
	filename := cat(filename, nu[r]):
	
	# Write the BZ polytope data in LattE format to the file filename.
	#
	writeto(filename):
	print(convert(new, symbol));
	writeto(terminal):
	close(filename):
od:

close(weights):

# Remove the file BZfromMaple that was created during this script's
# execution.
#
! rm BZfromMaple