// ================================================================
//  Exemple 6.2 (cas supersingulier) 
//  Parametres : k = F_5, E : y² + y = x³, l = 7, f = [2]
// ================================================================

// ================================================================
//  Etape 1 : Definition du corps de base et de la courbe
// ================================================================

p := 5;
k := GF(p);
k_x<X> := PolynomialRing(k);
a3 := 1;

E := EllipticCurve([k | 0,0,a3,0,0]);  // y^2 + y = x^3
printf "Courbe E : %o\n", E;
printf "Cardinalite de E(F_%o) : %o\n", p, #E;

// ================================================================
//  Etape 2 : Parametres principaux
// ================================================================

l := 7;
n := (l-1) div 2;
printf "l = %o, n = %o\n", l, n;

// ================================================================
//  Etape 3 : Calcul du polynome de division et selection du facteur psi_H
// ================================================================

psi_l := DivisionPolynomial(E, l);
printf "\nPolynome de division psi_%o (degre %o)\n", l, Degree(psi_l);

facteurs := Factorization(psi_l);
printf "\nFacteurs de psi_%o :\n", l;
for fa in facteurs do
    printf "  degre %o, multiplicite %o\n", Degree(fa[1]), fa[2];
end for;

// Selection du facteur de degre n dont la trace est non nulle
psiH := k_x!0;
selected := false;
for fa in facteurs do
    pol := fa[1];
    if Degree(pol) eq n and not selected then
        assert IsIrreducible(pol);
        Kc<wc> := SplittingField(pol);
        Kc_X<Xc> := PolynomialRing(Kc);
        racines_c := [rt[1] : rt in Roots(Kc_X!pol)];
        trace_c := &+racines_c;
        printf "\nTest du facteur de degre %o :\n", n;
        printf "  Trace = %o\n", trace_c;
        if trace_c ne Kc!0 then
            psiH := pol;
            selected := true;
            printf "  -> Facteur retenu (trace non nulle)\n";
        else
            printf "  -> Rejete (trace nulle)\n";
        end if;
    end if;
end for;

if psiH eq 0 then
    error "Aucun facteur de degre", n, "ne verifie la condition de trace non nulle.";
end if;

printf "\npsi_H(X) = %o\n", psiH;

// ================================================================
//  Etape 4 : Corps de decomposition et verification de la base normale
// ================================================================

printf "\n=== Verification de la condition de base normale ===\n";
K<w> := SplittingField(psiH);
K_X<XK> := PolynomialRing(K);
printf "Corps K = F_%o^%o\n", p, Degree(K, k);

racines_psiH := [rt[1] : rt in Roots(K_X!psiH)];
assert #racines_psiH eq n;

// Construction de la matrice de Hankel pour tester la base normale
mat := Matrix(K, n, n, []);
for i in [1..n] do
    for j in [1..n] do
        mat[i][j] := racines_psiH[i]^(p^(j-1));
    end for;
end for;

det := Determinant(mat);
printf "Determinant de la matrice de Hankel : %o\n", det;
if det ne 0 then
    printf "-> Les racines forment une base de K/k (base normale)\n";
else
    printf "-> Attention : les racines ne forment PAS une base normale\n";
end if;

// ================================================================
//  Etape 5 : Relevement d'un point P d'ordre l
// ================================================================

printf "\n=== Relevement d'un point P d'ordre %o ===\n", l;

x1 := racines_psiH[1];
rhs := x1^3;  // car y^2 + y = x^3

Y := PolynomialRing(K).1;
poly_y := Y^2 + Y - rhs;
racines_y := Roots(poly_y);

if #racines_y ge 1 then
    y1 := racines_y[1][1];
    L := K;
    printf "y1 = %o\n", y1;
else
    printf "Extension necessaire pour la racine carree\n";
    L<yext> := ext<K | poly_y>;
    y1 := yext;
end if;

E_L := BaseChange(E, L);
P := E_L ! [L!x1, L!y1];
printf "P = (%o, %o)\n", P[1], P[2];

