p-Adic ZZ_pX FM Element.

This file implements elements of Eisenstein and unramified extensions of \mathbb{Z}_p with fixed modulus precision.

For the parent class see padic_extension_leaves.pyx.

The underlying implementation is through NTL’s ZZ_pX class. Each element contains the following data:

  • value (ZZ_pX_c) – An ntl ZZ_pX storing the value. The variable x is the uniformizer in the case of Eisenstein extensions. This ZZ_pX is created with global ntl modulus determined by the parent’s precision cap and shared among all elements.
  • prime_pow (some subclass of PowComputer_ZZ_pX) – a class, identical among all elements with the same parent, holding common data.
    • prime_pow.deg – The degree of the extension
    • prime_pow.e – The ramification index
    • prime_pow.f – The inertia degree
    • prime_pow.prec_cap – the unramified precision cap. For Eisenstein extensions this is the smallest power of p that is zero.
    • prime_pow.ram_prec_cap – the ramified precision cap. For Eisenstein extensions this will be the smallest power of x that is indistinguishable from zero.
    • prime_pow.pow_ZZ_tmp, prime_pow.pow_mpz_t_tmp``, prime_pow.pow_Integer – functions for accessing powers of p. The first two return pointers. See sage/rings/padics/pow_computer_ext for examples and important warnings.
    • prime_pow.get_context, prime_pow.get_context_capdiv, prime_pow.get_top_context – obtain an ntl_ZZ_pContext_class corresponding to p^n. The capdiv version divides by prime_pow.e as appropriate. top_context corresponds to p^{prec_cap}.
    • prime_pow.restore_context, prime_pow.restore_context_capdiv, prime_pow.restore_top_context – restores the given context.
    • prime_pow.get_modulus, get_modulus_capdiv, get_top_modulus – Returns a ZZ_pX_Modulus_c* pointing to a polynomial modulus defined modulo p^n (appropriately divided by prime_pow.e in the capdiv case).

EXAMPLES:

An Eisenstein extension:

sage: R = ZpFM(5,5)
sage: S.<x> = R[]
sage: f = x^5 + 75*x^3 - 15*x^2 +125*x - 5
sage: W.<w> = R.ext(f); W
Eisenstein Extension of 5-adic Ring of fixed modulus 5^5 in w defined by (1 + O(5^5))*x^5 + (3*5^2 + O(5^5))*x^3 + (2*5 + 4*5^2 + 4*5^3 + 4*5^4 + O(5^5))*x^2 + (5^3 + O(5^5))*x + (4*5 + 4*5^2 + 4*5^3 + 4*5^4 + O(5^5))
sage: z = (1+w)^5; z
1 + w^5 + w^6 + 2*w^7 + 4*w^8 + 3*w^10 + w^12 + 4*w^13 + 4*w^14 + 4*w^15 + 4*w^16 + 4*w^17 + 4*w^20 + w^21 + 4*w^24 + O(w^25)
sage: y = z >> 1; y
w^4 + w^5 + 2*w^6 + 4*w^7 + 3*w^9 + w^11 + 4*w^12 + 4*w^13 + 4*w^14 + 4*w^15 + 4*w^16 + 4*w^19 + w^20 + 4*w^23 + 4*w^24 + O(w^25)
sage: y.valuation()
4
sage: y.precision_relative()
21
sage: y.precision_absolute()
25
sage: z - (y << 1)
1 + O(w^25)

An unramified extension:

sage: g = x^3 + 3*x + 3
sage: A.<a> = R.ext(g)
sage: z = (1+a)^5; z
(2*a^2 + 4*a) + (3*a^2 + 3*a + 1)*5 + (4*a^2 + 3*a + 4)*5^2 + (4*a^2 + 4*a + 4)*5^3 + (4*a^2 + 4*a + 4)*5^4 + O(5^5)
sage: z - 1 - 5*a - 10*a^2 - 10*a^3 - 5*a^4 - a^5
O(5^5)
sage: y = z >> 1; y
(3*a^2 + 3*a + 1) + (4*a^2 + 3*a + 4)*5 + (4*a^2 + 4*a + 4)*5^2 + (4*a^2 + 4*a + 4)*5^3 + O(5^5)
sage: 1/a
(3*a^2 + 4) + (a^2 + 4)*5 + (3*a^2 + 4)*5^2 + (a^2 + 4)*5^3 + (3*a^2 + 4)*5^4 + O(5^5)

