# Matthew Russell
# Experimental Math
# Homework 7
# I give permission to post this.

#####################################################################################

# 1. Write a procedure Gsomosk(Ini,k,n): that outputs the n-th term of the
# Somos-k recurrence, where Somos-k is the natural k-th order analog
# of the Somos-4 recurrence (first think what it is).

GSomosk:=proc(Ini,k,n)
option remember:

if n<=nops(Ini) then
	return Ini[n]:
else
 	return (add(GSomosk(Ini,k,n-j)*GSomosk(Ini,k,n-k+j),j=1..trunc(k/2)))/GSomosk(Ini,k,n-k):
fi:

end:

#####################################################################################

# 2. Find out, experimentally, for which k, Somos-k, with Ini=[1$k], seem to
# always produce integers.

# Somos-1 does not make sense.
# Somos-2 and Somos-3 are the constant all-ones sequence.
# Somos-4, Somos-5, Somos-6, and Somos-7 all appear to be integer sequences.
# Somos-k, for k>=8, appear to be not be integer sequences.

#####################################################################################

# 3. Find out, experimentally, for which k, Somos-k, with Ini=[x[1], ..., x[k]],
# seem to always produce Laurent polynomials (rational functions whose denominators are monomials)

# Somos-k for k<=7 seems to always produce Laurent polynomials.
# Somos-k for k>=8 seems to always not produce Laurent polynomials.

#####################################################################################

# 4. The total degree of a monomial x1^a1 x2^a2 x3^a3 ... is a1+a2+a3+...
# a homogenous polynomial of degree d is one for which each of its terms
# is a monomial of degree d. Adapt GenPol(X,d,a,co), to write a procedure
# GenHomPol(X,d,a,co) that yields a homog. polynomial of total degree d.

GenHomPol:=proc(X,d,a,counter) local x,m,counter1,i,P1,X1,P,var:

m:=nops(X):
counter1:=counter:
x:=X[1]:

if m=1 then
	return a[counter1]*x^d,{a[counter1]},counter1+1:
else
	P:=0:
	var:={}:
	X1:=[op(2..m,X)]:
	for i from 0 to d do
		P1:=GenHomPol(X1,i,a,counter1):
		P:=expand(P+P1[1]*x^(d-i)):
		var:=var union P1[2]:
		counter1:=P1[3]:
	od:
	return P,var,counter1:
fi:

end:

#####################################################################################

# 5. Using GenHomPol(X,d,a,co), write a program GenNonHomPol(X,d,a,co) of
# total degree d (in other words one that is the sum of a homog. polynomials of total degree i for i=0...d) 

GenNonHomPol:=proc(X,d,a,counter) local i,P,var,counter1,G:

P:=0:
var:={}:
counter1:=counter:
for i from 0 to d do
	G:=GenHomPol(X,i,a,counter1):
	P:=P+G[1]:
	var:=var union G[2]:
	counter1:=G[3]:
od:

return P,var,counter1:

end:

#####################################################################################