Anforderungen  |   Konzepte  |   Entwurf  |   Entwicklung  |   Qualitätssicherung  |   Lebenszyklus  |   Steuerung
 
 
 
 


Quelle  csetperm.gi   Sprache: unbekannt

 
Spracherkennung für: .gi vermutete Sprache: Unknown {[0] [0] [0]} [Methode: Schwerpunktbildung, einfache Gewichte, sechs Dimensionen]

#############################################################################
##
##  This file is part of GAP, a system for computational discrete algebra.
##  This file's authors include Alexander Hulpke, Heiko Theißen.
##
##  Copyright of GAP belongs to its developers, whose names are too numerous
##  to list here. Please refer to the COPYRIGHT file for details.
##
##  SPDX-License-Identifier: GPL-2.0-or-later
##
##  This file contains the operations for cosets of permutation groups
##

#############################################################################
##
#F  MinimizeExplicitTransversal( <U>, <maxmoved> )  . . . . . . . . . . local
##
InstallGlobalFunction( MinimizeExplicitTransversal, function( U, maxmoved )
    local   explicit,  lenflock,  flock,  lenblock,  index,  s;

    if     IsBound( U.explicit )
       and IsBound( U.stabilizer )  then
        explicit := U.explicit;
        lenflock := U.stabilizer.index * U.lenblock / Length( U.orbit );
        flock    := U.flock;
        lenblock := U.lenblock;
        index    := U.index;
        ChangeStabChain( U, [ 1 .. maxmoved ] );
        for s  in [ 1 .. Length( explicit ) ]  do
            explicit[ s ] := MinimalElementCosetStabChain( U, explicit[ s ] );
        od;
        Sort( explicit );
        U.explicit := explicit;
        U.lenflock := lenflock;
        U.flock    := flock;
        U.lenblock := lenblock;
        U.index    := index;
    fi;
end );

#############################################################################
##
#F  RightTransversalPermGroupConstructor( <filter>, <G>, <U> )  . constructor
##
MAX_SIZE_TRANSVERSAL := 100000;

# so far only orbits and perm groups -- TODO: Other deduced actions
InstallGlobalFunction(ActionRefinedSeries,function(G,U)
local o,A,ser,act,i;
  o:=List(Orbits(U,MovedPoints(G)),Set);
  SortBy(o,Length);
  A:=G;
  ser:=[A];
  act:=[0]; # dummy entry
  i:=1;
  while i<=Length(o) and Size(A)>Size(U) do
    A:=Stabilizer(A,o[i],OnSets);
    if Size(A)<Size(Last(ser)) then
      Add(ser,A);
      Add(act,[o[i],OnSets]);
    fi;
    i:=i+1;
  od;
  if Size(A)>Size(U) then
    Add(ser,U);
    Add(act,fail);
  fi;
  # refine large step?
  for i in [1..Length(ser)-1] do
    if IndexNC(ser[i],ser[i+1])>MAX_SIZE_TRANSVERSAL then
      A:=IntermediateGroup(ser[i],ser[i+1]:cheap);
      if A<>fail then
        # refine with action
        o:=ActionRefinedSeries(ser[i],A);
        ser:=Concatenation(ser{[1..i]},o[1]{[1..Length(o[1])-1]},
          ser{[i+1..Length(ser)]});
        act:=Concatenation(act{[1..i]},o[2]{[1..Length(o[2])-1]},
          act{[i+1..Length(act)]});
      else
        # no refinement, next step
        i:=i+1;
      fi;
    else
      # no refinement needed, next step
      i:=i+1;
    fi;
  od;
  # make ascending like AscendingSeries
  return [Reversed(ser),Reversed(act)];
end);