ord_P := Order(P);
printf "Ordre(P) = %o\n", ord_P;
assert ord_P eq l;

// ================================================================
//  Etape 6 : Construction explicite de la base (x_i) = x(iP)
// ================================================================

printf "\n=== Construction de la base (x_i) ===\n";
xcoord := function(Q) return Q[1]; end function;

xs := [xcoord(i*P) : i in [1..n]];
printf "x_i = x(iP) :\n";
for i in [1..n] do
    printf "  x_%o = %o\n", i, xs[i];
end for;

// Verification que chaque x_i est bien une racine de psi_H
for i in [1..n] do
    if not (xs[i] in racines_psiH) then
        error "Erreur : x_", i, "n'est pas une racine de psi_H";
    end if;
end for;
printf "Verification : les %o x_i sont bien des racines de psi_H\n", n;

// ================================================================
//  Etape 7 : Action de l'endomorphisme f = [2]
// ================================================================

printf "\n=== Action de f = [2] ===\n";
m := 2;

images := [xcoord(m*(i*P)) : i in [1..n]];
printf "Images directes :\n";
for i in [1..n] do
    printf "  f(x_%o) = x(%o*P) = %o\n", i, m*i mod l, images[i];
end for;

// Calcul de la permutation induite sur les indices
perm := [];
for i in [1..n] do
    cible := (m*i) mod l;
    if cible gt l - cible then
        cible := l - cible;
    end if;
    Append(~perm, cible);
end for;

printf "\nPermutation induite sur les indices :\n";
for i in [1..n] do
    printf "  %o -> %o\n", i, perm[i];
end for;

images_theoriques := [xs[perm[i]] : i in [1..n]];
if images eq images_theoriques then
    printf "\nVerification des images : OK\n";
else
    printf "\nAttention : les images ne correspondent pas\n";
end if;

sigma_perm := Sym(n) ! perm;
printf "Ordre de sigma_f : %o\n", Order(sigma_perm);
printf "Decomposition en cycles : %o\n", CycleDecomposition(sigma_perm);

// ================================================================
//  Etape 8 : Expression explicite de sigma_f
// ================================================================

printf "\n=== sigma_f explicite ===\n";
printf "Pour tout z = sum_{i=1}^{%o} a_i x_i :\n", n;
printf "sigma_f(z) = ";
for i in [1..n] do
    if i gt 1 then
        printf " + ";
    end if;
    printf "a_%o * x_%o", i, perm[i];
end for;
printf "\n";

// ================================================================
//  Etape 9 : Verification que sigma_f est un automorphisme de K/k
// ================================================================

printf "\n=== Verification que sigma_f est un automorphisme de K/k ===\n";

base_xs := [K | xs[i] : i in [1..n]];

// Matrice de passage de la base canonique vers la base (base_xs)
M := Matrix(K, n, n, []);
for i in [1..n] do
    coords := Eltseq(base_xs[i]);
    for j in [1..n] do
        M[j][i] := K!coords[j];
    end for;
end for;

detM := Determinant(M);
printf "Determinant de la matrice de passage M : %o\n", detM;
assert detM ne 0;

M_inv := M^-1;

// Fonction auxiliaire : coordonnees d'un element dans la base (x_i)
function coords_in_xs(z)
    coords := Eltseq(z);
    vec_z := Matrix(K, n, 1, coords);
    c_mat := M_inv * vec_z;
    return [c_mat[i][1] : i in [1..n]];
end function;

// Application de sigma_f a un element quelconque
function sigma_f_element(z)
    c := coords_in_xs(z);
    res := K!0;
    for i in [1..n] do
        res +:= c[i] * base_xs[perm[i]];
    end for;
    return res;
end function;

