# Note: In the previous version, HFBoundsAndBettis3-31-09, 
# the lwrbnd on ideals was the bound correpsonding to the lower bound on the HF of the scheme,
# but was in fact the upper bound on the HF of the ideal. Likewise the uppbnd on ideals
# was actually the lower bound, but corresponded to the upper bound on the scheme.
# This was confusing so has been changed.

# Also, the order of the entries in the the old version of the paper is the reverse of what
# we the paper now uses; this script now expects input using the new ordering.

# Enter reduction vector d = (d_1,...,d_n); get out first differences for
# both an upper and lower bound on Hilbert function of quotient ring
# and, if d is Bezout (i.e., if upper and lower Hilbert 
# functions coincide, also outputs bounds on the graded Betti numbers
# Run as:
# echo d_1 ... d_n | awk -f HFBoundsAndBettis10-17-2010
# Basic idea:
# This program computes what would be the Hilbert function 
# of the ideal for a given line count vector
# where all multiplicities are equal to 1 (this is all we need
# since (d_1,...d_n) is playing the role of a*m).
# It also prints out diag(d_1,...d_n).
# The lower bound comes from pretending the exact sequences 
# obtained by subtracting off a line at a time (starting with 
# the line corresponding to d_n) are exact on global sections
# (but just working formally, so there is no requirement
# that the mi are nondecreasing).
# The lower bound comes from assuming that each sequence
# is as far from exactness on the right as possible
# (but how far sequence i can be from exactness depends 
# on how far sequence i-1 is). If however the d_i are increasing,
# for example, our theorem implies that the lower and upper bounds
# coincide, and equal diag and we also get exact values for the
# graded Betti numbers. Indeed, the lower bound
# always equals diag, with no restrictions on (d_1,...,d_n)
# beyond having d_i > 0.

{tot=0
print "\n*******************************\n"
print "Input vector: "$0
inputcnt=NF
for(i=1;i<= inputcnt;i++) {
x[inputcnt-i+1]=$i
tot=tot+$i}

hhstr=""
hstr=""
z=0
sum=0
t=-1
alphalower=-1
while(sum<tot) {
t=t+1
Gamma=0
for(i=1;i<= inputcnt;i++) {
Gamma = Gamma + gg(t,i, inputcnt)}

# Gamma counts the number of rows in standard diagram for reduction vector such that
# the dots in the row are strictly to the left of line x+y=t+1.
# Hence, if n = max(t+1, number of rows), then diag(t) = n-Gamma

if(t< inputcnt) 
	y=t+1-Gamma
if(t>= inputcnt) 
	y= inputcnt-Gamma

diag[t+1]=y  # diag(t) = diag[t+1]; this just keeps diag[] an array whose first entry for t=0 is diag[1]

hhstr=hhstr" "y
z=z+y
h[t]=z

hstr=hstr" "z
hideal[t]=(t^2+3*t+2)/2-z

# print "hideal["t"]="hideal[t]

if(hideal[t]==0)
	alphalower = t+1

sum=sum+y
}

h[t+1]=tot
hideal[t+1]=((t+1)^2+3*t+5)/2-tot
# print "hideal["t+1"]="hideal[t+1]

diag[t+2]=0
diag[t+3]=0

sum=0
t=-1
prev=0
z=0
alphaupper=-1
Hstr=""
upperbndstr=""
while(sum<tot) {
t=t+1
tmp0=AAA(t, inputcnt)
tmp1=(t+2)*(t+1)/2-tmp0
# print t, tmp0
tmp2=tmp1-prev
upperbnd[t]=tmp2
upperbndstr= upperbndstr" "tmp2
z=z+tmp2
H[t]=z
Hstr=Hstr" "z
Hideal[t]=(t^2+3*t+2)/2-z
# print "Hideal["t"]="Hideal[t]
if(Hideal[t]==0)
	alphaupper = t+1
prev=tmp1
sum=sum+tmp2
}

H[t+1]=tot
Hideal[t+1]=((t+1)^2+3*t+5)/2-tot
# print "Hideal["t+1"]="Hideal[t+1]

i=0
difflag=0
while(H[i]<tot && h[i]<tot) {
	if(H[i]!=h[i])
		difflag=1
	i=i+1
}
if(H[i]!=h[i])
	difflag=1
if(difflag==1)
	print "Upper and lower Hilbert functions are NOT equal."
if(difflag==0)
	print "Upper and lower Hilbert functions ARE equal."
print "alpha lwrbnd= "alphaupper
tmpstr=""
i=alphaupper-1
while(H[i]<tot) {
	tmpstr=tmpstr" "Hideal[i]
	i=i+1}
tau=i

tmpstr=tmpstr" "Hideal[i]
tmpstr=tmpstr" "Hideal[i+1]
print "  lower bound on dim of I_t for t >= alpha-1:"tmpstr
print "alpha uppbnd= "alphalower
tmpstr=""
i=alphalower-1
while(h[i]<tot) {
	tmpstr=tmpstr" "hideal[i]
	i=i+1}
tmpstr=tmpstr" "hideal[i]
tmpstr=tmpstr" "hideal[i+1]
print "  upper bound on dim of I_t for t >= alpha-1:"tmpstr

print "Delta uprbnd= "upperbndstr
print "Delta lwrbnd= "hhstr

print "      uprbnd= "Hstr
print "      lwrbnd= "hstr
print " "

if(difflag ==0) {
	print "The input is Bezout, so we now give upper and lower bounds on the Betti numbers."


Bettiflag=0
for(i=1;i< inputcnt;i++) {
if(x[i]==x[i+1])
	Bettiflag= Bettiflag+1
}

if(Bettiflag==0)
	print "Note: the upper and lower bounds coincide."
if(Bettiflag==1) {
	for(i=1;i< inputcnt;i++) {
		if(x[i]==x[i+1])
			Bettiflag=i+1}
	newflag=0
	for(i=1;i< Bettiflag;i++) {
		if(x[i]!=i)
			newflag=1}
	if(Bettiflag+1<=inputcnt && x[Bettiflag+1]<2)
			newflag=1
	if(newflag==0)
		print "Note : the upper and lower bounds coincide."
}
print " "

	bettigenstr="nmbr of gens ="
	bettisyzstr="nmbr of syzs ="
	for(i=0;i<alphalower;i++) {
		bettigenstr = bettigenstr" 0"
		bettisyzstr = bettisyzstr" 0"}
	bettigenstr = bettigenstr" "Hideal[alphalower]
# print ".........", alphalower, Hideal[alphalower]
	bettisyzstr = bettisyzstr" 0"
	for(i=alphalower;i<=tau+2;i++) {
		neg=0
		jflag=0
		for(j=1;j<=inputcnt;j++) {
			if(i-(j-1)-x[inputcnt-j+1]>=0)
				 jflag = 1
			if(i-(j-1)-x[inputcnt-j+1]<0 && jflag == 0)
				neg=neg+1}
		tttmp1=neg - diag[i+2]
		ttttmp1=diag[i+1] - diag[i+2]
		bettigenstr = bettigenstr" "ttttmp1
		tttmp2=(diag[i+2]-diag[i+1])-(diag[i+1]-diag[i])+ttttmp1
		bettisyzstr = bettisyzstr" "tttmp2}
print "Upper bounds for the graded Betti numbers:"
print bettigenstr
print bettisyzstr

	bettigenstr="nmbr of gens ="
	bettisyzstr="nmbr of syzs ="
	for(i=0;i<alphalower;i++) {
		bettigenstr = bettigenstr" 0"
		bettisyzstr = bettisyzstr" 0"}
	bettigenstr = bettigenstr" "Hideal[alphalower]
	bettisyzstr = bettisyzstr" 0"

	for(i=alphalower;i<=tau+2;i++) {
		neg=0
		jflag=0
		for(j=1;j<=inputcnt;j++) {
			if(i-(j-1)-x[inputcnt-j+1]>=0)
				 jflag = 1
			if(i-(j-1)-x[inputcnt-j+1]<0 && jflag == 0)
				neg=neg+1}
		tttmp1=neg - diag[i+2]
		ttttmp1=2*diag[i+1]-neg-i
# print "deg="i,"neg="neg,"diag(i)="diag[i+1],"diag(i+1)="diag[i+2],tttmp1,ttttmp1
		if(ttttmp1<0)
			ttttmp1=0
		bettigenstr = bettigenstr" "tttmp1 + ttttmp1
		tttmp2=(diag[i+2]-diag[i+1])-(diag[i+1]-diag[i])+tttmp1 + ttttmp1
		bettisyzstr = bettisyzstr" "tttmp2}
print "Lower bounds for the graded Betti numbers:"
print bettigenstr
print bettisyzstr
}
		
}



