|
#############################################################################
##
## 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 );
[ zur Elbe Produktseite wechseln0.96Quellennavigators
Analyse erneut starten
]
|
