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


Quelle  padics.gi   Sprache: unbekannt

 
#############################################################################
##
##  This file is part of GAP, a system for computational discrete algebra.
##  This file's authors include Jens Hollmann.
##
##  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 implementation part of the padic numbers.
##


#############################################################################
##
#F  PrintPadicExpansion( <ppower>, <int>, <prime>, <precision> )
##
##  PrintPadicsExpansion prints   a pure p-adic  number x,  which is given as
##  <ppower> and <int> such  that  x =  <prime>^<ppower>*<int> as the  p-adic
##  expansion in the  pure p-adic numbers with   <precision> "digits".  <int>
##  may be  divisible by <prime>.  Each "digit"  ranges from  0 to <prime>-1.
##  If a "digit" has more  than one decimal digit  it  is embraced in  single
##  quotes.
##
##  For Example:
##      153.0(17) is  1*17^(-2) + 5*17^(-1) + 3*17^0 and
##    '15'3.0(17) is 15*17^(-1) + 3*17^0
##
BindGlobal( "PrintPadicExpansion", function( ppower, int, prime, precision )
    local   pos,  flag,  z,  k,  r;

    if int = 0  then
        Print( "0" );
    else

        # <int> might be divisible by <prime>
        while int mod prime = 0 do
            int := int / prime;
            ppower := ppower + 1;
        od;
        if ppower > 0 then

            # leading zeros
            for pos in [0..ppower-1] do
                if pos = 1 then
                    Print( "." );
                fi;
                Print( "0" );
            od;
        fi;
        pos := ppower;
        flag := false;

        # print the <int>
        z := int;
        k := 1;
        r := z mod prime;
        while (pos < 1) or ((k<=precision) and (z<>0)) do
            if pos = 1 then
                Print( "." );
            fi;
            if prime >= 10 then
                if flag then
                    Print( r, "'" );
                else
                    Print( "'" , r, "'" );
                fi;
                flag := true;
            else
                Print( r );
                flag := false;
            fi;
            z := (z - r) / prime;
            r := z mod prime;
            k := k + 1;
            pos := pos + 1;
        od;
    fi;
    Print( "(", prime, ")" );
end );


#############################################################################
##
#F  PadicExpansionByRat( <a>, <prime>, <precision>)
##
##  PadicExpansionByRat takes a rational <a> and returns a list [ppart, erg],
##  such  that  <a>  is <prime>^ppart*erg in   the pure  p-adic numbers  over
##  <prime> with <precision> "digits".
##
##  For Example:
##    PadicExpansionByRat(5, 3, 4)   -> [0,  5] = 12.00(3)
##    PadicExpansionByRat(1/2, 3, 4) -> [0, 41] = 2.111(3)
##
BindGlobal( "PadicExpansionByRat", function( a, prime, precision )
    local   c,  flag,  ppart,  z,  step,  erg,  ppot,  l,  digit;

    if a = 0 then
        c := [0,0];
    else

        # extract the p-part (num and den of rationals are coprime)
        flag := false;
        ppart := 0;
        if NumeratorRat(a) mod prime = 0 then
            z := NumeratorRat(a) / prime;
            step := 1;
            flag := true;
        fi;
        if DenominatorRat(a) mod prime = 0 then
            z := DenominatorRat(a) / prime;
            step := -1;
            flag := true;
        fi;
        if flag then

            # extract the <prime>-part
            ppart := step;
            while z mod prime = 0 do
                z := z / prime;
                ppart := ppart + step;
            od;
        fi;
        a := a / prime^ppart;
        erg := 0;
        ppot := 1;
        l := 1;
        while( l<= precision) and (a <> 0) do
            digit := a mod prime;
            erg := erg + digit * ppot;
            a := (a - digit) / prime;
            ppot := ppot * prime;
            l := l + 1;
        od;
        c := [ppart, erg];
    fi;
    return c;
end );


#############################################################################
##
#F  MultMatrixPadicNumbersByCoefficientsList( <list> )
##
##  MultMatrix...List  takes  the coeff.-list  <list>   of  a polynomial  and
##  returns  a (n  x   2*n-1)-matrix if n  is the   degree of  the polynomial
##  i.e. the length of <list> minus 1.  If you have an extension L of a field
##  K by a given polynomial f (i.e. L = K[X]/(f(X))) with coeff.-list <list>,
##  than the i-th row of the returned matrix gives the coeff.-presentation of
##  X^(i-1) in the basis {X^0, ..., X^n-1} of L.
##
##  For example:
##    Mult...List( [-1,5,-2,1] ) -> [1      ]
##      so the polynomial is        [   1   ]
##      X^3-2X^2+5X-1               [      1]
##                                  [1 -5  2]
##                                  [2 -9 -1]
##
##  So multiplying two elements of L (polynomials in X) is simply multiplying
##  the  polynomials and then multiplying the  resulting coeff.-list with the
##  matrix and get the  coeff.-presentation of the result  in the right basis
##  of L.
##
BindGlobal( "MultMatrixPadicNumbersByCoefficientsList", function ( list )
    local   n,  F,  zero,  one,  mat,  i,  j;

    n := Length(list) - 1;
    F := FamilyObj(list[1]);
    zero := Zero(F);
    one  := One(F);
    if n <= 1 then
        return [[one]];
    fi;
    # prepare a zero-matrix with ones on the main diagonal:
    mat := [  ];
    for i in [1..2*n-1] do
        mat[i] := [ ];
        for j in [1..n] do
            mat[i][j] := zero;
        od;
    od;
    for i in [1..n] do
        mat[i][i] := one;
    od;
    # the n+1-th row is simple:
    for j in [1..n] do
        mat[n+1][j] := - list[j];
    od;
    # the rest is a little more complex. Regard the fact, that
    # x^i is x*x^(i-1) and the coeff.presentation of x^(i-1) is
    # already known.
    for i in [n+2..2*n-1] do
        mat[i] := ShiftedCoeffs(mat[i-1],1) + mat[i-1][n] * mat[n+1];
    od;
    return mat;
end );