BindGlobal( "RightTransversalPermGroupConstructor", function( filter, G, U )
  local GC, UC, noyet, orbs, domain, GCC, UCC, ac, nc, bpt, enum, i,
    actions,nct;

    GC := CopyStabChain( StabChainImmutable( G ) );
    UC := CopyStabChain( StabChainImmutable( U ) );
    noyet:=ValueOption("noascendingchain")<>true;
    if not IsTrivial( G )  then
        orbs := ShallowCopy( OrbitsDomain( U, MovedPoints( G ) ) );
        Sort( orbs, function( o1, o2 )
            return Length( o1 ) < Length( o2 ); end );
        domain := Concatenation( orbs );
        GCC:=GC;
        UCC:=UC;
        while    Length( GCC.genlabels ) <> 0
              or Length( UCC.genlabels ) <> 0  do
#Print(SizeStabChain(GCC),"/",SizeStabChain(UCC),":",
#  SizeStabChain(GCC)/SizeStabChain(UCC),"\n");
          if noyet and (
          (SizeStabChain(GCC)/SizeStabChain(UCC) >MAX_SIZE_TRANSVERSAL) or
          (Length(UCC.genlabels)=0 and
            SizeStabChain(GCC)>MAX_SIZE_TRANSVERSAL)
            ) then

            # first get a factorization through actions
            ac:=ActionRefinedSeries(G,U);
            actions:=ac[2];
            ac:=ac[1];

            # go in biggish steps through the chain
            nc:=[ac[1]];
            nct:=[actions[1]];
            for i in [3..Length(ac)] do
              if Size(ac[i])/Size(Last(nc))>MAX_SIZE_TRANSVERSAL then
                Add(nc,ac[i-1]);
                Add(nct,actions[i-1]);
              fi;
            od;
            Add(nc,Last(ac));
            Add(nct,Last(actions));
            if Length(nc)>2 then
              ac:=[];
              for i in [Length(nc),Length(nc)-1..2] do
                Info(InfoCoset,4,"Recursive [",Size(nc[i]),",",Size(nc[i-1]));
                Add(ac,RightTransversal(nc[i],nc[i-1]
                      # do not try to factor again
                      :noascendingchain));
              od;
              return FactoredTransversal(G,U,ac);
            fi;
            noyet:=false;

          fi;
          bpt := First( domain, p -> not IsFixedStabilizer( GCC, p ) );
          ChangeStabChain( GCC, [ bpt ], true  );  GCC := GCC.stabilizer;
          ChangeStabChain( UCC, [ bpt ], false );  UCC := UCC.stabilizer;
        od;
    fi;

    AddCosetInfoStabChain(GC,UC,LargestMovedPoint(G));
    MinimizeExplicitTransversal(UC,LargestMovedPoint(G));

    enum := Objectify( NewType( FamilyObj( G ),
                           filter and IsList and IsDuplicateFreeList
                           and IsAttributeStoringRep ),
          rec( group := G,
            subgroup := U,
      stabChainGroup := GC,
   stabChainSubgroup := UC ) );

    return enum;
end );


#############################################################################
##
#R  IsRightTransversalPermGroupRep( <obj> ) . right transversal of perm group
##
DeclareRepresentation( "IsRightTransversalPermGroupRep",
    IsRightTransversalRep,
    [ "stabChainGroup", "stabChainSubgroup" ] );

InstallMethod( \[\],
    "for right transversal of perm. group, and pos. integer",
    true,
    [ IsList and IsRightTransversalPermGroupRep, IsPosInt ], 0,
    function( cs, num )
    return CosetNumber( cs!.stabChainGroup, cs!.stabChainSubgroup, num );
end );

InstallMethod( PositionCanonical,
    "for right transversal of perm. group, and permutation",
    IsCollsElms,
    [ IsList and IsRightTransversalPermGroupRep, IsPerm ], 0,
    function( cs, elm )
    return NumberCoset( cs!.stabChainGroup,
                        cs!.stabChainSubgroup,
                        elm );
end );

#############################################################################
##
#M  RightTransversalOp( <G>, <U> )  . . . . . . . . . . . . . for perm groups
##
InstallMethod( RightTransversalOp,
    "for two perm. groups",
    IsIdenticalObj,
    [ IsPermGroup, IsPermGroup ], 0,
    function( G, U )
    return RightTransversalPermGroupConstructor(
               IsRightTransversalPermGroupRep, G, U );
end );


#############################################################################
##
#F  AddCosetInfoStabChain( <G>, <U>, <maxmoved> ) . . . . . .  add coset info
##

