Quelle  rcwamap.gi   Sprache: unbekannt

#############################################################################
##
#W  rcwamap.gi                GAP4 Package `RCWA'                 Stefan Kohl
##
##  This file contains implementations of methods and functions for computing
##  with rcwa mappings of
##
##    - the ring Z of the integers, of
##    - the ring Z^2, of
##    - the semilocalizations Z_(pi) of the ring of integers, and of
##    - the polynomial rings GF(q)[x] in one variable over a finite field.
##
##  See the definitions given in the file rcwamap.gd.
##
#############################################################################

#############################################################################
##
#F  RCWAInfo . . . . . . . . . . . . . . . . . . set info level of `InfoRCWA'
##
InstallGlobalFunction( RCWAInfo,
                       function ( n ) SetInfoLevel( InfoRCWA, n ); end );

#############################################################################
##
#S  Shorthands for commonly used filters. ///////////////////////////////////
##
#############################################################################

BindGlobal( "IsRcwaMappingInStandardRep",
             IsRcwaMapping and IsRcwaMappingStandardRep );
BindGlobal( "IsRcwaMappingOfZInStandardRep",
             IsRcwaMappingOfZ and IsRcwaMappingStandardRep );
BindGlobal( "IsRcwaMappingOfZ_piInStandardRep",
             IsRcwaMappingOfZ_pi and IsRcwaMappingStandardRep );
BindGlobal( "IsRcwaMappingOfZOrZ_piInStandardRep",
             IsRcwaMappingOfZOrZ_pi and IsRcwaMappingStandardRep );
BindGlobal( "IsRcwaMappingOfZxZInStandardRep",
             IsRcwaMappingOfZxZ and IsRcwaMappingStandardRep );
BindGlobal( "IsRcwaMappingOfGFqxInStandardRep",
             IsRcwaMappingOfGFqx and IsRcwaMappingStandardRep );

BindGlobal( "IsRcwaMappingInSparseRep",
             IsRcwaMapping and IsRcwaMappingSparseRep );
BindGlobal( "IsRcwaMappingOfZInSparseRep",
             IsRcwaMappingOfZ and IsRcwaMappingSparseRep );
BindGlobal( "IsRcwaMappingOfZ_piInSparseRep",                 # unused so far
             IsRcwaMappingOfZ_pi and IsRcwaMappingSparseRep );
BindGlobal( "IsRcwaMappingOfZOrZ_piInSparseRep",
             IsRcwaMappingOfZOrZ_pi and IsRcwaMappingSparseRep );
BindGlobal( "IsRcwaMappingOfZxZInSparseRep",                  # unused so far
             IsRcwaMappingOfZxZ and IsRcwaMappingSparseRep );
BindGlobal( "IsRcwaMappingOfGFqxInSparseRep",                 # unused so far
             IsRcwaMappingOfGFqx and IsRcwaMappingSparseRep );

#############################################################################
##
#S  The families of rcwa mappings. //////////////////////////////////////////
##
#############################################################################

#############################################################################
##
#V  RcwaMappingsOfZFamily . . . . . . .  the family of all rcwa mappings of Z
##
BindGlobal( "RcwaMappingsOfZFamily",
            NewFamily( "RcwaMappingsFamily( Integers )",
                       IsRcwaMappingOfZ,
                       CanEasilySortElements, CanEasilySortElements ) );
SetFamilySource( RcwaMappingsOfZFamily, FamilyObj( 1 ) );
SetFamilyRange ( RcwaMappingsOfZFamily, FamilyObj( 1 ) );
SetUnderlyingRing( RcwaMappingsOfZFamily, Integers );

#############################################################################
##
#V  RcwaMappingsOfZxZFamily . . . . .  the family of all rcwa mappings of Z^2
##
BindGlobal( "RcwaMappingsOfZxZFamily",
            NewFamily( "RcwaMappingsFamily( Integers^2 )",
                       IsRcwaMappingOfZxZ,
                       CanEasilySortElements, CanEasilySortElements ) );
SetFamilySource( RcwaMappingsOfZxZFamily, FamilyObj( [ 1, 1 ] ) );
SetFamilyRange ( RcwaMappingsOfZxZFamily, FamilyObj( [ 1, 1 ] ) );
SetUnderlyingLeftModule( RcwaMappingsOfZxZFamily, Integers^2 );

## Internal variables storing the rcwa mapping families used in the
## current GAP session.

BindGlobal( "Z_PI_RCWAMAPPING_FAMILIES", [] );
BindGlobal( "GFQX_RCWAMAPPING_FAMILIES", [] );

#############################################################################
##
#F  RcwaMappingsOfZ_piFamily( <R> )
##
##  Returns the family of all rcwa mappings of a given semilocalization <R>
##  of the ring of integers.
##
InstallGlobalFunction( RcwaMappingsOfZ_piFamily,

  function ( R )

    local  fam, name;

    if   not IsZ_pi( R )
    then Error("usage: RcwaMappingsOfZ_piFamily( <R> )\n",
               "where <R> = Z_pi( <pi> ) for a set of primes <pi>.\n");
    fi;
    fam := First( Z_PI_RCWAMAPPING_FAMILIES,
                  fam -> UnderlyingRing( fam ) = R );
    if fam <> fail then return fam; fi;
    name := Concatenation( "RcwaMappingsFamily( ",
                           String( R ), " )" );
    fam := NewFamily( name, IsRcwaMappingOfZ_pi,
                      CanEasilySortElements, CanEasilySortElements );
    SetUnderlyingRing( fam, R );
    SetFamilySource( fam, FamilyObj( 1 ) );
    SetFamilyRange ( fam, FamilyObj( 1 ) );
    MakeReadWriteGlobal( "Z_PI_RCWAMAPPING_FAMILIES" );
    Add( Z_PI_RCWAMAPPING_FAMILIES, fam );
    MakeReadOnlyGlobal( "Z_PI_RCWAMAPPING_FAMILIES" );

    return fam;
  end );

#############################################################################
##
#F  RcwaMappingsOfGFqxFamily( <R> )
##
##  Returns the family of all rcwa mappings of a given polynomial ring <R>
##  in one variable over a finite field.
##
InstallGlobalFunction( RcwaMappingsOfGFqxFamily,

  function ( R )

    local  fam, x;

    if   not (     IsUnivariatePolynomialRing( R )
               and IsFiniteFieldPolynomialRing( R ) )
    then Error("usage: RcwaMappingsOfGFqxFamily( <R> ) for a ",
               "univariate polynomial ring <R> over a finite field.\n");
    fi;
    x := IndeterminatesOfPolynomialRing( R )[ 1 ];
    fam := First( GFQX_RCWAMAPPING_FAMILIES,
                  fam -> IsIdenticalObj( UnderlyingRing( fam ), R ) );
    if fam <> fail then return fam; fi;
    fam := NewFamily( Concatenation( "RcwaMappingsFamily( ",
                                      String( R ), " )" ),
                      IsRcwaMappingOfGFqx,
                      CanEasilySortElements, CanEasilySortElements );
    SetUnderlyingIndeterminate( fam, x );
    SetUnderlyingRing( fam, R );
    SetFamilySource( fam, FamilyObj( x ) );
    SetFamilyRange ( fam, FamilyObj( x ) );
    MakeReadWriteGlobal( "GFQX_RCWAMAPPING_FAMILIES" );
    Add( GFQX_RCWAMAPPING_FAMILIES, fam );
    MakeReadOnlyGlobal( "GFQX_RCWAMAPPING_FAMILIES" );

    return fam;
  end );

#############################################################################
##
#F  RcwaMappingsFamily( <R> ) . . . family of rcwa mappings over the ring <R>
##
InstallGlobalFunction( RcwaMappingsFamily,

  function ( R )

    if   IsIntegers( R ) then return RcwaMappingsOfZFamily;
    elif IsZxZ( R )      then return RcwaMappingsOfZxZFamily;
    elif IsZ_pi( R )     then return RcwaMappingsOfZ_piFamily( R );
    elif IsUnivariatePolynomialRing( R ) and IsFiniteFieldPolynomialRing( R )
    then return RcwaMappingsOfGFqxFamily( R );
    else Error("Sorry, rcwa mappings over ",R," are not yet implemented.\n");
    fi;
  end );

#############################################################################
##
#F  RcwaMappingsType( <R> ) . . . . . . . . . .  filter: rcwa mappings of <R>
##
InstallGlobalFunction( RcwaMappingsType,

  function ( R )
    if   IsIntegers( R ) then return IsRcwaMappingOfZ;
    elif IsZxZ( R )      then return IsRcwaMappingOfZxZ;
    elif IsZ_pi( R )     then return IsRcwaMappingOfZ_pi;
    elif IsUnivariatePolynomialRing( R ) and IsFiniteFieldPolynomialRing( R )
    then return IsRcwaMappingOfGFqx;
    else return fail; fi;
  end );

#############################################################################
##
#S  The methods for the general-purpose constructor for rcwa mappings. //////
##
#############################################################################