// Test de la preservation de la multiplication sur les elements de la base
printf "Test de la multiplication sur les elements de base :\n";
base_mul_ok := true;
for i in [1..n] do
    for j in [1..n] do
        a := base_xs[i];
        b := base_xs[j];
        if sigma_f_element(a*b) ne sigma_f_element(a)*sigma_f_element(b) then
            printf "  Echec pour i=%o, j=%o\n", i, j;
            base_mul_ok := false;
        end if;
    end for;
end for;
printf "Preservation de la multiplication sur la base : %o\n", base_mul_ok;

// Tests sur des elements aleatoires du corps
printf "\nTest sur des elements aleatoires :\n";
all_add_ok := true;
all_mul_ok := true;
for test in [1..10] do
    a := Random(K);
    b := Random(K);
    add_ok := sigma_f_element(a + b) eq sigma_f_element(a) + sigma_f_element(b);
    mul_ok := sigma_f_element(a * b) eq sigma_f_element(a) * sigma_f_element(b);
    printf "Test %o : add=%o, mul=%o\n", test, add_ok, mul_ok;
    if not add_ok then all_add_ok := false; end if;
    if not mul_ok then all_mul_ok := false; end if;
end for;
printf "Preservation de l'addition : %o\n", all_add_ok;
printf "Preservation de la multiplication : %o\n", all_mul_ok;
printf "sigma_f est un automorphisme de corps : %o\n", all_add_ok and all_mul_ok and base_mul_ok;

// ================================================================
//  Etape 10 : Identification de sigma_f avec le Frobenius (calcul direct)
// ================================================================

printf "\n=== Identification avec le Frobenius ===\n";

// Calcul de l'action du Frobenius pi_p sur la base
frob_indices := [];
for i in [1..n] do
    frob_xi := xs[i]^p;
    idx := Position(xs, frob_xi);
    if idx eq 0 then
        error "Erreur : l'image par le Frobenius n'est pas dans la base";
    end if;
    Append(~frob_indices, idx);
end for;
frob_perm := Sym(n) ! frob_indices;

printf "Frobenius pi_%o agit sur les racines par : %o\n", p, frob_perm;
printf "Ordre du Frobenius : %o\n", Order(frob_perm);
printf "Ordre attendu (degre de l'extension) : %o\n", n;

is_cyclic := Order(frob_perm) eq n;
printf "Gal(K/k) est cyclique d'ordre %o : %o\n", n, is_cyclic;

// Recherche de la puissance du Frobenius qui coincide avec sigma_f
found := false;
for k in [0..n-1] do
    if frob_perm^k eq sigma_perm then
        printf "sigma_f = pi_%o^%o\n", p, k;
        found := true;
        break;
    end if;
end for;
if not found then
    printf "sigma_f n'est pas une puissance du Frobenius dans ce groupe.\n";
end if;

// Affichage detaille de l'action du Frobenius
printf "\nAction du Frobenius sur les racines :\n";
for i in [1..n] do
    frob_xi := xs[i]^p;
    idx := Position(xs, frob_xi);
    printf "  pi_%o(x_%o) = x_%o^%o = %o", p, i, i, p, frob_xi;
    if idx ne 0 then
        printf " = x_%o\n", idx;
    else
        printf "\n";
    end if;
end for;

// ================================================================
//  Etape 11 : Resume final des resultats
// ================================================================

printf "\n=== RESUME FINAL ===\n";
printf "Extension K/k : F_%o^%o\n", p, n;
printf "Groupe de Galois Gal(K/k) ≃ (Z/%oZ)^*/{±1}\n", l;
printf "  -> ordre %o\n", (l-1) div 2;
printf "Frobenius (calcule directement) : pi_%o = %o\n", p, frob_perm;
printf "Element sigma_f (construit a partir de f=[2]) : %o\n", sigma_perm;
printf "sigma_f est un automorphisme de K/k : %o\n", all_add_ok and all_mul_ok and base_mul_ok;
printf "Le Frobenius engendre Gal(K/k) : %o\n", is_cyclic;

printf "\n=== FIN DE L'EXECUTION ===\n";