#############################################################################
##
#F  StructureConstantsPadicNumbers( <e>, <f> )
##
##  An  extended  p-adic field  is   given by  two polynomials.  One  for the
##  ramified  part g with  degree <e> and one for  the unramified part h with
##  degree <f>. The extended p-adic field is then the extension of Q_p i.e. L
##  = Q_p[x,y]/(g(y),h(x)).
##
##  L has a basis in x^i*y^j. This basis is ordered as follows:
##    {1, x, x^2, ..., y, x y, x^2 y, ..., y^2, x y^2, x^2 y^2, ...}
##
##  Let B_i be  the i-th basiselement.   StructureConstantsPadicNumbers takes
##  the two degrees <e> and <f> and returns a (n x n)-matrix (n = <e> * <f>).
##  In this matrix stands  at  position (i,j) the  index  of the row  of  the
##  multiplication-matrix that gives the coeff.-presentation of B_i * B_j The
##  mult.-matrix   of  an   extended   p-adic  field    is   given  as    the
##  Kronecker-product of the mult.-matrices  of the two  polynomials returned
##  by        MultMatrixPadicNumbersByCoefficientsList.        (see        in
##  PadicExtensionNumberFamily)
##
##  So I get the structure-constants m_ijk if
##    B_i*B_j = SUM(k=1,...,n) m_ijk B_k
##  simply by M[ B[i,j], k].
##
##  M is the mult.-matrix and B is the matrix returned by this function.
##
BindGlobal( "StructureConstantsPadicNumbers", function( e, f )
    local   mat,  i,  j,  a1,  a2,  b1,  b2;

    # there are <e>*<f> basis-elements and according to the above ordering
    # B_i is simply x^((i-1) mod f)*y^((i-1) div f)
    mat := [];
    for i in [1..e*f] do
        mat[i] := [];
        for j in [1..e*f] do
            a1 := (i-1) mod f;
            a2 := (j-1) mod f;
            b1 := QuoInt(i-1, f);
            b2 := QuoInt(j-1, f);
            mat[i][j] := (a1+a2) + (b1+b2)*(2*f-1) + 1;
        od;
    od;
    return mat;
end );


#############################################################################
##
#M  ShiftedPadicNumber( <padic>, <shift> )
##
##  ShiftedPadicNumber  takes a p-adic number <padic>  and an integer <shift>
##  and returns the  p-adic number   c, that is   <padic> *  p^<shift>.   The
##  <shift> is just added to the p-part.
##
InstallMethod( ShiftedPadicNumber,
    true,
    [ IsPadicNumber, IsInt ],
    0,

function( x, shift )
    local   fam,  c;

    fam := FamilyObj(x);
    c := Immutable( [ x![1]+shift, x![2] ] );
    return Objectify( fam!.defaultType, c );
end );


#############################################################################
##
#M  <padic> * <rat>
##
InstallMethod( \*,
     true,
     [ IsPadicNumber, IsRat ],
     0,

function( a, b )
    local   fam;

    fam := FamilyObj( a );
    return a * PadicNumber( fam, b );
end );


#############################################################################
##
#M  <rat> * <padic>
##
InstallMethod( \*,
     true,
     [ IsRat, IsPadicNumber ],
     0,

function( a, b )
    local   fam;

    fam := FamilyObj( b );
    return PadicNumber( fam, a ) * b;
end );


#############################################################################
##
#M  <padic-list> * <rat>
##
InstallMethod( \*,
     true,
     [ IsPadicNumberList, IsRat ],
     0,

function( a, b )
    b := PadicNumber( FamilyObj(a[1]), b );
    return List( a, x -> x * b );
end );


#############################################################################
##
#M  <rat> * <padic-list>
##
InstallMethod( \*,
     true,
     [ IsRat, IsPadicNumberList ],
     0,

function( a, b )
    a := PadicNumber( FamilyObj(b[1]), a );
    return List( b, x -> a * x );
end );


