Hot-keys on this page
r m x p toggle line displays
j k next/prev highlighted chunk
0 (zero) top of page
1 (one) first highlighted chunk
""" Floating point template for complete discrete valuation rings.
In order to use this template you need to write a linkage file and gluing file. For an example see mpz_linkage.pxi (linkage file) and padic_floating_point_element.pyx (gluing file).
The linkage file implements a common API that is then used in the class FPElement defined here. See sage/libs/linkages/padics/API.pxi for the functions needed.
The gluing file does the following:
- ctypedef's celement to be the appropriate type (e.g. mpz_t) - includes the linkage file - includes this template - defines a concrete class inheriting from FPElement, and implements any desired extra methods
AUTHORS:
- David Roe (2016-03-21) -- initial version """
#***************************************************************************** # Copyright (C) 2007-2016 David Roe <roed.math@gmail.com> # William Stein <wstein@gmail.com> # # Distributed under the terms of the GNU General Public License (GPL) # as published by the Free Software Foundation; either version 2 of # the License, or (at your option) any later version. # # http://www.gnu.org/licenses/ #*****************************************************************************
from sage.ext.stdsage cimport PY_NEW include "padic_template_element.pxi" from cpython.int cimport *
from sage.structure.element cimport Element from sage.rings.padics.common_conversion cimport comb_prec, _process_args_and_kwds from sage.rings.integer_ring import ZZ from sage.rings.rational_field import QQ from sage.categories.sets_cat import Sets from sage.categories.sets_with_partial_maps import SetsWithPartialMaps from sage.categories.homset import Hom
cdef inline bint overunderflow(long* ordp, celement unit, PowComputer_ prime_pow): """ Check for over and underflow. If detected, sets ordp and unit appropriately, and returns True. If not, returns False. """ ordp[0] = maxordp csetzero(unit, prime_pow) ordp[0] = minusmaxordp csetone(unit, prime_pow) else: return True
cdef inline bint overunderflow_mpz(long* ordp, mpz_t ordp_mpz, celement unit, PowComputer_ prime_pow): """ Check for over and underflow with an mpz_t ordp. If detected, sets ordp and unit appropriately, and returns True. If not, returns False. """ if mpz_sgn(ordp_mpz) > 0: ordp[0] = maxordp csetzero(unit, prime_pow) else: ordp[0] = minusmaxordp csetone(unit, prime_pow) return True
cdef inline bint very_pos_val(long ordp):
cdef inline bint very_neg_val(long ordp):
cdef inline bint huge_val(long ordp):
cdef class FPElement(pAdicTemplateElement): cdef int _set(self, x, long val, long xprec, absprec, relprec) except -1: """ Sets the value of this element from given defining data.
This function is intended for use in conversion, and should not be called on an element created with :meth:`_new_c`.
INPUT:
- ``x`` -- data defining a `p`-adic element: int, long, Integer, Rational, other `p`-adic element...
- ``val`` -- the valuation of the resulting element
- ``xprec -- an inherent precision of ``x``, if ``val`` is larger then the result will be zero.
- ``absprec`` -- an absolute precision cap for this element, if ``val`` is larger then the result will be zero.
- ``relprec`` -- a relative precision cap for this element (unused; for compatibility with other `p`-adic precision modes)
TESTS::
sage: R = ZpFP(5) sage: a = R(17,5); a #indirect doctest 2 + 3*5 sage: R(15) #indirect doctest 3*5
sage: R = ZpFP(5,5) sage: a = R(25/9); a #indirect doctest 4*5^2 + 2*5^3 + 5^5 + 2*5^6 sage: R(ZpCR(5)(25/9)) == a True sage: R(5) - R(5) 0
We check that :trac:`23966` is resolved::
sage: R = ZpFM(2) sage: K = R.fraction_field() sage: K(R.zero()) 0 """ self._set_infinity() else: ccopy(self.unit, (<FPElement>x).unit, self.prime_pow) else:
cdef int _set_exact_zero(self) except -1: """ Sets this element to zero.
TESTS::
sage: R = Zp(5); R(0) #indirect doctest 0 """
cdef int _set_infinity(self) except -1: """ Sets this element to zero.
TESTS::
sage: R = Zp(5); R(0) #indirect doctest 0 """
cdef FPElement _new_c(self): """ Creates a new element with the same basic info.
TESTS::
sage: R = ZpFP(5); R(6) * R(7) #indirect doctest 2 + 3*5 + 5^2 """
cdef pAdicTemplateElement _new_with_value(self, celement value, long absprec): """ Creates a new element with a given value and absolute precision.
Used by code that doesn't know the precision type. """
cdef int _get_unit(self, celement value) except -1: """ Sets ``value`` to the unit of this p-adic element. """
cdef int check_preccap(self) except -1: """ Check that the precision of this element does not exceed the precision cap. Does nothing for floating point elements.
TESTS::
sage: ZpFP(5)(1).lift_to_precision(30) # indirect doctest 1 """ pass
def __copy__(self): """ Return a copy of this element.
EXAMPLES::
sage: a = ZpFP(5,6)(17); b = copy(a) sage: a == b True sage: a is b False """
cdef int _normalize(self) except -1: """ Normalizes this element, so that ``self.ordp`` is correct.
TESTS::
sage: R = ZpFP(5) sage: R(6) + R(4) #indirect doctest 2*5 """ cdef long diff cdef bint is_zero self._set_exact_zero() self._set_infinity() else: else: self._set_exact_zero()
def __dealloc__(self): """ Deallocate the underlying data structure.
TESTS::
sage: R = ZpFP(5) sage: a = R(17) sage: del(a) """
def __reduce__(self): """ Return a tuple of a function and data that can be used to unpickle this element.
EXAMPLES::
sage: a = ZpFP(5)(-3) sage: type(a) <type 'sage.rings.padics.padic_floating_point_element.pAdicFloatingPointElement'> sage: loads(dumps(a)) == a True """
# def __richcmp__(self, right, int op): # """ # Compare this element to ``right`` using the comparison operator ``op``. # # TESTS:: # # sage: R = ZpFP(5) # sage: a = R(17) # sage: b = R(21) # sage: a == b # False # sage: a < b # True # """ # return (<Element>self)._richcmp(right, op)
cpdef _neg_(self): r""" Return the additive inverse of this element.
EXAMPLES::
sage: R = Zp(7, 4, 'floating-point', 'series') sage: -R(7) #indirect doctest 6*7 + 6*7^2 + 6*7^3 + 6*7^4 """ else:
cpdef _add_(self, _right): r""" Return the sum of this element and ``_right``.
EXAMPLES::
sage: R = Zp(7, 4, 'floating-point', 'series') sage: x = R(1721); x 6 + 5*7^3 sage: y = R(1373); y 1 + 4*7^3 sage: x + y #indirect doctest 7 + 2*7^3 """ cdef FPElement ans cdef long tmpL else: else: # Addition is commutative, swap so self.ordp < right.ordp ccopy(ans.unit, self.unit, ans.prime_pow) else:
cpdef _sub_(self, _right): r""" Return the difference of this element and ``_right``.
EXAMPLES::
sage: R = Zp(7, 4, 'floating-point', 'series') sage: x = R(1721); x 6 + 5*7^3 sage: y = R(1373); y 1 + 4*7^3 sage: x - y #indirect doctest 5 + 7^3 """ cdef FPElement ans cdef long tmpL else: ccopy(ans.unit, self.unit, ans.prime_pow) else: else: ccopy(ans.unit, self.unit, ans.prime_pow) else:
def __invert__(self): r""" Returns multiplicative inverse of this element.
EXAMPLES::
sage: R = Zp(7, 4, 'floating-point', 'series') sage: ~R(2) 4 + 3*7 + 3*7^2 + 3*7^3 sage: ~R(0) infinity sage: ~R(7) 7^-1 """ # Input should be normalized! csetone(ans.unit, ans.prime_pow) else:
cpdef _mul_(self, _right): r""" Return the product of this element and ``_right``.
EXAMPLES::
sage: R = Zp(7, 4, 'floating-point', 'series') sage: R(3) * R(2) #indirect doctest 6 sage: R(1/2) * R(2) 1 """ raise ZeroDivisionError("Cannot multiply 0 by infinity") return right return ans
cpdef _div_(self, _right): r""" Return the quotient of this element and ``right``.
EXAMPLES::
sage: R = Zp(7, 4, 'floating-point', 'series') sage: R(3) / R(2) #indirect doctest 5 + 3*7 + 3*7^2 + 3*7^3 sage: R(5) / R(0) infinity sage: R(7) / R(49) 7^-1 """ # Input should be normalized! if very_neg_val(self.ordp): raise ZeroDivisionError("Cannot divide infinity by infinity") ans._set_exact_zero() else: return ans
def __pow__(FPElement self, _right, dummy): # NOTE: dummy ignored, always use self.prime_pow.prec_cap """ Exponentiation by an integer
EXAMPLES::
sage: R = ZpFP(11, 5) sage: R(1/2)^5 10 + 7*11 + 11^2 + 5*11^3 + 4*11^4 sage: R(1/32) 10 + 7*11 + 11^2 + 5*11^3 + 4*11^4 sage: R(1/2)^5 == R(1/32) True sage: R(3)^1000 #indirect doctest 1 + 4*11^2 + 3*11^3 + 7*11^4 sage: R(11)^-1 11^-1 """ cdef long dummyL cdef mpz_t tmp cdef Integer right cdef FPElement base, pright, ans cdef bint exact_exp elif self.parent() is _right.parent(): ## For extension elements, we need to switch to the ## fraction field sometimes in highly ramified extensions. exact_exp = False pright = _right else: self, _right = canonical_coercion(self, _right) return self.__pow__(_right, dummy) # We may assume from above that right > 0 else: # log(0) and log(infinity) not defined raise ValueError("0^x and inf^x not defined for p-adic x") # exact_pow_helper is defined in padic_template_element.pxi return ans else: finally: else: # padic_pow_helper is defined in padic_template_element.pxi dummyL = padic_pow_helper(ans.unit, self.unit, self.ordp, self.prime_pow.prec_cap, pright.unit, pright.ordp, pright.prime_pow.prec_cap, self.prime_pow) ans.ordp = 0
cdef pAdicTemplateElement _lshift_c(self, long shift): """ Multiplies self by `\pi^{shift}`.
Negative shifts may truncate the result if the parent is not a field.
EXAMPLES:
We create a floating point ring::
sage: R = ZpFP(5, 20); a = R(1000); a 3*5^3 + 5^4
Shifting to the right is the same as dividing by a power of the uniformizer `\pi` of the `p`-adic ring.::
sage: a >> 1 3*5^2 + 5^3
Shifting to the left is the same as multiplying by a power of `\pi`::
sage: a << 2 3*5^5 + 5^6 sage: a*5^2 3*5^5 + 5^6
Shifting by a negative integer to the left is the same as right shifting by the absolute value::
sage: a << -3 3 + 5 sage: a >> 3 3 + 5 """ # check both in case of overflow in sum; this case also includes self.ordp = maxordp # need to check that we're not shifting infinity if very_neg_val(self.ordp): raise ZeroDivisionError("Cannot multiply zero by infinity") ans.ordp = maxordp csetzero(ans.unit, ans.prime_pow) else:
cdef pAdicTemplateElement _rshift_c(self, long shift): """ Divides by `\pi^{shift}`.
Positive shifts may truncate the result if the parent is not a field.
EXAMPLES::
sage: R = ZpFP(997, 7); a = R(123456878908); a 964*997 + 572*997^2 + 124*997^3
Shifting to the right divides by a power of `\pi`, but dropping terms with negative valuation::
sage: a >> 3 124
A negative shift multiplies by that power of `\pi`::
sage: a >> -3 964*997^4 + 572*997^5 + 124*997^6 """ raise ZeroDivisionError("Cannot divide zero by zero") if very_neg_val(shift): raise ZeroDivisionError("Cannot divide infinity by infinity") return self cdef long diff ans._set_infinity() ans._set_exact_zero() else: else: ans._set_exact_zero() else:
def _repr_(self, mode=None, do_latex=False): """ Returns a string representation of this element.
INPUT:
- ``mode`` -- allows one to override the default print mode of the parent (default: ``None``).
- ``do_latex`` -- whether to return a latex representation or a normal one.
EXAMPLES::
sage: ZpFP(5,5)(1/3) # indirect doctest 2 + 3*5 + 5^2 + 3*5^3 + 5^4 sage: ~QpFP(5,5)(0) infinity """
def add_bigoh(self, absprec): """ Returns a new element truncated modulo `\pi^{\mbox{absprec}}`.
INPUT:
- ``absprec`` -- an integer or infinity
OUTPUT:
- a new element truncated modulo `\pi^{\mbox{absprec}}`.
EXAMPLES::
sage: R = Zp(7,4,'floating-point','series'); a = R(8); a.add_bigoh(1) 1 """ cdef long aprec, newprec else: absprec = Integer(absprec) aprec = minusmaxordp else: else:
cpdef bint _is_exact_zero(self) except -1: """ Tests whether this element is exactly zero.
EXAMPLES::
sage: R = Zp(7,4,'floating-point','series'); a = R(8); a._is_exact_zero() False sage: b = R(0); b._is_exact_zero() True """
cpdef bint _is_inexact_zero(self) except -1: """ Returns True if self is indistinguishable from zero.
EXAMPLES::
sage: R = ZpFP(7, 5) sage: R(14)._is_inexact_zero() False sage: R(0)._is_inexact_zero() True """
def is_zero(self, absprec = None): r""" Returns whether self is zero modulo `\pi^{\mbox{absprec}}`.
INPUT:
- ``absprec`` -- an integer
EXAMPLES::
sage: R = ZpFP(17, 6) sage: R(0).is_zero() True sage: R(17^6).is_zero() False sage: R(17^2).is_zero(absprec=2) True """ return True return False return self.ordp >= absprec absprec = Integer(absprec)
def __nonzero__(self): """ Returns True if this element is distinguishable from zero.
For most applications, explicitly specifying the power of p modulo which the element is supposed to be nonzero is preferrable.
EXAMPLES::
sage: R = ZpFP(5); a = R(0); b = R(75) sage: bool(a), bool(b) # indirect doctest (False, True) """
def is_equal_to(self, _right, absprec=None): r""" Returns whether this element is equal to ``right`` modulo `p^{\mbox{absprec}}`.
If ``absprec`` is ``None``, determines whether self and right have the same value.
INPUT:
- ``right`` -- a p-adic element with the same parent - ``absprec`` -- a positive integer or ``None`` (default: ``None``)
EXAMPLES::
sage: R = ZpFP(2, 6) sage: R(13).is_equal_to(R(13)) True sage: R(13).is_equal_to(R(13+2^10)) True sage: R(13).is_equal_to(R(17), 2) True sage: R(13).is_equal_to(R(17), 5) False """ cdef FPElement right cdef long aprec, rprec else: right = self.parent().coerce(_right) return False return False absprec = Integer(absprec) return False return False return False return ccmp(self.unit, right.unit, self.prime_pow.prec_cap, False, False, self.prime_pow) == 0 rprec = self.prime_pow.prec_cap right.unit, rprec, rprec < self.prime_pow.prec_cap, rprec < right.prime_pow.prec_cap,
cdef int _cmp_units(self, pAdicGenericElement _right) except -2: """ Comparison of units, used in equality testing.
EXAMPLES::
sage: R = ZpFP(5) sage: a = R(17); b = R(0,3); c = R(85,7); d = R(2, 1) sage: any([a == b, a == c, b == c, b == d, c == d, a == d]) # indirect doctest False sage: all([a == a, b == b, c == c, d == d]) True """
cdef pAdicTemplateElement lift_to_precision_c(self, long absprec): """ Lifts this element to another with precision at least absprec.
Since floating point elements don't track precision, this function just returns the same element.
EXAMPLES::
sage: R = ZpFP(5); sage: a = R(77, 2); a 2 sage: a.lift_to_precision(17) # indirect doctest 2 """
def _teichmuller_set_unsafe(self): """ Sets this element to the Teichmuller representative with the same residue.
.. WARNING::
This function modifies the element, which is not safe. Elements are supposed to be immutable.
EXAMPLES::
sage: R = ZpFP(17,5); a = R(11) sage: a 11 sage: a._teichmuller_set_unsafe(); a 11 + 14*17 + 2*17^2 + 12*17^3 + 15*17^4 sage: E = a.expansion(lift_mode='teichmuller'); E 17-adic expansion of 11 + 14*17 + 2*17^2 + 12*17^3 + 15*17^4 (teichmuller) sage: list(E) [11 + 14*17 + 2*17^2 + 12*17^3 + 15*17^4, 0, 0, 0, 0]
Note that if you set an element which is congruent to 0 you get 0 to maximum precision::
sage: b = R(17*5); b 5*17 sage: b._teichmuller_set_unsafe(); b 0 """ else:
def polynomial(self, var='x'): """ Returns a polynomial over the base ring that yields this element when evaluated at the generator of the parent.
INPUT:
- ``var`` -- string, the variable name for the polynomial
EXAMPLES::
sage: K.<a> = QqFP(5^3) sage: a.polynomial() x sage: a.polynomial(var='y') y sage: (5*a^2 + K(25, 4)).polynomial() 5*x^2 + 5^2 """ return S([]) return S([~R(0)]) else:
def precision_absolute(self): """ The absolute precision of this element.
EXAMPLES::
sage: R = Zp(7,4,'floating-point'); a = R(7); a.precision_absolute() 5 sage: R(0).precision_absolute() +Infinity sage: (~R(0)).precision_absolute() -Infinity """
def precision_relative(self): r""" The relative precision of this element.
EXAMPLES::
sage: R = Zp(7,4,'floating-point'); a = R(7); a.precision_relative() 4 sage: R(0).precision_relative() 0 sage: (~R(0)).precision_relative() 0 """ else:
cpdef pAdicTemplateElement unit_part(FPElement self): r""" Returns the unit part of this element.
If the valuation of this element is positive, then the high digits of the result will be zero.
EXAMPLES::
sage: R = Zp(17, 4, 'floating-point') sage: R(5).unit_part() 5 sage: R(18*17).unit_part() 1 + 17 sage: R(0).unit_part() Traceback (most recent call last): ... ValueError: unit part of 0 and infinity not defined sage: type(R(5).unit_part()) <type 'sage.rings.padics.padic_floating_point_element.pAdicFloatingPointElement'> sage: R = ZpFP(5, 5); a = R(75); a.unit_part() 3 """
cdef long valuation_c(self): """ Returns the valuation of this element.
If this element is an exact zero, returns ``maxordp``, which is defined as ``(1L << (sizeof(long) * 8 - 2))-1``.
If this element is infinity, returns ``-maxordp``.
TESTS::
sage: R = ZpFP(5, 5); R(1).valuation() #indirect doctest 0 sage: R = Zp(17, 4,'floating-point') sage: a = R(2*17^2) sage: a.valuation() 2 sage: R = Zp(5, 4,'floating-point') sage: R(0).valuation() +Infinity sage: (~R(0)).valuation() -Infinity sage: R(1).valuation() 0 sage: R(2).valuation() 0 sage: R(5).valuation() 1 sage: R(10).valuation() 1 sage: R(25).valuation() 2 sage: R(50).valuation() 2 """
cpdef val_unit(self, p=None): """ Returns a 2-tuple, the first element set to the valuation of this element, and the second to the unit part.
If this element is either zero or infinity, raises an error.
EXAMPLES::
sage: R = ZpFP(5,5) sage: a = R(75); b = a - a sage: a.val_unit() (2, 3) sage: b.val_unit() Traceback (most recent call last): ... ValueError: unit part of 0 and infinity not defined """ raise ValueError('Ring (%s) residue field of the wrong characteristic.'%self.parent())
def __hash__(self): """ Hashing.
EXAMPLES::
sage: R = ZpFP(11, 5) sage: hash(R(3)) == hash(3) True """ return 314159
cdef class pAdicCoercion_ZZ_FP(RingHomomorphism): """ The canonical inclusion from the integer ring to a floating point ring.
EXAMPLES::
sage: f = ZpFP(5).coerce_map_from(ZZ); f Ring morphism: From: Integer Ring To: 5-adic Ring with floating precision 20
TESTS::
sage: TestSuite(f).run()
""" def __init__(self, R): """ Initialization.
EXAMPLES::
sage: f = ZpFP(5).coerce_map_from(ZZ); type(f) <type 'sage.rings.padics.padic_floating_point_element.pAdicCoercion_ZZ_FP'> """
cdef dict _extra_slots(self): """ Helper for copying and pickling.
EXAMPLES::
sage: f = ZpFP(5).coerce_map_from(ZZ) sage: g = copy(f) # indirect doctest sage: g == f True sage: g(6) 1 + 5 sage: g(6) == f(6) True """
cdef _update_slots(self, dict _slots): """ Helper for copying and pickling.
EXAMPLES::
sage: f = ZpFP(5).coerce_map_from(ZZ) sage: g = copy(f) # indirect doctest sage: g == f True sage: g(6) 1 + 5 sage: g(6) == f(6) True """
cpdef Element _call_(self, x): """ Evaluation.
EXAMPLES::
sage: f = ZpFP(5).coerce_map_from(ZZ) sage: f(0).parent() 5-adic Ring with floating precision 20 sage: f(5) 5 """
cpdef Element _call_with_args(self, x, args=(), kwds={}): """ This function is used when some precision cap is passed in (relative or absolute or both).
INPUT:
- ``x`` -- an Integer
- ``absprec``, or the first positional argument -- the maximum absolute precision (unused for floating point elements).
- ``relprec``, or the second positional argument -- the maximum relative precision (unused for floating point elements)
EXAMPLES::
sage: R = ZpFP(5,4) sage: type(R(10,2)) <type 'sage.rings.padics.padic_floating_point_element.pAdicFloatingPointElement'> sage: R(30,2) 5 sage: R(30,3,2) 5 + 5^2 sage: R(30,absprec=2) 5 sage: R(30,relprec=2) 5 + 5^2 sage: R(30,absprec=1) 0 sage: R(30,empty=True) 0 """ cdef long val, aprec, rprec
def section(self): """ Returns a map back to ZZ that approximates an element of this `p`-adic ring by an integer.
EXAMPLES::
sage: f = ZpFP(5).coerce_map_from(ZZ).section() sage: f(ZpFP(5)(-1)) - 5^20 -1 """ import copy self._section = copy.copy(self._section)
cdef class pAdicConvert_FP_ZZ(RingMap): """ The map from a floating point ring back to ZZ that returns the smallest non-negative integer approximation to its input which is accurate up to the precision.
If the input is not in the closure of the image of ZZ, raises a ValueError.
EXAMPLES::
sage: f = ZpFP(5).coerce_map_from(ZZ).section(); f Set-theoretic ring morphism: From: 5-adic Ring with floating precision 20 To: Integer Ring """ def __init__(self, R): """ Initialization.
EXAMPLES::
sage: f = ZpFP(5).coerce_map_from(ZZ).section(); type(f) <type 'sage.rings.padics.padic_floating_point_element.pAdicConvert_FP_ZZ'> sage: f.category() Category of homsets of sets """ else:
cpdef Element _call_(self, _x): """ Evaluation.
EXAMPLES::
sage: f = QpFP(5).coerce_map_from(ZZ).section() sage: f(QpFP(5)(-1)) - 5^20 -1 sage: f(QpFP(5)(5)) 5 sage: f(QpFP(5)(0)) 0 sage: f(~QpFP(5)(5)) Traceback (most recent call last): ... ValueError: negative valuation sage: f(~QpFP(5)(0)) Traceback (most recent call last): ... ValueError: Infinity cannot be converted to a rational """ else:
cdef class pAdicCoercion_QQ_FP(RingHomomorphism): """ The canonical inclusion from the rationals to a floating point field.
EXAMPLES::
sage: f = QpFP(5).coerce_map_from(QQ); f Ring morphism: From: Rational Field To: 5-adic Field with floating precision 20
TESTS::
sage: TestSuite(f).run()
""" def __init__(self, R): """ Initialization.
EXAMPLES::
sage: f = QpFP(5).coerce_map_from(QQ); type(f) <type 'sage.rings.padics.padic_floating_point_element.pAdicCoercion_QQ_FP'> """
cdef dict _extra_slots(self): """ Helper for copying and pickling.
EXAMPLES::
sage: f = QpFP(5).coerce_map_from(QQ) sage: g = copy(f) # indirect doctest sage: g Ring morphism: From: Rational Field To: 5-adic Field with floating precision 20 sage: g == f True sage: g is f False sage: g(6) 1 + 5 sage: g(6) == f(6) True """
cdef _update_slots(self, dict _slots): """ Helper for copying and pickling.
EXAMPLES::
sage: f = QpFP(5).coerce_map_from(QQ) sage: g = copy(f) # indirect doctest sage: g Ring morphism: From: Rational Field To: 5-adic Field with floating precision 20 sage: g == f True sage: g is f False sage: g(6) 1 + 5 sage: g(6) == f(6) True """
cpdef Element _call_(self, x): """ Evaluation.
EXAMPLES::
sage: f = QpFP(5).coerce_map_from(QQ) sage: f(0).parent() 5-adic Field with floating precision 20 sage: f(1/5) 5^-1 sage: f(1/4) 4 + 3*5 + 3*5^2 + 3*5^3 + 3*5^4 + 3*5^5 + 3*5^6 + 3*5^7 + 3*5^8 + 3*5^9 + 3*5^10 + 3*5^11 + 3*5^12 + 3*5^13 + 3*5^14 + 3*5^15 + 3*5^16 + 3*5^17 + 3*5^18 + 3*5^19 sage: f(1/4, 5) 4 + 3*5 + 3*5^2 + 3*5^3 + 3*5^4 """
cpdef Element _call_with_args(self, x, args=(), kwds={}): """ This function is used when some precision cap is passed in (relative or absolute or both).
See the documentation for :meth:`pAdicCappedRelativeElement.__init__` for more details.
EXAMPLES::
sage: R = QpFP(5,4) sage: type(R(10/3,2)) <type 'sage.rings.padics.padic_floating_point_element.pAdicFloatingPointElement'> sage: R(10/3,2) 4*5 sage: R(10/3,3,1) 4*5 sage: R(10/3,absprec=2) 4*5 sage: R(10/3,relprec=2) 4*5 + 5^2 sage: R(10/3,absprec=1) 0 sage: R(3/100,absprec=-1) 2*5^-2 """ cdef long val, aprec, rprec cdef FPElement ans
def section(self): """ Returns a map back to the rationals that approximates an element by a rational number.
EXAMPLES::
sage: f = QpFP(5).coerce_map_from(QQ).section() sage: f(QpFP(5)(1/4)) 1/4 sage: f(QpFP(5)(1/5)) 1/5 """ import copy self._section = copy.copy(self._section)
cdef class pAdicConvert_FP_QQ(RingMap): """ The map from the floating point ring back to the rationals that returns a rational approximation of its input.
EXAMPLES::
sage: f = QpFP(5).coerce_map_from(QQ).section(); f Set-theoretic ring morphism: From: 5-adic Field with floating precision 20 To: Rational Field """ def __init__(self, R): """ Initialization.
EXAMPLES::
sage: f = QpFP(5).coerce_map_from(QQ).section(); type(f) <type 'sage.rings.padics.padic_floating_point_element.pAdicConvert_FP_QQ'> sage: f.category() Category of homsets of sets with partial maps """
cpdef Element _call_(self, _x): """ Evaluation.
EXAMPLES::
sage: f = QpFP(5).coerce_map_from(QQ).section() sage: f(QpFP(5)(-1)) -1 sage: f(QpFP(5)(0)) 0 sage: f(QpFP(5)(1/5)) 1/5 """ raise ValueError("Infinity cannot be converted to a rational") else:
cdef class pAdicConvert_QQ_FP(Morphism): """ The inclusion map from QQ to a floating point ring.
EXAMPLES::
sage: f = ZpFP(5).convert_map_from(QQ); f Generic morphism: From: Rational Field To: 5-adic Ring with floating precision 20 """ def __init__(self, R): """ Initialization.
EXAMPLES::
sage: f = ZpFP(5).convert_map_from(QQ); type(f) <type 'sage.rings.padics.padic_floating_point_element.pAdicConvert_QQ_FP'> """
cdef dict _extra_slots(self): """ Helper for copying and pickling.
EXAMPLES::
sage: f = ZpFP(5).convert_map_from(QQ) sage: g = copy(f) # indirect doctest sage: g == f # todo: comparison not implemented True sage: g(1/6) 1 + 4*5 + 4*5^3 + 4*5^5 + 4*5^7 + 4*5^9 + 4*5^11 + 4*5^13 + 4*5^15 + 4*5^17 + 4*5^19 sage: g(1/6) == f(1/6) True """
cdef _update_slots(self, dict _slots): """ Helper for copying and pickling.
EXAMPLES::
sage: f = ZpFP(5).convert_map_from(QQ) sage: g = copy(f) # indirect doctest sage: g == f # todo: comparison not implemented True sage: g(1/6) 1 + 4*5 + 4*5^3 + 4*5^5 + 4*5^7 + 4*5^9 + 4*5^11 + 4*5^13 + 4*5^15 + 4*5^17 + 4*5^19 sage: g(1/6) == f(1/6) True """
cpdef Element _call_(self, x): """ Evaluation.
EXAMPLES::
sage: f = ZpFP(5,4).convert_map_from(QQ) sage: f(1/7) 3 + 3*5 + 2*5^3 sage: f(0/1) 0 """
cpdef Element _call_with_args(self, x, args=(), kwds={}): """ This function is used when some precision cap is passed in (relative or absolute or both).
INPUT:
- ``x`` -- a Rational
- ``absprec``, or the first positional argument -- the maximum absolute precision (unused for floating point elements).
- ``relprec``, or the second positional argument -- the maximum relative precision (unused for floating point elements)
EXAMPLES::
sage: R = ZpFP(5,4) sage: type(R(1/7,2)) <type 'sage.rings.padics.padic_floating_point_element.pAdicFloatingPointElement'> sage: R(1/7,2) 3 + 3*5 sage: R(1/7,3,2) 3 + 3*5 sage: R(1/7,absprec=2) 3 + 3*5 sage: R(5/7,relprec=2) 3*5 + 3*5^2 sage: R(1/7,absprec=1) 3 sage: R(1/7,empty=True) 0 """ cdef long val, aprec, rprec raise ValueError("p divides the denominator")
cdef class pAdicCoercion_FP_frac_field(RingHomomorphism): r""" The canonical inclusion of `\ZZ_q` into its fraction field.
EXAMPLES::
sage: R.<a> = ZqFP(27, implementation='FLINT') sage: K = R.fraction_field() sage: f = K.coerce_map_from(R); f Ring morphism: From: Unramified Extension in a defined by x^3 + 2*x + 1 with floating precision 20 over 3-adic Ring To: Unramified Extension in a defined by x^3 + 2*x + 1 with floating precision 20 over 3-adic Field
TESTS::
sage: TestSuite(f).run()
""" def __init__(self, R, K): r""" Initialization.
EXAMPLES::
sage: R.<a> = ZqFP(27, implementation='FLINT') sage: K = R.fraction_field() sage: f = K.coerce_map_from(R); type(f) <type 'sage.rings.padics.qadic_flint_FP.pAdicCoercion_FP_frac_field'> """
cpdef Element _call_(self, _x): r""" Evaluation.
EXAMPLES::
sage: R.<a> = ZqFP(27, implementation='FLINT') sage: K = R.fraction_field() sage: f = K.coerce_map_from(R) sage: f(a) a sage: f(R(0)) 0 """
cpdef Element _call_with_args(self, _x, args=(), kwds={}): r""" This function is used when some precision cap is passed in (relative or absolute or both).
See the documentation for :meth:`pAdicCappedAbsoluteElement.__init__` for more details.
EXAMPLES::
sage: R.<a> = ZqFP(27, implementation='FLINT') sage: K = R.fraction_field() sage: f = K.coerce_map_from(R) sage: f(a, 3) a sage: b = 9*a + 27 sage: f(b, 3) a*3^2 sage: f(b, 4, 1) a*3^2 sage: f(b, 4, 3) a*3^2 + 3^3 sage: f(b, absprec=4) a*3^2 + 3^3 sage: f(b, relprec=3) a*3^2 + 3^3 sage: f(b, absprec=1) 0 sage: f(R(0)) 0 """ cdef long aprec, rprec else: else:
def section(self): r""" Returns a map back to the ring that converts elements of non-negative valuation.
EXAMPLES::
sage: R.<a> = ZqFP(27, implementation='FLINT') sage: K = R.fraction_field() sage: f = K.coerce_map_from(R) sage: f(K.gen()) a """
cdef dict _extra_slots(self): r""" Helper for copying and pickling.
TESTS::
sage: R.<a> = ZqFP(27, implementation='FLINT') sage: K = R.fraction_field() sage: f = K.coerce_map_from(R) sage: g = copy(f) # indirect doctest sage: g Ring morphism: From: Unramified Extension in a defined by x^3 + 2*x + 1 with floating precision 20 over 3-adic Ring To: Unramified Extension in a defined by x^3 + 2*x + 1 with floating precision 20 over 3-adic Field sage: g == f True sage: g is f False sage: g(a) a sage: g(a) == f(a) True """
cdef _update_slots(self, dict _slots): r""" Helper for copying and pickling.
TESTS::
sage: R.<a> = ZqFP(9, implementation='FLINT') sage: K = R.fraction_field() sage: f = K.coerce_map_from(R) sage: g = copy(f) # indirect doctest sage: g Ring morphism: From: Unramified Extension in a defined by x^2 + 2*x + 2 with floating precision 20 over 3-adic Ring To: Unramified Extension in a defined by x^2 + 2*x + 2 with floating precision 20 over 3-adic Field sage: g == f True sage: g is f False sage: g(a) a sage: g(a) == f(a) True
"""
cdef class pAdicConvert_FP_frac_field(Morphism): r""" The section of the inclusion from `\ZZ_q` to its fraction field.
EXAMPLES::
sage: R.<a> = ZqFP(27, implementation='FLINT') sage: K = R.fraction_field() sage: f = R.convert_map_from(K); f Generic morphism: From: Unramified Extension in a defined by x^3 + 2*x + 1 with floating precision 20 over 3-adic Field To: Unramified Extension in a defined by x^3 + 2*x + 1 with floating precision 20 over 3-adic Ring """ def __init__(self, K, R): r""" Initialization.
EXAMPLES::
sage: R.<a> = ZqFP(27, implementation='FLINT') sage: K = R.fraction_field() sage: f = R.convert_map_from(K); type(f) <type 'sage.rings.padics.qadic_flint_FP.pAdicConvert_FP_frac_field'> """
cpdef Element _call_(self, _x): r""" Evaluation.
EXAMPLES::
sage: R.<a> = ZqFP(27, implementation='FLINT') sage: K = R.fraction_field() sage: f = R.convert_map_from(K) sage: f(K.gen()) a """
cpdef Element _call_with_args(self, _x, args=(), kwds={}): r""" This function is used when some precision cap is passed in (relative or absolute or both).
See the documentation for :meth:`pAdicCappedAbsoluteElement.__init__` for more details.
EXAMPLES::
sage: R.<a> = ZqFP(27, implementation='FLINT') sage: K = R.fraction_field() sage: f = R.convert_map_from(K); a = K(a) sage: f(a, 3) a sage: b = 9*a + 27 sage: f(b, 3) a*3^2 sage: f(b, 4, 1) a*3^2 sage: f(b, 4, 3) a*3^2 + 3^3 sage: f(b, absprec=4) a*3^2 + 3^3 sage: f(b, relprec=3) a*3^2 + 3^3 sage: f(b, absprec=1) 0 sage: f(K(0)) 0 """ cdef long aprec, rprec else: else:
cdef dict _extra_slots(self): r""" Helper for copying and pickling.
TESTS::
sage: R.<a> = ZqFP(27, implementation='FLINT') sage: K = R.fraction_field() sage: f = R.convert_map_from(K) sage: a = K(a) sage: g = copy(f) # indirect doctest sage: g Generic morphism: From: Unramified Extension in a defined by x^3 + 2*x + 1 with floating precision 20 over 3-adic Field To: Unramified Extension in a defined by x^3 + 2*x + 1 with floating precision 20 over 3-adic Ring sage: g == f True sage: g is f False sage: g(a) a sage: g(a) == f(a) True """
cdef _update_slots(self, dict _slots): r""" Helper for copying and pickling.
TESTS::
sage: R.<a> = ZqFP(9, implementation='FLINT') sage: K = R.fraction_field() sage: f = R.convert_map_from(K) sage: a = K(a) sage: g = copy(f) # indirect doctest sage: g Generic morphism: From: Unramified Extension in a defined by x^2 + 2*x + 2 with floating precision 20 over 3-adic Field To: Unramified Extension in a defined by x^2 + 2*x + 2 with floating precision 20 over 3-adic Ring sage: g == f True sage: g is f False sage: g(a) a sage: g(a) == f(a) True
"""
def unpickle_fpe_v2(cls, parent, unit, ordp): """ Unpickles a floating point element.
EXAMPLES::
sage: from sage.rings.padics.padic_floating_point_element import pAdicFloatingPointElement, unpickle_fpe_v2 sage: R = ZpFP(5) sage: a = unpickle_fpe_v2(pAdicFloatingPointElement, R, 17, 2); a 2*5^2 + 3*5^3 sage: a.parent() is R True """ |