InstallGlobalFunction( AddCosetInfoStabChain, function( G, U, maxmoved )
    local   orb,  pimg,  img,  vert,  s,  t,  index,
            block,  B,  blist,  pos,  sliced,  lenflock,  i,  j,
            ss,  tt,t1,t1lim,found,tl,vimg,
            sel,shortsel,gorpo,pisel,prepi,transinv;

    # iterated image
    vimg:=function(point,list)
    local i;
      for i in list do
        point:=point^i;
      od;
      return point;
    end;

    Info(InfoCoset,5,"AddCosetInfoStabChain [",
          SizeStabChain(G),",",SizeStabChain(U),"]");
    if IsEmpty( G.genlabels )  then
        U.index    := 1;
        U.explicit := [ U.identity ];
        U.lenflock := 1;
        U.flock    := U.explicit;
    else
        AddCosetInfoStabChain( G.stabilizer, U.stabilizer, maxmoved );

        # U.index := [G_1:U_1];
        U.index := U.stabilizer.index * Length( G.orbit ) / Length( U.orbit );
        Info(InfoCoset,5,"U.index=",U.index);

        # block := 1 ^ <U,G_1>; is a block for G.
        block := OrbitPerms( Concatenation( U.generators,
                 G.stabilizer.generators ), G.orbit[ 1 ] );
        U.lenblock := Length( block );
        lenflock := Length( G.orbit ) / U.lenblock;

        # For small indices,  permutations   are multiplied,  so  we  need  a
        # multiplied transversal.
        if     IsBound( U.stabilizer.explicit )
           and U.lenblock * maxmoved <= MAX_SIZE_TRANSVERSAL
           and U.index    * maxmoved <= MAX_SIZE_TRANSVERSAL * lenflock  then
            U.explicit := [  ];
            U.flock    := [ G.identity ];
            tt := [  ];  tt[ G.orbit[ 1 ] ] := G.identity;
            for t  in G.orbit  do
                tt[ t ] := tt[ t ^ G.transversal[ t ] ] /
                           G.transversal[ t ];
            od;
        fi;

        # flock := { G.transversal[ B[1] ] | B in block system };
        blist := BlistList( G.orbit, block );
        pos := Position( blist, false );
        while pos <> fail  do
            img := G.orbit[ pos ];
            B := block{ [ 1 .. U.lenblock ] };
            sliced := [  ];
            while img <> G.orbit[ 1 ]  do
                Add( sliced, G.transversal[ img ] );
                img := img ^ G.transversal[ img ];
            od;
            for i  in Reversed( [ 1 .. Length( sliced ) ] )  do
                for j  in [ 1 .. Length( B ) ]  do
                    B[ j ] := B[ j ] / sliced[ i ];
                od;
            od;
            Append( block, B );
            if IsBound( U.explicit )  then
                Add( U.flock, tt[ B[ 1 ] ] );
            fi;
            #UniteBlist( blist, BlistList( G.orbit, B ) );
            UniteBlistList(G.orbit, blist, B );
            pos := Position( blist, false, pos );
        od;
        G.orbit := block;

        # Let <s> loop over the transversal elements in the stabilizer.
        U.repsStab := List( [ 1 .. U.lenblock ], x ->
                           BlistList( [ 1 .. U.stabilizer.index ], [  ] ) );
        U.repsStab[ 1 ] := BlistList( [ 1 .. U.stabilizer.index ],
                                      [ 1 .. U.stabilizer.index ] );
        index := U.stabilizer.index * lenflock;
        s := 1;

        # For  large  indices, store only   the  numbers of  the  transversal
        # elements needed.
        if not IsBound( U.explicit )  then

            # If  the   stabilizer   is the   topmost  level   with  explicit
            # transversal, this must contain minimal coset representatives.
            MinimizeExplicitTransversal( U.stabilizer, maxmoved );

            # if there are over 200 points, do a cheap test first.
            t1lim:=Length(G.orbit);
            if t1lim>200 then
              t1lim:=50;
            fi;

            sel:=Filtered([1..Maximum(G.orbit)],x->IsBound(U.translabels[x]));
            shortsel:=Length(sel)<t1lim;

            if shortsel then
              # inverse transversal elements
              t1:=ShallowCopy(G.generators);
              Add(t1,G.identity);
              vert:=List(t1,x->x^-1);
              # inverses of transversal, stored compact
              # TODO: Instead of `Position`, use translabels entry
              #transinv:=List(G.transversal,x->vert[Position(t1,x)]);
              transinv:=[];
              for i in [1..Length(G.transversal)] do
                if IsBound(G.transversal[i]) then
                  t:=Position(t1,G.transversal[i]);
                  if t=fail then
                    Add(t1,G.transversal[i]);
                    Add(vert,G.transversal[i]^-1);
                    t:=Length(t1);
                  fi;
                  transinv[i]:=vert[t];
                fi;
              od;
              # Position in orbit
              gorpo:=[];
              for i in [1..Length(G.orbit)] do
                gorpo[G.orbit[i]]:=i;
              od;
            fi;

            orb := G.orbit{ [ 1 .. U.lenblock ] };
            pimg := [  ];
            while index < U.index  do

                pimg{ orb } := CosetNumber( G.stabilizer, U.stabilizer, s,
                                       orb );
                t := 2;

                if shortsel then

                  # test for the few wrong values, mapping backwards

                  pisel:=Filtered([1..Length(pimg)],
                    x->IsBound(pimg[x]) and pimg[x] in sel);

                  while t <= U.lenblock  and  index < U.index  do

                    # For this point  in the  block,  find the images  of the
                    # earlier points under the representative.
                    vert := G.orbit{ [ 1 .. t-1 ] };
                    img := G.orbit[ t ];
                    tl:=[];
                    while img <> G.orbit[ 1 ]  do
                      Add(tl,transinv[img]);
                      img  := img           ^ G.transversal[ img ];
                    od;
                    prepi:=pisel;
                    for t1 in [Length(tl),Length(tl)-1..1] do
                      prepi:=OnTuples(prepi,tl[t1]);
                    od;

                    # If $Ust = Us't'$ then $1t'/t/s in 1U$. Also if $1t'/t/s
                    # in 1U$ then $st/t' =  u.g_1$ with $u  in U, g_1 in G_1$
                    # and $g_1  =  u_1.s'$ with $u_1  in U_1,  s' in S_1$, so
                    # $Ust = Us't'$.
                    #if ForAll( [ 1 .. t-1 ], i -> not
                    #     vert[ i ] in prepi  )  then
                    if not ForAny(gorpo{prepi},x->x>=1 and x<=t-1) then
                      U.repsStab[ t ][ s ] := true;
                      index := index + lenflock;

                    fi;

                    t := t + 1;
                  od;

                else

                  while t <= U.lenblock  and  index < U.index  do

                    # do not test all points first if not necessary
                    # (test only at most t1lim points, if the test succeeds,
                    # test the rest)
                    # this gives a major speedup.
                    t1:=Minimum(t-1,t1lim);
                    # For this point  in the  block,  find the images  of the
                    # earlier points under the representative.
                    vert := G.orbit{ [ 1 .. t1 ] };
                    img := G.orbit[ t ];
                    while img <> G.orbit[ 1 ]  do
                      vert := OnTuples( vert, G.transversal[ img ] );
                      img  := img           ^ G.transversal[ img ];
                    od;

                    # If $Ust = Us't'$ then $1t'/t/s in 1U$. Also if $1t'/t/s
                    # in 1U$ then $st/t' =  u.g_1$ with $u  in U, g_1 in G_1$
                    # and $g_1  =  u_1.s'$ with $u_1  in U_1,  s' in S_1$, so
                    # $Ust = Us't'$.
                    if ForAll( [ 1 .. t1 ], i -> not IsBound
                        ( U.translabels[ pimg[ vert[ i ] ] ] ) )  then

                      # do all points
                      if t1<t-1 then
                        img := G.orbit[ t ];
                        if t<=10*t1lim then
                          vert := G.orbit{ [ 1 .. t - 1 ] };
                          while img <> G.orbit[ 1 ]  do
                            vert := OnTuples( vert, G.transversal[ img ] );
                            img  := img           ^ G.transversal[ img ];
                          od;
                          found:=ForAll( [ t1+1 .. t - 1 ], i -> not IsBound
                            ( U.translabels[ pimg[ vert[ i ] ] ] ) );
                        else
                          # avoid calculating tons of images of a long list
                          # instead calculate images on the fly
                          # this implicitly assumes that, if we get to so
                          # long a list, failure will happen quickly.
                          tl:=[];
                          while img <> G.orbit[ 1 ]  do
                            Add(tl,G.transversal[img]);
                            img  := img           ^ G.transversal[ img ];
                          od;
                          found:=ForAll( [ t1+1 .. t - 1 ], i -> not IsBound
                            ( U.translabels[ pimg[ vimg(G.orbit[ i ],tl) ] ] ) );
                        fi;

                        if found then
                          U.repsStab[ t ][ s ] := true;
                          index := index + lenflock;
                        fi;

                      else
                        U.repsStab[ t ][ s ] := true;
                        index := index + lenflock;
                      fi;
                    fi;

                    t := t + 1;
                  od;
                fi;

                s := s + 1;
            od;

        # For small indices, store a transversal explicitly.
        else
            for ss  in U.stabilizer.flock  do
                Append( U.explicit, U.stabilizer.explicit * ss );
            od;
            while index < U.index  do
                t := 2;
                while t <= U.lenblock  and  index < U.index  do
                    ss := U.explicit[ s ] * tt[ G.orbit[ t ] ];
                    if ForAll( [ 1 .. t - 1 ], i -> not IsBound
                           ( U.translabels[ G.orbit[ i ] / ss ] ) )  then
                        U.repsStab[ t ][ s ] := true;
                        Add( U.explicit, ss );
                        index := index + lenflock;
                    fi;
                    t := t + 1;
                od;
                s := s + 1;
            od;
            Unbind( U.stabilizer.explicit );
            Unbind( U.stabilizer.flock    );
        fi;

    fi;
end );