Different printing modes:

sage: R = ZpFM(5, print_mode='digits'); S.<x> = R[]; f = x^5 + 75*x^3 - 15*x^2 + 125*x -5; W.<w> = R.ext(f)
sage: z = (1+w)^5; repr(z)
'...4110403113210310442221311242000111011201102002023303214332011214403232013144001400444441030421100001'
sage: R = ZpFM(5, print_mode='bars'); S.<x> = R[]; g = x^3 + 3*x + 3; A.<a> = R.ext(g)
sage: z = (1+a)^5; repr(z)
'...[4, 4, 4]|[4, 4, 4]|[4, 4, 4]|[4, 4, 4]|[4, 4, 4]|[4, 4, 4]|[4, 4, 4]|[4, 4, 4]|[4, 4, 4]|[4, 4, 4]|[4, 4, 4]|[4, 4, 4]|[4, 4, 4]|[4, 4, 4]|[4, 4, 4]|[4, 4, 4]|[4, 4, 4]|[4, 3, 4]|[1, 3, 3]|[0, 4, 2]'
sage: R = ZpFM(5, print_mode='terse'); S.<x> = R[]; f = x^5 + 75*x^3 - 15*x^2 + 125*x -5; W.<w> = R.ext(f)
sage: z = (1+w)^5; z
6 + 95367431640505*w + 25*w^2 + 95367431640560*w^3 + 5*w^4 + O(w^100)
sage: R = ZpFM(5, print_mode='val-unit'); S.<x> = R[]; f = x^5 + 75*x^3 - 15*x^2 + 125*x -5; W.<w> = R.ext(f)
sage: y = (1+w)^5 - 1; y
w^5 * (2090041 + 95367431439401*w + 76293946571402*w^2 + 57220458985049*w^3 + 57220459001160*w^4) + O(w^100)

AUTHORS:

  • David Roe (2008-01-01) initial version
sage.rings.padics.padic_ZZ_pX_FM_element.make_ZZpXFMElement(parent, f)

Creates a new pAdicZZpXFMElement out of an ntl_ZZ_pX f, with parent parent. For use with pickling.

EXAMPLES:

sage: R = ZpFM(5,5)
sage: S.<x> = R[]
sage: f = x^5 + 75*x^3 - 15*x^2 +125*x - 5
sage: W.<w> = R.ext(f)
sage: z = (1 + w)^5 - 1
sage: loads(dumps(z)) == z # indirect doctest
True
class sage.rings.padics.padic_ZZ_pX_FM_element.pAdicZZpXFMElement

Bases: sage.rings.padics.padic_ZZ_pX_element.pAdicZZpXElement

is_equal_to(right, absprec=None)

Returns whether self is equal to right modulo self.uniformizer()^absprec.

If absprec is None, returns if self is equal to right modulo the precision cap.

EXAMPLES:

sage: R = Zp(5,5)
sage: S.<x> = R[]
sage: f = x^5 + 75*x^3 - 15*x^2 +125*x - 5
sage: W.<w> = R.ext(f)
sage: a = W(47); b = W(47 + 25)
sage: a.is_equal_to(b)
False
sage: a.is_equal_to(b, 7)
True
is_zero(absprec=None)

Returns whether the valuation of self is at least absprec. If absprec is None, returns whether self is indistinguishable from zero.

EXAMPLES:

sage: R = ZpFM(5,5)
sage: S.<x> = R[]
sage: f = x^5 + 75*x^3 - 15*x^2 +125*x - 5
sage: W.<w> = R.ext(f)
sage: O(w^189).is_zero()
True
sage: W(0).is_zero()
True
sage: a = W(675)
sage: a.is_zero()
False
sage: a.is_zero(7)
True
sage: a.is_zero(21)
False
lift_to_precision(absprec)

Returns self.

EXAMPLES:

sage: R = ZpFM(5,5)
sage: S.<x> = R[]
sage: f = x^5 + 75*x^3 - 15*x^2 +125*x - 5
sage: W.<w> = R.ext(f)
sage: w.lift_to_precision(10000)
w + O(w^25)
list(lift_mode='simple')

Returns a list giving a series representation of self.

  • If lift_mode == 'simple' or 'smallest', the returned list will consist of
    • integers (in the eisenstein case) or
    • lists of integers (in the unramified case).
  • self can be reconstructed as
    • a sum of elements of the list times powers of the uniformiser (in the eisenstein case), or
    • as a sum of powers of the p times polynomials in the generator (in the unramified case).
  • If lift_mode == 'simple', all integers will be in the range [0,p-1],
  • If lift_mode == 'smallest' they will be in the range [(1-p)/2,
p/2].
  • If lift_mode == 'teichmuller', returns a list of pAdicZZpXCRElements, all of which are Teichmuller representatives and such that self is the sum of that list times powers of the uniformizer.

EXAMPLES:

sage: R = ZpFM(5,5)
sage: S.<x> = R[]
sage: f = x^5 + 75*x^3 - 15*x^2 +125*x - 5
sage: W.<w> = R.ext(f)
sage: y = W(775); y
w^10 + 4*w^12 + 2*w^14 + w^15 + 2*w^16 + 4*w^17 + w^18 + w^20 + 2*w^21 + 3*w^22 + w^23 + w^24 + O(w^25)
sage: (y>>9).list()
[0, 1, 0, 4, 0, 2, 1, 2, 4, 1, 0, 1, 2, 3, 1, 1, 4, 1, 2, 4, 1, 0, 4, 3]
sage: (y>>9).list('smallest')
[0, 1, 0, -1, 0, 2, 1, 2, 0, 1, 2, 1, 1, -1, -1, 2, -2, 0, -2, -2, -2, 0, 2, -2, 2]
sage: w^10 - w^12 + 2*w^14 + w^15 + 2*w^16 + w^18 + 2*w^19 + w^20 + w^21 - w^22 - w^23 + 2*w^24
w^10 + 4*w^12 + 2*w^14 + w^15 + 2*w^16 + 4*w^17 + w^18 + w^20 + 2*w^21 + 3*w^22 + w^23 + w^24 + O(w^25)
sage: g = x^3 + 3*x + 3
sage: A.<a> = R.ext(g)
sage: y = 75 + 45*a + 1200*a^2; y
4*a*5 + (3*a^2 + a + 3)*5^2 + 4*a^2*5^3 + a^2*5^4 + O(5^5)
sage: y.list()
[[], [0, 4], [3, 1, 3], [0, 0, 4], [0, 0, 1]]
sage: y.list('smallest')
[[], [0, -1], [-2, 2, -2], [1], [0, 0, 2]]
sage: 5*((-2*5 + 25) + (-1 + 2*5)*a + (-2*5 + 2*125)*a^2)
4*a*5 + (3*a^2 + a + 3)*5^2 + 4*a^2*5^3 + a^2*5^4 + O(5^5)
sage: W(0).list()
[0]
sage: A(0,4).list()
[[]]
log(branch=None, same_ring=True)

Compute the p-adic logarithm of any unit.

See below for normalization.

INPUT:

  • branchpAdicZZpXFMElement (default None). The log of the uniformizer. If None, then an error is raised if self is not a unit.
  • same_ringbool or pAdicGeneric (default True). When e > p it is possible (even common) that the image of the log map is not contained in the ring of integers. If same_ring is True, then this function will return a value in self.parent() or raise an error if the answer would have negative valuation. If same_ring is False, then this function will raise an error (this behavior is for consistency with other p-adic types). If same_ring is a p-adic field into which this fixed mod ring can be successfully cast, then self is cast into that field and the log is taken there. Note that this casting will assume that self has the full precision possible.

OUTPUT:

The p-adic log of self.

