/**** Decomposition Group of A Triangular Ideal */

/*** Introduction 
Documented function and procedures above are the \textsc{Magma} implementation
of the algorithm presented in the paper \href{http://www.springerlink.com/link.asp?id=ked7r9667gdm798v}
{Computation of the decomposition group of a triangular
ideal}.

The function Decomposition_Group computes the decomposition group of a triangular ideal $I \in k[x_1,\ldots,x_n] $. The decomposition group of $I$ is the set-wise stabiliser of $I$ for the natural group action of $S_n$ (the symmetric group of degree $n$) on $k[x_1,\ldots ,x_n]$~: \\
$$S_n \times k[x_1,\ldots,x_n] \longrightarrow  k[x_1,\ldots,x_n] \qquad\;\;$$
$$(\sigma,P(x_1,\ldots,x_n)) \longmapsto P(x_{\sigma(1)},\ldots,x_{\sigma(n)}).$$
In other words, the decomposition group is defined by~: \\
$$Dec(I) = \{ \sigma \in S_n \, \mid \, \sigma.f \in I\, \}. $$
The input of this implementation is a triangular generating set $\{f_1,\ldots,f_n\}$ of $I$ such that $f_i$ depends only on the variables $x_1,\ldots,x_i$. As this algorithm uses ideal membership tests, in this implementation, the triangular set of generators $S$ is supposed to be a lexicographical Gröbner basis for the order induced by $x_1>\ldots>x_n$.\\

The output is a strong generating set of the decomposition group of the ideal generated by S. During the computation, the algorithm tests whether a backtrack appears (in this case, the Galois ideal $I$ is not pure) and at the end of the computation return the predicate

\begin{center}
Is I a pure Galois Ideal ? 
\end{center}

Authors : S.~Orange (LIP6, Paris 6 University), 
	  G.~Renault (LIP6, Paris 6 University).

Source of this file : \href{../Decomposition_group_aaecc.dmc}{Documented Magma Code file}.
*/



/*** Implementation 
*/

/** NewOrbits
  @param orbites the set of the orbits of the group the action on $\{1,...,n\}$ 
     of a sub-group $G$ of $S_n$.
  @param sigma a permutation of $S_n$
  @return the set of orbits of the group $\langle G, sigma\rangle$.
*/
/***************************************************************************/
function NewOrbits (orbites,sigma);