#############################################################################
##
#M  ZeroOp( <padic> )
##
InstallMethod( ZeroOp,
    "for a p-adic number",
    true,
    [ IsPadicNumber ],
    0,

function( padic )
    return Zero( FamilyObj( padic ) );
end );


#############################################################################
##
#M  OneOp( <padic> )
##
InstallMethod( OneOp,
    "for a p-adic number",
    true,
    [ IsPadicNumber ],
    0,

function( padic )
    return One( FamilyObj( padic ) );
end );


#############################################################################
##
#M  PurePadicNumberFamily( <p>, <precision> )
##
##  PurePadicNumberFamily returns the family of pure p-adic numbers over the
##  prime  <p> with  <precision>  "digits".  For the  representation  of pure
##  p-adic numbers see "PadicNumber" below.
##

BindGlobal("PADICS_FAMILIES", []);


InstallGlobalFunction( PurePadicNumberFamily, function( p, precision )
    local   str,  fam;

    if not IsPrimeInt( p ) then
        Error( "<p> must be a prime" );
    fi;
    if (not IsInt( precision )) or (precision < 0) then
        Error( "<precision> must be a positive integer" );
    fi;
    if not IsBound(PADICS_FAMILIES[p]) then
        PADICS_FAMILIES[p] := [];
    fi;
    if not IsBound(PADICS_FAMILIES[p][precision]) then
        str := "PurePadicNumberFamily(";
        Append( str, String(p) );
        Append( str, "," );
        Append( str, String(precision) );
        Append( str, ")" );
        fam := NewFamily( str, IsPurePadicNumber );
        fam!.prime:= p;
        fam!.precision:= precision;
        fam!.modulus:= p^precision;
        fam!.printPadicSeries:= true;
        fam!.defaultType := NewType( fam, IsPurePadicNumber and IsPositionalObjectRep );
        PADICS_FAMILIES[p][precision] := fam;
    fi;
    return PADICS_FAMILIES[p][precision];
end );


#############################################################################
##
#M  PadicNumber( <pure-padic-family>, <list> )
##
##  Make  a pure  p-adic number  out of  a list.    A pure  p-adic  number is
##  represented as a list of  length 2 such  that the  number is p^list[1]  *
##  list[2].  It is easily guaranteed that  list[2] is never divisible by the
##  prime p.  By that we have always maximum precision.
##
InstallMethod( PadicNumber, "for a pure p-adic family and a list",
    true,
    [ IsPurePadicNumberFamily,
      IsCyclotomicCollection ],
    0,

function( fam, list )
    if Length(list) <> 2 then
        Error( "<list> must have length 2" );
    elif not IsInt(list[1]) then
        Error( "<list>[1] must be an integer" );
    elif not IsInt(list[2]) or list[2] < 0 or list[2] >= fam!.modulus then
        Error( "<list>[2] must be an integer in {0..p^precision}" );
    fi;
    return Objectify( fam!.defaultType, list );
end );


#############################################################################
##
#M  PadicNumber( <pure-padic-family>, <rat> )
##
##  Make a pure p-adic number out of a rational.
##
InstallMethod( PadicNumber, "for a pure p-adic family and a rational",
    true,
    [ IsPurePadicNumberFamily,
      IsRat ],
    0,

function( fam, rat )
    local   c;

    c := Immutable( PadicExpansionByRat( rat, fam!.prime, fam!.precision ) );
    return Objectify( fam!.defaultType, c );
end );


#############################################################################
##
#M  PrintObj( <pure-padic> )
##
InstallMethod( PrintObj,
    true,
    [ IsPurePadicNumber ],
    0,

function ( x )
    local   fam;

    fam := FamilyObj(x);
    # printPadicSeries is just a boolean variable. Handy for checking
    # what REALLY happens if set to false.
    if fam!.printPadicSeries then
        PrintPadicExpansion( x![1], x![2], fam!.prime, fam!.precision );
    else
        Print( fam!.prime, "^", x![1], "*", x![2], "(" , fam!.prime , ")" );
    fi;
end );


#############################################################################
##
#M  Random( <pure-padic-family> )
##
##  This is just  something that actually returns  a pure p-adic number.  The
##  range of the  p-part is not totally covered  as it is  infinity.  But you
##  may get two pure p-adic  numbers that have no  "digit"  in common, so  by
##  adding them, one of the two vanishes.
##
InstallOtherMethodWithRandomSource( Random, "for a random source and a pure p-adic family",
    true,
    [ IsRandomSource, IsPurePadicNumberFamily ],
    0,

function ( rg, fam )
    local c;

    c := [];
    c[1] := Random( rg, -fam!.precision, fam!.precision );
    c[2] := Random( rg, 0, fam!.modulus-1 );
    while c[2] mod fam!.prime = 0  do
        c[1] := c[1] + 1;
        c[2] := c[2] / fam!.prime;
    od;
    return Objectify( fam!.defaultType, c );
end );


#############################################################################
##
#M  Zero( <pure-padic-family> )
##
InstallOtherMethod( Zero,
    true,
    [ IsPurePadicNumberFamily ],
    0,

function ( fam )
    return PadicNumber( fam, [0,0] );
end );


