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( , <precision> )
##
## PurePadicNumberFamily returns the family of pure p-adic numbers over the
## prime 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( " 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 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( , <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 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( " 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 );

Quelle padics.gi Sprache: unbekannt

[ Dauer der Verarbeitung: 0.40 Sekunden (vorverarbeitet) ]