
<<cox.m

(*------------------Salem 12.7 -- glue group 7x49, still not realizable -----------------*)
prime = 7
Print[" 			SALEM 12:  Twisting at ",prime]

(*--------------------------Ideal Theory ---------------*)
SetOptions[GroebnerBasis,CoefficientDomain->Integers];
GB[x_]:=(Print[""]; Print["GB for ",x,": ",GroebnerBasis[x]];GroebnerBasis[x])

Print[""]
Print["         Ideal Theory"]

basis1 = GB[{salem12, (x+1)^3}]
basis2 = GB[{u^7-1,2u^2+3u+2}]
basis3 = GB[{t^2+7,t^3}]

(* We have A1 = A3 under t = 3x+3 *)
GB[basis3/.t->3x+3];Print[" = ",basis1]

(* We have A1 = A2 under u = -x *)
GB[basis2/.u->-x];Print[" = ",basis1]

(* ------------------- Construct Isometries ---------------------- *)
(* Find matrix conjugating action on Z/7 x Z/49 to standard form *)
(* phi^-1 h phi = {{0,1},{-7,0}} and m = +/-{{7,0},{0,1}} *)

FindGlue[m_,g_]:=Block[{a,b,e1,e2,mod,phi,h},
		mod = {7,49};
		h = 3 g + 3 Id[g];
		If[Max[Mod[h.h+7Id[h],mod]] > 0,	
			Print[" * * * Incorrect module * * *"]];
		e2 = {0,1};
		e1 = h.e2;
		phi = {e1,e2};
		phi = Transpose[phi];
		a = 49 (Transpose[phi].m.phi)[[2,2]];
		b=Select[Range[0,48],Mod[#*#*a-1,49]==0 || Mod[#*#*a+1,49]==0 &][[1]];
		phi = Mod[b phi,mod];
		Return[phi]
]

(* Test action mod 49 *)
TestGlue[q_,f_]:=Block[{g,m,phi},
		g = GlueF[q,f,7];
		m = GlueQ[q,7];
		phi = FindGlue[m,g];
		m = Transpose[phi].m.phi;
		m = Mod[m,1];
		Return[49 m[[2,2]]]
]

(*-------------------Find viable extra twists -------------------*)
Extras[q_,f_]:=Block[{a,ps},
        a = f + Inverse[f];
        ps = (IntegerCube[30,4]-15).{1,x,x^2,x^3};
        ps = Select[ps,0<Abs[Det[PofM[#,a]]] <= 1000&];
        ps = Select[ps,Sig[QTwist[q,f,#]][[2]]<=2&];
        ps = Select[ps,EvenForm[QTwist[q,f,#]]&];
        ps = Select[ps,TestGlue[QTwist[q,f,#],f]!=1&];
        Return[ps]
]


(*--------------------------Salem part---------------*)
salem = salem12
qs = QEuler[salem]
fs = Mat[salem]
(* ps = AllTwists[salem] *)
twists = {-2 - x + 4*x^2 - x^4, -3*x + x^3, -1 + 5*x - 5*x^3 + x^5}

polys = twists[[1]]
prime = TwistN[fs,polys]

polys = FixSig[qs,fs,polys]
ps = polys[[2]]
Print["Fixing signature; twist is ",ps] 
qs = QTwist[qs,fs,ps]
charfs = GlueChar[qs,fs,prime]

ms = GlueQ[qs,prime];
gs = GlueF[qs,fs,prime];
phis = FindGlue[ms,gs];
vss = GlueVectors[qs,prime];

Print["\n 			Salem Factor"]
DiGlue[qs,fs]

(*--------------------------Cyclotomic part ---------------*)
fc = Mat[x^7-1]
qc = Id[fc]
polyc = 2x + 3
(* This is the same as 2f^2 + 3f + 2 *)
qc = QTwist[qc,fc,polyc]

Print["\n 			Cyclotomic Factor \n"]
DiGlue[qc,fc]

Print["For compatible gluing these should be opposite:  ",
	TestGlue[qs,fs]," != ",TestGlue[qc,-fc]]

(*--------------------------Extras ---------------*)
fc = Mat[x^7-1]
qc = Id[fc]
polyc = 2x + 3
(* This is the same as 2f^2 + 3f + 2 *)
qc = QTwist[qc,fc,polyc]
DiGlue[qc,fc]

(* Marvellously, this already has signature (5,2)! *)
(* Terribly, the form is not reversed *)
Print["-----------------------------EXTRAS--------------------"]

TryExtra[extra_]:=Block[{ans,c1,c2,f0,fc,fi,gc,mc,q0,qc,qi,s1,s2,vsc,vss,v1,v2},
	fc = Mat[x^7-1];
	qc = Id[fc];
	polyc = 2x + 3;
	qc = QTwist[qc,fc,polyc];
	qc = QTwist[qc,-fc,extra];
	mc = GlueQ[qc,prime];
	gc = GlueF[qc,fc,prime];
	phic = FindGlue[mc,-gc];
	vss = GlueVectors[qs,prime];
	vsc = GlueVectors[qc,prime];
	q0 = DirSum[qs,qc];
	f0 = DirSum[fs,-fc];
	{s1,s2} = Transpose[phis].vss;
	{c1,c2} = Transpose[phic].vsc;
	v1 = Join[s1,c1];
	v2 = Join[s2,c2];
	{qi,fi} = Glue[q0,f0,{v1,v2}];
	Print["Extra: ",extra," Sig ",Sig[qi]," Glue ",AbGp[qi]," Positive? ",ans=PosQ[qi,fi]];
	Return[ans]
]

(* extras = Extras[qc,-fc] *)
<<f12.7.extras.m
Print["Number of extras is ",Length[extras]]
i = Random[Integer,Length[extras]]
Print["We will sample a few"]
extras = Take[extras,2] 
ans = Map[TryExtra,extras]
Print["Some extra worked?? ",Or @@ ans]

(*--------------------------Cylotomic part ---------------*)
fc = Mat[x^7-1]
qc = Id[fc]
polyc = 2x + 3
(* This is the same as 2f^2 + 3f + 2 *)
qc = QTwist[qc,fc,polyc]

DiGlue[qc,fc]

extra = 2-x^2
Print["EXTRA is ",extra]

(* -----------Gluing Salem to Cyclo------------- *)
qc = QTwist[qc,-fc,extra];
mc = GlueQ[qc,prime];
gc = GlueF[qc,fc,prime];
phic = FindGlue[mc,-gc];   (* NOTE minus sign! *)

vsc = GlueVectors[qc,prime];
q0 = DirSum[qs,qc];
f0 = DirSum[fs,-fc];

{s1,s2} = Transpose[phis].vss;
{c1,c2} = Transpose[phic].vsc;
v1 = Join[s1,c1];
v2 = Join[s2,c2];
{qd,fd} = Glue[q0,f0,{v1,v2}];
DiGlue[qd,fd]

(*-------------Final jeopardy------------*)
qk =  {{2, -1, 0}, {-1, 2, -1}, {0, 1, 0}}
fk = Id[qk]

{q,f} = GlueAuto[qd,fd,qk,fk,2]

Print["\n 			Realized (but wrong signature)"]
DiGlue[q,f]
