/* ======================================================= */
/*             Computations for CM threefolds  	     	   */
/* ======================================================= */

load "MT.m";

/* 
   Initialise the list of possible Galois
   groups for the Galois closure of a
   sextic CM field 
*/
dimension := 3;
PossibleGaloisGroups := [];

// C3 x C2
G<sigma,i> := FPGroup<sigma,i | i^2=1, sigma*i=i*sigma, sigma^3=1 >;
Append(~PossibleGaloisGroups, G);

// A4 x C2
G := DirectProduct(AlternatingGroup(4), CyclicGroup(2));
G := FPGroup(G);
Append(~PossibleGaloisGroups, G);

// S4 x C2
G := DirectProduct(SymmetricGroup(4), CyclicGroup(2));
G := FPGroup(G);
Append(~PossibleGaloisGroups, G);

// D3 x C2
G<sigma,tau,i> := FPGroup<sigma,tau,i | i^2=1, sigma*i=i*sigma, tau*i=i*tau, sigma^3=1, tau^2=1, tau*sigma=sigma^(-1)*tau >;
Append(~PossibleGaloisGroups, G);

// Start the loop over possible Galois groups
for G in PossibleGaloisGroups do
	"-------------------------";
	"Group G", IdentifyGroup(G);
	pG, phi := PermutationGroup(G);
	psi := Inverse(phi);
	pGels := [g : g in pG];
	Gels := [psi(g) : g in pGels];
	Gdata := [* G, pG, phi, psi, pGels *];

	// Also loop over possible complex conjugations,
	// namely central involutions.
	// For our groups, there is only one choice
	invos := [psi(g) : g in Center(pG) | Order(g) eq 2 ];
	for iota in invos do
		"Involution", iota;

		// List subgroups that define a sextic CM field
		Subgps := AllSubgroups(phi(G));
		// Subgps := [ H`subgroup : H in Subgps | not phi(iota) in H`subgroup ];
		Subgps := [ H : H in Subgps | #G/#H eq 2*dimension ];
		Subgps := [ H : H in Subgps | #Core(pG,H) eq 1 ];
		Subgps := [ psi(H) : H in Subgps ];


		for H1 in Subgps do
			for H2 in Subgps do
				"*************************";
				// "Subgroup H1 =", {psi(h) : h in phi(H1)};
				// "Subgroup H2 =", {psi(h) : h in phi(H2)};

				"#Normaliser H1 modulo H1:", #(Normaliser(pG, phi(H1)) / phi(H1));
				"#Normaliser H2 modulo H2:", #(Normaliser(pG, phi(H2)) / phi(H2));

				CMTypes1 := ListSimpleCMTypes(H1, Gdata, iota);
				CMTypes2 := ListSimpleCMTypes(H2, Gdata, iota);

				"#Distinct CM types", #CMTypes1, #CMTypes2;
				for i1 in [1..#CMTypes1] do
				for i2 in [1..#CMTypes2] do			// only loop over distinct CM types
					"=========================";
					"CM type 1", CMTypes1[i1];
					"CM type 2", CMTypes2[i2];
					testIsogeny := ComputeData(H1, H2, CMTypes1[i1], CMTypes2[i2], Gdata);
					if testIsogeny then			// if the projections MT(A_1 x A_2) \to MT(A_i)
										// are isogenies, we check that the subgroups
										// H_1, H_2 are conjugated (hence define the
										// same CM field, up to isomorphism) and that
										// under this isomorphism the CM types are
										// (essentially) equal
						testConj, elConj := IsConjugate(phi(G), phi(H1), phi(H2));	// check and assert that H1, H2 are conjugate
						assert testConj;

						S2conj := [ psi(elConj * phi(s)) : s in CMTypes2[i2] ];			// conjugate S2 to a CM type for the fixed field of H1

						assert CheckCMType(H1, S2conj, Gdata, iota);				// check that it is indeed a CM type for H1\G
						assert AreEssentiallyEquivalent( H1, CMTypes1[i1], S2conj, Gdata );	// check that the CM types are essentially equivalent

					end if;
				end for;
				end for;
			end for;
		end for;
	end for;

end for;