#############################################################################
##
#F  RCWAMAPPING_COMPRESS_COEFFICIENT_LIST( <coeffs> ) . . . . . . . . utility
##
##  This function takes care of that equal coefficient triples are always
##  also identical, in order to save memory.
##
BindGlobal( "RCWAMAPPING_COMPRESS_COEFFICIENT_LIST",

  function ( coeffs )

    local  cset, i;

    cset := Set(coeffs);
    for i in [1..Length(coeffs)] do
      coeffs[i] := cset[PositionSorted(cset,coeffs[i])];
    od;
  end );

#############################################################################
##
#M  RcwaMapping( <R>, <modulus>, <coeffs> ) . . . .  method (a) in the manual
##
InstallMethod( RcwaMapping,
               "rcwa mapping by ring, modulus and coefficients (RCWA)",
               ReturnTrue, [ IsRing, IsRingElement, IsList ], 0,

  function ( R, modulus, coeffs )

    if not modulus in R then TryNextMethod(); fi;
    if   IsIntegers(R) or IsZ_pi(R)
    then return RcwaMapping(R,coeffs);
    elif IsPolynomialRing(R)
    then return RcwaMapping(Size(LeftActingDomain(R)),modulus,coeffs);
    else TryNextMethod(); fi;
  end );

#############################################################################
##
#M  RcwaMapping( <R>, <modulus>, <coeffs> ) . . . .  method (a) in the manual
##
InstallMethod( RcwaMapping,
               "rcwa mapping by ring = Z^2, modulus and coefficients (RCWA)",
               ReturnTrue, [ IsRowModule, IsMatrix, IsList ], 0,

  function ( R, modulus, coeffs )

    local  residues, errormessage, i;

    errormessage := Concatenation("construction of an rcwa mapping of Z^2:",
                                  "\nmathematically incorrect arguments.\n");

    if   not IsZxZ(R) or DimensionsMat(modulus) <> [2,2]
      or not ForAll(Flat(modulus),IsInt) or DeterminantMat(modulus) = 0
      or Length(coeffs) <> DeterminantMat(modulus)
      or not ForAll(coeffs,IsList)
      or not Set(List(coeffs,Length)) in [[2],[3]]
      or not (    Length(coeffs[1])=2
              and ForAll( coeffs, c ->    c[1] in R and IsList(c[2])
                                      and Length(c[2])=3
                                      and IsMatrix(c[2][1])
                                      and ForAll(Flat(c[2][1]),IsInt)
                                      and c[2][2] in R
                                      and IsInt(c[2][3]) and c[2][3] <> 0 )
            or    Length(coeffs[1])=3
              and ForAll( coeffs, c ->    IsMatrix(c[1])
                                      and ForAll(Flat(c[1]),IsInt)
                                      and c[2] in R
                                      and IsInt(c[3]) and c[3] <> 0 ) )
    then Error(errormessage); return fail; fi;

    modulus  := HermiteNormalFormIntegerMat(modulus);
    residues := AllResidues(R,modulus);

    if Length(coeffs[1]) = 2 then
      for i in [1..Length(coeffs)] do
        coeffs[i][1] := coeffs[i][1] mod modulus;
      od;
      Sort( coeffs, function ( c1, c2 ) return c1[1] < c2[1]; end );
      if   List(coeffs,c->c[1]) <> residues
      then Error(errormessage); return fail; fi;
      coeffs := List(coeffs,c->c[2]);
    fi;

    if   not ForAll( [1..Length(residues)],
                     i ->   ( residues[i]*coeffs[i][1] + coeffs[i][2] )
                          mod coeffs[i][3] = [ 0, 0 ] )
    then Error(errormessage); return fail; fi;

    return RcwaMappingNC(R,modulus,coeffs);
  end );

#############################################################################
##
#M  RcwaMappingNC( <R>, <modulus>, <coeffs> ) . . NC-method (a) in the manual
##
InstallMethod( RcwaMappingNC,
               "rcwa mapping by ring, modulus and coefficients (RCWA)",
               ReturnTrue, [ IsRing, IsRingElement, IsList ], 0,

  function ( R, modulus, coeffs )

    if not modulus in R then TryNextMethod(); fi;
    if   IsIntegers(R) or IsZ_pi(R)
    then return RcwaMappingNC(R,coeffs);
    elif IsPolynomialRing(R)
    then return RcwaMappingNC(Size(LeftActingDomain(R)),modulus,coeffs);
    else TryNextMethod(); fi;
  end );

#############################################################################
##
#M  RcwaMappingNC( <R>, <modulus>, <coeffs> ) . . NC-method (a) in the manual
##
InstallMethod( RcwaMappingNC,
               "rcwa mapping by ring = Z^2, modulus and coefficients (RCWA)",
               ReturnTrue, [ IsRowModule, IsMatrix, IsList ], 0,

  function ( R, modulus, coeffs )

    local  ReduceRcwaMappingOfZxZ, result;

    ReduceRcwaMappingOfZxZ := function ( f )

      local  m, c, d, divs, res, resRed, mRed, cRed,
             nraffs, identres, pos, i;

      m := f!.modulus; c := f!.coeffs;
      for i in [1..Length(c)] do
        c[i] := c[i]/Gcd(Flat(c[i]));
        if c[i][3] < 0 then c[i] := -c[i]; fi;
      od;
      nraffs := Length(Set(c));
      res    := AllResidues(R,m);
      divs   := Superlattices(m);
      mRed := m; cRed := c;
      for d in divs do
        if DeterminantMat(d) < nraffs then continue; fi;
        resRed   := List(res,r->r mod d);
        identres := EquivalenceClasses([1..Length(c)],i->resRed[i]);
        if ForAll(identres,res->Length(Set(c{res}))=1) then
          mRed   := d;
          pos    := List(identres,cl->cl[1]);
          resRed := res{pos};
          cRed   := Permuted(c{pos},SortingPerm(resRed));
          break;
        fi;
      od;
      RCWAMAPPING_COMPRESS_COEFFICIENT_LIST(cRed);
      f!.modulus := Immutable(mRed);
      f!.coeffs  := Immutable(cRed);
    end;

    if not IsZxZ( R ) then TryNextMethod( ); fi;

    modulus := HermiteNormalFormIntegerMat( modulus );

    result := Objectify( NewType( RcwaMappingsOfZxZFamily,
                                  IsRcwaMappingOfZxZInStandardRep ),
                         rec( modulus := modulus,
                              coeffs  := coeffs ) );
    SetSource( result, R );
    SetRange ( result, R );

    ReduceRcwaMappingOfZxZ( result );

    return result;
  end );

#############################################################################
##
#M  RcwaMapping( <R>, <coeffs> ) . . . . . . . . . . method (b) in the manual
##
InstallMethod( RcwaMapping,
               "rcwa mapping by ring and coefficients (RCWA)",
               ReturnTrue, [ IsRing, IsList ], 0,

  function ( R, coeffs )

    if   IsIntegers(R)
    then return RcwaMapping(coeffs);
    elif IsZ_pi(R)
    then return RcwaMapping(NoninvertiblePrimes(R),coeffs);
    else TryNextMethod(); fi;
  end );

#############################################################################
##
#M  RcwaMappingNC( <R>, <coeffs> ) . . . . . . .  NC-method (b) in the manual
##
InstallMethod( RcwaMappingNC,
               "rcwa mapping by ring and coefficients (RCWA)",
               ReturnTrue, [ IsRing, IsList ], 0,

  function ( R, coeffs )

    if   IsIntegers(R)
    then return RcwaMappingNC(coeffs);
    elif IsZ_pi(R)
    then return RcwaMappingNC(NoninvertiblePrimes(R),coeffs);
    else TryNextMethod(); fi;
  end );

#############################################################################
##
#M  RcwaMapping( <coeffs> ) . . . . . . . . . . . .  method (c) in the manual
##
InstallMethod( RcwaMapping,
               "rcwa mapping of Z by coefficients (RCWA)",
               true, [ IsList ], 10,

  function ( coeffs )

    local  quiet;

    if not IsList( coeffs[1] ) or not IsInt( coeffs[1][1] )
      or Length( coeffs[1] ) = 5
    then TryNextMethod( ); fi;
    quiet := ValueOption("BeQuiet") = true;
    if not (     ForAll(Flat(coeffs),IsInt)
             and ForAll(coeffs, IsList)
             and ForAll(coeffs, c -> Length(c) = 3)
             and ForAll([0..Length(coeffs) - 1],
                        n -> coeffs[n + 1][3] <> 0 and
                             (n * coeffs[n + 1][1] + coeffs[n + 1][2])
                             mod coeffs[n + 1][3] = 0 and
                             (  (n + Length(coeffs)) * coeffs[n + 1][1] 
                              +  coeffs[n + 1][2])
                             mod coeffs[n + 1][3] = 0))
    then if quiet then return fail; fi;
         Error("the coefficients ",coeffs," do not define a proper ",
               "rcwa mapping of Z.\n");
    fi;
    return RcwaMappingNC( coeffs );
  end );

