/* vim: set syntax=magma : */

/*

    Magma implementation for
        "Polarizations of abelian varieties over finite fields via canonical liftings", 
         by Jonas Bergström, Valentijn Karemaker and Stefano Marseglia.

    Please send bug reports and comments to 
        Stefano Marseglia, stefano.marseglia89@gmail.com

    Here we include the code to reporduce the examples.
    Last tested on Magma V2.25-8

*/

/////////////////////////////////
// Almost Ordinary Example 
// (not in the paper)
/////////////////////////////////

    AttachSpec("~/packages_github/AbVarFq/packages.spec");
    Attach("~/packages_github/PolsAbVarFpCanLift/ResRefCond.m");
    load "~/packages_github/PolsAbVarFpCanLift/AlmOrd.m"; 

    _<x>:=PolynomialRing(Integers());
    f:=x^6 - 11*x^5 + 65*x^4 - 240*x^3 + 585*x^2 - 891*x + 729;
    If:=IsogenyClass(f:Check:=false);
    If;
    q:=FiniteField(If);
    R:=ZFVOrder(If);
    L:=Algebra(R);
    oo:=FindOverOrders(R);
    oo_max_at_ss,Lss:=overorders_maximal_at_ss(If);
    printf "In total there are %o overorder of Z[F,V], but only %o are maximal at the supersingular part 
            and hence are Endomrophism Rings of some abelian variety in the isogene class\n",#oo,#oo_max_at_ss;
    O:=MaximalOrder(L);

    iso:=[ ComputeIsomorphismClassesWithEndomorphismRing(If,S) : S in oo_max_at_ss ];
    printf "There are %o isomorphism classes of AVs\n",&+[ #I : I in iso];
    "The isogeny class is ramified",IsRamified(Lss); 
    "Hence there is a 1-1 map between AVs and ideals";

    all_cm:=AllCMTypes(If);
    CM_ST:=[ PHI : PHI in all_cm | ShimuraTaniyama(If,PHI)];
    assert #CM_ST eq 2;
    "There are 2 CM-types that satisfy the Shimura Tanyiama formula. Both will give the same output (in terms of polarizations).";
    PHI:=CM_ST[1];
    for isoS in iso do
        "new endomorphism ring";
        for I in isoS do
            "\tnew isomorphism class of AVs...";
            _,pps:=IsPrincipallyPolarized(I,PHI);
            printf "\t... with %o principal polarizations, up to isomorphism\n",#pps;
        end for;
    end for;

/////////////////////////////////
// Example 5.10
/////////////////////////////////

    AttachSpec("~/packages_github/AbVarFq/packages.spec");
    Attach("~/packages_github/PolsAbVarFpCanLift/ResRefCond.m");
    load "~/packages_github/PolsAbVarFpCanLift/AlmOrd.m"; 

    _<x>:=PolynomialRing(Integers());
    f:=(x^2 + 9)*(x^4 - 6*x^3 + 19*x^2 - 54*x + 81);
    If:=IsogenyClass(f:Check:=false);
    If;
    q:=FiniteField(If);
    R:=ZFVOrder(If);
    L:=Algebra(R);
    oo:=FindOverOrders(R);
    oo_max_at_ss,Lss:=overorders_maximal_at_ss(If);
    printf "In total there are %o overorder of Z[F,V], but only %o are maximal at the supersingular part 
            and hence are Endomrophism Rings of some abelian variety in the isogene class\n",#oo,#oo_max_at_ss;
    O:=MaximalOrder(L);

    iso:=[ ComputeIsomorphismClassesWithEndomorphismRing(If,S) : S in oo_max_at_ss ];
    printf "There are %o isomorphism classes of AVs\n",&+[ #I : I in iso];
    "The isogeny class not is ramified : ",IsRamified(Lss); 
    "Hence there is a 2-1 map between AVs and ideals";

    all_cm:=AllCMTypes(If);
    CM_ST:=[ PHI : PHI in all_cm | ShimuraTaniyama(If,PHI)];
    assert #CM_ST eq 2;
    "There are 2 CM-types that satisfy the Shimura Tanyiama formula.\n
    Each fractional ideal give rise to two AVs, one per CM-Type.\n
    For each of these two CM Types and each ideal, we compute the principal polarizations.";
    for isoS in iso do
        "new endomorphism ring";
        for I in isoS do
            "\tnew ideal";
            for PHI in CM_ST do
                "\t\tnew isomorphism class ...";
                _,pps:=IsPrincipallyPolarized(I,PHI);
                printf "\t\t... with %o principal polarizations, up to isomorphism\n",#pps;
            end for;
        end for;
    end for;

/////////////////////////////////
// Example 5.11
/////////////////////////////////

    AttachSpec("~/packages_github/AbVarFq/packages.spec");
    Attach("~/packages_github/PolsAbVarFpCanLift/ResRefCond.m");
   
    _<x>:=PolynomialRing(Integers());
    h1:=x^2+3;
    h2:=x^4+9;
    AVh1:=IsogenyClass(h1);
    "CM types of AVh1 satisfying RRC (sorted as in AllCMTypes)=\n",
            [ PHI in ResidualReflexCondition(AVh1) : PHI in AllCMTypes(AVh1) ];
    AVh2:=IsogenyClass(h2);
    "CM types of AVh2 satisfying RRC (sorted as in AllCMTypes)=\n",
            [ PHI in ResidualReflexCondition(AVh2) : PHI in AllCMTypes(AVh2) ];
    AVh:=IsogenyClass(h1*h2);
    "CM types of AVh satisfying RRC (sorted as in AllCMTypes)=\n",
            [ PHI in ResidualReflexCondition(AVh) : PHI in AllCMTypes(AVh) ];
    "One can observe that the Shimura-Taniyama formula is satisfied by all cm-types of AVh1, AVH2 and AVh ... ";
            [ ShimuraTaniyama(AVh1,PHI) : PHI in AllCMTypes(AVh1) ];
            [ ShimuraTaniyama(AVh2,PHI) : PHI in AllCMTypes(AVh2) ];
            [ ShimuraTaniyama(AVh,PHI) : PHI in AllCMTypes(AVh) ];
    "... so the failure of RRC is in the embeddability of the residue field of the reflex field.";
            [ IsResidueReflexFieldEmbeddable(AVh1,PHI) : PHI in AllCMTypes(AVh1) ];
            [ IsResidueReflexFieldEmbeddable(AVh2,PHI) : PHI in AllCMTypes(AVh2) ];
            [ IsResidueReflexFieldEmbeddable(AVh,PHI) : PHI in AllCMTypes(AVh) ];

/////////////////////////////////
// Example 5.12
/////////////////////////////////

    AttachSpec("~/packages_github/AbVarFq/packages.spec");
    Attach("~/packages_github/PolsAbVarFpCanLift/ResRefCond.m");
   
    _<x>:=PolynomialRing(Integers());
    f:=x^6 - x^5 + 6*x^4 - 6*x^3 + 18*x^2 - 9*x + 27;
    If:=IsogenyClass(f);
    If;
    q:=FiniteField(If);
    pRank(If);
    R:=ZFVOrder(If);
    L:=Algebra(R);
    oo:=FindOverOrders(R);
    O:=MaximalOrder(L);
    "there is only one non-maximal order of index coprime to p"; 
    S:=[T : T in oo | IsCoprime(Index(O,T),q) ][1];
    US,uS:=UnitGroup2(S);
    USr:=TotallyRealUnitGroup(S);
    gens_USr:=[USr.i@uS : i in [1..Ngens(USr)]]; 
    for T in oo do
        "index =",Index(O,T);
        "IsGorenstein =",IsGorenstein(T);
        "# classes with End=T =",#PicardGroup(T);
        UT,uT:=UnitGroup2(T);
        UTr:=TotallyRealUnitGroup(T); 
        "does Theorem 5.5.1 hold?",USr subset UTr;
        Tau:=[ u@uT : u in Transversal(UT,sub<UT|[ UT.i*(ComplexConjugate(UT.i@uT))@@uT : i in [1..Ngens(UT)] ]>)];
        "does Theorem 5.5.2 hold?",forall{ g : g in gens_USr | exists{ u : u in Tau | IsTotallyRealPositive(g/u) }};
        "\n";
    end for;
    all_cm:=AllCMTypes(If);
    "number of CM Types satisfying RRC=",#ResidualReflexCondition(If);
    "CM types satisfying RRC (sorted as in AllCMTypes)=\n",[ PHI in ResidualReflexCondition(If) : PHI in all_cm ];
    for T in oo do
        "\nindex of End =",Index(O,T);
        for A in ComputeIsomorphismClassesWithEndomorphismRing(If,T) do
            "new isomorphism class";
            printf "princ. pols. for each CM-type : ";
            for PHI in all_cm do
                _,pp:=IsPrincipallyPolarized(A,PHI);
                printf "%o ",#pp;
            end for;
            printf "\n";
        end for;
    end for;

/////////////////////////////////
// Example 5.13
/////////////////////////////////

    AttachSpec("~/packages_github/AbVarFq/packages.spec");
    Attach("~/packages_github/PolsAbVarFpCanLift/ResRefCond.m");
   
    _<x>:=PolynomialRing(Integers());
    h:=x^8+16;
    Ih:=IsogenyClass(h);
    Ih;
    q:=FiniteField(Ih);
    pRank(Ih);
    R:=ZFVOrder(Ih);
    L:=Algebra(R);
    O:=MaximalOrder(L);
    oo:=FindOverOrders(R);
    "there are 19 endomorphism rings, all with index in the maximal order a power of 2";
    [ Index(O,S) : S in oo ];
    all_cm:=AllCMTypes(Ih);
    "all cm-types satisfy RRC:";
    forall{ PHI : PHI in all_cm | ResidualReflexCondition(Ih,PHI)};

    "by looking at the indices we deduce the only order in \mathcal{S} is O, that is, we can only ensure that we can canonically lift an abelin variety with maximal endomorphism ring.";
    S:=O;
    US,uS:=UnitGroup2(S);
    USr:=TotallyRealUnitGroup(S);
    gens_USr:=[USr.i@uS : i in [1..Ngens(USr)]]; 
    for T in Reverse(Setseq(oo)) do
        UT,uT:=UnitGroup2(T);
        conj:=T eq ComplexConjugate(T);
        if conj then 
            Tau:=[ u@uT : u in Transversal(UT,sub<UT|[ UT.i*(ComplexConjugate(UT.i@uT))@@uT : i in [1..Ngens(UT)] ]>)];
            test_cond2:=forall{ g : g in gens_USr | exists{ u : u in Tau | IsTotallyRealPositive(g/u) }};
        else
            test_cond2:="N/A";
        end if;
        printf "\nindex = %o\tconj stable =%o\tThm 5.5.2 =%o",Index(O,T),conj,test_cond2;
    end for;

    "we verify that each overorder T is generated equals ZZ+(T:O) and that the conductor (T:O) is a principal O-ideal"
    for T in Reverse(Setseq(oo)) do 
        t1:=OneIdeal(T) eq ideal<T|ZBasis(Conductor(T)) cat [1]>;
        t2,g:=IsPrincipal(O!Conductor(T));
        printf "\nindex = %o\tT=ZZ+(T:O) ? %o\tis principal (T:O) ? %o",Index(O,T),t1,t2;
    end for;

    "we compute the isomorphism classes for each over-order T: this computation takes quite a while ...";
    for T in Reverse(Setseq(oo)) do
        "\nindex of End =",Index(O,T);
        for A in ComputeIsomorphismClassesWithEndomorphismRing(Ih,T) do
            "new isomorphism class";
            printf "princ. pols. for each CM-type : ";
            for PHI in all_cm do
                _,pp:=IsPrincipallyPolarized(A,PHI);
                printf "%o ",#pp;
            end for;
            printf "\n";
        end for;
    end for;