#############################################################################
##
#F  NumberCoset( <G>, <U>, <r> )  . . . . . . . . . . . . . . coset to number
##
InstallGlobalFunction( NumberCoset, function( G, U, r )
    local   num,  b,  t,  u,  g1,  pnt,  bpt;

    if IsEmpty( G.genlabels )  or  U.index = 1  then
        return 1;
    fi;

    # Find the block number of $r$.
    bpt := G.orbit[ 1 ];
    b := QuoInt( Position( G.orbit, bpt ^ r ) - 1, U.lenblock );

    # For small indices, look at the explicit transversal.
    if IsBound( U.explicit )  then
        return b * U.lenflock + Position( U.explicit,
               MinimalElementCosetStabChain( U, r / U.flock[ b + 1 ] ) );
    fi;

    pnt := G.orbit[ b * U.lenblock + 1 ];
    while pnt <> bpt  do
        r   := r   * G.transversal[ pnt ];
        pnt := pnt ^ G.transversal[ pnt ];
    od;

    # Now $r$ stabilises the block. Find the first $t in G/G_1$ such that $Ur
    # = Ust$ for $s in G_1$. In this code, G.orbit[ <t> ] = bpt ^ $t$.
    num := b * U.stabilizer.index * U.lenblock / Length( U.orbit );
             # \_________This is [<U,G_1>:U] = U.lenflock_________/
    t := 1;
    pnt := G.orbit[ t ] / r;
    while not IsBound( U.translabels[ pnt ] )  do
        num := num + SizeBlist( U.repsStab[ t ] );
        t := t + 1;
        pnt := G.orbit[ t ] / r;
    od;

    # $r/t = u.g_1$ with $u in U, g_1 in G_1$, hence $t/r.u = g_1^-1$.
    u := U.identity;
    while pnt ^ u <> bpt  do
        u := u * U.transversal[ pnt ^ u ];
    od;
    g1 := LeftQuotient( u, r );  # Now <g1> = $g_1.t = u mod r$.
    while bpt ^ g1 <> bpt  do
        g1 := g1 * G.transversal[ bpt ^ g1 ];
    od;

    # The number of $r$  is the number of $g_1$  plus an offset <num> for
    # the earlier values of $t$.
    return num + SizeBlist( U.repsStab[ t ]{ [ 1 ..
                   NumberCoset( G.stabilizer, U.stabilizer, g1 ) ] } );

end );