#############################################################################
##
#M  IsZero( <pure-padic> )
##
InstallMethod( IsZero,
    true,
    [ IsPurePadicNumber ],
    0,

function ( x )
    return x![2] = 0;
end );


########################################################################
##
#M  One( <pure-padic-family> )
##
InstallOtherMethod( One,
    true,
    [ IsPurePadicNumberFamily ],
    0,

function ( fam )
    return PadicNumber( fam, [0,1] );
end );


#############################################################################
##
#M  Valuation( <pure-padic>
##
##  The Valuation is the p-part of the p-adic number.
##
InstallMethod( Valuation,
    true,
    [ IsPurePadicNumber ],
    0,

function( x )
    if IsZero(x) then
        return infinity;
    fi;
    return x![1];
end );


#############################################################################
##
#M  AdditiveInverseOp( <pure-padic> )
##
InstallMethod( AdditiveInverseOp,
     true,
     [ IsPurePadicNumber ],
     0,

function( a )
    local   fam,  c;

    fam := FamilyObj( a );
    c := [ a![1], -a![2] mod fam!.modulus];
    return Objectify( fam!.defaultType, c );
end );


#############################################################################
##
#M  InverseOp( <pure-padic> )
##
InstallMethod( InverseOp,
     true,
     [ IsPurePadicNumber ],
     0,

function( a )
    local   fam,  c;

    if IsZero(a)  then
        Error("division by zero");
    fi;
    fam:= FamilyObj( a );
    c:= [ -a![1] , 1/a![2] mod fam!.modulus ];
    return Objectify( fam!.defaultType, c );
end );


#############################################################################
##
#M  <pure-padic> + <pure-padic>
##
InstallMethod( \+,
     IsIdenticalObj,
     [ IsPurePadicNumber, IsPurePadicNumber ],
     0,

function( a, b )
    local   fam,  c,  r;

    # if <a> or <b> is zero, return the other one
    if IsZero(a) then
        return b;
    fi;
    if IsZero(b) then
        return a;
    fi;
    fam:= FamilyObj( a );

    # different valuation: c[2] is NOT divisible by p!
    if a![1] < b![1] then
        c := [ a![1],
               (a![2]+fam!.prime^(b![1]-a![1])*b![2]) mod fam!.modulus ];

    # equal valuation: c[2] MAY BE divisible by p! So check that.
    elif a![1] = b![1] then
        c := [];
        c[1] := a![1];
        c[2] := (a![2] + b![2]) mod fam!.modulus;

        # c[2] might be divisible by p
        r := c[2] mod fam!.prime;
        while (r=0) and (c[2]>1) do
            c[1] := c[1] + 1;
            c[2] := c[2] / fam!.prime;
            r := c[2] mod fam!.prime;
        od;

    # different valuation: again c[2] is NOT divisible by p!
    else
        c:= [ b![1],
              (fam!.prime^(a![1]-b![1])*a![2]+b![2]) mod fam!.modulus ];
    fi;
    return Objectify( fam!.defaultType, c );
end );


#############################################################################
##
#M  <pure-padic> * <pure-padic>
##
InstallMethod( \*,
     IsIdenticalObj,
     [ IsPurePadicNumber, IsPurePadicNumber ],
     0,

function( a, b )
    local   fam,  c;

    fam:= FamilyObj( a );
    if IsZero(a) then
        c:= [0, 0];
    elif IsZero(b) then
        c:= [0, 0];
    else
        c:= [ a![1]+b![1] , a![2]*b![2] mod fam!.modulus ];
    fi;
    return Objectify( fam!.defaultType, c );
end );


#############################################################################
##
#M  <pure-padic> / <pure-padic>
##
InstallMethod( \/,
     IsIdenticalObj,
     [ IsPurePadicNumber, IsPurePadicNumber ],
     0,

function( a, b )
    local   fam,  c;

    if IsZero(b) then
        Error("division by zero");
    fi;
    fam:= FamilyObj( a );
    c:= [ a![1]-b![1] , a![2]/b![2] mod fam!.modulus ];
    return Objectify( fam!.defaultType, c );
end );


#############################################################################
##
#M  <pure-padic> = <pure-padic>
##
InstallMethod( \=,
     IsIdenticalObj,
     [ IsPurePadicNumber, IsPurePadicNumber ],
     0,

function( a, b )
    return (a![1] = b![1]) and (a![2] = b![2]);
end );


#############################################################################
##
#M  <pure-padic> < <pure-padic>
##
##  This is just something to keep GAP quiet
##
InstallMethod( \<,
     IsIdenticalObj,
     [ IsPurePadicNumber, IsPurePadicNumber ],
     0,

function( a, b )
    if a![1] = b![1]  then
        return a![2] < b![2];
    else
        return a![1] < b![1];
    fi;
end );