function gg(t,i,n,       a, b) 
{
a=t-(n-i)-x[i]+1
b=t-(n-i)-x[i]
if(a<0)
	a=0
if(b<0)
	b=0
return a-b
}

function dd(t,i,n,         a) 
{
if(t<i-1)
	a=0
if(t>=i-1)
	a=1
if(t>=x[n-i+1]+i-1)
	a=0
return a
}

# This next function call is no longer used
function AA(t,n,  s,i,p,q)
{
A=(t-n+2)*(t-n+1)/2
if(t-n+2<2)
	A=0
B=0
# print t, 0, A, B
# if(t==10)
#	print i" .. "A, B" .. "t-x[1]-(n-1)" t="t" ci="x[1]
p=0
for(i=1;i<=n;i++) {
s=t-x[i]-(n-i)+1-B
p=p+x[i]
q=(t-n+i+2)*(t-n+i+1)/2
if(t-n+i<2)
	q=0
q=q-p
if(s<0)
	s=0
A=A+s
if(A<q)
	A=q
B=A-q
# print t, i, A, B
# if(t==10)
#	print i" .. "A, B" .. "t-x[i]-(n-i)" t="t" ci="x[i]
}
return A
}

function AAA(t,n,  s,i,p,q,A,B)
{
A=(t-n+2)*(t-n+1)/2
if(t-n+2<2)
	A=0
B=0
p=0
for(i=1;i<=n;i++) {
s=t-x[i]-(n-i)+1-B
p=p+x[i]
q=(t-n+i+2)*(t-n+i+1)/2
if(t-n+i<0)
	q=0
q=q-p
if(s<0)
	s=0
A=A+s
B=A-q
# print "  t="t," i="i, " l="A, " l-l'="q, " l'="B, " r="t-n+i-x[i]+1

}
return A
}