#############################################################################
##
#F  CosetNumber( <arg> )  . . . . . . . . . . . . . . . . . . number to coset
##
InstallGlobalFunction( CosetNumber, function( arg )
    local   G,  U,  num,  tup,  b,  t,  rep,  pnt,  bpt,  index,  len;

    # Get the arguments.
    G := arg[ 1 ];  U := arg[ 2 ];  num := arg[ 3 ];
    if Length( arg ) > 3  then  tup := arg[ 4 ];
                          else  tup := false;     fi;

    if num = 1  then
        if tup = false  then  return G.identity;
                        else  return tup;         fi;
    fi;

    # Find the block $b$ addressed by <num>.
    if IsBound( U.explicit )  then
        index := U.lenflock;
    else
        index := U.stabilizer.index * U.lenblock / Length( U.orbit );
               # \_________This is [<U,G_1>:U] = U.lenflock_________/
    fi;
    b := QuoInt( num - 1, index );
    num := ( num - 1 ) mod index + 1;

    # For small indices, look at the explicit transversal.
    if IsBound( U.explicit )  then
        if tup = false  then
            return U.explicit[ num ] * U.flock[ b + 1 ];
        else
            return List( tup, t -> t / U.flock[ b + 1 ] / U.explicit[ num ] );
        fi;
    fi;

    # Otherwise, find the point $t$ addressed by <num>.
    t := 1;
    len := SizeBlist( U.repsStab[ t ] );
    while num > len  do
        num := num - len;
        t := t + 1;
        len := SizeBlist( U.repsStab[ t ] );
    od;
    if len < U.stabilizer.index  then
        num := PositionNthTrueBlist( U.repsStab[ t ], num );
    fi;

    # Find the representative $s$ in   the stabilizer addressed by <num>  and
    # return $st$.
    rep := G.identity;
    bpt := G.orbit[ 1 ];
    if tup = false  then
        pnt := G.orbit[ b * U.lenblock + 1 ];
        while pnt <> bpt  do
            rep := rep * G.transversal[ pnt ];
            pnt := pnt ^ G.transversal[ pnt ];
        od;
        pnt := G.orbit[ t ];
        while pnt <> bpt  do
            rep := rep * G.transversal[ pnt ];
            pnt := pnt ^ G.transversal[ pnt ];
        od;
        return CosetNumber( G.stabilizer, U.stabilizer, num ) / rep;
    else
        pnt := G.orbit[ b * U.lenblock + 1 ];
        while pnt <> bpt  do
            tup := OnTuples( tup, G.transversal[ pnt ] );
            pnt := pnt ^ G.transversal[ pnt ];
        od;
        pnt := G.orbit[ t ];
        while pnt <> bpt  do
            tup := OnTuples( tup, G.transversal[ pnt ] );
            pnt := pnt ^ G.transversal[ pnt ];
        od;
        return CosetNumber( G.stabilizer, U.stabilizer, num, tup );
    fi;
end );

