
(* Study uniqueness of Salem module for f8 *)
<<cox.m

(*-------------------Realize lambda_8 -------------------*)
salem = salem8
real = CxToReal[salem]

p = 3
(*-------------------Salem Part-------------------*)
qs = QEuler[salem];
fs = Mat[salem];

polys=  (2 + x)^2*(-2*x - x^2 + x^3);
qs = QTwist[qs,fs,polys];

(*  This unit can also be applied.... is the result really different?? *)
u = FixSig[qs,fs,1][[2]]

unit = PofM[u,fs + Inverse[fs]]

qu = QTwist[qs,fs,u]

Anti[n_]:=Block[{i,j},Table[If[i+j==n+1,1,0],{i,1,n},{j,1,n}]]

Gal[f_]:=Block[{n},
		n = Length[f];
		MatrixPower[f,1-n].Anti[n]
]

(*-------------------Check other twists for positivity -------------------*)
qe = QEuler[salem];
fe = Mat[salem];

(* {p1,p2} =  PrimeTwists[qe,fe,3] *)
{p1,p2} =  {-1 - 3*x + x^2 + x^3, -4*x + x^3}

(* This one has much too big glue group *)
q1 = QTwist[qe,fe,p1] 
(* This one s Z/3 x Z/9 *)
q2 = QTwist[qe,fe,p2] 

{s1,s2} = FixSig[qs,fs,p2]
q1 = QTwist[qe,fe,s1] 
q2 = QTwist[qe,fe,s2] 

(* Now look for obstructions! *)
Print["Obstructions for q1:  ",Length[CubeOb[q1,fs]]]
Print["Obstructions for q2:  ",Length[CubeOb[q2,fs]]]
Print["So no extension of (fs,qs) will work!"]

(*-----------Study twists by units-----------------*)
as = fs + Inverse[fs]
 
q8 = {{-2, 1, 1, 0, 0, 0, 1, 0}, {1, -2, 1, 1, 0, 0, 0, 1}, 
     {1, 1, -2, 1, 1, 0, 0, 0}, {0, 1, 1, -2, 1, 1, 0, 0}, 
     {0, 0, 1, 1, -2, 1, 1, 0}, {0, 0, 0, 1, 1, -2, 1, 1}, 
     {1, 0, 0, 0, 1, 1, -2, 1}, {0, 1, 0, 0, 0, 1, 1, -2}}
 
f8 = {{0, 0, 0, 0, 0, 0, 0, -1}, {1, 0, 0, 0, 0, 0, 0, 0}, 
     {0, 1, 0, 0, 0, 0, 0, 0}, {0, 0, 1, 0, 0, 0, 0, 1}, 
     {0, 0, 0, 1, 0, 0, 0, 1}, {0, 0, 0, 0, 1, 0, 0, 1}, 
     {0, 0, 0, 0, 0, 1, 0, 0}, {0, 0, 0, 0, 0, 0, 1, 0}}

Units[r_,n_]:=Block[{gs,i,us},
		us=Units[r];
		gs=Join[us,{1}];
		us = gs; 
		For[i=0, i<n, i++,
			us = Flatten[Outer[Times,us,gs],1];
			us = Map[PolynomialRemainder[#,r,x]&,us];
			us = Union[us]
		];
		Return[us]
]
		

(*
depth = 25
ru = Units[real,depth];
Print["Depth ",depth,"  Number of units: ",Length[ru]]

list = Map[QTwist[qs,fs,#]&,ru];
Print[{Position[list,q8],Position[list,qs],Position[list,qu]}]
*)

Size[q_]:=Plus @@ (Flatten[q]^2)
Lower[q_,f_]:=Block[{i,r,ss,qs,us},
		r = CxToReal[CharPoly[f,x]];
		us = Units[r,20];
		qs = Map[QTwist[q,f,#]&,us];
		ss = Map[Size,qs];
		i = Position[ss,Min[ss]][[1,1]];
		Return[qs[[i]]]
]