#############################################################################
##
#M  RcwaMappingNC( <coeffs> ) . . . . . . . . . . NC-method (c) in the manual
##
InstallMethod( RcwaMappingNC,
               "rcwa mapping of Z by coefficients (RCWA)",
               true, [ IsList ], 10,

  function ( coeffs )

    local  ReduceRcwaMappingOfZ, Result;

    ReduceRcwaMappingOfZ := function ( f )

      local  c, m, fact, p, cRed, cRedBuf, n;

      c := f!.coeffs; m := f!.modulus;
      for n in [1..Length(c)] do
        c[n] := c[n]/Gcd(c[n]);
        if c[n][3] < 0 then c[n] := -c[n]; fi;
      od;
      fact := Set(FactorsInt(m)); cRed := c;
      for p in fact do
        repeat
          cRedBuf := StructuralCopy(cRed);
          cRed := List([1..p], i -> cRedBuf{[(i - 1) * m/p + 1 .. i * m/p]});
          if   Length(Set(cRed)) = 1
          then cRed := cRed[1]; m := m/p; else cRed := cRedBuf; fi;
        until cRed = cRedBuf or m mod p <> 0;
      od;
      RCWAMAPPING_COMPRESS_COEFFICIENT_LIST(cRed);
      f!.coeffs  := Immutable(cRed);
      f!.modulus := Length(cRed);
    end;

    if not IsList( coeffs[1] ) or not IsInt( coeffs[1][1] )
      or Length( coeffs[1] ) = 5
    then TryNextMethod( ); fi;
    Result := Objectify( NewType(    RcwaMappingsOfZFamily,
                                     IsRcwaMappingOfZInStandardRep ),
                         rec( coeffs  := coeffs,
                              modulus := Length(coeffs) ) );
    SetSource(Result, Integers);
    SetRange (Result, Integers);
    ReduceRcwaMappingOfZ(Result);
    return Result;
  end );

#############################################################################
##
#M  RcwaMapping( <perm>, <range> ) . . . . . . . . . method (d) in the manual
##
InstallMethod( RcwaMapping,
               "rcwa mapping of Z by permutation and range (RCWA)",
               true, [ IsPerm, IsRange ], 0,

  function ( perm, range )

    local  quiet;

    quiet := ValueOption("BeQuiet") = true;
    if   Permutation(perm,range) = fail
    then if quiet then return fail; fi;
         Error("the permutation ",perm," does not act on the range ",
               range,".\n");
    fi;
    return RcwaMappingNC( perm, range );
  end );

#############################################################################
##
#M  RcwaMappingNC( <perm>, <range> ) . . . . . .  NC-method (d) in the manual
##
InstallMethod( RcwaMappingNC,
               "rcwa mapping of Z by permutation and range (RCWA)",
               true, [ IsPerm, IsRange ], 0,

  function ( perm, range )

    local  result, coeffs, min, max, m, n, r;

    min := Minimum(range); max := Maximum(range);
    m := max - min + 1; coeffs := [];
    for n in [min..max] do
      r := n mod m + 1;
      coeffs[r] := [1, n^perm - n, 1];
    od;
    result := RcwaMappingNC( coeffs );
    SetIsBijective(result,true);
    SetIsTame(result,true); SetIsIntegral(result,true);
    SetOrder(result,Order(RestrictedPerm(perm,range)));
    return result;
  end );

#############################################################################
##
#M  RcwaMapping( <modulus>, <values> ) . . . . . . . method (e) in the manual
##
InstallMethod( RcwaMapping,
               "rcwa mapping of Z by modulus and values (RCWA)",
               true, [ IsInt, IsList ], 0,

  function ( modulus, values )

    local  f, coeffs, pts, r, quiet;

    quiet := ValueOption("BeQuiet") = true;
    coeffs := [];
    for r in [1..modulus] do
      pts := Filtered(values, pt -> pt[1] mod modulus = r - 1);
      if   Length(pts) < 2
      then if quiet then return fail; fi;
           Error("the mapping is not given at at least 2 points <n> ",
                 "with <n> mod ",modulus," = ",r - 1,".\n");
      fi;
    od;
    f := RcwaMappingNC( modulus, values );
    if Mod(f) mod Div(f) <> 0 or not ForAll(values,t -> t[1]^f = t[2])
    then if quiet then return fail; fi;
         Error("the values ",values," do not define a proper ",
               "rcwa mapping of Z.\n"); 
    fi;
    return f;
  end );

#############################################################################
##
#M  RcwaMappingNC( <modulus>, <values> ) . . . .  NC-method (e) in the manual
##
InstallMethod( RcwaMappingNC,
               "rcwa mapping of Z by modulus and values (RCWA)",
               true, [ IsInt, IsList ], 0,

  function ( modulus, values )

    local  coeffs, pts, r;

    coeffs := [];
    for r in [1..modulus] do
      pts := Filtered(values, pt -> pt[1] mod modulus = r - 1);
      coeffs[r] := [  pts[1][2] - pts[2][2],
                      pts[1][2] * (pts[1][1] - pts[2][1])
                    - pts[1][1] * (pts[1][2] - pts[2][2]),
                      pts[1][1] - pts[2][1]];
    od;
    return RcwaMappingNC( coeffs );
  end );

#############################################################################
##
#M  RcwaMapping( <pi>, <coeffs> ) . . . . . . . . .  method (f) in the manual
##
InstallMethod( RcwaMapping,
               "rcwa mapping by noninvertible primes and coeff's (RCWA)",
               true, [ IsObject, IsList ], 0,

  function ( pi, coeffs )

    local  R, quiet;

    quiet := ValueOption("BeQuiet") = true;
    if IsInt(pi) then pi := [pi]; fi; R := Z_pi(pi);
    if not (     IsList(pi) and ForAll(pi,IsInt)
             and IsSubset(Union(pi,[1]),Set(Factors(Length(coeffs))))
             and ForAll(Flat(coeffs), x -> IsRat(x) and Intersection(pi,
                                        Set(Factors(DenominatorRat(x))))=[])
             and ForAll(coeffs, IsList)
             and ForAll(coeffs, c -> Length(c) = 3)
             and ForAll([0..Length(coeffs) - 1],
                        n -> coeffs[n + 1][3] <> 0 and
                             NumeratorRat(n * coeffs[n + 1][1]
                                            + coeffs[n + 1][2])
                             mod StandardAssociate(R,coeffs[n + 1][3]) = 0
                         and NumeratorRat(  (n + Length(coeffs))
                                           * coeffs[n + 1][1]
                                           + coeffs[n + 1][2])
                             mod StandardAssociate(R,coeffs[n + 1][3]) = 0))
    then if quiet then return fail; fi;
         Error("the coefficients ",coeffs," do not define a proper ",
               "rcwa mapping of Z_(",pi,").\n");
    fi;
    return RcwaMappingNC(pi,coeffs);
  end );

#############################################################################
##
#M  RcwaMappingNC( <pi>, <coeffs> ) . . . . . . . NC-method (f) in the manual
##
InstallMethod( RcwaMappingNC,
               "rcwa mapping by noninvertible primes and coeff's (RCWA)",
               true, [ IsObject, IsList ], 0,

  function ( pi, coeffs )

    local  ReduceRcwaMappingOfZ_pi, f, R, fam;

    ReduceRcwaMappingOfZ_pi := function ( f )

      local  c, m, pi, d_pi, d_piprime, divs, d, cRed, n, i;

      c := f!.coeffs; m := f!.modulus;
      pi := NoninvertiblePrimes(Source(f));
      for n in [1..Length(c)] do
        if c[n][3] < 0 then c[n] := -c[n]; fi;
        d_pi := Gcd(Product(Filtered(Factors(Gcd(NumeratorRat(c[n][1]),
                                                 NumeratorRat(c[n][2]))),
                                     p -> p in pi or p = 0)),
                    NumeratorRat(c[n][3]));
        d_piprime := c[n][3]/Product(Filtered(Factors(NumeratorRat(c[n][3])),
                                              p -> p in pi));
        c[n] := c[n] / (d_pi * d_piprime);
      od;
      divs := DivisorsInt(m); i := 1;
      repeat
        d := divs[i]; i := i + 1;
        cRed := List([1..m/d], i -> c{[(i - 1) * d + 1 .. i * d]});
      until Length(Set(cRed)) = 1;
      cRed := cRed[1];
      RCWAMAPPING_COMPRESS_COEFFICIENT_LIST(cRed);
      f!.coeffs  := Immutable(cRed);
      f!.modulus := Length(cRed);
    end;

    if IsInt(pi) then pi := [pi]; fi;
    if   not IsList(pi) or not ForAll(pi,IsInt) or not ForAll(coeffs,IsList)
    then TryNextMethod(); fi;
    R := Z_pi(pi); fam := RcwaMappingsFamily( R );
    f := Objectify( NewType( fam, IsRcwaMappingOfZ_piInStandardRep ),
                    rec( coeffs  := coeffs,
                         modulus := Length(coeffs) ) );
    SetSource(f,R); SetRange(f,R);
    ReduceRcwaMappingOfZ_pi(f);
    return f;
  end );