#############################################################################
##
#M  AscendingChainOp(<G>,<pnt>) . . . approximation of
##
InstallMethod( AscendingChainOp, "PermGroup", IsIdenticalObj,
  [IsPermGroup,IsPermGroup],0,
function(G,U)
local s,c,mp,o,i,step,a;
  s:=G;
  c:=[G];
  repeat
    mp:=MovedPoints(s);
    o:=ShallowCopy(OrbitsDomain(s,mp));
    SortBy(o,Length);
    i:=1;
    step:=false;
    while i<=Length(o) and step=false do
      if not IsTransitive(U,o[i]) then
        Info(InfoCoset,2,"AC: orbit");
        o:=ShallowCopy(OrbitsDomain(U,o[i]));
        SortBy(o,Length);
        # union of same length -- smaller index
        a:=Union(Filtered(o,x->Length(x)=Length(o[1])));
        if Length(a)=Sum(o,Length) then
          a:=Set(o[1]);
        fi;
        s:=Stabilizer(s,a,OnSets);
        step:=true;
      elif Index(G,U)>NrMovedPoints(U)
          and IsPrimitive(s,o[i]) and not IsPrimitive(U,o[i]) then
        Info(InfoCoset,2,"AC: blocks");
        s:=Stabilizer(s,Set(MaximalBlocks(U,o[i]),Set),
                      OnSetsDisjointSets);
        step:=true;
      else
        i:=i+1;
      fi;
    od;
    if step then
      Add(c,s);
    fi;
  until step=false or Index(s,U)=1; # we could not refine better
  if Index(s,U)>1 then
    Add(c,U);
  fi;
  Info(InfoCoset,2,"Indices",List([1..Length(c)-1],i->Index(c[i],c[i+1])));
  return RefinedChain(G,Reversed(c));
end);

InstallMethod(CanonicalRightCosetElement,"Perm",IsCollsElms,
  [IsPermGroup,IsPerm],0,
function(U,e)
  return MinimalElementCosetStabChain(MinimalStabChain(U),e);
end);

InstallMethod(\<,"RightCosets of perm group",IsIdenticalObj,
  [IsRightCoset and IsPermCollection,IsRightCoset and IsPermCollection],0,
function(a,b)
  # for permutation groups the canonical rep is the smallest element of the
  # coset
  if ActingDomain(a)<>ActingDomain(b) then
    return ActingDomain(a)<ActingDomain(b);
  fi;
  return CanonicalRepresentativeOfExternalSet(a)
         <CanonicalRepresentativeOfExternalSet(b);
end);



