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
""" Lazy real and complex numbers
These classes are very lazy, in the sense that it doesn't really do anything but simply sits between exact rings of characteristic 0 and the real numbers. The values are actually computed when they are cast into a field of fixed precision.
The main purpose of these classes is to provide a place for exact rings (e.g. number fields) to embed for the coercion model (as only one embedding can be specified in the forward direction). """
#***************************************************************************** # Copyright (C) 2008 Robert Bradshaw <robertwb@math.washington.edu> # # 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 __future__ import absolute_import, division, print_function
import math
cdef add, sub, mul, truediv, pow, neg, inv from operator import add, sub, mul, pow, neg, inv, truediv
cdef canonical_coercion from sage.structure.element import canonical_coercion from sage.structure.element cimport parent from sage.structure.richcmp cimport richcmp
import sage.categories.map from sage.categories.morphism cimport Morphism from sage.rings.ring cimport Field import sage.rings.infinity from sage.rings.integer import Integer
cdef QQ, RR, CC, RealField, ComplexField from sage.rings.rational_field import QQ from sage.rings.real_mpfr import RR, RealField from sage.rings.complex_field import ComplexField CC = ComplexField(53)
cdef _QQx = None
cdef QQx(): global _QQx _QQx = QQ['x']
cdef named_unops = [ 'sqrt', 'erf', 'gamma', 'abs', 'floor', 'ceil', 'trunc', 'exp', 'log', 'log10', 'log2', 'sin', 'cos', 'tan', 'arcsin', 'arccos', 'arctan', 'csc', 'sec', 'cot', 'sinh', 'cosh', 'tanh', 'arcsinh', 'arccosh', 'arctanh' ]
cdef named_constants = [ 'pi', 'e', 'euler_constant', 'catalan_constant' ]
cdef class LazyField(Field): """ The base class for lazy real fields.
.. WARNING::
LazyField uses :meth:`__getattr__`, to implement::
sage: CLF.pi 3.141592653589794?
I (NT, 20/04/2012) did not manage to have ``__getattr__`` call :meth:`Parent.__getattr__` in case of failure; hence we can't use this ``__getattr__`` trick for extension types to recover the methods from categories. Therefore, at this point, no concrete subclass of this class should be an extension type (which is probably just fine)::
sage: RLF.__class__ <class 'sage.rings.real_lazy.RealLazyField_class_with_category'> sage: CLF.__class__ <class 'sage.rings.real_lazy.ComplexLazyField_class_with_category'> """ def __init__(self, base=None, names=None, normalize=True, category=None): """ Initialize ``self``.
EXAMPLES::
sage: RLF # indirect doctest Real Lazy Field
""" Field.__init__(self,base or self, names=names, normalize=normalize, category=category)
Element = LazyWrapper
def __getattr__(self, name): """ Simulates a list of methods found on the real/complex rings.
EXAMPLES::
sage: a = CLF.pi() * CLF.I(); a 3.141592653589794?*I sage: CDF(a) 3.141592653589793*I """ else:
cpdef _coerce_map_from_(self, R): r""" The only things that coerce into this ring are exact rings that embed into `\RR` or `\CC` (depending on whether this field is real or complex), that is, exact rings that coerce into all rings into which this ring coerces.
.. NOTE::
The rings into which this ring coerces are currently the corresponding floating-point fields (RealField(p) or ComplexField(p)), machine-precision floating-point fields (RDF, CDF), and interval fields (RealIntervalField(p), ComplexIntervalField(p)). This method should be updated if a new parent is added that declares a coercion from RLF/CLF but not from one of these, otherwise coercions of elements of type LazyWrapper into the new ring might fail.
EXAMPLES::
sage: RLF.has_coerce_map_from(ZZ) # indirect doctest True sage: RLF.has_coerce_map_from(QQ) True sage: RLF.has_coerce_map_from(AA) True sage: CLF.has_coerce_map_from(QQbar) True sage: RLF.has_coerce_map_from(RDF) False
sage: CLF.has_coerce_map_from(QQ) True sage: CLF.has_coerce_map_from(QQbar) True sage: CLF.has_coerce_map_from(CC) False sage: CLF.has_coerce_map_from(RLF) True """ # Indirect coercions might lead to loops both in the coercion # discovery algorithm and when trying to convert LazyWrappers, # so we only consider direct coercions. # We can skip the test for a coercion to RDF/CDF since RR/CC # already coerce into it.
def algebraic_closure(self): """ Returns the algebraic closure of ``self``, i.e., the complex lazy field.
EXAMPLES::
sage: RLF.algebraic_closure() Complex Lazy Field
sage: CLF.algebraic_closure() Complex Lazy Field """
cpdef interval_field(self, prec=None): """ Abstract method to create the corresponding interval field.
TESTS::
sage: RLF.interval_field() # indirect doctest Real Interval Field with 53 bits of precision """ raise NotImplementedError("subclasses must override this method")
class RealLazyField_class(LazyField): r""" This class represents the set of real numbers to unspecified precision. For the most part it simply wraps exact elements and defers evaluation until a specified precision is requested.
Its primary use is to connect the exact rings (such as number fields) to fixed precision real numbers. For example, to specify an embedding of a number field `K` into `\RR` one can map into this field and the coercion will then be able to carry the mapping to real fields of any precision.
EXAMPLES::
sage: a = RLF(1/3) sage: a 0.3333333333333334? sage: a + 1/5 0.5333333333333334? sage: a = RLF(1/3) sage: a 0.3333333333333334? sage: a + 5 5.333333333333334? sage: RealField(100)(a+5) 5.3333333333333333333333333333
::
sage: CC.0 + RLF(1/3) 0.333333333333333 + 1.00000000000000*I sage: ComplexField(200).0 + RLF(1/3) 0.33333333333333333333333333333333333333333333333333333333333 + 1.0000000000000000000000000000000000000000000000000000000000*I
TESTS::
sage: TestSuite(RLF).run() """ def interval_field(self, prec=None): """ Returns the interval field that represents the same mathematical field as ``self``.
EXAMPLES::
sage: RLF.interval_field() Real Interval Field with 53 bits of precision sage: RLF.interval_field(200) Real Interval Field with 200 bits of precision """ else:
def construction(self): """ Returns the functorial construction of ``self``, namely, the completion of the rationals at infinity to infinite precision.
EXAMPLES::
sage: c, S = RLF.construction(); S Rational Field sage: RLF == c(S) True """
def _latex_(self): r""" Return a latex representation of ``self``.
EXAMPLES::
sage: latex(RLF) # indirect doctest \Bold{R} """
def gen(self, i=0): """ Return the `i`-th generator of ``self``.
EXAMPLES::
sage: RLF.gen() 1 """ else: raise ValueError("RLF has only one generator.")
def _repr_(self): """ Return a string representation of ``self``.
EXAMPLES::
sage: RealLazyField() Real Lazy Field """
def __hash__(self): """ Return the hash of ``self``.
EXAMPLES::
sage: hash(RLF) % 2^32 == hash(str(RLF)) % 2^32 True """
def __reduce__(self): """ For pickling.
TESTS::
sage: RLF == loads(dumps(RLF)) True sage: RLF is loads(dumps(RLF)) True """
RLF = RealLazyField_class()
def RealLazyField(): """ Return the lazy real field.
EXAMPLES:
There is only one lazy real field::
sage: RealLazyField() is RealLazyField() True """
class ComplexLazyField_class(LazyField): """ This class represents the set of complex numbers to unspecified precision. For the most part it simply wraps exact elements and defers evaluation until a specified precision is requested.
For more information, see the documentation of the :class:`RLF <sage.rings.real_lazy.RealLazyField_class>`.
EXAMPLES::
sage: a = CLF(-1).sqrt() sage: a 1*I sage: CDF(a) 1.0*I sage: ComplexField(200)(a) 1.0000000000000000000000000000000000000000000000000000000000*I
TESTS::
sage: TestSuite(CLF).run(skip=["_test_prod"])
.. NOTE::
The following ``TestSuite`` failure::
sage: CLF._test_prod() Traceback (most recent call last): ... AssertionError: False is not true
is due to (acceptable?) numerical noise::
sage: x = CLF.I sage: x*x == x^2 False sage: x*x -1 sage: x^2 -0.9999999999999999? + 0.?e-15*I """ def __init__(self): """ This lazy field doesn't evaluate its elements until they are cast into a field of fixed precision.
EXAMPLES: sage: a = RLF(1/3); a 0.3333333333333334? sage: Reals(200)(a) 0.33333333333333333333333333333333333333333333333333333333333 """ LazyField.__init__(self, base=RLF) self._populate_coercion_lists_(coerce_list=[LazyWrapperMorphism(RLF, self)])
def interval_field(self, prec=None): """ Returns the interval field that represents the same mathematical field as ``self``.
EXAMPLES::
sage: CLF.interval_field() Complex Interval Field with 53 bits of precision sage: CLF.interval_field(333) Complex Interval Field with 333 bits of precision sage: CLF.interval_field() is CIF True """ else:
def gen(self, i=0): """ Return the `i`-th generator of ``self``.
EXAMPLES::
sage: CLF.gen() 1*I sage: ComplexField(100)(CLF.gen()) 1.0000000000000000000000000000*I """ else: raise ValueError("CLF has only one generator.")
def construction(self): """ Returns the functorial construction of ``self``, namely, algebraic closure of the real lazy field.
EXAMPLES::
sage: c, S = CLF.construction(); S Real Lazy Field sage: CLF == c(S) True """
def _latex_(self): r""" Return a latex representation of ``self``.
EXAMPLES::
sage: latex(CLF) # indirect doctest \Bold{C} """
def _repr_(self): """ Return a string representation of ``self``.
EXAMPLES::
sage: CLF Complex Lazy Field """
def __hash__(self): """ Return the hash of ``self``.
EXAMPLES::
sage: hash(CLF) % 2^32 == hash(str(CLF)) % 2^32 True """
def __reduce__(self): """ For pickling.
TESTS::
sage: CLF == loads(dumps(CLF)) True sage: CLF is loads(dumps(CLF)) True """
CLF = ComplexLazyField_class()
def ComplexLazyField(): """ Returns the lazy complex field.
EXAMPLES:
There is only one lazy complex field::
sage: ComplexLazyField() is ComplexLazyField() True """
cdef int get_new_prec(R, int depth) except -1: """ There are depth operations, so we want at least that many more digits of precision.
Field creation may be expensive, so we want to avoid incrementing by 1 so that it is more likely for cached fields to be used. """ cdef int needed_prec = R.prec() needed_prec += depth if needed_prec % 10 != 0: needed_prec += 10 - needed_prec % 10 return needed_prec
cdef class LazyFieldElement(FieldElement):
cpdef _add_(left, right): """ Add ``left`` with ``right``.
EXAMPLES::
sage: RLF(5) + RLF(1/2) # indirect doctest 5.5000000000000000? """ except TypeError: pass
cpdef _sub_(left, right): """ Subtract ``right`` from ``left``.
EXAMPLES::
sage: CLF(5) - 2 # indirect doctest 3 """ except TypeError: pass return left._new_binop(left, right, sub)
cpdef _mul_(left, right): """ Mutliply ``left`` with ``right``.
EXAMPLES::
sage: CLF(10) * RLF(5) # indirect doctest 50 """ except TypeError: pass
cpdef _div_(left, right): """ Divide ``left`` by ``right``.
EXAMPLES::
sage: a = RLF(1) / RLF(6); a # indirect doctest 0.1666666666666667? sage: Reals(300)(a) 0.166666666666666666666666666666666666666666666666666666666666666666666666666666666666666667 """ except TypeError: pass
def __pow__(left, right, dummy): """ Raise ``left`` to the ``right`` power.
EXAMPLES::
sage: a = RLF(2) ^ (1/2); a 1.414213562373095? sage: Reals(300)(a) 1.41421356237309504880168872420969807856967187537694807317667973799073247846210703885038753 """ try: return left._new_wrapper((<LazyWrapper>left)._value ** (<LazyWrapper>right)._value) except TypeError: pass left = (<LazyFieldElement>right)._new_wrapper(left)
def __neg__(self): """ Return the negation of ``self``.
EXAMPLES::
sage: -RLF(7) -7 """
def __invert__(self): """ Take the reciprocal of ``self``.
EXAMPLES::
sage: a = ~RLF(6); a 0.1666666666666667? sage: Reals(90)(a) 0.16666666666666666666666667 """
cpdef _richcmp_(self, other, int op): """ If things are being wrapped, tries to compare values. That failing, it tries to compare intervals, which may return a false negative.
EXAMPLES::
sage: RLF(3) == RLF(9/3) True sage: RLF(3) == RLF(4) False sage: RLF(3) < RLF(5/3) False
TESTS::
sage: from sage.rings.real_lazy import LazyBinop sage: RLF(3) < LazyBinop(RLF, 5, 3, operator.truediv) False sage: from sage.rings.real_lazy import LazyWrapper sage: LazyWrapper(RLF, 3) < LazyWrapper(RLF, 5/3) False sage: from sage.rings.real_lazy import LazyUnop sage: RLF(3) < LazyUnop(RLF, 2, sqrt) False sage: from sage.rings.real_lazy import LazyNamedUnop sage: RLF(3) < LazyNamedUnop(RLF, 0, 'sin') False sage: from sage.rings.real_lazy import LazyConstant sage: RLF(3) < LazyConstant(RLF, 'e') False """ pass
def __hash__(self): """ Return the hash value of ``self``.
EXAMPLES::
sage: a = RLF(3) sage: hash(a) 3 """
cdef LazyFieldElement _new_wrapper(self, value):
cdef LazyFieldElement _new_binop(self, LazyFieldElement left, LazyFieldElement right, op):
cdef LazyFieldElement _new_unop(self, LazyFieldElement arg, op):
def _repr_(self): """ The string representation of ``self`` is an interval in which ``self`` is contained.
EXAMPLES::
sage: RLF(3) # indirect doctest 3 sage: RLF(1/3) 0.3333333333333334? """
def approx(self): """ Returns ``self`` as an element of an interval field.
EXAMPLES::
sage: CLF(1/6).approx() 0.1666666666666667? sage: CLF(1/6).approx().parent() Complex Interval Field with 53 bits of precision
When the absolute value is involved, the result might be real::
sage: z = exp(CLF(1 + I/2)); z 2.38551673095914? + 1.303213729686996?*I sage: r = z.abs(); r 2.71828182845905? sage: parent(z.approx()) Complex Interval Field with 53 bits of precision sage: parent(r.approx()) Real Interval Field with 53 bits of precision """
def _real_double_(self, R): """ Return ``self`` as a real double.
EXAMPLES::
sage: a = RLF(3) sage: RDF(a) # indirect doctest 3.0 """
def _complex_double_(self, R): """ Return ``self`` as a complex double.
EXAMPLES::
sage: a = RLF(5) sage: CDF(a) # indirect doctest 5.0 sage: a = CLF(-1)^(1/4) sage: CDF(a) 0.7071067811865476 + 0.7071067811865475*I """ return self.eval(R)
def _generic_(self, R): """ Return ``self`` in a generic ring ``R``.
EXAMPLES::
sage: a = RLF(2/3) sage: RR(a) # indirect doctest 0.666666666666667 sage: RR(a^2) 0.444444444444444 """
_real_mpfi_ = _complex_mpfi_ =_mpfr_ = _complex_mpfr_field_ = _generic_
def __complex__(self): """ Return ``self`` as a complex.
EXAMPLES::
sage: complex(CLF(-1)^(1/4)) (0.707106781186547...+0.707106781186547...j) """
cpdef eval(self, R): """ Abstract method for converting ``self`` into an element of ``R``.
EXAMPLES::
sage: a = RLF(12) sage: a.eval(ZZ) 12 """ raise NotImplementedError("Subclasses must override this method.")
cpdef int depth(self): """ Abstract method for returning the depth of ``self`` as an arithmetic expression.
This is the maximum number of dependent intermediate expressions when evaluating ``self``, and is used to determine the precision needed to get the final result to the desired number of bits.
It is equal to the maximum of the right and left depths, plus one.
EXAMPLES::
sage: from sage.rings.real_lazy import LazyBinop sage: a = LazyBinop(RLF, 6, 8, operator.mul) sage: a.depth() 1 """ raise NotImplementedError("Subclasses must override this method.")
def __dir__(self): """ Adds the named_unops to ``__dir__`` so that tab completion works.
TESTS::
sage: "log" in RLF(sqrt(8)).__dir__() True
"""
def __getattribute__(self, name): """ Simulates a list of methods found on the real/complex mpfr classes.
EXAMPLES::
sage: a = RLF(3) sage: a.sqrt() 1.732050807568878? sage: sin(a) 0.1411200080598673? sage: RealField(160)(tanh(RLF(3))) 0.99505475368673045133188018525548847509781385470 """ else:
def continued_fraction(self): r""" Return the continued fraction of self.
EXAMPLES::
sage: a = RLF(sqrt(2)) + RLF(sqrt(3)) sage: cf = a.continued_fraction() sage: cf [3; 6, 1, 5, 7, 1, 1, 4, 1, 38, 43, 1, 3, 2, 1, 1, 1, 1, 2, 4, ...] sage: cf.convergent(100) 444927297812646558239761867973501208151173610180916865469/141414466649174973335183571854340329919207428365474086063 """ from sage.rings.continued_fraction import ContinuedFraction_real return ContinuedFraction_real(self)
def make_element(parent, *args): """ Create an element of ``parent``.
EXAMPLES::
sage: a = RLF(pi) + RLF(sqrt(1/2)) # indirect doctest sage: bool(loads(dumps(a)) == a) True """
cdef class LazyWrapper(LazyFieldElement):
cpdef int depth(self): """ Returns the depth of ``self`` as an expression, which is always 0.
EXAMPLES::
sage: RLF(4).depth() 0 """
def __init__(self, LazyField parent, value, check=True): """ A lazy element that simply wraps an element of another ring.
EXAMPLES::
sage: from sage.rings.real_lazy import LazyWrapper sage: a = LazyWrapper(RLF, 3) sage: a._value 3 """
def __neg__(self): """ Return the negation of ``self``.
EXAMPLES::
sage: from sage.rings.real_lazy import LazyWrapper sage: a = LazyWrapper(RLF, 3) sage: (-a)._value -3 """
def __invert__(self): """ Return the reciprocal of ``self``.
EXAMPLES::
sage: from sage.rings.real_lazy import LazyWrapper sage: a = LazyWrapper(RLF, 23) sage: ~a 0.04347826086956522? sage: (~a)._value 1/23 """
def __float__(self): """ EXAMPLES::
sage: from sage.rings.real_lazy import LazyWrapper sage: a = LazyWrapper(CLF, 19) sage: float(a) 19.0 """
def __nonzero__(self): """ Check to see if ``self`` is not zero.
EXAMPLES::
sage: from sage.rings.real_lazy import LazyWrapper sage: not LazyWrapper(RLF, 1) False sage: not LazyWrapper(RLF, 0) True """
def __hash__(self): """ Return the hash value of ``self``.
EXAMPLES::
sage: hash(CLF(-1)) -2 sage: hash(RLF(9/4)) == hash(9/4) True """
cpdef eval(self, R): """ Convert ``self`` into an element of ``R``.
EXAMPLES::
sage: a = RLF(12) sage: a.eval(ZZ) 12 sage: a.eval(ZZ).parent() Integer Ring """ else: raise TypeError("unable to convert {!r} to an element of {}".format(self._value, R))
def __reduce__(self): """ For pickling.
TESTS::
sage: a = RLF(2) sage: loads(dumps(a)) == a True """
def continued_fraction(self): r""" Return the continued fraction of self.
EXAMPLES::
sage: a = RLF(sqrt(2)) sage: a.continued_fraction() [1; 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, ...] """
cdef class LazyBinop(LazyFieldElement):
def __init__(self, LazyField parent, left, right, op): """ A lazy element representing a binary (usually arithmetic) operation between two other lazy elements.
EXAMPLES::
sage: from sage.rings.real_lazy import LazyBinop sage: a = LazyBinop(RLF, 2, 1/3, operator.add) sage: a 2.333333333333334? sage: Reals(200)(a) 2.3333333333333333333333333333333333333333333333333333333333 """
cpdef int depth(self): """ Return the depth of ``self`` as an arithmetic expression.
This is the maximum number of dependent intermediate expressions when evaluating ``self``, and is used to determine the precision needed to get the final result to the desired number of bits.
It is equal to the maximum of the right and left depths, plus one.
EXAMPLES::
sage: from sage.rings.real_lazy import LazyBinop sage: a = LazyBinop(RLF, 6, 8, operator.mul) sage: a.depth() 1 sage: b = LazyBinop(RLF, 2, a, operator.sub) sage: b.depth() 2 """
cpdef eval(self, R): """ Convert the operands to elements of ``R``, then perform the operation on them.
EXAMPLES::
sage: from sage.rings.real_lazy import LazyBinop sage: a = LazyBinop(RLF, 6, 8, operator.add) sage: a.eval(RR) 14.0000000000000
A bit absurd::
sage: a.eval(str) '68' """ return left - right else: # We only do a call after testing the above because it is a python call. return self._op(left, right)
def __float__(self): """ EXAMPLES::
sage: from sage.rings.real_lazy import LazyBinop sage: a = LazyBinop(RLF, 3, 1/2, operator.sub) sage: float(a) 2.5 sage: type(float(a)) <... 'float'> """ elif self._op is truediv: return left / right elif self._op is pow: return left ** right else: # We only do a call here because it is a python call. return self._op(left, right)
def __hash__(self): """ Return the hash value of ``self``.
EXAMPLES::
sage: from sage.rings.real_lazy import LazyBinop sage: a = LazyBinop(RLF, 5, 1/2, operator.sub) sage: hash(a) -1607638785 # 32-bit -7461864723258187521 # 64-bit """
def __reduce__(self): """ For pickling.
TESTS::
sage: from sage.rings.real_lazy import LazyBinop sage: a = LazyBinop(CLF, 3, 2, operator.truediv) sage: loads(dumps(a)) == a True """
cdef class LazyUnop(LazyFieldElement):
def __init__(self, LazyField parent, arg, op): """ Represents a unevaluated single function of one variable.
EXAMPLES::
sage: from sage.rings.real_lazy import LazyUnop sage: a = LazyUnop(RLF, 3, sqrt); a 1.732050807568878? sage: a._arg 3 sage: a._op <function sqrt at ...> sage: Reals(100)(a) 1.7320508075688772935274463415 sage: Reals(100)(a)^2 3.0000000000000000000000000000 """
cpdef int depth(self): """ Return the depth of ``self`` as an arithmetic expression.
This is the maximum number of dependent intermediate expressions when evaluating ``self``, and is used to determine the precision needed to get the final result to the desired number of bits.
It is equal to one more than the depth of its operand.
EXAMPLES::
sage: from sage.rings.real_lazy import LazyUnop sage: a = LazyUnop(RLF, 3, sqrt) sage: a.depth() 1 sage: b = LazyUnop(RLF, a, sin) sage: b.depth() 2 """
cpdef eval(self, R): """ Convert ``self`` into an element of ``R``.
EXAMPLES::
sage: from sage.rings.real_lazy import LazyUnop sage: a = LazyUnop(RLF, 3, sqrt) sage: a.eval(ZZ) sqrt(3) """
def __hash__(self): """ Return the hash value of ``self``.
EXAMPLES::
sage: hash(RLF(sin(1))) == hash(RLF(sin(1))) True """ return hash(self._op(hash(self._arg)))
def __float__(self): """ Convert ``self`` into a floating point.
EXAMPLES::
sage: from sage.rings.real_lazy import LazyUnop sage: a = LazyUnop(RLF, 3, sqrt) sage: float(a) 1.7320508075688772 """
def __reduce__(self): """ For pickling.
TESTS::
sage: from sage.rings.real_lazy import LazyUnop sage: a = LazyUnop(RLF, 7, sqrt) sage: float(loads(dumps(a))) == float(a) True """
cdef class LazyNamedUnop(LazyUnop):
def __init__(self, LazyField parent, arg, op, extra_args=None): """ This class is used to represent the many named methods attached to real numbers, and is instantiated by the ``__getattr__`` method of :class:`LazyElements`.
EXAMPLES::
sage: from sage.rings.real_lazy import LazyNamedUnop sage: a = LazyNamedUnop(RLF, 1, 'arcsin') sage: RR(a) 1.57079632679490 sage: a = LazyNamedUnop(RLF, 9, 'log', extra_args=(3,)) sage: RR(a) 2.00000000000000 """ raise TypeError("extra args must be a tuple")
cpdef eval(self, R): """ Convert ``self`` into an element of ``R``.
TESTS::
sage: from sage.rings.real_lazy import LazyNamedUnop sage: a = LazyNamedUnop(RLF, 4, 'sqrt') sage: RR(a) # indirect doctest 2.00000000000000 sage: a.sqrt() 1.414213562373095? sage: RealField(212)(a) 2.00000000000000000000000000000000000000000000000000000000000000 sage: float(a) 2.0
Now for some extra arguments::
sage: a = RLF(100) sage: a.log(10) 2 sage: float(a.log(10)) 2.0 """ else: else: else:
def approx(self): """ Does something reasonable with functions that are not defined on the interval fields.
TESTS::
sage: from sage.rings.real_lazy import LazyNamedUnop sage: LazyNamedUnop(RLF, 8, 'sqrt') # indirect doctest 2.828427124746190? """ except AttributeError: # not everything defined on interval fields # this is less info though, but mostly just want to print it interval_field = self._parent.interval_field() return self.eval(interval_field.middle_field())
def __hash__(self): """ Return the hash value of ``self``.
EXAMPLES::
sage: from sage.rings.real_lazy import LazyNamedUnop sage: a = LazyNamedUnop(RLF, 1, 'sin') sage: hash(a) 2110729788 """
def __float__(self): """ TESTS::
sage: from sage.rings.real_lazy import LazyNamedUnop sage: a = LazyNamedUnop(RLF, 1, 'sin') sage: float(a) 0.8414709848078965 """
def __call__(self, *args): """ TESTS::
sage: a = RLF(32) sage: a.log(2) 5 sage: float(a.log(2)) 5.0
What is going on here in the background is::
sage: from sage.rings.real_lazy import LazyNamedUnop sage: b = LazyNamedUnop(RLF, a, 'log') sage: b(2) 5 sage: b(2)._extra_args (2,) """
def __reduce__(self): """ TESTS::
sage: from sage.rings.real_lazy import LazyNamedUnop sage: a = LazyNamedUnop(RLF, 1, 'sin') sage: float(loads(dumps(a))) == float(a) True """
cdef class LazyConstant(LazyFieldElement):
cdef readonly _name cdef readonly _extra_args cdef readonly bint _is_special
def __init__(self, LazyField parent, name, extra_args=None): """ This class represents a real or complex constant (such as ``pi`` or ``I``).
TESTS::
sage: a = RLF.pi(); a 3.141592653589794? sage: RealField(300)(a) 3.14159265358979323846264338327950288419716939937510582097494459230781640628620899862803482
sage: from sage.rings.real_lazy import LazyConstant sage: a = LazyConstant(RLF, 'euler_constant') sage: RealField(200)(a) 0.57721566490153286060651209008240243104215933593992359880577 """
cpdef eval(self, R): """ Convert ``self`` into an element of ``R``.
EXAMPLES::
sage: from sage.rings.real_lazy import LazyConstant sage: a = LazyConstant(RLF, 'e') sage: RDF(a) # indirect doctest 2.718281828459045 sage: a = LazyConstant(CLF, 'I') sage: CC(a) 1.00000000000000*I """ else: raise TypeError("The complex constant I is not in this real field.") else:
def __call__(self, *args): """ TESTS::
sage: CLF.I() 1*I sage: CDF(CLF.I()) 1.0*I """
def __hash__(self): """ Return the hash value of ``self``.
TESTS::
sage: from sage.rings.real_lazy import LazyConstant sage: a = LazyConstant(RLF, 'e') sage: hash(a) 2141977644 """
def __float__(self): """ TESTS::
sage: from sage.rings.real_lazy import LazyConstant sage: a = LazyConstant(RLF, 'pi') sage: float(a) 3.141592653589793 """
def __reduce__(self): """ TESTS::
sage: from sage.rings.real_lazy import LazyConstant sage: a = LazyConstant(RLF, 'pi') sage: float(loads(dumps(a))) == float(a) True """
cdef class LazyAlgebraic(LazyFieldElement):
cdef readonly _poly cdef readonly _root_approx cdef readonly int _prec cdef readonly _quadratic_disc cdef readonly _root
def __init__(self, parent, poly, approx, int prec=0): r""" This represents an algebraic number, specified by a polynomial over `\QQ` and a real or complex approximation.
EXAMPLES::
sage: x = polygen(QQ) sage: from sage.rings.real_lazy import LazyAlgebraic sage: a = LazyAlgebraic(RLF, x^2-2, 1.5) sage: a 1.414213562373095? """ prec = approx.parent().prec() else:
cpdef eval(self, R): """ Convert ``self`` into an element of ``R``.
EXAMPLES::
sage: from sage.rings.real_lazy import LazyAlgebraic sage: a = LazyAlgebraic(CLF, QQ['x'].cyclotomic_polynomial(7), 0.6+0.8*CC.0) sage: a 0.6234898018587335? + 0.7818314824680299?*I sage: ComplexField(150)(a) # indirect doctest 0.62348980185873353052500488400423981063227473 + 0.78183148246802980870844452667405775023233452*I
sage: a = LazyAlgebraic(CLF, QQ['x'].0^2-7, -2.0) sage: RR(a) -2.64575131106459 sage: RR(a)^2 7.00000000000000 """ # Carl Witty said: # Quadratic equation faster and more accurate than roots(), # but the current code doesn't do the right thing with interval # arithmetic (it returns a point interval) so it's being disabled # for now # if self._quadratic_disc is not None: # c, b, a = self._poly.list() # if self._root_approx.real() < -b/2*a: # z = (-b - R(self._quadratic_disc).sqrt()) / (2*a) # else: # z = (-b + R(self._quadratic_disc).sqrt()) / (2*a) # if z.parent() is not R: # z = R(z) # self._root_approx = z # from sage.rings.complex_interval_field import is_IntervalField # if is_IntervalField(R): # self._root_approx = (self._root_approx.upper() + self._root_approx.lower()) / 2 # self._prec = R.prec() # return R(self._root_approx) # This could be done much more efficiently with Newton iteration, # but will require some care to make sure we get the right root, and # to the correct precision.
def __float__(self): """ TESTS::
sage: x = polygen(QQ) sage: from sage.rings.real_lazy import LazyAlgebraic sage: a = LazyAlgebraic(RLF, x^3-10, 1.5) sage: float(a) 2.154434690031883... """
def __reduce__(self): """ TESTS::
sage: from sage.rings.real_lazy import LazyAlgebraic sage: a = LazyAlgebraic(RLF, x^2-2, 1.5) sage: float(loads(dumps(a))) == float(a) True """
cdef class LazyWrapperMorphism(Morphism):
def __init__(self, domain, LazyField codomain): """ This morphism coerces elements from anywhere into lazy rings by creating a wrapper element (as fast as possible).
EXAMPLES::
sage: from sage.rings.real_lazy import LazyWrapperMorphism sage: f = LazyWrapperMorphism(QQ, RLF) sage: a = f(3); a 3 sage: type(a) <type 'sage.rings.real_lazy.LazyWrapper'> sage: a._value 3 sage: a._value.parent() Rational Field """
cpdef Element _call_(self, x): """ EXAMPLES::
sage: from sage.rings.real_lazy import LazyWrapperMorphism sage: f = LazyWrapperMorphism(QQ, CLF) sage: a = f(1/3); a # indirect doctest 0.3333333333333334? sage: type(a) <type 'sage.rings.real_lazy.LazyWrapper'> sage: Reals(100)(a) 0.33333333333333333333333333333
Note that it doesn't double-wrap lazy elements::
sage: f = LazyWrapperMorphism(RLF, CLF) sage: x = RLF(20) sage: f(x) 20 sage: f(x)._value 20 """ else:
|