#############################################################################
##
#M  RcwaMapping( <q>, <modulus>, <coeffs> ) . . . .  method (g) in the manual
##
InstallMethod( RcwaMapping,
               Concatenation("rcwa mapping by finite field size, ",
                             "modulus and coefficients (RCWA)"),
               true, [ IsInt, IsPolynomial, IsList ], 0,

  function ( q, modulus, coeffs )

    local  d, x, P, p, quiet;

    quiet := ValueOption("BeQuiet") = true;
    if not (    IsPosInt(q) and IsPrimePowerInt(q) 
            and ForAll(coeffs, IsList)
            and ForAll(coeffs, c -> Length(c) = 3) 
            and ForAll(Flat(coeffs), IsPolynomial)
            and Length(Set(List(Flat(coeffs),
                                IndeterminateNumberOfLaurentPolynomial)))=1)
    then if quiet then return fail; fi;
         Error("see RCWA manual for information on how to construct\n",
               "an rcwa mapping of a polynomial ring.\n");
    fi;
    d := DegreeOfLaurentPolynomial(modulus);
    x := IndeterminateOfLaurentPolynomial(coeffs[1][1]);
    P := AllGFqPolynomialsModDegree(q,d,x);
    if not ForAll([1..Length(P)],
                  i -> IsZero(   (coeffs[i][1]*P[i] + coeffs[i][2])
                              mod coeffs[i][3]))
    then Error("the coefficients ",coeffs," do not define a proper ",
               "rcwa mapping.\n");
    fi;
    return RcwaMappingNC( q, modulus, coeffs );
  end );

#############################################################################
##
#M  RcwaMappingNC( <q>, <modulus>, <coeffs> ) . . NC-method (g) in the manual
##
InstallMethod( RcwaMappingNC,
               Concatenation("rcwa mapping by finite field size, ",
                             "modulus and coefficients (RCWA)"),
               true, [ IsInt, IsPolynomial, IsList ], 0,

  function ( q, modulus, coeffs )

    local  ReduceRcwaMappingOfGFqx, f, R, fam, ind;

    ReduceRcwaMappingOfGFqx := function ( f )

      local  c, m, F, q, x, deg, r, fact, d, degd,
             sigma, csorted, numresred, numresd, mred, rred,
             n, l, i;

      c := f!.coeffs; m := f!.modulus;
      for n in [1..Length(c)] do
        d := Gcd(c[n]);
        c[n] := c[n]/(d * LeadingCoefficient(c[n][3]));
      od;
      deg := DegreeOfLaurentPolynomial(m);
      F := CoefficientsRing(UnderlyingRing(FamilyObj(f)));
      q := Size(F);
      x := UnderlyingIndeterminate(FamilyObj(f));
      r := AllGFqPolynomialsModDegree(q,deg,x);
      fact := Difference(Factors(m),[One(m)]);
      for d in fact do 
        degd := DegreeOfLaurentPolynomial(d);
        repeat
          numresd := q^degd; numresred := q^(deg-degd);
          mred  := m/d;
          rred  := List(r, P -> P mod mred);
          sigma := SortingPerm(rred);
          csorted := Permuted(c,sigma);
          if ForAll([1..numresred],
                    i->Length(Set(csorted{[(i-1)*numresd+1..i*numresd]}))=1)
          then m   := mred;
               deg := deg - degd;
               r := AllGFqPolynomialsModDegree(q,deg,x);
               c := csorted{[1, 1 + numresd .. 1 + (numresred-1) * numresd]};
          fi;
        until m <> mred or not IsZero(m mod d);
      od;
      RCWAMAPPING_COMPRESS_COEFFICIENT_LIST(c);
      f!.coeffs  := Immutable(c);
      f!.modulus := m;
    end;

    ind := IndeterminateNumberOfLaurentPolynomial( coeffs[1][1] );
    R   := PolynomialRing( GF( q ), [ ind ] );
    fam := RcwaMappingsFamily( R );
    f   := Objectify( NewType( fam, IsRcwaMappingOfGFqxInStandardRep ),
                      rec( coeffs  := coeffs,
                           modulus := modulus ) );
    SetUnderlyingField( f, CoefficientsRing( R ) );
    SetSource( f, R ); SetRange( f, R );
    ReduceRcwaMappingOfGFqx( f );
    return f;
  end );

#############################################################################
##
#M  RcwaMapping( <P1>, <P2> ) . . . . . . . . . . .  method (h) in the manual
##
InstallMethod( RcwaMapping,
               "rcwa mapping by two class partitions (RCWA)",
               true, [ IsList, IsList ], 0,

  function ( P1, P2 )

    local  result;

    if not (     ForAll(Concatenation(P1,P2),IsResidueClass)
             and Length(P1) = Length(P2)
             and Sum(List(P1,Density)) = 1
             and Union(P1) = UnderlyingRing(FamilyObj(P1[1])))
    then TryNextMethod(); fi;
    result := RcwaMappingNC(P1,P2);
    IsBijective(result);
    return result;
  end );

#############################################################################
##
#M  RcwaMappingNC( <P1>, <P2> ) . . . . . . . . . NC-method (h) in the manual
##
InstallMethod( RcwaMappingNC,
               "rcwa mapping by two class partitions (RCWA)",
               true, [ IsList, IsList ], 0,

  function ( P1, P2 )

    local  R, coeffs, m, res, r1, m1, r2, m2, i, j;

    if not IsResidueClassUnion(P1[1]) then TryNextMethod(); fi;
    R := UnderlyingRing(FamilyObj(P1[1]));
    m := Lcm(R,List(P1,Modulus)); res := AllResidues(R,m);
    coeffs := List(res,r->[1,0,1]*One(R));
    for i in [1..Length(P1)] do
      r1 := Residue(P1[i]); m1 := Modulus(P1[i]);
      r2 := Residue(P2[i]); m2 := Modulus(P2[i]);
      for j in Filtered([1..Length(res)],j->res[j] mod m1 = r1) do
        coeffs[j] := [m2,m1*r2-m2*r1,m1];
      od;
    od;
    return RcwaMappingNC(R,m,coeffs);
  end );

#############################################################################
##
#M  RcwaMappingNC( <P1>, <P2> ) . . . .  NC-method (h) in the manual, for Z^2
##
InstallMethod( RcwaMappingNC,
               "rcwa mapping by two class partitions of Z^2 (RCWA)",
               true, [ IsList, IsList ], 0,

  function ( P1, P2 )

    local  R, coeffs, m, res, affectedpos, t, r1, m1, r2, m2, i, j;

    if not IsResidueClassUnionOfZxZ(P1[1]) then TryNextMethod(); fi;
    R := UnderlyingRing(FamilyObj(P1[1]));
    m := Lcm(List(P1,Modulus)); res := AllResidues(R,m);
    coeffs := List(res,r->[[[1,0],[0,1]],[0,0],1]);
    for i in [1..Length(P1)] do
      r1 := Residue(P1[i]); m1 := Modulus(P1[i]);
      r2 := Residue(P2[i]); m2 := Modulus(P2[i]);
      affectedpos := Filtered([1..Length(res)],j->res[j] mod m1 = r1);
      t := [m1^-1*m2,r2-r1*m1^-1*m2,1];
      t := t * Lcm(List(Flat(t),DenominatorRat));
      for i in affectedpos do coeffs[i] := t; od;
    od;
    return RcwaMappingNC(R,m,coeffs);
  end );

#############################################################################
##
#M  RcwaMapping( <cycles> ) . . . . . . . . . . . .  method (i) in the manual
##
InstallMethod( RcwaMapping,
               "rcwa mapping by class cycles (RCWA)", true, [ IsList ], 0,

  function ( cycles )

    local  CheckClassCycles, R;

    CheckClassCycles := function ( R, cycles )

      if not (    ForAll(cycles,IsList)
              and ForAll(Flat(cycles),S->IsResidueClass(S)
              and IsSubset(R,S)))
         or  ForAny(Combinations(Flat(cycles),2),s->Intersection(s) <> [])
      then Error("there is no rcwa mapping of ",R," having the class ",
                 "cycles ",cycles,".\n"); 
      fi;
    end;

    if   not IsList(cycles[1]) or not IsResidueClass(cycles[1][1])
    then TryNextMethod(); fi;
    R := cycles[1][1];
    if   not (IsRing(R) or IsZxZ(R))
    then R := UnderlyingRing(FamilyObj(R)); fi;
    CheckClassCycles(R,cycles);
    return RcwaMappingNC(cycles);
  end );

#############################################################################
##
#M  RcwaMapping( <cycles> ) .  method (i), variation for rc. with fixed rep's
##
InstallMethod( RcwaMapping,
               "rcwa mapping by class cycles (fixed rep's) (RCWA)",
               true, [ IsList ], 0,

  function ( cycles )

    local  CheckClassCycles, R;

    CheckClassCycles := function ( R, cycles )

      if not (    ForAll(cycles,IsList)
              and ForAll(Flat(cycles),S->IsResidueClassWithFixedRep(S)
              and UnderlyingRing(FamilyObj(S)) = R))
         or  ForAny(Combinations(Flat(cycles),2),
                    s->Intersection(List([s[1],s[2]],
                                    AsOrdinaryUnionOfResidueClasses)) <> [])
      then Error("there is no rcwa mapping of ",R," having the class ",
                 "cycles ",cycles,".\n"); 
      fi;
    end;

    if   not IsList(cycles[1])
      or not IsResidueClassWithFixedRepresentative(cycles[1][1])
    then TryNextMethod(); fi;
    R := UnderlyingRing(FamilyObj(cycles[1][1]));
    CheckClassCycles(R,cycles);
    return RcwaMappingNC(cycles);
  end );