#############################################################################
##
#M  PadicExtensionNumberFamily( <p>, <precision>, <unram>, <ram> )
##
##  An   extended p-adic field  L  is given by two   polynomials h and g with
##  coeff.-lists   <unram> (for  the  unramified  part)  and <ram>  (for  the
##  ramified part). Then L  is Q_p[x,y]/(h(x),g(y)).  This function takes the
##  prime number  <p> and the two coeff.-lists  <unram> and <ram> for the two
##  polynomials.  It   is  not checked  BUT <unram>   should be  a cyclotomic
##  polynomial and <ram> should be  an Eisenstein-polynomial or [1,1].  Every
##  number  out   of  L is  represented   as   a coeff.-list   for  the basis
##  {1,x,x^2,...,y,xy,x^2y,...} of L.   <precision> is the number of "digits"
##  that all the coeff. have.
##
InstallGlobalFunction( PadicExtensionNumberFamily,
    function( p, precision, unram, ram )
    local   str,  fam,  yem1;

    if not IsPrimeInt( p ) then
        Error( "<p> must be a prime" );
    fi;
    if (not IsInt( precision )) or (precision < 0) then
        Error( "<precision> must be a positive integer" );
    fi;
    str := "PadicExtensionNumberFamily(";
    Append( str, String(p) );
    Append( str, "," );
    Append( str, String(precision) );
    Append( str, ",...)" );
    fam := NewFamily( str, IsPadicExtensionNumber );
    fam!.defaultType := NewType( fam, IsPadicExtensionNumber and IsPositionalObjectRep );
    fam!.prime       := p;
    fam!.precision   := precision;
    fam!.modulus     := p^precision;
    fam!.unramified  := unram;
    fam!.f           := Length(unram)-1;
    fam!.ramified    := ram;
    fam!.e           := Length(ram)-1;
    fam!.n           := fam!.e * fam!.f;
    fam!.M           := KroneckerProduct(
                         MultMatrixPadicNumbersByCoefficientsList(ram),
                         MultMatrixPadicNumbersByCoefficientsList(unram) );
    fam!.B           := StructureConstantsPadicNumbers(fam!.e, fam!.f);

    yem1 := List( [1..fam!.n], i->0 );
    yem1[ (fam!.e-1) * fam!.f + 1 ] := 1;
    fam!.yem1 := Objectify( fam!.defaultType, [0, yem1] );

    fam!.printPadicSeries := true;

    return fam;
end );


#############################################################################
##
##  General  comment:    In  PadicExtensionNumberFamily  you  give   the two
##  polynomials,  that define  the  extension of  Q_p.  You have  to care for
##  yourself, that these polynomials  are  really irreducible over Q_p!   Try
##  PadicExtensionNumberFamily(3, 4, [1,1,1], [1,1]) for example.  You think
##  this is ok? It is not, because  x^2+x+1 is NOT  irreducible over Q_p. The
##  result being,  that you get non-invertible extended  p-adic numbers.  So,
##  if that happens, check your polynomials!
##


#############################################################################
##
#M  PadicNumber( <extended-padic-family>, <list> )
##
##  Make an extended p-adic number out of  a list.  An extended p-adic number
##  is represented as a list L of length 2.
##
##  L[2] is  the list of coeff. for  the  Basis {1,..,x^(f-1)*y^(e-1)} of the
##  extended p-adic field.
##
##  L[1] is a common p-part of all the coeff.
##
##  It is  NOT  guaranteed that all  or  at least one   of the coeff.  is not
##  divisible by the prime p.
##
##  For example: in PadicExtensionNumberFamily(3, 5, [1,1,1], [1,1])
##    the number (1.2000, 0.1210)(3) may be
##      [ 0, [ 1.2000, 0.1210 ] ]   or
##      [-1, [ 12.000, 1.2100 ] ]  here the coeff. have to be multiplied
##                                 by p^(-1)
##
##    so there may be a number (1.2, 2.2)(3) and you may ask: "Where are my 5
##    digits? There  are only two! Where  is the complain  department!"  But
##    the number is  intern: [-3, [ 0.0012, 0.0022  ] ]  and  so has in  fact
##    maximum precision.
##
##  So watch it!
##
InstallMethod( PadicNumber, "for a p-adic extension family and a list",
    true,
    [ IsPadicExtensionNumberFamily,
      IsList ],
    0,

function( fam, list )
    local   range;

    range := [ 0 .. fam!.modulus-1 ];
    if not IsInt(list[1]) then
        Error( "<list>[1] must be an integer" );
    elif not IsList(list[2]) or Length(list[2]) <> fam!.n  then
        Error( "<list>[2] must be a list of length ", fam!.n );
    elif not ForAll( list[2], x -> x in range )  then
        Error( "<list>[2] must be a list of integers in ", range );
    fi;
    return Objectify( fam!.defaultType, list );
end );


#############################################################################
##
#M  PadicNumber( <extended-padic-family>, <rat> )
##
##  Make an extended p-adic number  out of a rational.   That means take  the
##  result of PadicExpansionByRat  and put it at  the  first position of  the
##  coeff.list.
##
InstallMethod( PadicNumber, "for a p-adic extension family and a rational",
    true,
    [ IsPadicExtensionNumberFamily,
      IsRat ],
    0,

function( fam, rat )
    local   c,  erg;

    c := PadicExpansionByRat( rat, fam!.prime, fam!.precision );
    erg    := List( [1..fam!.n], i->0 );
    erg[1] := c[2];
    c      := [ c[1], erg ];
    return Objectify( fam!.defaultType, c );
end );


