/////////////////////////////////////////////////////////////////////////
//example2.m
//
//Nils Bruin, 11/09/2002
//
//Accompanying MAGMA script for article
//
//Nils Bruin, Visualising Sha[2] in Abelian surfaces
//
//This script verifies the computations referred to in Section 7, Example 2
//The script has been tested on Magma V2.9-12

//First load the package for computing Selmer groups
AttachSpec("Algaem/Algaem.spec");
_<x>:=PolynomialRing(Rationals());
//Oddly enough, we need the *number field* Q.
Q:=NumberField(x-1:DoLinearExtension);
K:=Q;
Kx<xK>:=PolynomialRing(K);

f:=x^3-22*x^2+21*x+1;

delta:=-11;

//The elliptic curves ...
E1:=PolToEllipticCurve(Kx!f);
E2:=PolToEllipticCurve(Kx!f*delta);

//Getting the 2-Selmer groups ...
two1:=MultiplicationByMMap(E1,2);
two2:=MultiplicationByMMap(E2,2);
mu1:=IsogenyMu(two1);
A<theta>:=Codomain(mu1);
S1small,mp1:=SelmerGroup(two1);
Aa,toAa1:=AbsoluteAlgebra(Codomain(mu1));
mu2:=IsogenyMu(two2:Fields:={Component(Aa,i):i in [1..NumberOfComponents(Aa)]});
A2:=Codomain(mu2);
A2toA:=hom<A2->A|[theta*delta]>;
S2small,mp2:=SelmerGroup(two2);

//And compute a structure in which both Selmer groups fit.
S:={LookupPrime(p[1]): p in Factorisation(2*Disc*IntegerRing(BaseRing(E2)))}
  where Disc:=Discriminant(E2);
H1,toH1:=pSelmerGroup(Codomain(mu1),2,S);

S1:=sub<H1|[(g@@mp1)@toH1:g in OrderedGenerators(S1small)]>;
S2:=sub<H1|[(g@@mp2)@A2toA@toH1:g in OrderedGenerators(S2small)]>;


Slower:=sub<H1|OrderedGenerators(S1) cat OrderedGenerators(S2)>;
Supper:=sub<H1|OrderedGenerators(S1 meet S2)>;

//Some stuff necessary to be able to compute the Norm map later on
Lrel:=NumberField(xK^2-delta);
L:=NumberField(x^2-delta);
LisLrel,LtoLrel:=IsIsomorphic(L,Lrel);
assert LisLrel;

EL:=BaseChange(E1,L);
twoL:=MultiplicationByMMap(EL,2);
muL:=IsogenyMu(twoL);
AL<thetaL>:=Codomain(muL);
AaL:=AbsoluteAlgebra(AL);
//We're lazy with the classgroup computation here. If you want, you can remove
//the "Bound:=300" and still arrive at the same (but now proven-to-be-correct)
//answer.
[ClassGroup(AaL[i]:Bound:=300):
    i in [1..NumberOfComponents(AaL)]];

SL,mpL:=SelmerGroup(twoL);

//Unfortunately, @LtoLrel@@LtoLrel isn't always the identity.
//Conjugation of L/K is easy enough for Lrel.
Lrelcon:=hom<Lrel->Lrel|[-Lrel.1]>;
//To compute what this does to L.1, we consider two options and pick the
//nontrivial one
L1con:=[a:a in [LtoLrel(L.1)@@LtoLrel,Lrelcon(LtoLrel(L.1))@@LtoLrel]| a ne L.1];
//We should really end up with only one option
assert #L1con eq 1;
//That supposedly defines conjugation
Lcon:=hom<L->L|[L1con[1]]>;
//But better check that it really does, i.e., leaves K invariant.
assert k eq Lcon(k) where k:=(Lrel!K.1)@@LtoLrel;
//And, just to be sure, that conjugation is of order 2 (it is not the identity
//by picking L1con different from L.1)
assert Lcon(Lcon(L.1)) eq L.1;

ALcon:=hom<AL->AL|a:->AL![Lcon(b):b in Eltseq(a)]>;
ALnorm:=hom<AL->A|a:->A![Eltseq(LtoLrel(b))[1]:b in Eltseq(a*ALcon(a))]>;
NS:=sub<H1|[ALnorm(g@@mpL)@toH1:g in OrderedGenerators(SL)]>;

///////////////////////////////////////////////////////////////////////////
//Given the computations above, we now have computed the objects with which we
//can check the claims in Example 2 in Section 7 of the article.

//First if S1=S^(2)(E1/Q) is what is claimed in the article
sub<H1|[toH1(g):g in [-theta,1-theta,1-2*theta,1-8*theta]]> eq S1;
   
//And that rational points on E2 generate S2.
sub<H1|[toH1(g):g in [(1-theta)*(1-2*theta)]]> eq
   sub<H1|[toH1(A2toA(mu2(P))):P in [E2![4,113]]]>;
sub<H1|[toH1(g):g in [(1-theta)*(1-2*theta)]]> eq S2;

      
//And that the Norm of the Selmer group of E1 over Q(sqrt(-11)) is equal
//to the intersection of S1 and S2 and therefore yields no new information.
NS eq S1 meet S2;