#############################################################################
##
#M  RcwaMappingNC( <cycles> ) . . . . . . . . . . NC-method (i) in the manual
##
InstallMethod( RcwaMappingNC,
               "rcwa mapping by class cycles (RCWA)", true, [ IsList ], 0,

  function ( cycles )

    local  result, R, coeffs, m, res, cyc, pre, im, affectedpos,
           r1, r2, m1, m2, pos, i;

    if    not IsResidueClass(cycles[1][1])
      and not IsResidueClassWithFixedRepresentative(cycles[1][1])
    then TryNextMethod(); fi;

    R := cycles[1][1];
    if   not (IsRing(R) or IsZxZ(R))
    then R := UnderlyingRing(FamilyObj(R)); fi;
    if not IsIntegers(R) and not IsZ_pi(R)
      and not (     IsUnivariatePolynomialRing(R)
                and IsFiniteFieldPolynomialRing(R))
    then TryNextMethod(); fi;

    m      := Lcm(List(Union(cycles),Modulus));
    res    := AllResidues(R,m);
    coeffs := List(res,r->[1,0,1]*One(R));
    for cyc in cycles do
      if Length(cyc) <= 1 then continue; fi;
      for pos in [1..Length(cyc)] do
        pre := cyc[pos]; im := cyc[pos mod Length(cyc) + 1];
        r1 := Residue(pre); m1 := Modulus(pre);
        r2 := Residue(im);  m2 := Modulus(im);
        affectedpos := Filtered([1..Length(res)],
                                i->res[i] mod m1 = r1 mod m1);
        for i in affectedpos do coeffs[i] := [m2,m1*r2-m2*r1,m1]; od;
      od;
    od;
    if   IsIntegers(R)
    then result := RcwaMappingNC(coeffs);
    elif IsZ_pi(R)
    then result := RcwaMappingNC(R,coeffs);
    elif IsPolynomialRing(R)
    then result := RcwaMappingNC(R,Lcm(List(Flat(cycles),Modulus)),coeffs);
    fi;
    Assert(4,Order(result)=Lcm(List(cycles,Length)));
    SetIsBijective(result,true); SetIsTame(result,true);
    SetOrder(result,Lcm(List(cycles,Length)));
    return result;
  end );

#############################################################################
##
#M  RcwaMappingNC( <cycles> ) . . . . .  NC-method (i) in the manual, for Z^2
##
InstallMethod( RcwaMappingNC,
               "rcwa mapping of Z^2 by class cycles (RCWA)", true,
               [ IsList ], 0,

  function ( cycles )

    local  result, R, coeffs, m, res, cyc, pre, im, affectedpos, t,
           r1, r2, m1, m2, pos, i;

    if   not IsResidueClass(cycles[1][1])
      or not IsResidueClassUnionOfZxZ(cycles[1][1])
    then TryNextMethod(); fi;

    R      := UnderlyingRing(FamilyObj(cycles[1][1]));
    m      := Lcm(List(Union(cycles),Modulus));
    res    := AllResidues(R,m);
    coeffs := List(res,r->[[[1,0],[0,1]],[0,0],1]);
    for cyc in cycles do
      if Length(cyc) <= 1 then continue; fi;
      for pos in [1..Length(cyc)] do
        pre := cyc[pos]; im := cyc[pos mod Length(cyc) + 1];
        r1 := Residue(pre); m1 := Modulus(pre);
        r2 := Residue(im);  m2 := Modulus(im);
        affectedpos := Filtered([1..Length(res)],i->res[i] mod m1 = r1);
        t := [m1^-1*m2,r2-r1*m1^-1*m2,1];
        t := t * Lcm(List(Flat(t),DenominatorRat));
        for i in affectedpos do coeffs[i] := t; od;
      od;
    od;
    result := RcwaMapping(R,m,coeffs);
    Assert(4,Order(result)=Lcm(List(cycles,Length)));
    SetIsBijective(result,true); SetIsTame(result,true);
    SetOrder(result,Lcm(List(cycles,Length)));
    return result;
  end );

#############################################################################
##
#M  RcwaMapping( <expression> ) . . . . . . . . . .  method (j) in the manual
##
InstallMethod( RcwaMapping,
               "rcwa mapping of Z by expression, given as a string (RCWA)",
               true, [ IsString ], 0,

  function ( expression )
    if   IsSubset( "0123456789-,()crs^*/", expression )
    then return RcwaMappingNC( expression );
    else TryNextMethod(); fi;
  end );

#############################################################################
##
#M  RcwaMappingNC( <expression> ) . . . . . . . . NC-method (j) in the manual
##
InstallMethod( RcwaMappingNC,
               "rcwa mapping of Z by expression, given as a string (RCWA)",
               true, [ IsString ], 0,

  function ( expression )

    local  ValueExpression, ValueElementaryExpression;

    ValueElementaryExpression := function ( exp )

      local  ints, cycle, i;

      if IsSubset("0123456789-",exp) then return Int(exp); fi;
      ints := List(Filtered(SplitString(exp,"()[],crs"),
                            s->s<>""),Int);
      if   Length(ints) = 2 then
        if   's' in exp
        then return ClassShift(ints);
        else return ClassReflection(ints); fi;
      elif Length(ints) = 4 then
        return ClassTransposition(ints);
      elif Length(ints) mod 2 = 0 then
        cycle := [];
        for i in [1,3..Length(ints)-1] do
          Add(cycle,ResidueClass(ints[i],ints[i+1]));
        od;
        return RcwaMapping([cycle]);
      else Error("unknown type of rcwa permutation\n"); fi;
    end;

    ValueExpression := function ( exp )

      local  brackets, parts, part, operators,
             values, valuesexp, value, i, j;

      if   IsSubset("0123456789-,()crs",exp)
      then return ValueElementaryExpression(exp); fi;

      brackets := 0; parts := []; operators := []; part := "";
      for i in [1..Length(exp)] do
        Add(part,exp[i]);
        if   exp[i] = '(' then
          brackets := brackets + 1;
        elif exp[i] = ')' then
          brackets := brackets - 1;
          if brackets = 0 then
            Add(parts,part); part := "";
          fi;
        elif brackets = 0 and exp[i] in "*/^" then
          Add(operators,exp[i]);
          Add(parts,part{[1..Length(part)-1]}); part := "";
        fi;
      od;
      Add(parts,part);
      parts := Filtered(parts,part->Intersection(part,"0123456789")<>"");
      for i in [1..Length(parts)] do
        if   parts[i][1] = '('
        then parts[i] := parts[i]{[2..Length(parts[i])-1]}; fi;
      od;
      values    := List(parts,ValueExpression);
      valuesexp := ShallowCopy(values);
      for i in [1..Length(operators)] do
        if operators[i] = '^' then
          valuesexp[i] := valuesexp[i]^valuesexp[i+1];
          valuesexp[i+1] := fail;
        fi;
      od;
      valuesexp := Filtered(valuesexp,val->val<>fail);

      operators := Filtered(operators,op->op<>'^');
      value     := valuesexp[1];
      for i in [1..Length(operators)] do
        if   operators[i] = '*'
        then value := value * valuesexp[i+1];
        elif operators[i] = '/'
        then value := value / valuesexp[i+1];
        else Error("RcwaMapping: unknown operator: ",operators[i],"\n"); fi; 
      od;

      return value;
    end;

    return ValueExpression( BlankFreeString( expression ) );
  end );

#############################################################################
##
#M  RcwaMapping( <R>, <f>, <g> ) . rcwa mapping of Z^2 by two rcwa map's of Z
#M  RcwaMapping( <f>, <g> )
##
InstallMethod( RcwaMapping,
               "rcwa mapping of Z^2 by two rcwa mappings of Z (RCWA)",
               ReturnTrue,
               [ IsRowModule, IsRcwaMappingOfZ, IsRcwaMappingOfZ ], 0,
               function ( R, f, g )
                 if IsZxZ(R) then return RcwaMappingNC(f,g);
                             else TryNextMethod(); fi;
               end );

InstallMethod( RcwaMapping,
               "rcwa mapping of Z^2 by two rcwa mappings of Z (RCWA)",
               IsIdenticalObj, [ IsRcwaMappingOfZ, IsRcwaMappingOfZ ], 0,
               function ( f, g ) return RcwaMappingNC(f,g); end );