InstallMethod(Intersection2, "perm cosets", IsIdenticalObj,
  [IsRightCoset and IsPermCollection,IsRightCoset and IsPermCollection],0,
function(cos1,cos2)
    local H1, H2, x1, x2, shift, sigma, listMoved_H1, listMoved_H2,
          listMoved_sigma, U, repr, H1_sigma, H2_sigma, H12, swap, rho, diff;
    # We set cosInt = cos1 \cap cos2 = H1 x1 \cap H2 x2
    H1:=ActingDomain(cos1);
    H2:=ActingDomain(cos2);
    x1:=Representative(cos1);
    x2:=Representative(cos2);
    if H1=H2 then
        if cos1=cos2 then
            return cos1;
        else
            return [];
        fi;
    fi;
    # We are using that
    #    H1*x1 \cap H2*x2 = (H1 \cap H2*x2/x1)*x1 = (H1 \cap H2*sigma)*shift,
    # where shift and sigma are defined as below:
    shift:=x1;
    sigma:=x2 / x1;

    # Reducing as much as possible in advance by using various relatively
    # cheap to compute criteria
    while true do
        listMoved_H1:=MovedPoints(H1);
        listMoved_H2:=MovedPoints(H2);
        listMoved_sigma:=MovedPoints(sigma);

        # If the coset intersection is non-empty, then there is h1 \in H1 and
        # h2 \in H2 such that h1 = h2*sigma. Therefore sigma is contained in
        # the group generated by H1 and H2. A necessary condition for this is
        # that the points moved by sigma are a subset of the points moved by
        # H1 and H2.
        if not IsSubset(Union(listMoved_H1, listMoved_H2), listMoved_sigma) then
            return [];
        fi;

        # Suppose x is an element of H1 \cap H2*sigma. Then for any positive
        # integer n, we know that n^x is contained in n^H1 but also in
        # (n^H2)^\sigma. Thus if the intersection of n^H1 and (n^H2)^\sigma is
        # empty, then the intersection of the cosets is also empty. Clearly
        # the orbit intersection contains n whenever n is fixed by sigma, so
        # we only have to consider this for n moved by sigma.
        if ForAny(listMoved_sigma, n -> IsEmpty(Intersection(Orbit(H1,n), OnTuples(Orbit(H2,n),sigma)))) then
            return [];
        fi;

        # If there are points that are moved by sigma and by H2 but not by H1,
        # then there must be an element in H2 which matches the action of sigma
        # on these points, or else the intersection is empty
        diff:=Difference(Intersection(listMoved_H2, listMoved_sigma), listMoved_H1);
        if Length(diff) > 0 then
            repr:=RepresentativeAction(H2, diff, OnTuples(diff, Inverse(sigma)), OnTuples);
            if repr=fail then
                return [];
            fi;
            # Since repr is in H2, we can replace sigma by repr*sigma
            # without changing the coset H2*sigma. The new sigma then fixes
            # all points in diff, as does H1. Hence replacing H2 by the
            # stabilizer in H2 of diff does not change H1 \cap H2 sigma.
            sigma:=repr * sigma;
            H2:=Stabilizer(H2, diff, OnTuples);
            continue;
        fi;

        # Mirror to the previous check:
        # If there are points that are moved by sigma and by H1 but not by H2,
        # then there must be an element in H1 which matches the action of sigma
        # on these points, or else the intersection is empty
        diff:=Difference(Intersection(listMoved_H1, listMoved_sigma), listMoved_H2);
        if Length(diff) > 0 then
            repr:=RepresentativeAction(H1, diff, OnTuples(diff, sigma), OnTuples);
            if repr=fail then
                return [];
            fi;
            # Again this is similar to the case before, except that we are adjusting
            # H1 now and thus also need to take `shift` into account. The situation
            # is as follows:
            #
            # cosInt = (H1 \cap H2 sigma) shift
            #        = (H1 sigma^{-1} \cap H2) sigma shift
            #        = (Stab_{H1}(diff) repr sigma^{-1} \cap H2) sigma shift
            #        = (Stab_{H1}(diff) \cap H2 sigma repr^{-1} ) repr shift
            H1:=Stabilizer(H1, diff, OnTuples);
            sigma:=sigma / repr;
            shift:=repr * shift;
            continue;
        fi;

        # easy termination criterion: reduction to group intersection
        if sigma in H2 then
            U:=Intersection(H1, H2);
            return RightCoset(U, shift);
        fi;

        # easy termination criterion: reduction to group intersection
        if sigma in H1 then
            # cosInt = (H1 \cap H2 sigma) shift
            #        = (H1 sigma^{-1} \cap H2) sigma shift
            U:=Intersection(H1, H2);
            return RightCoset(U, sigma * shift);
        fi;

        # any element of H1 which moves points not moved by sigma or anything
        # in H2 can not be in the intersection, so we may as well remove them
        # by a stabilizer computation
        diff:=Difference(listMoved_H1, Union(listMoved_H2, listMoved_sigma));
        if Length(diff) > 0 then
            H1:=Stabilizer(H1, diff, OnTuples);
            continue;
        fi;

        # the same but with the roles of H1 and H2 reversed
        diff:=Difference(listMoved_H2, Union(listMoved_H1, listMoved_sigma));
        if Length(diff) > 0 then
            H2:=Stabilizer(H2, diff, OnTuples);
            continue;
        fi;

        # More general but more expensive than previous check
        H1_sigma:=ClosureGroup(H1, sigma);
        if not IsSubgroup(H1_sigma, H2) then
            H2:=Intersection(H1_sigma, H2);
            continue;
        fi;
        H2_sigma:=ClosureGroup(H2, sigma);
        if not IsSubgroup(H2_sigma, H1) then
            H1:=Intersection(H2_sigma, H1);
            continue;
        fi;
        # No more reduction tricks available
        break;
    od;

    # A final termination criterion
    H12:=ClosureGroup(H1, H2);
    if not sigma in H12 then
        return [];
    fi;

    # We are now inspired by the algorithm from
    # Lazlo Babai, Coset Intersection in Moderately Exponential Time
    #
    # We use the algorithm from Page 10 of coset analysis and we reformulate
    # it here in order to avoid errors:
    # --- The naive algorithm for computing H1 \cap H2 sigma is to iterate
    # over elements of H1 and testing if one belongs to H2 sigma. If we find
    # one such z then the result is the coset RightCoset(U, z). If not
    # then it is empty.
    # --- Since the result is independent of the cosets U, what we can
    # do is iterate over the RightCosets(H1, U). The algorithm is the
    # one of Proposition 3.2
    # for r in RightCosets(H1, U) do
    #   if r in H1*sigma then
    #     return RightCoset(U, r * shift)
    #   fi;
    # od;
    # --- (TODO for future work): The question is how to make it faster.
    # One idea is to use an ascending chain between U and H1.
    # Section 3.4 of above paper gives statement related to that but not a
    # useful algorithm. The question deserves further exploration.
    #
    # We select the smallest group for that computation in order to have
    # as few cosets as possible
    if Order(H2) < Order(H1) then
        # cosInt = (H1 \cap H2 sigma) shift
        #        = (H1 sigma^{-1} \cap H2) sigma shift
        swap:=H1;
        H1:=H2;
        H2:=swap;
        shift:=sigma * shift;
        sigma:=Inverse(sigma);
    fi;
    # So now Order(H1) <= Order(H2)
    U:=Intersection(H1, H2);
    for rho in RightTransversal(H1, U) do
        if rho / sigma in H2 then
            return RightCoset(U, rho * shift);
        fi;
    od;
    return [];
end);