nvorbites := {};
while not(IsEmpty(orbites)) do
    e := Random(orbites) ;
    p := e join Image(sigma,e) ;
    orbites := Exclude(orbites,e) ;
    while not(p eq e) do
        for oprime in orbites do
            if not(#(p meet oprime) eq 0) then
                e := e join oprime;
                orbites := Exclude(orbites,oprime) ;
            end if;
        end for;
        p := e join Image(sigma,e) ;
    end while ;
    nvorbites := Include(nvorbites,e) ;
end while;
return nvorbites;

end function;
/***************************************************************************/

/** OnePermutation
  @param l a sequence of integers representing the prefix.
  @param S a sequence of polynomials representing the triangular basis
     of an ideal.
  @param SP the sequence S modulo a compatible prime (for the modular pre-test).
  @return if there exists a permutation of prefix \texttt{l} in $Dec(<S>)$, one of 
  such a permutation is returned else $Id(S_n)$ is returned.
*/
/***************************************************************************/
function OnePermutation(l,S,SP)
local r,Sn,n,boucle1,boucle2,a,ens,s,ll;
n := #S;
Sn := Sym(n);
r := n - #l;

boucle1 := true;
while (r gt 0) and boucle1 do
    ens:=Reverse(Sort(SetToSequence({1..n} diff Set(l))));
    boucle2 := true;
    while (r gt 0) and boucle2 do
          if #ens eq 0 then //A backtrack appears
                boucle1:=false;
                boucle2:=false;
          else
          a := ens[1];
          ens := ens[2 .. #ens];
          ll := [a] cat l;
          s:= Sn!(SetToSequence({1..n} diff Set(ll)) cat ll);  
          if NormalForm((SP[r])^s ,SP) eq 0 then
                if NormalForm((S[r])^s,S) eq 0 then
                   r:=r-1;
                   boucle2 := false;
                   l := ll;
                end if;
          end if;
          end if;
    end while;
end while;

if r eq 0 then
   return s;
else //A backtrack appears
   return Id(Sn);
end if;

end function;
/***************************************************************************/


/** deGkaGkm
  This procedure determines a set of generators $G$ of the set-wise stabiliser of $\{k-1,\ldots,n\}$ 
  of the decomposition group of a triangular ideal $\langle S \rangle$ from the set-wise stabiliser of $\{k-1,\ldots,n\}$ 
 

  In order to improve efficiency, a preliminary modular test is performed. 

  @param k an integer
  @param G sequence of generators of the successive generators
  @param S triangular basis.
  @param SP the image of \texttt{S} modulo a prime number $p$.
  @param n the length of \texttt{S}.
  @param sn the symmetric group.
  @param orbites the set of orbits of $\{1,...,n\}$ under the action of <G>
  @param IsPureGaloisIdeal a boolean indicating if the case where we compute the
  predicate \textit{IsPureGaloisIdeal?}.
  @return G the set of the generators of the computed set-wise stabiliser.
  @warn when IsPureGaloisIdeal is true TRUE the procedure stops and returns FALSE.
*/
/***************************************************************************/
procedure deGkaGkm(k,~G,~S,~SP,~n,~sn,~orbites,~IsPureGaloisIdeal)
elt:={1..(k-1)} meet {Max(o): o in orbites};
while (not(IsEmpty(elt))) do
    ak := Max(elt);
    elt := elt diff {ak};
    ll:= SetToSequence({1..k} diff {ak}) cat [ak] cat [(k+1)..n] ;
    if NormalForm((SP[k])^(sn! ll) ,SP) eq 0 then
        if NormalForm((S[k])^(sn! ll) ,S) eq 0 then
            sigma2:= OnePermutation([ak] cat [(k+1)..n],S,SP);
            if not (sigma2 eq Id(sn)) then
                G:=G cat [sigma2];
                orbites:=NewOrbits(orbites,sigma2);
                elt:=elt meet {Max(o): o in orbites};
            else  
                if IsPureGaloisIdeal then
                    IsPureGaloisIdeal := false;//A backtrack appears;
                    break;
                end if;
            end if;
        end if;
    end if;
end while;
end procedure;
/***************************************************************************/


/** Decomposition_group

This is the central function. It determines the decomposition group of a
triangular ideal and permits the computation of the predicate
\textit{IsPureGaloisIdeal?}.
 
 
  @param S a triangular basis given as a sequence.
  @param GalId a boolean. 
  @warn If \texttt{GalId} is set to TRUE, then \texttt{S} must be
  separable and the computation stops as soon as the ideal $\langle S
  \rangle$ is detected to not be a pure Galois ideal.  
  @return the group of decomposition of the ideal <S>.
  @return the value of \textit{IsPureGaloisIdeal?}.
  @warn when \texttt{GalId} is set to TRUE then the second output is
  correct but the former is not guaranteed, and vice-versa.  
*/
/***************************************************************************/
function Decomposition_group(S,GalId)

S:=Reverse(Sort(S));
n:=Rank(Parent (S[1]));
sn:=Sym(n);
G:=[];
orbites := {{i}: i in {1..n}} ;

// We search for a prime p which is compatible.

p:=2;
rep:=false;
while rep eq false do
    polmod:=PolynomialRing(FiniteField(p),n);
    k:=1;
    rep:=true;
    while (rep eq true) and (k le n) do
        rep:=IsCoercible(polmod,S[k]);
        k:=k+1;
    end while;
    if rep eq false then p:=NextPrime(p); end if;
end while;

polmod:=PolynomialRing(FiniteField(p),n);
SP := [(polmod ! S[i]) : i in [1..n]];

// Computation of the group of decomposition.

if GalId then
for k:= 2 to n do
    deGkaGkm(k,~G,~S,~SP,~n,~sn,~orbites,~GalId);
    if not(GalId) then
       break;
    end if;
end for;
else
for k:= 2 to n do
    deGkaGkm(k,~G,~S,~SP,~n,~sn,~orbites,~GalId);
end for;
end if;

// Computation of Card(V(<S>)).
// Here S must be separable.

card := &*[LeadingTotalDegree(f) : f in S];

// Output

G:=PermutationGroup<n|G>;
return G, (card eq #G) and GalId;
end function;
/***************************************************************************/

/*** Examples of uses in \textsc{Magma} of the function \texttt{Decomposition_group}
*/

/** Computation of the decomposition group of a triangular ideal and of the predicate {\it Is Pure Galois Ideal ?}

The ideal $I$ of $\mathbb{Q}[x_1,\ldots,x_n]$ generated by the triangular Gröbner basis (with respect to the lexicographical order induced by $x_8<\ldots<x_1$)

$ S = \{ x_1 - x_4x_6x_8^7 + x_4x_6x_8^3,\\
    x2^2 + x2*x3 + x2*x4 + x3^2 + x3*x4 + x4^2,\\
    x3^3 + x3^2*x4 + x3*x4^2 + x4^3,\\
    x4^4 + x8^4 - 1,\\
    x_5 + x_6,\\
    x_6^2 + x_8^2,\\
    x_7 + x_8,\\
    x_8^8 - x_8^4 - 1\} \\
$
is an Galois ideal of $ f(x) = x_8^8 - x_8^4 - 1 \in \mathbb{Q}[x] $. 

The \textsc{Magma} function Decomposition_Group return the decomposition group of $I$~:

\texttt{
> PR8<x1,x2,x3,x4,x5,x6,x7,x8>:= PolynomialRing(RationalField(), 8);\\
> \\
> S:=[\\  
>     x1 + x2 + x3 + x4,\\
>     x2^2 + x2*x3 + x2*x4 + x3^2 + x3*x4 + x4^2,\\
>     x3^3 + x3^2*x4 + x3*x4^2 + x4^3,\\
>     x4^4 + x8^4 - 1,\\
>     x5 + x6,\\
>     x6^2 + x8^2,\\
>     x7 + x8,\\
>     x8^8 - x8^4 - 1\\
> ];\\
> \\
}
The decomposition group of the ideal $I$ and the predicate "Is $I$ a pure Galois Ideal ?" are computed by the function call~:

\texttt{
> Decomposition_group(S,true); \\
Permutation group acting on a set of cardinality 8\\
Order = 192 = 2^6 * 3\\
    (1, 2)\\
    (2, 3)\\
    (3, 4)\\
    (5, 6)\\
    (7, 8)\\
    (5, 7)(6, 8)\\
false             \\
> \\
}

The second ouput shows that this Galois ideal is not a pure Galois ideal.
*/
/** Computation of the decomposition group

The ideal $I$ of $\mathbb{Q}[x_1,\ldots,x_n]$ generated by the triangular Gröbner basis (with respect to the lexicographical order induced by $x_8<\ldots<x_1$)

$ S = \{ x_1 - x_4x_6x_8^7 + x_4x_6x_8^3,\\
    x_2 + x_4x_6x_8^7 - x_4x_6x_8^3,\\
    x_3 + x_4,\\
    x_4^2 - x_6x_8^5 + x_6x_8,\\
    x_5 + x_6,\\
    x_6^2 + x_8^2,\\
    x_7 + x_8,\\
    x_8^8 - x_8^4 - 1\} \\
$
is a relations ideal of $ f(x) = x_8^8 - x_8^4 - 1 \in \mathbb{Q}[x] $. 

The \textsc{Magma} function Decomposition_Group return the decomposition group of $I$~:

\texttt{
> PR8<x1,x2,x3,x4,x5,x6,x7,x8>:= PolynomialRing(RationalField(), 8);\\
> \\
> S:=[\\  
>     x1 + x2 + x3 + x4,\\
>     x_2 + x_4x_6x_8^7 - x_4x_6x_8^3,\\
>     x_3 + x_4,\\
>     x_4^2 - x_6x_8^5 + x_6x_8,\\
>     x5 + x6,\\
>     x6^2 + x8^2,\\
>     x7 + x8,\\
>     x8^8 - x8^4 - 1\\
> ];\\
> \\
> Decomposition_group(S,true);                   \\         
Permutation group acting on a set of cardinality 8\\
Order = 32 = 2^5\\
    (1, 2)(3, 4)\\
    (1, 3)(2, 4)(5, 6)\\
    (1, 3)(2, 4)(7, 8)\\
    (1, 2)(5, 7)(6, 8)\\
    (1, 8, 4, 6, 2, 7, 3, 5)\\
true\\
> \\
}
 
The second parameter of \textit{Decomposition_group} is set to
\texttt{true} and, in this case, the predicate is calculated. The second output shows that the ideal $\langle S \rangle$ is a pure Galois ideal.
This Galois ideal is actually a relations ideal of $f$ and therefore the decomposition group of $I$ is the the symetric representation of the Galois group of $f$ corresponding to $I$.
  
*/