#############################################################################
##
#M  RcwaMappingNC( <f>, <g> ) . rcwa mapping of Z^2 by two rcwa mappings of Z
##
InstallMethod( RcwaMappingNC,
               "rcwa mapping of Z^2 by two rcwa mappings of Z (RCWA)",
               IsIdenticalObj, [ IsRcwaMappingOfZ, IsRcwaMappingOfZ ], 0,

  function ( f, g )

    local  result, m, mf, mg, c, cf, cg, res, r, t, d, d1, d2;

    mf := Modulus(f);      mg  := Modulus(g);
    m  := [[mf,0],[0,mg]]; res := AllResidues(Integers^2,m);
    cf := Coefficients(f); cg  := Coefficients(g); c := [];
    for r in res do
      t := [cf[r[1]+1],cg[r[2]+1]];
      d := Lcm(t[1][3],t[2][3]); d1 := d/t[1][3]; d2 := d/t[2][3];
      Add(c,[[[t[1][1]*d1,0],[0,t[2][1]*d2]],[t[1][2]*d1,t[2][2]*d2],d]);
    od;
    result := RcwaMapping(Integers^2,m,c);
    if   ForAny([f,g],HasIsClassTransposition)
    then IsClassTransposition(result); fi;
    return result;
  end );

#############################################################################
##
#M  RcwaMapping( <coeffs> ) . . .  rcwa mapping of Z in sparse representation
##
InstallMethod( RcwaMapping,
               "rcwa mapping of Z in sparse representation (RCWA)",
               true, [ IsList ], 5,

  function ( coeffs )

    local  result;

    if   not ForAll(coeffs,c->IsList(c) and Length(c)=5 and ForAll(c,IsInt)) 
    then TryNextMethod(); fi;

    if ForAny(coeffs,c->c[2] = 0 or c[5] = 0) then
      Error("zero in modulus or denominator not allowed.\n");
      return fail;
    fi;
    if ForAny(Combinations([1..Length(coeffs)],2),
              i->(    coeffs[i[1]][1]-coeffs[i[2]][1])
                  mod Gcd(coeffs[i[1]][2],coeffs[i[2]][2]) = 0)
    then Error("rcwa mapping must be single-valued.\n"); return fail; fi;
    if   Sum(List(coeffs,c->1/c[2])) <> 1
    then Error("rcwa mapping must be total.\n"); return fail; fi;

    return RcwaMappingNC(coeffs);
  end );

#############################################################################
##
#M  RcwaMappingNC( <coeffs> ) . .  rcwa mapping of Z in sparse representation
##
InstallMethod( RcwaMappingNC,
               "rcwa mapping of Z in sparse representation (RCWA)",
               true, [ IsList ], 5,

  function ( coeffs )

    local  result, modulus, coeffset, coeffcoll, cls, redcls, cl,
           d, d_red, res, res_d, div, m, r, r_red, range, i, j;

    if   not ForAll(coeffs,c->IsList(c) and Length(c)=5)
      or not ForAll(coeffs[1],IsInt)
    then TryNextMethod(); fi;

    for i in [1..Length(coeffs)] do
      coeffs[i][2] := AbsInt(coeffs[i][2]);
      coeffs[i][1] := coeffs[i][1] mod coeffs[i][2];
      coeffs[i]{[3..5]} := coeffs[i]{[3..5]}/Gcd(coeffs[i]{[3..5]});
      if coeffs[i][5] < 0 then coeffs[i]{[3..5]} := -coeffs[i]{[3..5]}; fi;
    od;

    coeffset  := Set(List(coeffs,c->c{[3..5]}));
    coeffcoll := List(coeffset,c->[]);
    for i in [1..Length(coeffs)] do
      j := PositionSorted(coeffset,coeffs[i]{[3..5]});
      Add(coeffcoll[j],coeffs[i]);
    od;

    coeffs := [];
    for i in [1..Length(coeffset)] do

      cls := Set(coeffcoll[i],c->c{[1,2]});
      if Length(cls) > 1 then

        d_red  := Gcd(DifferencesList(List(cls,cl->cl[1])));
        d_red  := Gcd(d_red,Gcd(List(cls,cl->cl[2])));
        r_red  := cls[1][1] mod d_red;
        redcls := List(cls,cl->[(cl[1]-r_red)/d_red,cl[2]/d_red]);
        m      := Lcm(List(redcls,cl->cl[2]));
        res    := [];
        for cl in redcls do
          for j in [0..m/cl[2]-1] do Add(res,cl[1]+j*cl[2]); od;
        od;
        res := Set(res);

        redcls := [];
        div    := DivisorsInt(m);
        for d in div do
          res_d := Set(res mod d);
          for r in res_d do
            range := List([0..m/d-1],j->r+j*d);
            if IsSubset(res,range) then
              Add(redcls,[r,d]);
              res := Difference(res,range);
              if res = [] then break; fi;
              if Length(res) < m/d then break; fi;
            fi;
          od;
          if res = [] then break; fi;
        od;

        cls := List(redcls,cl->[cl[1]*d_red+r_red,cl[2]*d_red]);

      # Erroneous reduction process: doesn't cope with cases like
      # 0(3) U 1(6) U 5(6) = 1(2) U 0(6).
      # 
      # repeat
      #   oldcls := ShallowCopy(cls);
      #   mods   := Set(List(oldcls,cl->cl[2]));
      #   for m in mods do
      #     resm := List(Set(Filtered(cls,cl->cl[2]=m)),c->c[1]);
      #     if Length(resm) >= 2 then
      #       for p in Set(Factors(m)) do
      #         for r in Filtered(resm,res->res<m/p) do
      #           # range := [r,r+m/p..m-m/p+r]; -> range restriction (<2^28)
      #           range := List([0..p-1],j->j*(m/p)+r);
      #           if IsSubset(resm,range) then
      #             cls := Difference(cls,List(range,r->[r,m]));
      #             Add(cls,[r,m/p]);
      #           fi;
      #         od;
      #       od;
      #     fi;
      #   od;
      #   cls := Set(cls);
      # until cls = oldcls;

      fi;

      for cl in cls do
        Add(coeffs,Concatenation(cl,coeffset[i]));
      od;

    od;

    # SortParallel(List(coeffs,c->[c[2],c[1]]),coeffs); -> problematic?
    coeffs  := Set(coeffs);
    modulus := Lcm(List(coeffs,c->c[2]));

    MakeImmutable(coeffs);
    result := Objectify( NewType( RcwaMappingsOfZFamily,
                                  IsRcwaMappingOfZ
                              and IsRcwaMappingSparseRep ),
                         rec( modulus := modulus,
                              coeffs  := coeffs ) );
    SetSource(result,Integers);
    SetRange (result,Integers);

    return result;
  end );

#############################################################################
##
#S  Conversion of rcwa mappings between standard and sparse representation. /
##
#############################################################################

#############################################################################
##
#M  SparseRepresentation( <f> )
##
InstallMethod( SparseRepresentation,
               "for rcwa mappings in standard representation (RCWA)",
               true, [ IsRcwaMappingInStandardRep ], 0,

  function ( f )

    local  result, attrs, attrsetters, props, propsetters, attr, prop,
           R, coeffs, src, sparse, cl, r, m, c, i;

    R      := Source(f);
    coeffs := f!.coeffs;

    src := Set(Flat(List(LargestSourcesOfAffineMappings(f),
                         AsUnionOfFewClasses)));
    sparse := [];
    for cl in src do
      r := Residue(cl); m := Modulus(cl);
      if   IsRcwaMappingOfZ(f)
      then c := coeffs[r+1];
      else c := coeffs[Position(AllResidues(R,m),r)]; fi;
      Add(sparse,[r,m,c[1],c[2],c[3]]);
    od;

    sparse := Set(sparse);

    result := Objectify(NewType(FamilyObj(f),RcwaMappingsType(R)
                                         and IsRcwaMappingSparseRep),
                        rec(modulus := f!.modulus, coeffs := sparse));

    # Copy over representation-independent attributes and properties.

    attrs := List( Intersection( RCWA_REP_INDEPENDENT_ATTRIBUTES,
                                 KnownAttributesOfObject( f ) ),
                   ValueGlobal );
    attrsetters := List(attrs,Setter);
    for i in [1..Length(attrs)] do attrsetters[i](result,attrs[i](f)); od;

    props := List( Intersection( RCWA_REP_INDEPENDENT_PROPERTIES,
                                 KnownPropertiesOfObject( f ) ),
                   ValueGlobal );
    propsetters := List(props,Setter);
    for i in [1..Length(props)] do propsetters[i](result,props[i](f)); od;

    return result;
  end );

#############################################################################
##
#M  SparseRepresentation( <f> )
##
InstallMethod( SparseRepresentation,
               "for rcwa mappings in sparse representation (RCWA)",
               true, [ IsRcwaMappingInSparseRep ], 0, f -> f );