Let K be the parent of self, \pi be a uniformizer of K
and w be a generator for the group of roots of unity in K. The usual power series for log with values in the additive group of K only converges for 1-units (units congruent to 1 modulo \pi). However, there is a unique extension of log to a homomorphism defined on all the units. If u = a \cdot v is a unit with v = 1
\pmod{p}, then we define \log(u) = \log(v). This is the correct extension because the units U of K split as a product U = V \times \langle w \rangle, where V is the subgroup of 1-units. The \langle w \rangle factor is torsion, so must go to 0 under any homomorphism to the torsion free group (K, +).

NOTES:

What some other systems do with regard to non-1-units:

  • PARI: Seems to define log the same way as we do.
  • MAGMA: Gives an error when unit is not a 1-unit.

In addition, if branch is specified, then the log map will work on non-units

..math

\log(\pi^k \cdot u) = k \cdot branch + \log(u)

ALGORITHM:

Input: Some unit u.

  1. Check that the input is really a unit (i.e., valuation 0), or that branch is specified.
  2. Let 1-x = u^{q-1}, which is a 1-unit, where q is the order of the residue field of K.
  3. Use the series expansion

..math

\log(1-x) = F(x) = -x - 1/2*x^2 - 1/3*x^3 - 1/4*x^4 - 1/5*x^5 - ...

to compute the logarithm log(u^{q-1}).

Add on terms until x^k is zero modulo the precision cap, and then determine if there are further terms that contribute to the sum (those where k is slightly above the precision cap but divisible by p).

  1. Then

..math

\log(u) = \log(u^{q-1})/(q-1) = F(1-u^{q-1})/(q-1).

EXAMPLES:

First, the Eisenstein case.:

sage: R = ZpFM(5,5)
sage: S.<x> = R[]
sage: f = x^4 + 15*x^2 + 625*x - 5
sage: W.<w> = R.ext(f)
sage: z = 1 + w^2 + 4*w^7; z
1 + w^2 + 4*w^7 + O(w^20)
sage: z.log()
4*w^2 + 3*w^4 + w^6 + w^7 + w^8 + 4*w^9 + 3*w^10 + w^12 + w^13 + 3*w^14 + w^15 + 4*w^16 + 4*w^17 + 3*w^18 + 3*w^19 + O(w^20)

Check that log is multiplicative:

sage: y = 1 + 3*w^4 + w^5
sage: y.log() + z.log() - (y*z).log()
O(w^20)

Now an unramified example.:

sage: g = x^3 + 3*x + 3
sage: A.<a> = R.ext(g)
sage: b = 1 + 5*(1 + a^2) + 5^3*(3 + 2*a)
sage: b.log()
(4*a^2 + 4)*5 + (a^2 + a + 2)*5^2 + (a^2 + 2*a + 4)*5^3 + (a^2 + 2*a + 2)*5^4 + O(5^5)

Check that log is multiplicative:

sage: c = 3 + 5^2*(2 + 4*a)
sage: b.log() + c.log() - (b*c).log()
O(5^5)

AUTHORS:

  • David Roe: initial version

TODO:

matrix_mod_pn()

Returns the matrix of right multiplication by the element on the power basis 1, x, x^2, \ldots, x^{d-1} for this extension field. Thus the emph{rows} of this matrix give the images of each of the x^i. The entries of the matrices are IntegerMod elements, defined modulo p^(self.absprec() / e).

Raises an error if self has negative valuation.

EXAMPLES:

sage: R = ZpFM(5,5)
sage: S.<x> = R[]
sage: f = x^5 + 75*x^3 - 15*x^2 +125*x - 5
sage: W.<w> = R.ext(f)
sage: a = (3+w)^7
sage: a.matrix_mod_pn()
[2757  333 1068  725 2510]
[  50 1507  483  318  725]
[ 500   50 3007 2358  318]
[1590 1375 1695 1032 2358]
[2415  590 2370 2970 1032]
norm(base=None)

Return the absolute or relative norm of this element.

NOTE! This is not the p-adic absolute value. This is a field theoretic norm down to a ground ring.

If you want the p-adic absolute value, use the abs() function instead.

If K is given then K must be a subfield of the parent L of self, in which case the norm is the relative norm from L to K. In all other cases, the norm is the absolute norm down to \mathbb{Q}_p or \mathbb{Z}_p.

EXAMPLES:

sage: R = ZpCR(5,5)
sage: S.<x> = R[]
sage: f = x^5 + 75*x^3 - 15*x^2 +125*x - 5
sage: W.<w> = R.ext(f)
sage: ((1+2*w)^5).norm()
1 + 5^2 + O(5^5)
sage: ((1+2*w)).norm()^5
1 + 5^2 + O(5^5)
precision_absolute()

Returns the absolute precision of self, ie the precision cap of self.parent().

EXAMPLES:

sage: R = ZpFM(5,5)
sage: S.<x> = R[]
sage: f = x^5 + 75*x^3 - 15*x^2 +125*x - 5
sage: W.<w> = R.ext(f)
sage: a = W(75); a
3*w^10 + 2*w^12 + w^14 + w^16 + w^17 + 3*w^18 + 3*w^19 + 2*w^21 + 3*w^22 + 3*w^23 + O(w^25)
sage: a.valuation()
10
sage: a.precision_absolute()
25
sage: a.precision_relative()
15
sage: a.unit_part()
3 + 2*w^2 + w^4 + w^6 + w^7 + 3*w^8 + 3*w^9 + 2*w^11 + 3*w^12 + 3*w^13 + w^15 + 4*w^16 + 2*w^17 + w^18 + w^22 + 3*w^24 + O(w^25)
precision_relative()

Returns the relative precision of self, ie the precision cap of self.parent() minus the valuation of self.

EXAMPLES:

sage: R = ZpFM(5,5)
sage: S.<x> = R[]
sage: f = x^5 + 75*x^3 - 15*x^2 +125*x - 5
sage: W.<w> = R.ext(f)
sage: a = W(75); a
3*w^10 + 2*w^12 + w^14 + w^16 + w^17 + 3*w^18 + 3*w^19 + 2*w^21 + 3*w^22 + 3*w^23 + O(w^25)
sage: a.valuation()
10
sage: a.precision_absolute()
25
sage: a.precision_relative()
15
sage: a.unit_part()
3 + 2*w^2 + w^4 + w^6 + w^7 + 3*w^8 + 3*w^9 + 2*w^11 + 3*w^12 + 3*w^13 + w^15 + 4*w^16 + 2*w^17 + w^18 + w^22 + 3*w^24 + O(w^25)
trace(base=None)

Return the absolute or relative trace of this element.

If K is given then K must be a subfield of the parent L of self, in which case the norm is the relative norm from L to K. In all other cases, the norm is the absolute norm down to \mathbb{Q}_p or \mathbb{Z}_p.

EXAMPLES:

sage: R = ZpCR(5,5)
sage: S.<x> = R[]
sage: f = x^5 + 75*x^3 - 15*x^2 +125*x - 5
sage: W.<w> = R.ext(f)
sage: a = (2+3*w)^7
sage: b = (6+w^3)^5
sage: a.trace()
3*5 + 2*5^2 + 3*5^3 + 2*5^4 + O(5^5)
sage: a.trace() + b.trace()
4*5 + 5^2 + 5^3 + 2*5^4 + O(5^5)
sage: (a+b).trace()
4*5 + 5^2 + 5^3 + 2*5^4 + O(5^5)
unit_part()

Returns the unit part of self, ie self / uniformizer^(self.valuation())

EXAMPLES:

sage: R = ZpFM(5,5)
sage: S.<x> = R[]
sage: f = x^5 + 75*x^3 - 15*x^2 +125*x - 5
sage: W.<w> = R.ext(f)
sage: a = W(75); a
3*w^10 + 2*w^12 + w^14 + w^16 + w^17 + 3*w^18 + 3*w^19 + 2*w^21 + 3*w^22 + 3*w^23 + O(w^25)
sage: a.valuation()
10
sage: a.precision_absolute()
25
sage: a.precision_relative()
15
sage: a.unit_part()
3 + 2*w^2 + w^4 + w^6 + w^7 + 3*w^8 + 3*w^9 + 2*w^11 + 3*w^12 + 3*w^13 + w^15 + 4*w^16 + 2*w^17 + w^18 + w^22 + 3*w^24 + O(w^25)

Previous topic

p-Adic ZZ_pX CA Element.

Next topic

PowComputer.

This Page