#############################################################################
##
#M  PrintObj( <extended-padic> )
##
InstallMethod( PrintObj,
    true,
    [ IsPadicExtensionNumber ],
    0,

function( x )
    local   fam,  l;

    fam := FamilyObj(x);
    if fam!.printPadicSeries then
        Print( "padic(" );
        for l in [1..fam!.n-1] do
            PrintPadicExpansion( x![1], x![2][l], fam!.prime, fam!.precision );
            Print( "," );
        od;
        PrintPadicExpansion( x![1], x![2][fam!.n], fam!.prime, fam!.precision );
        Print( ")" );
    else
        Print( "padic(", fam!.prime, "^", x![1], "*", x![2], ")" );
    fi;
end );


#############################################################################
##
#M  Random( <extended-padic-family>
##
##  Again this is just something that returns an extended p-adic number.  The
##  p-part is not totally covered (just ranges from -precision to precision).
##
InstallOtherMethodWithRandomSource( Random, "for a random source and p-adic extension family",
    true,
    [ IsRandomSource, IsPadicExtensionNumberFamily ],
    0,

function ( rg, fam )
    local   c,  l;

    c := [];
    c[1] := Random( rg, -fam!.precision, fam!.precision );
    c[2] := [];
    for l  in [ 1 .. fam!.n ] do
        c[2][l] := Random( rg, 0, fam!.modulus-1 );
    od;
    while ForAll( c[2], x-> x mod fam!.prime = 0 ) do
        c[1] := c[1] + 1;
        c[2] := c[2] / fam!.prime;
    od;
    return Objectify( fam!.defaultType, c );
end );


#############################################################################
##
#M  Zero( <extended-padic-family> )
##
InstallOtherMethod( Zero,
    true,
    [ IsPadicExtensionNumberFamily ],
    0,

function( fam )
    return Objectify( fam!.defaultType, [0, List( [1..fam!.n], i->0 )] );
end );


########################################################################
##
#M  IsZero( <extended-padic> )
##
InstallMethod( IsZero,
    true,
    [ IsPadicExtensionNumber ],
    0,

function ( x )
    if PositionNonZero( x![2] ) > Length( x![2] ) then
        return true;
    fi;
    return false;
end );


#############################################################################
##
#M  One( <extended-padic-family> )
##
InstallOtherMethod( One,
    true,
    [ IsPadicExtensionNumberFamily ],
    0,

function( fam )
    local   c;

    c := List( [ 1 .. fam!.n ], i -> 0 );
    c[1] := 1;
    return Objectify( fam!.defaultType, [0, c] );
end );


#############################################################################
##
#M  Valuation( <extended-padic>
##
##  In an  extended p-adic field the  prime p has  valuation e (the degree of
##  the totally ramified part) and y has valuation 1. y^e  is p so this makes
##  sense.
##
##  The valuation of a sum is (in this case) the minimum of the valuations of
##  the summands. As an extended p-adic number is
##
##    SUM(i=0..(f-1)) SUM(j=0..(e-1)) a_ij x^i y^j
##
##  the valuation nu of that number is the minimum over i and j of
##
##    nu(a_ij x^i y^j)
##
##  The valuation of a product is the sum of the single valuations, so
##
##    nu(a_ij x^i y^j) = nu(a_ij) + nu(x^i) + nu(y^j)
##                     = nu(a_ij) + j
##
InstallMethod( Valuation,
    true,
    [ IsPadicExtensionNumber ],
    0,

function ( x )
    local   fam,  min,  j,  wert;

    fam := FamilyObj(x);
    min := Minimum( List([1..fam!.f], i -> Valuation(x![2][i],fam!.prime)) );
    if min <> infinity then
        min := fam!.e * min;
    fi;
    for j in [1..fam!.e-1] do
        wert := Minimum( List( [1..fam!.f], i ->
                    Valuation(x![2][i+j*fam!.f], fam!.prime) ) );
        if wert <> infinity then
            wert := fam!.e * wert + j;
        fi;
        if wert < min then
            min := wert;
        fi;
    od;
    if min <> infinity then
        min := min + x![1] * fam!.e;
    fi;
    return min;
end );


#############################################################################
##
#M  AdditiveInverseOp( <extended-padic> )
##
InstallMethod( AdditiveInverseOp,
     true,
     [ IsPadicExtensionNumber ],
     0,

function( a )
    local   fam,  c;

    fam := FamilyObj( a );
    c := [ a![1], List( a![2], x -> (-x) mod fam!.modulus ) ];
    return Objectify( fam!.defaultType, c );
end );