#############################################################################
##
#M  StandardRepresentation( <f> )
##
InstallMethod( StandardRepresentation,
               "for rcwa mappings in sparse representation (RCWA)",
               true, [ IsRcwaMappingInSparseRep ], 0,

  function ( f )

    local  result, attrs, attrsetters, props, propsetters,
           R, modulus, coeffs, src, sparse, res, cl, r, m, c, n, i;

    R       := Source(f);
    modulus := f!.modulus;
    sparse  := f!.coeffs;

    if IsRing(R) then
      coeffs := List([1..NumberOfResidues(R,modulus)],r->[1,0,1]*One(R));
    elif IsZxZ(R) then
      coeffs := List([1..NumberOfResidues(R,modulus)],
                     r->[[[1,0],[0,1]],[0,0],1]);
    fi;

    res := AllResidues(R,modulus);
    for i in [1..Length(sparse)] do
      c := sparse[i];
      r := c[1]; m := c[2];
      if IsRcwaMappingOfZ(f) then
        for n in [r+1,r+m+1..modulus-m+r+1] do
          coeffs[n] := c{[3..5]};
        od;
      else
        for n in PositionsProperty(res,el->el mod m = r) do
          coeffs[n] := c{[3..5]};
        od;
      fi;
    od;

    result := Objectify(NewType(FamilyObj(f),RcwaMappingsType(R)
                                         and IsRcwaMappingStandardRep),
                        rec( modulus := modulus, coeffs := coeffs ));

    # Copy over representation-independent attributes and properties.

    attrs := List( Intersection( RCWA_REP_INDEPENDENT_ATTRIBUTES,
                                 KnownAttributesOfObject( f ) ),
                   ValueGlobal );
    attrsetters := List(attrs,Setter);
    for i in [1..Length(attrs)] do attrsetters[i](result,attrs[i](f)); od;

    props := List( Intersection( RCWA_REP_INDEPENDENT_PROPERTIES,
                                 KnownPropertiesOfObject( f ) ),
                   ValueGlobal );
    propsetters := List(props,Setter);
    for i in [1..Length(props)] do propsetters[i](result,props[i](f)); od;

    return result;
  end );

#############################################################################
##
#M  StandardRepresentation( <f> )
##
InstallMethod( StandardRepresentation,
               "for rcwa mappings in standard representation (RCWA)",
               true, [ IsRcwaMappingInStandardRep ], 0, f -> f );

#############################################################################
##
#S  ExtRepOfObj / ObjByExtRep for rcwa mappings. ////////////////////////////
##
#############################################################################

#############################################################################
##
#M  ExtRepOfObj( <f> ) . . . . . . . . . . . . . . . . . .  for rcwa mappings
##
InstallMethod( ExtRepOfObj,
               "for rcwa mappings (RCWA)", true, [ IsRcwaMapping ], 0,
               f -> [ Modulus( f ), ShallowCopy( Coefficients( f ) ) ] );

#############################################################################
##
#M  ExtRepOfObj( <ct> ) . . . . . . . . . . . . . .  for class transpositions
##
InstallMethod( ExtRepOfObj,
               "for class transpositions (RCWA)", true,
               [ IsRcwaMapping and IsClassTransposition ], 0,

  function ( ct )

    local  cls;

    cls := TransposedClasses(ct);
    return [ Residue(cls[1]), Mod(cls[1]), Residue(cls[2]), Mod(cls[2]) ];
  end );

#############################################################################
##
#M  ObjByExtRep( <fam>, <l> ) . rcwa mapping, by list [ <modulus>, <coeffs> ]
##
InstallMethod( ObjByExtRep,
               "rcwa mapping, by list [ <modulus>, <coefficients> ] (RCWA)",
               ReturnTrue, [ IsFamily, IsList ], 0,

  function ( fam, l )

    local  R;

    if not HasUnderlyingRing(fam) or Length(l) <> 2 then TryNextMethod(); fi;
    R := UnderlyingRing(fam);
    if fam <> RcwaMappingsFamily(R) then TryNextMethod(); fi;
    if   ForAll(l[2],c->IsList(c) and Length(c)=5)
    then return RcwaMappingNC(l[2]); # sparse rep., of Z
    else return RcwaMappingNC(R,l[1],l[2]); fi;
  end );

#############################################################################
##
#S  Creating new rcwa mappings from rcwa mappings of the same ring. /////////
##
#############################################################################

#############################################################################
##
#M  PiecewiseMapping( <sources>, <maps> ) . . . . . .  for rcwa mappings of Z
##
InstallMethod( PiecewiseMapping,
               "for rcwa mappings of Z (RCWA)",
               ReturnTrue, [ IsList, IsList ], 0,

  function ( sources, maps )

    local  map, coeffs, f, c, S, t, cls, cl, i;

    if  Length(sources) <> Length(maps)
      or not ForAll(sources,IsListOrCollection)
      or not ForAll(maps,IsMapping)
    then return fail; fi;
    if   not ForAll(sources,IsResidueClassUnionOfZ)
      or not ForAll(maps,IsRcwaMappingOfZ)
      or not IsIntegers(Union(sources))
      or Sum(List(sources,Density)) <> 1
    then TryNextMethod(); fi;

    maps := List(maps,SparseRep);
    coeffs := [];
    for i in [1..Length(maps)] do
      S := sources[i];
      f := maps[i];
      c := Coefficients(f);
      for t in c do
        cls := AsUnionOfFewClasses(Intersection(ResidueClass(t[1],t[2]),S));
        for cl in cls do
          Add(coeffs,[Residue(cl),Modulus(cl),t[3],t[4],t[5]]);
        od;
      od;
    od;
    map := RcwaMapping(coeffs);
    return SparseRep(map);
  end );

#############################################################################
##
#S  Creating new rcwa mappings from rcwa mappings of different rings. ///////
##
#############################################################################

#############################################################################
##
#F  LocalizedRcwaMapping( <f>, <p> )
##
InstallGlobalFunction( LocalizedRcwaMapping,

  function ( f, p )
    if   not IsRcwaMappingOfZ(f) or not IsInt(p) or not IsPrimeInt(p)
    then Error("usage: see ?LocalizedRcwaMapping( f, p )\n"); fi;
    return SemilocalizedRcwaMapping( f, [ p ] );
  end );

#############################################################################
##
#F  SemilocalizedRcwaMapping( <f>, <pi> )
##
InstallGlobalFunction( SemilocalizedRcwaMapping,

  function ( f, pi )
    if    IsRcwaMappingOfZ(f) and IsList(pi) and ForAll(pi,IsPosInt)
      and ForAll(pi,IsPrimeInt) and IsSubset(pi,Factors(Modulus(f)))
    then return RcwaMapping(Z_pi(pi),ShallowCopy(Coefficients(f)));
    else Error("usage: see ?SemilocalizedRcwaMapping( f, pi )\n"); fi;
  end );

#############################################################################
##
#M  ProjectionsToCoordinates( <f> ) . . . . . . . .  for rcwa mappings of Z^2
##
InstallMethod( ProjectionsToCoordinates,
               "rcwa mapping of Z^2 to two rcwa mappings of Z (RCWA)", true,
               [ IsRcwaMappingOfZxZ ], 0,

  function ( f )

    local  m, mf, mg, c, cf, cg, res, t, t1, t2, r, i;

    m := Modulus(f); c := Coefficients(f);
    res := AllResidues(Integers^2,m);
    mf := m[1][1]; mg := m[2][2]; cf := []; cg := [];
    for i in [1..Length(res)] do
      t := c[i]; r := res[i];
      t1 := [t[1][1][1],t[2][1],t[3]]; t1 := t1/Gcd(t1);
      t2 := [t[1][2][2],t[2][2],t[3]]; t2 := t2/Gcd(t2);
      if   not IsBound(cf[r[1]+1]) then cf[r[1]+1] := t1;
      elif cf[r[1]+1] <> t1        then return fail; fi;
      if   not IsBound(cg[r[2]+1]) then cg[r[2]+1] := t2;
      elif cg[r[2]+1] <> t2        then return fail; fi;
      if t[1][1][2] <> 0 or t[1][2][1] <> 0 then return fail; fi;
    od;
    return [ RcwaMapping(cf), RcwaMapping(cg) ];
  end );

#############################################################################
##
#M  Projection( <f>, <coord> ) .  proj. of an rcwa mapping of Z^2 to 1 coord.
##
InstallOtherMethod( Projection,
                    "for rcwa mappings of Z^2 (RCWA)",
                    ReturnTrue, [ IsRcwaMappingOfZxZ, IsPosInt ], 0,

  function ( f, coord )

    local  m, c, c_proj, proj;

    if not coord in [1,2] then return fail; fi; # there are only 2 coord's

    proj := ProjectionsToCoordinates(f); # maybe even both projections exist
    if proj <> fail then return proj[coord]; fi;

    m := Modulus(f); # check whether the choice of the affine partial mapping
                     # depends only on the specified coordinate:
    if   m[1][2] <> 0 or m[2][1] <> 0 or m[3-coord][3-coord] <> 1
    then return fail; fi;
 
    c := Coefficients(f); # check for dependency on other coordinate:
    if not ForAll(c,t->t[1][3-coord][coord]=0) then return fail; fi;

    # build coefficient list in one dimension:
    c_proj := List(c,t->[t[1][coord][coord],t[2][coord],t[3]]);

    return RcwaMapping(c_proj);
  end );

#############################################################################
##
#S  The automorphism switching action on negative and nonnegative integers. /
##
#############################################################################

