# This script determines a reduction vector for a given
# fat point scheme Z = m_1p_1+...+m_rp_r with associated matrix M.

# NOTE: the reduction vector you get may not always be the best one.
# It's chosen using a greedy algorithm.

# M is any 0,1 matrix with r columns such that the dot products 
# of any two distinct rows is at most 1 (since different lines can't
# both pass through the same two points), the dot product of any
# two distinct columns is also at most 1 (since different points
# can't both be on the same two lines), and every column has a 1 in it.
# Rows indicate which subsets of points are to be taken
# as being collinear: for each maximal subset S
# of collinear points there is a row in M having 1's
# where the ith column of the row is 1 iff point p_i is in S.

# FOR THIS SCRIPT TO WORK, EVERY COLUMN OF M MUST HAVE A 1 IN IT AND
# AT LEAST ONE OF THE MULTIPLICITIES m_i MUST BE POSITIVE.

# Let File1 contain the matrix M; i.e., the rows of File1 are
# just the rows of M.

# Let each row of File2 be the multiplicities m_1 m_2 ... m_r
# of a fat point scheme whose reduction vector with respect to M
# you want to compute. Be sure that the m_i specify the multiplicities 
# of the points where the points are taken in the same order as the
# corresponding columns in M. So column 3 of M and multiplicity m_3
# should both correspond to the same point, point p_3.

# The reduction vector d is just the sequence d_1 d_2 ... d_n
# where d_1 is the max entry of Mm^T, where m=(m_1 m_2 ... m_r).
# If the max entry is in row i, then replace m by (m-row(i))_+
# where row(i) is the ith row of M, and for any vector v,
# v_+ denotes zeroing out any negative entries. 
# Then d_2 is the max entry of Mm^T with respect to this new m.
# Repeat until the max entry of Mm^T is 0 (i.e., d_i = 0).

# Note: the ordering of the entries of d now is the reverse of what it was
# in earlier versions of this script.

# Run as:
# awk -f reductionvectorscript10-17-10 File1 File2
# The output is d (or, if m has more than one row, the output is,
# line by line, the d you get for each row of m).

# A typical use is in conjunction with the script HFBoundsAndBettis10-17-2010, to
# compute our upper and lower bounds:
# awk -f reductionvectorscript10-17-10 File1 File2 | awk -f HFBoundsAndBettis10-17-2010


FILENAME == ARGV[1] {
matroid[FNR]=$0
matroidn=FNR}

FILENAME == ARGV[2] {
fatptmults[FNR]=$0
fatptmultsn=FNR}

END{
for(i=1;i<= fatptmultsn;i++) {
	reductionvector=""
	tmpmultvector = fatptmults[i]
	while(totweight(tmpmultvector)>0) {
		maxdot=0
		maxdotrow=0
		for(j=1;j<=matroidn;j++) {
			tmpdot=dot(tmpmultvector, matroid[j])
			if(maxdot<tmpdot) {
				maxdot=tmpdot
				maxdotrow=j
			}
		}
		if(maxdot>0) {
			reductionvector= reductionvector" "maxdot
			tmpmultvector=subtract(tmpmultvector, matroid[maxdotrow])
			# DEBUGGING CODE: print tmpmultvector" .. "reductionvector
		}
		
	}
	print reductionvector
}}


function totweight(F,   D,n,s,i)
{
n=split(F,D," ")
s=0
for(i=1;i<=n;i++)
	s=s+D[i]
return s
}

function dot(F,G,    n,D,i,E,m,s)
{
n=split(F,D," ")
m=split(G,E," ")
if(n>m)
  n=m
s=0
for(i=1;i<=n;i++)
  s=s+D[i]*E[i]
return s
}


function subtract(F,G,  C,D,m,n,i,x,s,t)
{n=split(F,C," ")
m=split(G,D," ")
if(n<m) {
	for(i=1;i<=n;i++)
		x[i]=C[i]-D[i]
	for(i=n+1;i<=m;i++)
		x[i]=C[i]
}
if(m<=n) {
	for(i=1;i<=m;i++)
		x[i]=C[i]-D[i]
	for(i=m+1;i<=n;i++)
		x[i]=-D[i]
}
if(n<m)
	n=m
s=x[1]
if(s<0)
	s=0
for(i=2;i<=n;i++) {
	t=x[i]
	if(t<0)
		t=0
	s=s" "t}

return s
}