#############################################################################
##
#M  InverseOp( <extended-padic> )
##
InstallMethod( InverseOp,
    true,
    [ IsPadicExtensionNumber ],
    0,

function(x)
    local   fam,  val,  coeffppart,  coeffypart,  ppart,  z,  L,  k,  j,
            Lp,  E,  Beta,  ppot,  Beta_k,  c,  addppart;

    if IsZero(x)  then
        Error( "<x> must be non-zero" );
    fi;
    fam := FamilyObj(x);
    # need a copy of x later:
    z := Objectify( fam!.defaultType, [x![1], ShallowCopy(x![2])] );

    # if x = [ppart, [x_1,...,x_n]] then
    # Valuation(x) = ppart*fam!.e + coeffppart*fam!.e + coeffypart
    val := Valuation(x);
    if fam!.e > 1 then
        coeffypart := val mod fam!.e;
    else
        coeffypart := 0;
    fi;
    ppart := x![1];
    coeffppart := (val - ppart*fam!.e - coeffypart) / fam!.e;
    # so x = p^(ppart + coeffppart) * y^coeffypart * z
    # and z is divisible neither by y nor by p, so it has Valuation 0
    # and can be inverted.
    # We don't have y^(-1) but y^(e-1) which is p*y^(-1)
    # so z = x * y^(e-1)^coeffypart * p^(-ppart-coeffppart-coeffypart).
    # at least there is the coeffppart in z![2]:
    if coeffppart > 0 then
        z![1] := z![1] + coeffppart;
        z![2] := z![2] / fam!.prime^coeffppart;
    fi;
    addppart := 0;
    if coeffypart > 0 then
        z := z * fam!.yem1^coeffypart;
        # BUT by multiplying with y^(e-1) one may get additional p-parts
        # in the coeffs
        while ForAll( z![2], y -> y mod fam!.prime = 0 ) do
            addppart := addppart + 1;
            z![1] := z![1] + 1;
            z![2] := z![2] / fam!.prime;
        od;
    fi;
    # z![1] := z![1] - ppart - coeffppart - coeffypart - addppart;

    # NOW z![2] has an entry that is not divisible by p and L should be
    # invertible
    L := [];
    for k in [1..fam!.n] do
        L[k] := [];
        for j in [1..fam!.n] do
            L[k][j]:=Sum(List([1..fam!.n],
                             i -> z![2][i] * fam!.M[ fam!.B[i][j] ][ k ] ));
        od;
    od;

    Lp := InverseMatMod(L, fam!.prime);

    # E is the right side (1,0,...,0) and Beta is just (0,...,0)
    E := List([1..fam!.n], i->0);
    Beta := ShallowCopy(E);
    E[1] := 1;
    ppot := 1;

    # now solve L * Beta = E mod p:
    for k in [0..fam!.precision-1-coeffppart] do
        Beta_k := List( Lp * E, x -> x mod fam!.prime );
        Beta := Beta + Beta_k * ppot;
        E := (E - L * Beta_k) / fam!.prime;
        ppot := ppot * fam!.prime;
    od;
    c := [];
    c[1] := -z![1];
    c[2] := Beta;

    # as z^(-1) is now calculated, the formula above gives
    # x^(-1) = z^(-1) * y^(e-1)^m * p^(-ppart-coeffppart-coeffypart)
    c[1] := c[1] - coeffppart - addppart;
    c[2] := c[2] * fam!.prime^(coeffppart + addppart);
    c[2] := List( c[2], x -> x mod fam!.modulus );
    Objectify( fam!.defaultType, c );
    return ( c * fam!.yem1^coeffypart );
end );


#############################################################################
##
#M  <extended-padic> + <extended-padic>
##
InstallMethod( \+,
     IsIdenticalObj,
     [ IsPadicExtensionNumber,
       IsPadicExtensionNumber ],
     0,

function( x, y )
    local   fam,  ppot,  c;

    if IsZero(y) then
        return x;
    elif IsZero(x) then
        return y;
    fi;
    fam:= FamilyObj( x );
    ppot := Minimum( x![1], y![1] );
    c := [];
    c[1] := ppot;
    c[2] := fam!.prime^(x![1]-ppot)*x![2] + fam!.prime^(y![1]-ppot)*y![2];
    c[2] := List( c[2], x -> x mod fam!.modulus );
    return Objectify( fam!.defaultType, c );
end );


#############################################################################
##
#M  <extended-padic> - <extended-padic>
##
InstallMethod( \-,
     IsIdenticalObj,
     [ IsPadicExtensionNumber,
       IsPadicExtensionNumber ],
     0,

function( x, y )
    local   fam,  ppot,  c;

    if IsZero(y) then
        return x;
    elif IsZero(x) then
        return y;
    fi;
    fam:= FamilyObj( x );
    ppot := Minimum( x![1], y![1] );
    c := [];
    c[1] := ppot;
    c[2] := fam!.prime^(x![1]-ppot)*x![2] - fam!.prime^(y![1]-ppot)*y![2];
    c[2] := List( c[2], x -> x mod fam!.modulus );
    return Objectify( fam!.defaultType, c );
end );


