\\ ---------------  GP code  ---------------------------------------
\\
\\ Time-stamp: <Wed, Jan 15, 2003 - 11:29:31 - apacetti>
\\
\\ Description: Routines for operating with quadratic forms.
\\
\\ File: /home/villegas/gp/qforms.gp
\\
\\ Original author:	  Ariel Pacetti
\\			  apacetti@math.utexas.edu
\\			  University of Texas at Austin
\\
\\ Created:             Wed Feb  19  2000
\\

read(polynomials)
read(vectors)

\\======================================================================
\\ Determines is a symmetric form is definite positive.

matisposdef(A)=
	{local(s);
	if(matsize(A)[1]!=matsize(A)[2],error("not a square matrix"));
	if(A~!=A,error("not a symmetric matrix"));
	qfsign(A)[1]==matsize(A)[1]
}

\\======================================================================
\\ Determines if the symmetrix matrix A is integral

matisintegral(A)=
	{local(aux);
	if(A~!=A,error("not a symmetric matrix"));
	if(abs(content(A))>=1,if(prod(k=1,matsize(A)[1],(A[k,k]%2)+1)==1,1,0),0)
}

\\======================================================================
\\ Computes the level of a quadratic form Q

qflevel(Q)=
	{local(aux);
	if(matisposdef(A),,error("not definite positive"));
	aux=abs(content(A^(-1)));
	if(aux<1,aux=aux^(-1),1);
	if(matisintegral(aux*A^(-1)),,aux=aux*2);
	aux
}
\\======================================================================
\\ Given a form, computes it's adjoint form

qfadjoint(Q)=qflevel(Q)*Q^(-1)

\\======================================================================
\\ Given a quadratic polynomial in several variables mod q, and a
\\ prime p, diagonalize the quadratic form mod p (if p!=2).

qdiagmod(P,q)=
	{local(aux,aux2,ans,n,var);
	P=polcoefmod(P,q);
	var=variables(P);
	n=length(var);
	ans=P;
	if(q==2,error("case p=2 not implemented yet"));
	for(k=1,n,
		aux=polvar(ans,var[k]);
		if(poldegree(aux,var[k])>=1,
			aux2=Mod(pollead(aux,var[k]),q);
			aux=polcoefmod(aux*lift(aux2^(-1)),q);
			ans=simplify(subst(ans,var[k],2*var[k]-deriv(aux,var[k])*lift(Mod(2,q)^(-1))));
			ans=polcoefmod(ans,q)));
	ans
}

\\======================================================================
\\ Given a quadratic form, answers if there's some nontrivial solution 
\\ modulo q (a prime different from 2) of Q(x_0)=0, where r is some
\\ integer. Q(x) is given by a polynomial, and need not to be in
\\ reduced form. 

qhassolution(pol,q)=
	{local(aux,var,aux2);
	var=variables(pol);
	if(q==2,
		while(polcoefmod(pol,2)==0,
			pol=pol/2);
		for(k=0,8^(length(var))-1,
		aux=vectorkill(base(k+8^(length(var)),8),1);
		aux2=pol;
		for(j=1,length(var),aux2=subst(aux2,var[j],aux[j]));
		if(aux2%8==0,return(1)));
		return(-1));
	while(polcoefmod(pol,q)==0,pol=pol/q);
	aux=simplify(qdiagmod(pol,q));
	var=variables(aux);
	if(length(var)>=3,return(1));
	if(length(var)==1,return(kronecker(-subst(aux,var[1],0)*polleadint(aux),q)));
	if(length(var)==0,return(-1));
	1
	}

{helpqforms()=
print("-------------------------------------------------------------------");
print("      		Computational number Theory                        "); 
print("-------------------------------------------------------------------");
print("");
print("Description: Routines for working with quadratic forms in
several variables");
print("");
print("Original Authors:     Fernando Rodriguez-Villegas");
print("                     villegas@math.utexas.edu");
print("                     University of Texas at Austin");
print("");
print("                      Ariel Pacetti");
print("                      apacetti@math.utexas.edu");
print("                      University of Texas at Austin");
print("-------------------------------------------------------------------");
print("");
print("Help functions");
print("");
print("helpmatisposdef	helpmatisintegral	helpqflevel	helpqfadjoint");
print("helpqdiagmod	helpqhassolution");
print("");
}

addhelp(matisposdef,"Given a symmetric form, answers if the quadratic form associated to it is positive definite")

addhelp(matisintegral,"Given a symmetric integer matrix, answers if the quadratic form associated to it is integral")

addhelp(qflevel,"Given a symmetric matrix, computes the level of the quadratic forma ssociated to it")

addhelp(qfadjoint,"Given a symmetric matrix, answers the matrix of the adjoint quadratic form assocaited to it")

addhelp(qdiagmod,"Given a quadratic form, and a prime p, diagonalizes the form mod p. Observe that the prime must be different from 2")

addhelp(qhassolution,"Given a quadratic form, and a prime number, answers if the form has solution module the prime number.")