#############################################################################
##
#F  FactorCosetAction( <G>, <U>, [<N>] )  operation on the right cosets Ug
##                                        with possibility to indicate kernel
##
BindGlobal("DoFactorCosetActionPerm",function(arg)
local G,u,op,h,N,rt,ac,actions,hom,i,q;
  G:=arg[1];
  u:=arg[2];
  if Length(arg)>2 then
    N:=arg[3];
  else
    N:=false;
  fi;
  if IsList(u) and Length(u)=0 then
    u:=G;
    Error("only trivial operation ?  I Set u:=G;");
  fi;
  if IsSubset(u, G) then
    return DoFactorCosetAction(G, u, G);
  fi;
  if N=false then
    N:=Core(G,u);
  fi;
  ac:=ActionRefinedSeries(G,u);
  actions:=ac[2];
  ac:=ac[1];
  hom:=false;
  for i in [2..Length(ac)] do
    if actions[i-1]<>fail
      # allow 2GB memory use for writing down orbit
      and SIZE_OBJ(actions[i-1][1])*IndexNC(ac[i],ac[i-1])<2*10^9 then

      op:=rec();
      h:=Orbit(ac[i],actions[i-1][1],actions[i-1][2]:permutations:=op);
      if IsBound(op.permutations) then
        rt:=List(op.permutations,PermList);
        q:=Group(rt);
        SetSize(q,IndexNC(G,N));
        h:=GroupHomomorphismByImagesNC(ac[i],Group(rt),
          op.generators,rt);
      else
        h:=ActionHomomorphism(ac[i],h,actions[i-1][2],"surjective");
      fi;
    else
      rt:=RightTransversal(ac[i],ac[i-1]);
      if not IsRightTransversalRep(rt) then
        # the right transversal has no special `PositionCanonical' method.
        rt:=List(rt,i->RightCoset(ac[i-1],i));
      fi;
      h:=ActionHomomorphism(ac[i],rt,OnRight,"surjective");

    fi;
    Unbind(op);
    Unbind(rt);
    if i=2 then
      hom:=h;
    else
      hom:=KuKGenerators(ac[i],h,hom);;
      q:=Group(hom);
      StabChainOptions(q).limit:=Size(ac[i]);
      hom:=GroupHomomorphismByImagesNC(ac[i],q,GeneratorsOfGroup(ac[i]),hom);;
    fi;
  od;

  op:=Image(hom,G);
  SetSize(op,IndexNC(G,N));

  # and note our knowledge
  SetKernelOfMultiplicativeGeneralMapping(hom,N);
  AddNaturalHomomorphismsPool(G,N,hom);
  return hom;
end);

InstallMethod(FactorCosetAction,"by right transversal operation",
  IsIdenticalObj,[IsPermGroup,IsPermGroup],0,
function(G,U)
  return DoFactorCosetActionPerm(G,U);
end);

InstallOtherMethod(FactorCosetAction,
  "by right transversal operation, given kernel",IsFamFamFam,
  [IsPermGroup,IsPermGroup,IsPermGroup],0,
function(G,U,N)
  return DoFactorCosetActionPerm(G,U,N);
end);


[ Dauer der Verarbeitung: 0.59 Sekunden  ]

                                                                                                                                                                                                                                                                                                                                                                                                     


Neuigkeiten

     Aktuelles
     Motto des Tages

Software

     Produkte
     Quellcodebibliothek

Aktivitäten

     Artikel über Sicherheit
     Anleitung zur Aktivierung von SSL

Muße

     Gedichte
     Musik
     Bilder

Jenseits des Üblichen ....

Besucherstatistik

Besucherstatistik

Monitoring

Montastic status badge