#############################################################################
##
#M  <pure-padic> * <extended-padic>
##
InstallMethod( \*,
     true,
     [ IsPurePadicNumber,
       IsPadicExtensionNumber ],
     0,

function( a, x )
    local   Qpxy,  Qp,  c;

    Qpxy := FamilyObj(x);
    Qp   := FamilyObj(a);
    if Qpxy!.prime <> Qp!.prime then
        Error( "different primes" );
    fi;
    if Qpxy!.precision <> Qp!.precision  then
        Error( "different precision" );
    fi;
    c := [ a![1] + x![1] ];
    c[2] := List( a![2] * x![2], x -> x mod Qpxy!.modulus );
    return Objectify( Qpxy!.defaultType, c );
end );


#############################################################################
##
#M  <extended-padic> * <pure-padic>
##
InstallMethod( \*,
     true,
     [ IsPadicExtensionNumber,
       IsPurePadicNumber ],
     0,

function( x, a )
    local   Qpxy,  Qp,  c;

    Qpxy := FamilyObj(x);
    Qp   := FamilyObj(a);
    if Qpxy!.prime <> Qp!.prime then
        Error( "different primes" );
    fi;
    if Qpxy!.precision <> Qp!.precision  then
        Error( "different precision" );
    fi;
    c := [ a![1] + x![1] ];
    c[2] := List( x![2] * a![2], x -> x mod Qpxy!.modulus );
    return Objectify( Qpxy!.defaultType, c );
end );


#############################################################################
##
#M  <extended-padic> * <extended-padic>
##
InstallMethod( \*,
     IsIdenticalObj,
     [ IsPadicExtensionNumber,
       IsPadicExtensionNumber ],
     0,

function (a, b)
    local   fam,  vec,  addvec,  bj,  bi,  aj,  ai,  c;

    fam := FamilyObj( a );

    ## zwei Nullvektoren der Laenge (2f-1)(2e-1):
    vec := List( [1..(2*fam!.f-1)*(2*fam!.e-1)], i->0 );
    addvec := List( [1..(2*fam!.f-1)*(2*fam!.e-1)], i->0 );

    ## Koeff. von b eintragen
    for bj in [1..fam!.e] do
        for bi in [1..fam!.f] do
            vec[ (bj-1)*(2*fam!.f-1) + bi ] := b![2][ (bj-1)*fam!.f + bi ];
        od;
    od;

    ## Das eigentliche Multiplizieren:
    for aj in [1..fam!.e] do
        for ai in [1..fam!.f] do
            addvec := addvec + vec * a![2][ (aj-1)*fam!.f + ai ];
            vec := ShiftedCoeffs( vec, 1 );
        od;
        vec := ShiftedCoeffs( vec, fam!.f-1 );
    od;
    c := [ a![1] + b![1], List(addvec * fam!.M, x -> x mod fam!.modulus) ];
    return Objectify( fam!.defaultType, c );
end );


#############################################################################
##
#M  <extended-padic> = <extended-padic>
##
InstallMethod( \=,
     IsIdenticalObj,
     [ IsPadicExtensionNumber,
       IsPadicExtensionNumber ],
     0,

function( a, b )
    local   fam;

    if IsZero(a)  then
        return IsZero(b);
    elif IsZero(b)  then
        return false;
    fi;
    # A little work is needed, because p^1 * 10 is equal to one
    fam := FamilyObj(a);
    a := [ a![1], ShallowCopy(a![2]) ];
    b := [ b![1], ShallowCopy(b![2]) ];
    while ForAll( a[2], z -> z mod fam!.prime = 0 ) do
        a[2] := a[2] / fam!.prime;
        a[1] := a[1] + 1;
    od;
    while ForAll( b[2], z -> z mod fam!.prime = 0 ) do
        b[2] := b[2] / fam!.prime;
        b[1] := b[1] + 1;
    od;
    return (a[1] = b[1]) and (a[2] = b[2]);
end );


#############################################################################
##
#M  <extended-padic> < <extended-padic>
##
##  Again just something to have it.
##
InstallMethod( \<,
     IsIdenticalObj,
     [ IsPadicExtensionNumber,
       IsPadicExtensionNumber ],
     0,

function( a, b )
    local   fam;

    if IsZero(a)  then
        return not IsZero(b);
    elif IsZero(b)  then
        return false;
    fi;
    fam := FamilyObj(a);
    a := [ a![1], ShallowCopy(a![2]) ];
    b := [ b![1], ShallowCopy(b![2]) ];
    while ForAll( a[2], z -> z mod fam!.prime = 0 ) do
        a[2] := a[2] / fam!.prime;
        a[1] := a[1] + 1;
    od;
    while ForAll( b[2], z -> z mod fam!.prime = 0 ) do
        b[2] := b[2] / fam!.prime;
        b[1] := b[1] + 1;
    od;
    if a[1] = b[1]  then
        return a[2] < b[2];
    else
        return a[1] < b[1];
    fi;
end );

[ Dauer der Verarbeitung: 0.5 Sekunden  (vorverarbeitet)  ]

                                                                                                                                                                                                                                                                                                                                                                                                     


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