#############################################################################
##
#M  Mirrored( <f> ) . . . . . . . . . for rcwa mappings of Z in standard rep.
#M  Mirrored( <f> ) . . . . . . . . . for rcwa mappings of Z in sparse rep.
##
InstallOtherMethod( Mirrored,
                    "for rcwa mappings of Z in standard rep. (RCWA)",
                    true, [ IsRcwaMappingOfZInStandardRep ], 0,
                    f -> f^RcwaMapping( [ [ -1, -1, 1 ] ] ) );
InstallOtherMethod( Mirrored,
                    "for rcwa mappings of Z in sparse rep. (RCWA)",
                    true, [ IsRcwaMappingOfZInSparseRep ], 0,
                    f -> f^RcwaMapping( [ [ 0, 1, -1, -1, 1 ] ] ) );

#############################################################################
##
#S  Constructors and basic methods for special types of rcwa permutations. //
##
#############################################################################

#############################################################################
##
#F  ClassShift( <R>, <r>, <m> ) . . . . . . . . . . . . . class shift nu_r(m)
#F  ClassShift( <r>, <m> )  . . . . . . . . . . . . . . . . . . . . .  (dito)
#F  ClassShift( <R>, <cl> ) . . . . . .  class shift nu_r(m), where cl = r(m)
#F  ClassShift( <cl> )  . . . . . . . . . . . . . . . . . . . . . . .  (dito)
#F  ClassShift( <R> ) . . . . . . . . . . . . .  class shift nu_R: n -> n + 1
##
##  (Enclosing the argument list in list brackets is permitted.)
##
InstallGlobalFunction( ClassShift,

  function ( arg )

    local  result, R, coeff, idcoeff, res, pos, r, m, latex;

    if Length(arg) = 1 and IsList(arg[1]) then arg := arg[1]; fi;

    if   IsZxZ(arg[1])
      or IsRowVector(arg[1]) and Length(arg[1]) = 2 and ForAll(arg[1],IsInt)
      or IsResidueClassOfZxZ(arg[1])
    then return CallFuncList(ClassShiftOfZxZ,arg); fi;

    if not Length(arg) in [1..3]
      or     Length(arg) = 1 and not IsResidueClass(arg[1])
      or     Length(arg) = 2
         and not (   ForAll(arg,IsRingElement)
                  or     IsRing(arg[1])
                     and IsResidueClass(arg[2])
                     and arg[1] = UnderlyingRing(FamilyObj(arg[2])))
      or     Length(arg) = 3
         and not (    IsRing(arg[1])
                  and IsSubset(arg[1],arg{[2,3]}))
    then Error("usage: see ?ClassShift( r, m )\n"); fi;

    if IsRing(arg[1]) then R := arg[1]; arg := arg{[2..Length(arg)]}; fi;
    if   IsBound(R) and IsEmpty(arg)
    then arg := [0,1] * One(R);
    elif IsResidueClass(arg[1])
    then if not IsBound(R) then R := UnderlyingRing(FamilyObj(arg[1])); fi;
         arg := [Residue(arg[1]),Modulus(arg[1])] * One(R);
    elif not IsBound(R) then R := DefaultRing(arg[2]); fi;
    arg := arg * One(R); # Now we know R, and we have arg = [r,m].

    m          := StandardAssociate(R,arg[2]);
    r          := arg[1] mod m;
    res        := AllResidues(R,m);
    idcoeff    := [1,0,1]*One(R);
    coeff      := List(res,r->idcoeff);
    pos        := PositionSorted(res,r);
    coeff[pos] := [1,m,1]*One(R);
    result     := RcwaMapping(R,m,coeff);
    SetIsClassShift(result,true); SetIsBijective(result,true);
    if   Characteristic(R) = 0
    then SetOrder(result,infinity);
    else SetOrder(result,Characteristic(R)); fi;
    SetIsTame(result,true);
    SetBaseRoot(result,result); SetPowerOverBaseRoot(result,1);
    if IsIntegers(R) then
      SetSmallestRoot(result,result); SetPowerOverSmallestRoot(result,1);
    fi;
    SetFactorizationIntoCSCRCT(result,[result]);

    latex := ValueOption("LaTeXString");
    if latex = fail then
      if IsOne(m) then
        SetLaTeXString(result,"\\nu");
      else
        SetLaTeXString(result,Concatenation("\\nu_{",String(r),"(",
                                                     String(m),")}"));
      fi;
    elif not IsEmpty(latex) then SetLaTeXString(result,latex); fi;

    return result;
  end );

#############################################################################
##
#F  ClassShiftOfZxZ( <R>, <r>, <m>, <coord> )  class shift nu_r(m),c; c=coord
#F  ClassShiftOfZxZ( <r>, <m>, <coord> )  . . . . . . . . . . . . . .  (dito)
#F  ClassShiftOfZxZ( <R>, <cl>, <coord> ) . .  class shift nu_r(m),c; cl=r(m)
#F  ClassShiftOfZxZ( <cl>, <coord> )  . . . . . . . . . . . . . . . .  (dito)
#F  ClassShiftOfZxZ( <R>, <coord> ) . . . . . . . . . . .  class shift nu_R_c
##
##  This function is called by `ClassShift' if the first argument is either
##  Integers^2, a row vector of length 2 with integer entries or a residue
##  class of Integers^2. Enclosing the argument list in list brackets is
##  permitted.
##
InstallGlobalFunction( ClassShiftOfZxZ,

  function ( arg )

    local  result, R, M, r, m, coord, cl, coeff, idcoeff, res, pos, latex;

    R := Integers^2; M := FullMatrixAlgebra(Integers,2);

    if Length(arg) = 1 and IsList(arg[1]) then arg := arg[1]; fi;

    coord := Remove(arg);
    if IsZxZ(arg[1]) then arg := arg{[2..Length(arg)]}; fi;

    if   not coord in [1,2]
      or not Length(arg) in [0..2]
      or Length(arg) = 1 and not IsResidueClassOfZxZ(arg[1])
      or Length(arg) = 2 and not (arg[1] in R and arg[2] in M)
    then Error("usage: see ?ClassShift( r, m, coord )\n"); fi;

    if   arg = [] then cl := R;
    elif Length(arg) = 1 then cl := arg[1];
    else cl := ResidueClass(arg[1],arg[2]); fi;

    r := Residue(cl); m := Modulus(cl);

    res        := AllResidues(R,m);
    idcoeff    := [[[1,0],[0,1]],[0,0],1];
    coeff      := ListWithIdenticalEntries(Length(res),idcoeff);
    pos        := PositionSorted(res,r);
    coeff[pos] := [[[1,0],[0,1]],m[coord],1];
    result     := RcwaMapping(R,m,coeff);
    SetIsClassShift(result,true); SetIsBijective(result,true);
    SetOrder(result,infinity);
    SetIsTame(result,true);
    SetBaseRoot(result,result); SetPowerOverBaseRoot(result,1);
    SetFactorizationIntoCSCRCT(result,[result]);

    latex := ValueOption("LaTeXString");
    if latex = fail then
      latex := Concatenation("\\nu_{",ViewString(cl),",",String(coord),"}");
      latex := ReplacedString(latex,"Z","\\mathbb{Z}");
      SetLaTeXString(result,latex);
    elif not IsEmpty(latex) then SetLaTeXString(result,latex); fi;

    return result;
  end );

#############################################################################
##
#M  IsClassShift( <sigma> ) . . . . . . . . . . . . . . . . for rcwa mappings
##
InstallMethod( IsClassShift,
               "for rcwa mappings (RCWA)", true, [ IsRcwaMapping ], 0,
               sigma -> IsResidueClass(Support(sigma))
                        and sigma = ClassShift(Support(sigma)) );

#############################################################################
##
#M  IsPowerOfClassShift( <sigma> ) . . . . . . . . . . for rcwa mappings of Z
##
InstallMethod( IsPowerOfClassShift, "for rcwa mappings of Z (RCWA)", true,
               [ IsRcwaMappingOfZ ], 0,

  function ( sigma )

    local  cl, c;

    if not IsClassWiseOrderPreserving(sigma) then return false; fi;
    if IsSignPreserving(sigma) then return IsOne(sigma); fi;
    if not IsBijective(sigma) then return false; fi;
    cl := Support(sigma);
    if not IsResidueClass(cl) then return false; fi;
    if Mod(sigma) <> Mod(cl) then return false; fi;
    return true;
  end );

#############################################################################
##
#M  String( <cs> ) . . . . . . . . . . . . . . . . . . . . . for class shifts
#M  ViewString( <cs> ) . . . . . . . . . . . . . . . . . . . for class shifts
#M  PrintObj( <cs> ) . . . . . . . . . . . . . . . . . . . . for class shifts
#M  ViewObj( <cs> )  . . . . . . . . . . . . . . . . . . . . for class shifts
##
InstallMethod( String, "for class shifts (RCWA)", true,
               [ IsRcwaMapping and IsClassShift ], SUM_FLAGS,

  function ( cs )

    local  str;

    if IsRcwaMappingOfZ(cs) then
      str := Concatenation(List(["ClassShift(",Residue(Support(cs)),",",
                                 Modulus(Support(cs))],String));
    else
      str := Concatenation(List(["ClassShift(",Source(cs),",",
                                 Residue(Support(cs)),",",
--> --------------------

--> maximum size reached

--> --------------------