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
# -*- coding: utf-8 r""" Arbitrary precision complex balls using Arb
This is a binding to the `Arb library <http://fredrikj.net/arb/>`_; it may be useful to refer to its documentation for more details.
Parts of the documentation for this module are copied or adapted from Arb's own documentation, licenced under the GNU General Public License version 2, or later.
.. SEEALSO::
- :mod:`Real balls using Arb <sage.rings.real_arb>` - :mod:`Complex interval field (using MPFI) <sage.rings.complex_interval_field>` - :mod:`Complex intervals (using MPFI) <sage.rings.complex_interval>`
Data Structure ==============
A :class:`ComplexBall` represents a complex number with error bounds. It wraps an Arb object of type ``acb_t``, which consists of a pair of real number balls representing the real and imaginary part with separate error bounds. (See the documentation of :mod:`sage.rings.real_arb` for more information.)
A :class:`ComplexBall` thus represents a rectangle `[m_1-r_1, m_1+r_1] + [m_2-r_2, m_2+r_2] i` in the complex plane. This is used in Arb instead of a disk or square representation (consisting of a complex floating-point midpoint with a single radius), since it allows implementing many operations more conveniently by splitting into ball operations on the real and imaginary parts. It also allows tracking when complex numbers have an exact (for example exactly zero) real part and an inexact imaginary part, or vice versa.
The parents of complex balls are instances of :class:`ComplexBallField`. The name ``CBF`` is bound to the complex ball field with the default precision of 53 bits::
sage: CBF is ComplexBallField() is ComplexBallField(53) True
Comparison ==========
.. WARNING::
In accordance with the semantics of Arb, identical :class:`ComplexBall` objects are understood to give permission for algebraic simplification. This assumption is made to improve performance. For example, setting ``z = x*x`` sets `z` to a ball enclosing the set `\{t^2 : t \in x\}` and not the (generally larger) set `\{tu : t \in x, u \in x\}`.
Two elements are equal if and only if they are exact and equal (in spite of the above warning, inexact balls are not considered equal to themselves)::
sage: a = CBF(1, 2) sage: b = CBF(1, 2) sage: a is b False sage: a == a True sage: a == b True
::
sage: a = CBF(1/3, 1/5) sage: b = CBF(1/3, 1/5) sage: a.is_exact() False sage: b.is_exact() False sage: a is b False sage: a == a False sage: a == b False
A ball is non-zero in the sense of usual comparison if and only if it does not contain zero::
sage: a = CBF(RIF(-0.5, 0.5)) sage: a != 0 False sage: b = CBF(1/3, 1/5) sage: b != 0 True
However, ``bool(b)`` returns ``False`` for a ball ``b`` only if ``b`` is exactly zero::
sage: bool(a) True sage: bool(b) True sage: bool(CBF.zero()) False
Coercion ========
Automatic coercions work as expected::
sage: bpol = 1/3*CBF(i) + AA(sqrt(2)) + (polygen(RealBallField(20), 'x') + QQbar(i)) sage: bpol x + [1.41421 +/- 5.09e-6] + [1.33333 +/- 3.97e-6]*I sage: bpol.parent() Univariate Polynomial Ring in x over Complex ball field with 20 bits of precision sage: bpol/3 ([0.333333 +/- 4.93e-7])*x + [0.47140 +/- 5.39e-6] + [0.44444 +/- 4.98e-6]*I
TESTS::
sage: polygen(CBF, 'x')^3 x^3
::
sage: SR.coerce(CBF(0.42 + 3.33*I)) [0.4200000000000000 +/- 1.56e-17] + [3.330000000000000 +/- 7.11e-17]*I
Check that :trac:`19839` is fixed::
sage: log(SR(CBF(0.42))).pyobject().parent() Complex ball field with 53 bits of precision
:trac:`24621`::
sage: CBF(NumberField(polygen(QQ, 'y')^3 + 20, 'a', embedding=CC(1.35,2.35)).gen()) [1.357208808297453 +/- 4.96e-16] + [2.350754612451197 +/- 7.67e-16]*I
Classes and Methods =================== """ #***************************************************************************** # Copyright (C) 2014 Clemens Heuberger <clemens.heuberger@aau.at> # 2017 Vincent Delecroix <20100.delecroix@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 __future__ import absolute_import
import operator, sys from cysignals.signals cimport sig_on, sig_str, sig_off, sig_error
import sage.categories.fields import sage.rings.number_field.number_field as number_field
cimport sage.rings.integer cimport sage.rings.rational
from cpython.float cimport PyFloat_AS_DOUBLE from cpython.int cimport PyInt_AS_LONG from cpython.object cimport Py_LT, Py_LE, Py_EQ, Py_NE, Py_GT, Py_GE from cpython.complex cimport PyComplex_FromDoubles
from sage.ext.stdsage cimport PY_NEW
from sage.libs.mpfr cimport MPFR_RNDU, MPFR_RNDD, mpfr_get_d_2exp from sage.libs.arb.types cimport ARF_RND_NEAR from sage.libs.arb.arb cimport * from sage.libs.arb.acb cimport * from sage.libs.arb.acb_calc cimport * from sage.libs.arb.acb_hypgeom cimport * from sage.libs.arb.acb_elliptic cimport * from sage.libs.arb.acb_modular cimport * from sage.libs.arb.arf cimport arf_init, arf_get_d, arf_get_mpfr, arf_set_mpfr, arf_clear, arf_set_mag, arf_set, arf_is_nan from sage.libs.arb.mag cimport (mag_init, mag_clear, mag_add, mag_set_d, MAG_BITS, mag_is_inf, mag_is_finite, mag_zero, mag_set_ui_2exp_si, mag_mul_2exp_si) from sage.libs.flint.fmpz cimport fmpz_t, fmpz_init, fmpz_get_mpz, fmpz_set_mpz, fmpz_clear, fmpz_abs from sage.libs.flint.fmpq cimport fmpq_t, fmpq_init, fmpq_set_mpq, fmpq_clear from sage.libs.gmp.mpz cimport mpz_fits_ulong_p, mpz_fits_slong_p, mpz_get_ui, mpz_get_si, mpz_sgn from sage.libs.gsl.complex cimport gsl_complex_rect
from sage.rings.real_double cimport RealDoubleElement from sage.rings.complex_double cimport ComplexDoubleElement from sage.rings.complex_field import ComplexField from sage.rings.complex_interval_field import ComplexIntervalField from sage.rings.integer_ring import ZZ from sage.rings.real_arb cimport mpfi_to_arb, arb_to_mpfi from sage.rings.real_arb import RealBallField from sage.rings.real_mpfr cimport RealField_class, RealField, RealNumber from sage.rings.ring import Field from sage.structure.element cimport Element, ModuleElement from sage.structure.parent cimport Parent from sage.structure.unique_representation import UniqueRepresentation
cdef void ComplexIntervalFieldElement_to_acb( acb_t target, ComplexIntervalFieldElement source): """ Convert a :class:`ComplexIntervalFieldElement` to an ``acb``.
INPUT:
- ``target`` -- an ``acb_t``
- ``source`` -- a :class:`ComplexIntervalFieldElement`
OUTPUT:
None. """ cdef long precision
cdef int acb_to_ComplexIntervalFieldElement( ComplexIntervalFieldElement target, const acb_t source) except -1: """ Convert an ``acb`` to a :class:`ComplexIntervalFieldElement`.
INPUT:
- ``target`` -- a :class:`ComplexIntervalFieldElement`
- ``source`` -- an ``acb_t``
OUTPUT:
A :class:`ComplexIntervalFieldElement`. """
cdef class IntegrationContext: r""" Used to wrap the integrand and hold some context information during numerical integration. """ cdef object f cdef object parent cdef object exn_type cdef object exn_obj cdef object exn_tb
cdef int acb_calc_func_callback(acb_ptr out, const acb_t inp, void * param, long order, long prec): r""" Callback used for numerical integration
TESTS::
sage: CBF.integral(lambda x, flag: 24, 0, 2) 48.00000000000000
sage: CBF.integral(lambda x, flag: "a", 0, 1) Traceback (most recent call last): ... TypeError: no canonical coercion ... to Complex ball field with 53 bits of precision
sage: def foo(*args): ....: raise RuntimeError sage: CBF.integral(foo, 0, 2) Traceback (most recent call last): ... RuntimeError
sage: points = [] sage: def foo(x, flag): ....: points.append(x) ....: return x sage: CBF.integral(foo, 0, 1) [0.50000000000000...] sage: points [[+/- 1.01], ..., [0.788...], [0.211...]] """ cdef IntegrationContext ctx cdef ComplexBall x finally: sig_on()
class ComplexBallField(UniqueRepresentation, Field): r""" An approximation of the field of complex numbers using pairs of mid-rad intervals.
INPUT:
- ``precision`` -- an integer `\ge 2`.
EXAMPLES::
sage: CBF(1) 1.000000000000000
TESTS::
sage: ComplexBallField(0) Traceback (most recent call last): ... ValueError: precision must be at least 2 sage: ComplexBallField(1) Traceback (most recent call last): ... ValueError: precision must be at least 2 """ Element = ComplexBall
@staticmethod def __classcall__(cls, long precision=53): r""" Normalize the arguments for caching.
TESTS::
sage: ComplexBallField(53) is ComplexBallField() True """
def __init__(self, long precision=53): r""" Initialize the complex ball field.
INPUT:
- ``precision`` -- an integer `\ge 2`.
EXAMPLES::
sage: CBF(1) 1.000000000000000
TESTS::
sage: CBF.base() Real ball field with 53 bits of precision sage: CBF.base_ring() Real ball field with 53 bits of precision
There are direct coercions from ZZ and QQ (for which arb provides construction functions)::
sage: CBF.coerce_map_from(ZZ) Coercion map: From: Integer Ring To: Complex ball field with 53 bits of precision sage: CBF.coerce_map_from(QQ) Coercion map: From: Rational Field To: Complex ball field with 53 bits of precision
Various other coercions are available through real ball fields or CLF::
sage: CBF.has_coerce_map_from(RLF) True sage: CBF.has_coerce_map_from(AA) True sage: CBF.has_coerce_map_from(QuadraticField(-1)) True sage: CBF.has_coerce_map_from(QQbar) True sage: CBF.has_coerce_map_from(CLF) True """
def _real_field(self): """ TESTS::
sage: CBF._real_field() is RBF True """
def _repr_(self): r""" String representation of ``self``.
EXAMPLES::
sage: ComplexBallField() Complex ball field with 53 bits of precision sage: ComplexBallField(106) Complex ball field with 106 bits of precision """
def construction(self): """ Return the construction of a complex ball field as the algebraic closure of the real ball field with the same precision.
EXAMPLES::
sage: functor, base = CBF.construction() sage: functor, base (AlgebraicClosureFunctor, Real ball field with 53 bits of precision) sage: functor(base) is CBF True """
def complex_field(self): """ Return the complex ball field with the same precision, i.e. ``self``
EXAMPLES::
sage: CBF.complex_field() is CBF True """
def ngens(self): r""" Return 1 as the only generator is the imaginary unit.
EXAMPLES::
sage: CBF.ngens() 1 """
def gen(self, i): r""" For i = 0, return the imaginary unit in this complex ball field.
EXAMPLES::
sage: CBF.0 1.000000000000000*I sage: CBF.gen(1) Traceback (most recent call last): ... ValueError: only one generator """ else:
def gens(self): r""" Return the tuple of generators of this complex ball field, i.e. ``(i,)``.
EXAMPLES::
sage: CBF.gens() (1.000000000000000*I,) sage: CBF.gens_dict() {'1.000000000000000*I': 1.000000000000000*I} """
def _coerce_map_from_(self, other): r""" Parents that canonically coerce into complex ball fields include:
- anything that coerces into the corresponding real ball field;
- real and complex ball fields with a larger precision;
- various exact or lazy parents representing subsets of the complex numbers, such as ``QQbar``, ``CLF``, and number fields equipped with complex embeddings.
TESTS::
sage: CBF.coerce_map_from(CBF) Identity endomorphism of Complex ball field with 53 bits of precision sage: CBF.coerce_map_from(ComplexBallField(100)) Coercion map: From: Complex ball field with 100 bits of precision To: Complex ball field with 53 bits of precision sage: CBF.has_coerce_map_from(ComplexBallField(42)) False sage: CBF.has_coerce_map_from(RealBallField(54)) True sage: CBF.has_coerce_map_from(RealBallField(52)) False sage: CBF.has_coerce_map_from(QuadraticField(-2)) True sage: CBF.has_coerce_map_from(QuadraticField(2, embedding=None)) False
Check that there are no coercions from interval or floating-point parents::
sage: CBF.has_coerce_map_from(RIF) False sage: CBF.has_coerce_map_from(CIF) False sage: CBF.has_coerce_map_from(RR) False sage: CBF.has_coerce_map_from(CC) False
Check that the map goes through the ``_acb_`` method::
sage: CBF.coerce_map_from(QuadraticField(-2, embedding=AA(-2).sqrt())) Conversion via _acb_ method map: ... sage: CBF.convert_map_from(QuadraticField(-2)) Conversion via _acb_ method map: ... """
def _element_constructor_(self, x=None, y=None): r""" Convert (x, y) to an element of this complex ball field, perhaps non-canonically.
INPUT:
- ``x``, ``y`` (optional) -- either a complex number, interval or ball, or two real ones (see examples below for more information on accepted number types).
EXAMPLES::
sage: CBF() 0 sage: CBF(1) # indirect doctest 1.000000000000000 sage: CBF(1, 1) 1.000000000000000 + 1.000000000000000*I sage: CBF(pi, sqrt(2)) [3.141592653589793 +/- 5.61e-16] + [1.414213562373095 +/- 4.10e-16]*I sage: CBF(I) 1.000000000000000*I sage: CBF(pi+I/3) [3.141592653589793 +/- 5.61e-16] + [0.3333333333333333 +/- 7.04e-17]*I sage: CBF(QQbar(i/7)) [0.1428571428571428 +/- 9.09e-17]*I sage: CBF(AA(sqrt(2))) [1.414213562373095 +/- 4.10e-16] sage: CBF(CIF(0, 1)) 1.000000000000000*I sage: CBF(RBF(1/3)) [0.3333333333333333 +/- 7.04e-17] sage: CBF(RBF(1/3), RBF(1/6)) [0.3333333333333333 +/- 7.04e-17] + [0.1666666666666667 +/- 7.04e-17]*I sage: CBF(1/3) [0.3333333333333333 +/- 7.04e-17] sage: CBF(1/3, 1/6) [0.3333333333333333 +/- 7.04e-17] + [0.1666666666666667 +/- 7.04e-17]*I sage: ComplexBallField(106)(1/3, 1/6) [0.33333333333333333333333333333333 +/- 6.94e-33] + [0.16666666666666666666666666666666 +/- 7.70e-33]*I sage: NF.<a> = QuadraticField(-2) sage: CBF(1/5 + a/2) [0.2000000000000000 +/- 4.45e-17] + [0.707106781186547 +/- 5.73e-16]*I sage: CBF(infinity, NaN) [+/- inf] + nan*I sage: CBF(x) Traceback (most recent call last): ... TypeError: unable to convert x to a ComplexBall
.. SEEALSO::
:meth:`sage.rings.real_arb.RealBallField._element_constructor_`
TESTS::
sage: CBF(1+I, 2) Traceback (most recent call last): ... TypeError: unable to convert I + 1 to a RealBall """ pass
pass pass pass else:
def _an_element_(self): r""" Construct an element.
EXAMPLES::
sage: CBF.an_element() # indirect doctest [0.3333333333333333 +/- 1.49e-17] - [0.1666666666666667 +/- 4.26e-17]*I """
def precision(self): """ Return the bit precision used for operations on elements of this field.
EXAMPLES::
sage: ComplexBallField().precision() 53 """
def is_exact(self): """ Complex ball fields are not exact.
EXAMPLES::
sage: ComplexBallField().is_exact() False """
def is_finite(self): """ Complex ball fields are infinite.
They already specify it via their category, but we currently need to re-implement this method due to the legacy implementation in :class:`sage.rings.ring.Ring`.
EXAMPLES::
sage: ComplexBallField().is_finite() False """
def characteristic(self): """ Complex ball fields have characteristic zero.
EXAMPLES::
sage: ComplexBallField().characteristic() 0 """
def some_elements(self): """ Complex ball fields contain elements with exact, inexact, infinite, or undefined real and imaginary parts.
EXAMPLES::
sage: CBF.some_elements() [1.000000000000000, -0.5000000000000000*I, 1.000000000000000 + [0.3333333333333333 +/- 1.49e-17]*I, [-0.3333333333333333 +/- 1.49e-17] + 0.2500000000000000*I, [-2.175556475109056e+181961467118333366510562 +/- 1.29e+181961467118333366510545], [+/- inf], [0.3333333333333333 +/- 1.49e-17] + [+/- inf]*I, [+/- inf] + [+/- inf]*I, nan, nan + nan*I, [+/- inf] + nan*I] """
def _sum_of_products(self, terms): r""" Compute a sum of product of complex balls without creating temporary Python objects
The input objects should be complex balls, but need not belong to this parent. The computation is performed at the precision of this parent.
EXAMPLES::
sage: Pol.<x> = ComplexBallField(1000)[] sage: pol = (x + 1/3)^100 sage: CBF._sum_of_products((c, c) for c in pol) [6.3308767660842e+23 +/- 4.59e+9]
TESTS::
sage: CBF._sum_of_products([]) 0 sage: CBF._sum_of_products([[]]) 1.000000000000000 sage: CBF._sum_of_products([["a"]]) Traceback (most recent call last): ... TypeError: Cannot convert str to sage.rings.complex_arb.ComplexBall """ cdef ComplexBall factor cdef acb_t tmp finally:
# Constants
def pi(self): r""" Return a ball enclosing `\pi`.
EXAMPLES::
sage: CBF.pi() [3.141592653589793 +/- 5.61e-16] sage: ComplexBallField(128).pi() [3.1415926535897932384626433832795028842 +/- 1.65e-38]
sage: CBF.pi().parent() Complex ball field with 53 bits of precision """
def integral(self, func, a, b, params=None, rel_tol=None, abs_tol=None, deg_limit=None, eval_limit=None, depth_limit=None, use_heap=None, verbose=None): r""" Compute a rigorous enclosure of the integral of ``func`` on the interval [``a``, ``b``].
INPUT:
- ``func`` -- a callable object accepting two parameters, a complex ball ``x`` and a boolean flag ``analytic``, and returning an element of this ball field (or some value that coerces into this ball field), such that:
- ``func(x, False)`` evaluates the integrand `f` on the ball ``x``. There are no restrictions on the behavior of `f` on ``x``; in particular, it can be discontinuous.
- ``func(x, True)`` evaluates `f(x)` if `f` is analytic on the whole ``x``, and returns some non-finite ball (e.g., ``self(NaN)``) otherwise.
(The ``analytic`` flag only needs to be checked for integrands that are non-analytic but bounded in some regions, typically complex functions with branch cuts, like `\sqrt{z}`. In particular, it can be ignored for meromorphic functions.)
- ``a``, ``b`` -- integration bounds. The bounds can be real or complex balls, or elements of any parent that coerces into this ball field, e.g. rational or algebraic numbers.
- ``rel_tol`` (optional, default `2^{-p}` where `p` is the precision of the ball field) -- relative accuracy goal
- ``abs_tol`` (optional, default `2^{-p}` where `p` is the precision of the ball field) -- absolute accuracy goal
Additionally, the following optional parameters can be used to control the integration algorithm. See the `Arb documentation <http://arblib.org/acb_calc.html>`_ for more information.
- ``deg_limit`` -- maximum quadrature degree for each subinterval
- ``eval_limit`` -- maximum number of function evaluations
- ``depth_limit`` -- maximum search depth for adaptive subdivision
- ``use_heap`` (boolean, default ``False``) -- if ``True``, use a priority queue instead of a stack to manage subintervals. This sometimes gives better results for integrals with slow convergence but may require more memory and increasing ``depth_limit``.
- ``verbose`` (integer, default 0) -- If set to 1, some information about the overall integration process is printed to standard output. If set to 2, information about each subinterval is printed.
EXAMPLES:
Some analytic integrands::
sage: CBF.integral(lambda x, _: x, 0, 1) [0.500000000000000 +/- 2.09e-16]
sage: CBF.integral(lambda x, _: x.gamma(), 1 - CBF(i), 1 + CBF(i)) [+/- 3.95e-15] + [1.5723926694981 +/- 4.53e-14]*I
sage: C = ComplexBallField(100) sage: C.integral(lambda x, _: x.cos() * x.sin(), 0, 1) [0.35403670913678559674939205737 +/- 8.89e-30]
sage: CBF.integral(lambda x, _: (x + x.exp()).sin(), 0, 8) [0.34740017266 +/- 6.36e-12]
sage: C = ComplexBallField(2000) sage: C.integral(lambda x, _: (x + x.exp()).sin(), 0, 8) # long time [0.34740017...55347713 +/- 6.72e-598]
Here the integration path crosses the branch cut of the square root::
sage: def my_sqrt(z, analytic): ....: if (analytic and not z.real() > 0 ....: and z.imag().contains_zero()): ....: return CBF(NaN) ....: else: ....: return z.sqrt() sage: CBF.integral(my_sqrt, -1 + CBF(i), -1 - CBF(i)) [+/- 1.14e-14] + [-0.4752076627926 +/- 5.18e-14]*I
Note, though, that proper handling of the ``analytic`` flag is required even when the path does not touch the branch cut::
sage: correct = CBF.integral(my_sqrt, 1, 2); correct [1.21895141649746 +/- 3.73e-15] sage: RBF(integral(sqrt(x), x, 1, 2)) [1.21895141649746 +/- 1.79e-15] sage: wrong = CBF.integral(lambda z, _: z.sqrt(), 1, 2) # WRONG! sage: correct - wrong [-5.640636259e-5 +/- 6.80e-15]
We can integrate the real absolute value function by defining a piecewise holomorphic extension::
sage: def real_abs(z, analytic): ....: if z.real().contains_zero(): ....: if analytic: ....: return z.parent()(NaN) ....: else: ....: return z.union(-z) ....: elif z.real() > 0: ....: return z ....: else: ....: return -z sage: CBF.integral(real_abs, -1, 1) [1.00000000000...] sage: CBF.integral(lambda z, analytic: real_abs(z.sin(), analytic), 0, 2*CBF.pi()) [4.00000000000...]
Here the integrand has a pole on or very close to the integration path, but there is no need to explicitly handle the ``analytic`` flag since the integrand is unbounded::
sage: CBF.integral(lambda x, _: 1/x, -1, 1) [+/- inf] + [+/- inf]*I sage: CBF.integral(lambda x, _: 1/x, 10^-1000, 1) [+/- inf] + [+/- inf]*I sage: CBF.integral(lambda x, _: 1/x, 10^-1000, 1, abs_tol=1e-10) [2302.5850930 +/- 1.26e-8]
Tolerances::
sage: CBF.integral(lambda x, _: x.exp(), -1020, -1010) [+/- 2.31e-438] sage: CBF.integral(lambda x, _: x.exp(), -1020, -1010, abs_tol=1e-450) [2.304377150950e-439 +/- 9.74e-452] sage: CBF.integral(lambda x, _: x.exp(), -1020, -1010, abs_tol=0) [2.304377150949e-439 +/- 7.53e-452] sage: CBF.integral(lambda x, _: x.exp(), -1020, -1010, rel_tol=1e-4, abs_tol=0) [2.30438e-439 +/- 3.90e-445]
sage: CBF.integral(lambda x, _: x*(1/x).sin(), 0, 1) [+/- 0.644] sage: CBF.integral(lambda x, _: x*(1/x).sin(), 0, 1, use_heap=True) [0.3785300 +/- 4.32e-8]
ALGORITHM:
Uses the `acb_calc <http://arblib.org/acb_calc.html>`_ module of the Arb library.
TESTS::
sage: CBF.integral(lambda x, _: x, 0, 10, rel_tol=1e-10, ....: abs_tol=1e-10, deg_limit=1, eval_limit=20, depth_limit=4, ....: use_heap=False) [50.00000000 +/- 2.20e-9]
sage: i = QuadraticField(-1).gen() sage: CBF.integral(lambda x, _: (1 + i*x).gamma(), -1, 1) [1.5723926694981 +/- 4.53e-14] + [+/- 3.95e-15]*I
sage: ComplexBallField(10000).integral(lambda x, _: x.sin(), 0, 1, rel_tol=1e-400) [0.459... +/- ...e-4...] sage: CBF.integral(lambda x, _: x.sin(), 0, 100, rel_tol=10) [+/- 7.61]
sage: ComplexBallField(10000).integral(lambda x, _: x.sin(), 0, 1, abs_tol=1e-400) [0.459697... +/- ...e-4...] sage: CBF.integral(lambda x, _: x.sin(), 0, 1, abs_tol=10) [+/- 0.980]
sage: ComplexBallField(100).integral(lambda x, _: sin(x), RBF(0), RBF(1)) [0.4596976941318602825990633926 +/- 5.09e-29] """ cdef acb_calc_integrate_opt_t arb_opts cdef long cgoal, expo cdef mag_t ctol cdef RealNumber tmp cdef ComplexBall ca, cb
else: else:
arb_opts.verbose = verbose
else:
else:
res.value, <acb_calc_func_t> acb_calc_func_callback, <void *> ctx, ca.value, cb.value, finally:
cdef inline bint _do_sig(long prec): """ Whether signal handlers should be installed for calls to arb. """
cdef inline long prec(ComplexBall ball):
cdef inline real_ball_field(ComplexBall ball):
cdef class ComplexBall(RingElement): """ Hold one ``acb_t`` of the `Arb library <http://fredrikj.net/arb/>`_
EXAMPLES::
sage: a = ComplexBallField()(1, 1) sage: a 1.000000000000000 + 1.000000000000000*I """ def __cinit__(self): """ Allocate memory for the encapsulated value.
EXAMPLES::
sage: ComplexBallField(2)(0) # indirect doctest 0 """
def __dealloc__(self): """ Deallocate memory of the encapsulated value.
EXAMPLES::
sage: a = ComplexBallField(2)(0) # indirect doctest sage: del a """
def __init__(self, parent, x=None, y=None): """ Initialize the :class:`ComplexBall`.
INPUT:
- ``parent`` -- a :class:`ComplexBallField`.
- ``x``, ``y`` (optional) -- either a complex number, interval or ball, or two real ones.
.. SEEALSO:: :meth:`ComplexBallField._element_constructor_`
TESTS::
sage: from sage.rings.complex_arb import ComplexBall sage: CBF53, CBF100 = ComplexBallField(53), ComplexBallField(100) sage: ComplexBall(CBF100) 0 sage: ComplexBall(CBF100, ComplexBall(CBF53, ComplexBall(CBF100, 1/3))) [0.333333333333333333333333333333 +/- 4.65e-31] sage: ComplexBall(CBF100, RBF(pi)) [3.141592653589793 +/- 5.61e-16]
sage: ComplexBall(CBF100, -3r) Traceback (most recent call last): ... TypeError: unsupported initializer sage: CBF100(-3r) -3.000000000000000000000000000000
sage: ComplexBall(CBF100, 10^100) 1.000000000000000000000000000000e+100 sage: ComplexBall(CBF100, CIF(1, 2)) 1.000000000000000000000000000000 + 2.000000000000000000000000000000*I sage: ComplexBall(CBF100, RBF(1/3), RBF(1)) [0.3333333333333333 +/- 7.04e-17] + 1.000000000000000000000000000000*I sage: NF.<a> = QuadraticField(-1, embedding=CC(0, -1)) sage: CBF(a) -1.000000000000000*I
sage: NF.<a> = QuadraticField(-1, embedding=None) sage: CBF(a) 1.000000000000000*I sage: CBF.coerce(a) Traceback (most recent call last): ... TypeError: no canonical coercion ...
sage: NF.<a> = QuadraticField(-2) sage: CBF(1/3 + a).real() [0.3333333333333333 +/- 7.04e-17]
sage: ComplexBall(CBF, 1, 1/2) 1.000000000000000 + 0.5000000000000000*I sage: ComplexBall(CBF, 1, 1) 1.000000000000000 + 1.000000000000000*I sage: ComplexBall(CBF, 1, 1/2) 1.000000000000000 + 0.5000000000000000*I sage: ComplexBall(CBF, 1/2, 1) 0.5000000000000000 + 1.000000000000000*I sage: ComplexBall(CBF, 1/2, 1/2) 0.5000000000000000 + 0.5000000000000000*I sage: ComplexBall(CBF, 1/2, 'a') Traceback (most recent call last): ... TypeError: unsupported initializer sage: ComplexBall(CBF, 'a') Traceback (most recent call last): ... TypeError: unsupported initializer
""" cdef fmpz_t tmpz cdef fmpq_t tmpq cdef long myprec
<ComplexIntervalFieldElement> x) else:
raise TypeError("unsupported initializer") else:
def __hash__(self): """ TESTS::
sage: hash(CBF(1/3)) == hash(RBF(1/3)) True sage: hash(CBF(1/3 + 2*i)) != hash(CBF(1/3 + i)) True """ else:
def _repr_(self): """ Return a string representation of ``self``.
OUTPUT:
A string.
EXAMPLES::
sage: CBF(1/3) [0.3333333333333333 +/- 7.04e-17] sage: CBF(0, 1/3) [0.3333333333333333 +/- 7.04e-17]*I sage: CBF(1/3, 1/6) [0.3333333333333333 +/- 7.04e-17] + [0.1666666666666667 +/- 7.04e-17]*I
TESTS::
sage: CBF(1-I/2) 1.000000000000000 - 0.5000000000000000*I """ else:
def _is_atomic(self): r""" Declare that complex balls print atomically in some cases.
TESTS::
sage: CBF(-1/3)._is_atomic() True
This method should in principle ensure that ``CBF['x']([1, -1/3])`` is printed as::
sage: CBF['x']([1, -1/3]) # todo - not tested [-0.3333333333333333 +/- 7.04e-17]*x + 1.000000000000000
However, this facility is not really used in Sage at this point, and we still get::
sage: CBF['x']([1, -1/3]) ([-0.3333333333333333 +/- 7.04e-17])*x + 1.000000000000000 """
# Conversions
cpdef ComplexIntervalFieldElement _complex_mpfi_(self, parent): """ Return :class:`ComplexIntervalFieldElement` of the same value.
EXAMPLES::
sage: CIF(CBF(1/3, 1/3)) # indirect doctest 0.3333333333333333? + 0.3333333333333333?*I """
def _integer_(self, _): """ Check that this ball contains a single integer and return that integer.
EXAMPLES::
sage: ZZ(CBF(-42, RBF(.1, rad=.2))) # indirect doctest -42 sage: ZZ(CBF(i)) Traceback (most recent call last): ... ValueError: 1.000000000000000*I does not contain a unique integer """ cdef sage.rings.integer.Integer res cdef fmpz_t tmp else: finally:
def _rational_(self): """ Check that this ball contains a single rational number and return that number.
EXAMPLES::
sage: QQ(CBF(12345/2^5)) 12345/32 sage: QQ(CBF(i)) Traceback (most recent call last): ... ValueError: 1.000000000000000*I does not contain a unique rational number """ else:
def _complex_mpfr_field_(self, parent): r""" Convert this complex ball to a complex number.
INPUT:
- ``parent`` - :class:`~sage.rings.complex_field.ComplexField_class`, target parent.
EXAMPLES::
sage: CC(CBF(1/3, 1/3)) 0.333333333333333 + 0.333333333333333*I sage: ComplexField(100)(CBF(1/3, 1/3)) 0.33333333333333331482961625625 + 0.33333333333333331482961625625*I
Check nan and inf::
sage: CC(CBF('nan', 1/3)) NaN + 0.333333333333333*I sage: CC(CBF('+inf')) +infinity """
def _real_mpfi_(self, parent): r""" Try to convert this complex ball to a real interval.
Fail if the imaginary part is not exactly zero.
INPUT:
- ``parent`` - :class:`~sage.rings.real_mpfi.RealIntervalField_class`, target parent.
EXAMPLES::
sage: RIF(CBF(RBF(1/3, rad=1e-5))) 0.3334? sage: RIF(CBF(RBF(1/3, rad=1e-5), 1e-10)) Traceback (most recent call last): ... ValueError: nonzero imaginary part """ else:
def _mpfr_(self, parent): r""" Try to convert this complex ball to a real number.
Fail if the imaginary part is not exactly zero.
INPUT:
- ``parent`` - :class:`~sage.rings.real_mpfr.RealField_class`, target parent.
EXAMPLES::
sage: RR(CBF(1/3)) 0.333333333333333 sage: RR(CBF(1, 1/3) - CBF(0, 1/3)) Traceback (most recent call last): ... ValueError: nonzero imaginary part """ else:
def __float__(self): """ Convert ``self`` to a ``float``.
EXAMPLES::
sage: float(CBF(1)) 1.0 sage: float(CBF(1/3)) 0.3333333333333333 sage: float(CBF(1,1)) Traceback (most recent call last): ... TypeError: can't convert complex ball to float """
def __complex__(self): """ Convert ``self`` to a ``complex``.
EXAMPLES::
sage: complex(CBF(1)) (1+0j) sage: complex(CBF(1,1)) (1+1j)
Check nan and inf::
sage: complex(CBF(0, 'nan')) nanj sage: complex(CBF('+inf', '-inf')) (inf-infj) """ arf_get_d(arb_midref(acb_realref(self.value)), ARF_RND_NEAR),
def _real_double_(self, parent): r""" Convert this complex ball to a real double.
EXAMPLES::
sage: RDF(CBF(3)) 3.0 sage: RDF(CBF(3/7)) 0.42857142857142855 sage: RDF(CBF(1 + I)) Traceback (most recent call last): ... TypeError: can't convert complex ball to float
sage: RDF(CBF(3)).parent() Real Double Field
Check nan and inf::
sage: RDF(CBF('nan')) NaN sage: RDF(CBF('nan')).parent() Real Double Field sage: RDF(CBF('+inf')) +infinity sage: RDF(CBF('+inf')).parent() Real Double Field
TESTS:
Check that conversions go through this method::
sage: RDF.convert_map_from(CBF) Conversion via _real_double_ method map: From: Complex ball field with 53 bits of precision To: Real Double Field """ cdef RealDoubleElement x
def _complex_double_(self, parent): r""" Convert this complex ball to a complex double.
EXAMPLES::
sage: CDF(CBF(1+I)) 1.0 + 1.0*I sage: CDF(CBF(3/7, -21/13)) 0.42857142857142855 - 1.6153846153846152*I
Check nan and inf::
sage: CDF(CBF('nan', 'nan')) NaN + NaN*I sage: CDF(CBF('+inf', 'nan')) +infinity + NaN*I sage: CDF(CBF('+inf', '-inf')) +infinity - +infinity*I
TESTS:
The conversion map should go through this method. However, it is manually called in the constructor of complex double elements::
sage: CDF.convert_map_from(CBF) # bug Conversion map: From: Complex ball field with 53 bits of precision To: Complex Double Field """ arf_get_d(arb_midref(acb_realref(self.value)), ARF_RND_NEAR), arf_get_d(arb_midref(acb_imagref(self.value)), ARF_RND_NEAR))
# Real and imaginary part, midpoint, radius
cpdef RealBall real(self): """ Return the real part of this ball.
OUTPUT:
A :class:`~sage.rings.real_arb.RealBall`.
EXAMPLES::
sage: a = CBF(1/3, 1/5) sage: a.real() [0.3333333333333333 +/- 7.04e-17] sage: a.real().parent() Real ball field with 53 bits of precision """
cpdef RealBall imag(self): """ Return the imaginary part of this ball.
OUTPUT:
A :class:`~sage.rings.real_arb.RealBall`.
EXAMPLES::
sage: a = CBF(1/3, 1/5) sage: a.imag() [0.2000000000000000 +/- 4.45e-17] sage: a.imag().parent() Real ball field with 53 bits of precision """
def __abs__(self): """ Return the absolute value of this complex ball.
EXAMPLES::
sage: CBF(1 + i).abs() # indirect doctest [1.414213562373095 +/- 2.99e-16] sage: abs(CBF(i)) 1.000000000000000
sage: CBF(1 + i).abs().parent() Real ball field with 53 bits of precision """
def below_abs(self, test_zero=False): """ Return a lower bound for the absolute value of this complex ball.
INPUT:
- ``test_zero`` (boolean, default ``False``) -- if ``True``, make sure that the returned lower bound is positive, raising an error if the ball contains zero.
OUTPUT:
A ball with zero radius
EXAMPLES::
sage: b = ComplexBallField(8)(1+i).below_abs() sage: b [1.4 +/- 0.0141] sage: b.is_exact() True sage: QQ(b)*128 181 sage: (CBF(1/3) - 1/3).below_abs() 0 sage: (CBF(1/3) - 1/3).below_abs(test_zero=True) Traceback (most recent call last): ... ValueError: ball contains zero
.. SEEALSO:: :meth:`above_abs` """
def above_abs(self): """ Return an upper bound for the absolute value of this complex ball.
OUTPUT:
A ball with zero radius
EXAMPLES::
sage: b = ComplexBallField(8)(1+i).above_abs() sage: b [1.4 +/- 0.0219] sage: b.is_exact() True sage: QQ(b)*128 182
.. SEEALSO:: :meth:`below_abs` """
def arg(self): """ Return the argument of this complex ball.
EXAMPLES::
sage: CBF(1 + i).arg() [0.785398163397448 +/- 3.91e-16] sage: CBF(-1).arg() [3.141592653589793 +/- 5.61e-16] sage: CBF(-1).arg().parent() Real ball field with 53 bits of precision """
def mid(self): """ Return the midpoint of this ball.
OUTPUT:
:class:`~sage.rings.complex_number.ComplexNumber`, floating-point complex number formed by the centers of the real and imaginary parts of this ball.
EXAMPLES::
sage: CBF(1/3, 1).mid() 0.333333333333333 + 1.00000000000000*I sage: CBF(1/3, 1).mid().parent() Complex Field with 53 bits of precision sage: CBF('inf', 'nan').mid() +infinity + NaN*I sage: CBF('nan', 'inf').mid() NaN + +infinity*I sage: CBF('nan').mid() NaN sage: CBF('inf').mid() +infinity sage: CBF(0, 'inf').mid() +infinity*I
.. SEEALSO:: :meth:`squash` """
def squash(self): """ Return an exact ball with the same midpoint as this ball.
OUTPUT:
A :class:`ComplexBall`.
EXAMPLES::
sage: mid = CBF(1/3, 1/10).squash() sage: mid [0.3333333333333333 +/- 1.49e-17] + [0.09999999999999999 +/- 1.68e-18]*I sage: mid.parent() Complex ball field with 53 bits of precision sage: mid.is_exact() True
.. SEEALSO:: :meth:`mid` """
def rad(self): """ Return an upper bound for the error radius of this ball.
OUTPUT:
A :class:`~sage.rings.real_mpfr.RealNumber` of the same precision as the radii of real balls.
.. WARNING::
Unlike a :class:`~sage.rings.real_arb.RealBall`, a :class:`ComplexBall` is *not* defined by its midpoint and radius. (Instances of :class:`ComplexBall` are actually rectangles, not balls.)
EXAMPLES::
sage: CBF(1 + i).rad() 0.00000000 sage: CBF(i/3).rad() 1.1102230e-16 sage: CBF(i/3).rad().parent() Real Field with 30 bits of precision
.. SEEALSO:: :meth:`diameter`, :meth:`mid`
TESTS::
sage: (CBF(0, 1/3) << (2^64)).rad() Traceback (most recent call last): ... RuntimeError: unable to convert the radius to MPFR (exponent out of range?) """ # Should we return a real number with rounding towards +∞ (or away from # zero if/when implemented)? cdef arf_t tmp sig_error()
# Should we implement rad_as_ball? If we do, should it return an enclosure # of the radius (which radius?), or an upper bound?
def diameter(self): r""" Return the diameter of this ball.
EXAMPLES::
sage: CBF(1 + i).diameter() 0.00000000 sage: CBF(i/3).diameter() 2.2204460e-16 sage: CBF(i/3).diameter().parent() Real Field with 30 bits of precision sage: CBF(CIF(RIF(1.02, 1.04), RIF(2.1, 2.2))).diameter() 0.20000000
.. SEEALSO:: :meth:`rad`, :meth:`mid` """
def union(self, other): r""" Return a ball containing the convex hull of ``self`` and ``other``.
EXAMPLES::
sage: b = CBF(1 + i).union(0) sage: b.real().endpoints() (-9.31322574615479e-10, 1.00000000093133) """
# Precision and accuracy
def round(self): """ Return a copy of this ball rounded to the precision of the parent.
EXAMPLES:
It is possible to create balls whose midpoint is more precise that their parent's nominal precision (see :mod:`~sage.rings.real_arb` for more information)::
sage: b = CBF(exp(I*pi/3).n(100)) sage: b.mid() 0.50000000000000000000000000000 + 0.86602540378443864676372317075*I
The ``round()`` method rounds such a ball to its parent's precision::
sage: b.round().mid() 0.500000000000000 + 0.866025403784439*I
.. SEEALSO:: :meth:`trim` """
def accuracy(self): """ Return the effective relative accuracy of this ball measured in bits.
This is computed as if calling :meth:`~sage.rings.real_arb.RealBall.accuracy()` on the real ball whose midpoint is the larger out of the real and imaginary midpoints of this complex ball, and whose radius is the larger out of the real and imaginary radii of this complex ball.
EXAMPLES::
sage: CBF(exp(I*pi/3)).accuracy() 52 sage: CBF(I/2).accuracy() == CBF.base().maximal_accuracy() True sage: CBF('nan', 'inf').accuracy() == -CBF.base().maximal_accuracy() True
.. SEEALSO::
:meth:`~sage.rings.real_arb.RealBallField.maximal_accuracy` """
def trim(self): """ Return a trimmed copy of this ball.
Return a copy of this ball with both the real and imaginary parts trimmed (see :meth:`~sage.rings.real_arb.RealBall.trim()`).
EXAMPLES::
sage: b = CBF(1/3, RBF(1/3, rad=.01)) sage: b.mid() 0.333333333333333 + 0.333333333333333*I sage: b.trim().mid() 0.333333333333333 + 0.333333015441895*I
.. SEEALSO:: :meth:`round` """
def add_error(self, ampl): """ Increase the radii of the real and imaginary parts by (an upper bound on) ``ampl``.
If ``ampl`` is negative, the radii remain unchanged.
INPUT:
- ``ampl`` - A **real** ball (or an object that can be coerced to a real ball).
OUTPUT:
A new complex ball.
EXAMPLES::
sage: CBF(1+i).add_error(10^-16) [1.000000000000000 +/- 1.01e-16] + [1.000000000000000 +/- 1.01e-16]*I """
# Comparisons and predicates
def is_NaN(self): """ Return ``True`` iff either the real or the imaginary part is not-a-number.
EXAMPLES::
sage: CBF(NaN).is_NaN() True sage: CBF(-5).gamma().is_NaN() True sage: CBF(oo).is_NaN() False sage: CBF(42+I).is_NaN() False """
def is_zero(self): """ Return ``True`` iff the midpoint and radius of this ball are both zero.
EXAMPLES::
sage: CBF(0).is_zero() True sage: CBF(RIF(-0.5, 0.5)).is_zero() False
.. SEEALSO:: :meth:`is_nonzero` """
def is_nonzero(self): """ Return ``True`` iff zero is not contained in the interval represented by this ball.
.. NOTE::
This method is not the negation of :meth:`is_zero`: it only returns ``True`` if zero is known not to be contained in the ball.
Use ``bool(b)`` (or, equivalently, ``not b.is_zero()``) to check if a ball ``b`` **may** represent a nonzero number (for instance, to determine the “degree” of a polynomial with ball coefficients).
EXAMPLES::
sage: CBF(pi, 1/3).is_nonzero() True sage: CBF(RIF(-0.5, 0.5), 1/3).is_nonzero() True sage: CBF(1/3, RIF(-0.5, 0.5)).is_nonzero() True sage: CBF(RIF(-0.5, 0.5), RIF(-0.5, 0.5)).is_nonzero() False
.. SEEALSO:: :meth:`is_zero` """
def __nonzero__(self): """ Return ``True`` iff this complex ball is not the zero ball, i.e. if the midpoint and radius of its real and imaginary parts are not all zero.
This is the preferred way, for instance, to determine the “degree” of a polynomial with ball coefficients.
.. WARNING::
A “nonzero” ball in the sense of this method may represent the value zero. Use :meth:`is_nonzero` to check that a complex number represented by a ``ComplexBall`` object is known to be nonzero.
EXAMPLES::
sage: bool(CBF(0)) # indirect doctest False sage: bool(CBF(i)) True sage: bool(CBF(RIF(-0.5, 0.5))) True """
def is_exact(self): """ Return ``True`` iff the radius of this ball is zero.
EXAMPLES::
sage: CBF(1).is_exact() True sage: CBF(1/3, 1/3).is_exact() False """
def is_real(self): """ Return ``True`` iff the imaginary part of this ball is exactly zero.
EXAMPLES::
sage: CBF(1/3, 0).is_real() True sage: (CBF(i/3) - CBF(1, 1/3)).is_real() False sage: CBF('inf').is_real() True """
cpdef _richcmp_(left, right, int op): """ Compare ``left`` and ``right``.
For more information, see :mod:`sage.rings.complex_arb`.
EXAMPLES::
sage: a = CBF(1) sage: b = CBF(1) sage: a is b False sage: a == b True sage: a = CBF(1/3) sage: a.is_exact() False sage: b = CBF(1/3) sage: b.is_exact() False sage: a == b False sage: a = CBF(1, 2) sage: b = CBF(1, 2) sage: a is b False sage: a == b True
TESTS:
Balls whose intersection consists of one point::
sage: a = CBF(RIF(1, 2), RIF(1, 2)) sage: b = CBF(RIF(2, 4), RIF(2, 4)) sage: a < b Traceback (most recent call last): ... TypeError: No order is defined for ComplexBalls. sage: a > b Traceback (most recent call last): ... TypeError: No order is defined for ComplexBalls. sage: a <= b Traceback (most recent call last): ... TypeError: No order is defined for ComplexBalls. sage: a >= b Traceback (most recent call last): ... TypeError: No order is defined for ComplexBalls. sage: a == b False sage: a != b False
Balls with non-trivial intersection::
sage: a = CBF(RIF(1, 4), RIF(1, 4)) sage: a = CBF(RIF(2, 5), RIF(2, 5)) sage: a == b False sage: a != b False
One ball contained in another::
sage: a = CBF(RIF(1, 4), RIF(1, 4)) sage: b = CBF(RIF(2, 3), RIF(2, 3)) sage: a == b False sage: a != b False
Disjoint balls::
sage: a = CBF(1/3, 1/3) sage: b = CBF(1/5, 1/5) sage: a == b False sage: a != b True
Exact elements::
sage: a = CBF(2, 2) sage: b = CBF(2, 2) sage: a.is_exact() True sage: b.is_exact() True sage: a == b True sage: a != b False """ cdef ComplexBall lt, rt cdef acb_t difference
def identical(self, ComplexBall other): """ Return whether ``self`` and ``other`` represent the same ball.
INPUT:
- ``other`` -- a :class:`ComplexBall`.
OUTPUT:
Return True iff ``self`` and ``other`` are equal as sets, i.e. if their real and imaginary parts each have the same midpoint and radius.
Note that this is not the same thing as testing whether both ``self`` and ``other`` certainly represent the complex real number, unless either ``self`` or ``other`` is exact (and neither contains NaN). To test whether both operands might represent the same mathematical quantity, use :meth:`overlaps` or ``in``, depending on the circumstance.
EXAMPLES::
sage: CBF(1, 1/3).identical(1 + CBF(0, 1)/3) True sage: CBF(1, 1).identical(1 + CBF(0, 1/3)*3) False """
def overlaps(self, ComplexBall other): """ Return True iff ``self`` and ``other`` have some point in common.
INPUT:
- ``other`` -- a :class:`ComplexBall`.
EXAMPLES::
sage: CBF(1, 1).overlaps(1 + CBF(0, 1/3)*3) True sage: CBF(1, 1).overlaps(CBF(1, 'nan')) True sage: CBF(1, 1).overlaps(CBF(0, 'nan')) False """
def contains_exact(self, other): """ Return ``True`` *iff* ``other`` is contained in ``self``.
Use ``other in self`` for a test that works for a wider range of inputs but may return false negatives.
INPUT:
- ``other`` -- :class:`ComplexBall`, :class:`~sage.rings.integer.Integer`, or :class:`~sage.rings.rational.Rational`
EXAMPLES::
sage: CBF(RealBallField(100)(1/3), 0).contains_exact(1/3) True sage: CBF(1).contains_exact(1) True sage: CBF(1).contains_exact(CBF(1)) True
sage: CBF(sqrt(2)).contains_exact(sqrt(2)) Traceback (most recent call last): ... TypeError: unsupported type: <type 'sage.symbolic.expression.Expression'> """ cdef fmpz_t tmpz cdef fmpq_t tmpq else: finally:
def __contains__(self, other): """ Return True if ``other`` can be verified to be contained in ``self``.
Depending on the type of ``other``, the test may use interval arithmetic with a precision determined by the parent of ``self`` and may return false negatives.
EXAMPLES::
sage: 1/3*i in CBF(0, 1/3) True
A false negative::
sage: RLF(1/3) in CBF(RealBallField(100)(1/3), 0) False
.. SEEALSO:: :meth:`contains_exact` """
def contains_zero(self): """ Return ``True`` iff this ball contains zero.
EXAMPLES::
sage: CBF(0).contains_zero() True sage: CBF(RIF(-1,1)).contains_zero() True sage: CBF(i).contains_zero() False """
def contains_integer(self): """ Return ``True`` iff this ball contains any integer.
EXAMPLES::
sage: CBF(3, RBF(0.1)).contains_integer() False sage: CBF(3, RBF(0.1,0.1)).contains_integer() True """
# Arithmetic
def __neg__(self): """ Return the opposite of this ball.
EXAMPLES::
sage: -CBF(1/3 + I) [-0.3333333333333333 +/- 7.04e-17] - 1.000000000000000*I """
def conjugate(self): """ Return the complex conjugate of this ball.
EXAMPLES::
sage: CBF(-2 + I/3).conjugate() -2.000000000000000 + [-0.3333333333333333 +/- 7.04e-17]*I """
cpdef _add_(self, other): """ Return the sum of two balls, rounded to the ambient field's precision.
The resulting ball is guaranteed to contain the sums of any two points of the respective input balls.
EXAMPLES::
sage: CBF(1) + CBF(I) 1.000000000000000 + 1.000000000000000*I """
cpdef _sub_(self, other): """ Return the difference of two balls, rounded to the ambient field's precision.
The resulting ball is guaranteed to contain the differences of any two points of the respective input balls.
EXAMPLES::
sage: CBF(1) - CBF(I) 1.000000000000000 - 1.000000000000000*I """
def __invert__(self): """ Return the inverse of this ball.
The result is guaranteed to contain the inverse of any point of the input ball.
EXAMPLES::
sage: ~CBF(i/3) [-3.00000000000000 +/- 9.44e-16]*I sage: ~CBF(0) [+/- inf] sage: ~CBF(RIF(10,11)) [0.1 +/- 9.53e-3] """
cpdef _mul_(self, other): """ Return the product of two balls, rounded to the ambient field's precision.
The resulting ball is guaranteed to contain the products of any two points of the respective input balls.
EXAMPLES::
sage: CBF(-2, 1)*CBF(1, 1/3) [-2.333333333333333 +/- 5.37e-16] + [0.333333333333333 +/- 4.82e-16]*I """
def __lshift__(val, shift): r""" If ``val`` is a ``ComplexBall`` and ``shift`` is an integer, return the ball obtained by shifting the center and radius of ``val`` to the left by ``shift`` bits.
INPUT:
- ``shift`` -- integer, may be negative.
EXAMPLES::
sage: CBF(i/3) << 2 [1.333333333333333 +/- 4.82e-16]*I sage: CBF(i) << -2 0.2500000000000000*I
TESTS::
sage: CBF(i) << (2^65) [3.636549880934858e+11106046577046714264 +/- 1.91e+11106046577046714248]*I sage: 'a' << CBF(1/3) Traceback (most recent call last): ... TypeError: unsupported operand type(s) for <<: 'str' and 'ComplexBall' sage: CBF(1) << 1/2 Traceback (most recent call last): ... TypeError: shift should be an integer """ cdef fmpz_t tmpz # the ComplexBall might be shift, not val acb_mul_2exp_si(res.value, self.value, PyInt_AS_LONG(shift)) else:
def __rshift__(val, shift): r""" If ``val`` is a ``ComplexBall`` and ``shift`` is an integer, return the ball obtained by shifting the center and radius of ``val`` to the right by ``shift`` bits.
INPUT:
- ``shift`` -- integer, may be negative.
EXAMPLES::
sage: CBF(1+I) >> 2 0.2500000000000000 + 0.2500000000000000*I
TESTS::
sage: 'a' >> CBF(1/3) Traceback (most recent call last): ... TypeError: unsupported operand type(s) for >>: 'str' and 'ComplexBall' """ # the ComplexBall might be shift, not val else:
cpdef _div_(self, other): """ Return the quotient of two balls, rounded to the ambient field's precision.
The resulting ball is guaranteed to contain the quotients of any two points of the respective input balls.
EXAMPLES::
sage: CBF(-2, 1)/CBF(1, 1/3) [-1.500000000000000 +/- 8.83e-16] + [1.500000000000000 +/- 5.64e-16]*I sage: CBF(2+I)/CBF(0) [+/- inf] + [+/- inf]*I sage: CBF(1)/CBF(0) [+/- inf] sage: CBF(1)/CBF(RBF(0, 1.)) nan """
def __pow__(base, expo, _): """ EXAMPLES::
sage: CBF(-1)**(1/2) 1.000000000000000*I sage: CBF(e)**CBF(i*pi) [-1.0000000000000 +/- 1.98e-15] + [+/- 2.32e-15]*I sage: CBF(0, 1)**AA(2)**(1/2) [-0.60569986707881 +/- 4.36e-15] + [0.79569320156748 +/- 2.53e-15]*I
sage: CBF(i)**RBF(2**1000) [+/- 2.51] + [+/- 2.87]*I sage: CBF(i)**(2**1000) 1.000000000000000
sage: CBF(0)^(1/3) 0 sage: CBF(0)^(-1) [+/- inf] sage: CBF(0)^(-2) [+/- inf] + [+/- inf]*I
TESTS::
sage: (CBF(e)**CBF(i))**RBF(pi) [-1.0000000000000 +/- 5.48e-15] + [+/- 4.14e-15]*I sage: CBF(2*i)**10r -1024.000000000000 sage: CBF(1,1) ^ -1r 0.5000000000000000 - 0.5000000000000000*I
""" cdef fmpz_t tmpz return sage.structure.element.bin_op(base, expo, operator.pow) else:
def sqrt(self): """ Return the square root of this ball.
If either the real or imaginary part is exactly zero, only a single real square root is needed.
EXAMPLES::
sage: CBF(-2).sqrt() [1.414213562373095 +/- 2.99e-16]*I """
def rsqrt(self): """ Return the reciprocal square root of ``self``.
If either the real or imaginary part is exactly zero, only a single real reciprocal square root is needed.
EXAMPLES::
sage: CBF(-2).rsqrt() [-0.707106781186547 +/- 5.73e-16]*I sage: CBF(0, 1/2).rsqrt() 1.000000000000000 - 1.000000000000000*I sage: CBF(0).rsqrt() nan """
def cube(self): """ Return the cube of this ball.
The result is computed efficiently using two real squarings, two real multiplications, and scalar operations.
EXAMPLES::
sage: CBF(1, 1).cube() -2.000000000000000 + 2.000000000000000*I """
def rising_factorial(self, n): """ Return the ``n``-th rising factorial of this ball.
The `n`-th rising factorial of `x` is equal to `x (x+1) \cdots (x+n-1)`.
For complex `n`, it is a quotient of gamma functions.
EXAMPLES::
sage: CBF(1).rising_factorial(5) 120.0000000000000 sage: CBF(1/3, 1/2).rising_factorial(300) [-3.87949484514e+612 +/- 5.23e+600] + [-3.52042209763e+612 +/- 5.55e+600]*I
sage: CBF(1).rising_factorial(-1) nan sage: CBF(1).rising_factorial(2**64) [+/- 2.30e+347382171305201370464] sage: ComplexBallField(128)(1).rising_factorial(2**64) [2.343691126796861348e+347382171305201285713 +/- 4.71e+347382171305201285694] sage: CBF(1/2).rising_factorial(CBF(2,3)) [-0.123060451458124 +/- 4.43e-16] + [0.040641263167655 +/- 3.72e-16]*I
"""
# Elementary functions
def log(self, base=None): """ General logarithm (principal branch).
INPUT:
- ``base`` (optional, complex ball or number) -- if ``None``, return the principal branch of the natural logarithm ``ln(self)``, otherwise, return the general logarithm ``ln(self)/ln(base)``
EXAMPLES::
sage: CBF(2*i).log() [0.6931471805599453 +/- 4.16e-17] + [1.570796326794897 +/- 6.65e-16]*I sage: CBF(-1).log() [3.141592653589793 +/- 5.61e-16]*I
sage: CBF(2*i).log(2) [1.000000000000000 +/- 8.01e-17] + [2.26618007091360 +/- 4.23e-15]*I sage: CBF(2*i).log(CBF(i)) [1.000000000000000 +/- 2.83e-16] + [-0.441271200305303 +/- 2.82e-16]*I
sage: CBF('inf').log() nan + nan*I sage: CBF(2).log(0) nan + nan*I """ cdef ComplexBall cst
def log1p(self): """ Return ``log(1 + self)``, computed accurately when ``self`` is close to zero.
EXAMPLES::
sage: eps = RBF(1e-50) sage: CBF(1+eps, eps).log() [+/- 2.23e-16] + [1.000000000000000e-50 +/- 2.30e-66]*I sage: CBF(eps, eps).log1p() [1.000000000000000e-50 +/- 7.63e-68] + [1.00000000000000e-50 +/- 2.30e-66]*I """
def exp(self): """ Return the exponential of this ball.
.. SEEALSO:: :meth:`exppii`
EXAMPLES::
sage: CBF(i*pi).exp() [-1.00000000000000 +/- 6.67e-16] + [+/- 5.68e-16]*I """
def exppii(self): """ Return ``exp(pi*i*self)``.
EXAMPLES::
sage: CBF(1/2).exppii() 1.000000000000000*I sage: CBF(0, -1/pi).exppii() [2.71828182845904 +/- 6.20e-15] """
def sin(self): """ Return the sine of this ball.
EXAMPLES::
sage: CBF(i*pi).sin() [11.5487393572577 +/- 5.34e-14]*I """
def cos(self): """ Return the cosine of this ball.
EXAMPLES::
sage: CBF(i*pi).cos() [11.59195327552152 +/- 8.38e-15] """
def tan(self): """ Return the tangent of this ball.
EXAMPLES::
sage: CBF(pi/2, 1/10).tan() [+/- 2.87e-14] + [10.0333111322540 +/- 2.36e-14]*I sage: CBF(pi/2).tan() [+/- inf] """
def cot(self): """ Return the cotangent of this ball.
EXAMPLES::
sage: CBF(pi, 1/10).cot() [+/- 5.74e-14] + [-10.0333111322540 +/- 2.81e-14]*I sage: CBF(pi).cot() [+/- inf] """
def arcsin(self): """ Return the arcsine of this ball.
EXAMPLES::
sage: CBF(1+i).arcsin() [0.66623943249252 +/- 5.40e-15] + [1.06127506190504 +/- 5.04e-15]*I """
def arccos(self): """ Return the arccosine of this ball.
EXAMPLES::
sage: CBF(1+i).arccos() [0.90455689430238 +/- 2.18e-15] + [-1.06127506190504 +/- 5.04e-15]*I """
def arctan(self): """ Return the arctangent of this ball.
EXAMPLES::
sage: CBF(1+i).arctan() [1.017221967897851 +/- 4.93e-16] + [0.4023594781085251 +/- 8.52e-17]*I sage: CBF(i).arctan() nan + nan*I """
def arcsinh(self): """ Return the hyperbolic arcsine of this ball.
EXAMPLES::
sage: CBF(1+i).arcsinh() [1.06127506190504 +/- 5.04e-15] + [0.66623943249252 +/- 5.40e-15]*I """
def arccosh(self): """ Return the hyperbolic arccosine of this ball.
EXAMPLES::
sage: CBF(1+i).arccosh() [1.061275061905035 +/- 8.44e-16] + [0.904556894302381 +/- 8.22e-16]*I """
def arctanh(self): """ Return the hyperbolic arctangent of this ball.
EXAMPLES::
sage: CBF(1+i).arctanh() [0.4023594781085251 +/- 8.52e-17] + [1.017221967897851 +/- 4.93e-16]*I """
# Special functions
def gamma(self, z=None): """ Return the image of this ball by the Euler Gamma function (if ``z = None``) or the incomplete Gamma function (otherwise).
EXAMPLES::
sage: CBF(1, 1).gamma() [0.498015668118356 +/- 9.16e-16] + [-0.154949828301811 +/- 7.08e-16]*I sage: CBF(-1).gamma() nan sage: CBF(1, 1).gamma(0) [0.498015668118356 +/- 9.16e-16] + [-0.154949828301811 +/- 7.08e-16]*I sage: CBF(1, 1).gamma(100) [-3.6143867454139e-45 +/- 6.88e-59] + [-3.7022961377791e-44 +/- 4.41e-58]*I sage: CBF(1, 1).gamma(CLF(i)) [0.32886684193500 +/- 5.04e-15] + [-0.18974945045621 +/- 1.26e-15]*I """ cdef ComplexBall my_z else:
def log_gamma(self): r""" Return the image of this ball by the logarithmic Gamma function.
The branch cut of the logarithmic gamma function is placed on the negative half-axis, which means that ``log_gamma(z) + log z = log_gamma(z+1)`` holds for all `z`, whereas ``log_gamma(z) != log(gamma(z))`` in general.
EXAMPLES::
sage: CBF(1000, 1000).log_gamma() [5466.22252162990 +/- 3.05e-12] + [7039.33429191119 +/- 3.81e-12]*I sage: CBF(-1/2).log_gamma() [1.265512123484645 +/- 8.82e-16] + [-3.141592653589793 +/- 5.68e-16]*I sage: CBF(-1).log_gamma() nan + [-3.141592653589793 +/- 5.68e-16]*I """
def rgamma(self): """ Compute the reciprocal gamma function with argument ``self``.
EXAMPLES::
sage: CBF(6).rgamma() [0.00833333333333333 +/- 4.96e-18] sage: CBF(-1).rgamma() 0 """
def psi(self, n=None): """ Compute the digamma function with argument ``self``.
If ``n`` is provided, compute the polygamma function of order ``n`` and argument ``self``.
EXAMPLES::
sage: CBF(1, 1).psi() [0.0946503206224770 +/- 7.74e-17] + [1.076674047468581 +/- 2.58e-16]*I sage: CBF(-1).psi() nan sage: CBF(1,1).psi(10) [56514.8269344249 +/- 4.70e-11] + [56215.1218005823 +/- 5.70e-11]*I
""" cdef ComplexBall my_n else:
def zeta(self, a=None): """ Return the image of this ball by the Hurwitz zeta function.
For ``a = None``, this computes the Riemann zeta function.
EXAMPLES::
sage: CBF(1, 1).zeta() [0.5821580597520036 +/- 5.26e-17] + [-0.9268485643308071 +/- 2.80e-17]*I sage: CBF(1, 1).zeta(1) [0.5821580597520036 +/- 5.26e-17] + [-0.9268485643308071 +/- 2.80e-17]*I sage: CBF(1, 1).zeta(1/2) [1.497919876084167 +/- 2.91e-16] + [0.2448655353684164 +/- 4.22e-17]*I sage: CBF(1, 1).zeta(CBF(1, 1)) [-0.3593983122202835 +/- 3.01e-17] + [-2.875283329756940 +/- 4.52e-16]*I sage: CBF(1, 1).zeta(-1) nan + nan*I """ cdef ComplexBall a_ball else:
def polylog(self, s): """ Return the polylogarithm `\operatorname{Li}_s(\mathrm{self})`.
EXAMPLES::
sage: CBF(2).polylog(1) [+/- 4.65e-15] + [-3.14159265358979 +/- 8.15e-15]*I sage: CBF(1, 1).polylog(CBF(1, 1)) [0.3708160030469 +/- 2.38e-14] + [2.7238016577979 +/- 4.20e-14]*I
TESTS::
sage: CBF(2).polylog(1r) [+/- 4.65e-15] + [-3.14159265358979 +/- 8.15e-15]*I """ cdef ComplexBall s_as_ball cdef sage.rings.integer.Integer s_as_Integer pass
def barnes_g(self): """ Return the Barnes G-function of ``self``.
EXAMPLES::
sage: CBF(-4).barnes_g() 0 sage: CBF(8).barnes_g() 24883200.00000000 sage: CBF(500,10).barnes_g() [4.54078781e+254873 +/- 5.43e+254864] + [8.65835455e+254873 +/- 7.28e+254864]*I """
def log_barnes_g(self): """ Return the logarithmic Barnes G-function of ``self``.
EXAMPLES::
sage: CBF(10^100).log_barnes_g() [1.14379254649702e+202 +/- 4.09e+187] sage: CBF(0,1000).log_barnes_g() [-2702305.04929258 +/- 2.60e-9] + [-790386.325561423 +/- 9.72e-10]*I """
def agm1(self): """ Return the arithmetic-geometric mean of 1 and ``self``.
The arithmetic-geometric mean is defined such that the function is continuous in the complex plane except for a branch cut along the negative half axis (where it is continuous from above). This corresponds to always choosing an "optimal" branch for the square root in the arithmetic-geometric mean iteration.
EXAMPLES::
sage: CBF(0, -1).agm1() [0.599070117367796 +/- 3.9...e-16] + [-0.599070117367796 +/- 5.5...e-16]*I """
def hypergeometric(self, a, b, bint regularized=False): r""" Return the generalized hypergeometric function of ``self``.
INPUT:
- ``a`` -- upper parameters, list of complex numbers that coerce into this ball's parent;
- ``b`` -- lower parameters, list of complex numbers that coerce into this ball's parent.
- ``regularized`` -- if True, the regularized generalized hypergeometric function is computed.
OUTPUT:
The generalized hypergeometric function defined by
.. MATH::
{}_pF_q(a_1,\ldots,a_p;b_1,\ldots,b_q;z) = \sum_{k=0}^\infty \frac{(a_1)_k\dots(a_p)_k}{(b_1)_k\dots(b_q)_k} \frac {z^k} {k!}
extended using analytic continuation or regularization when the sum does not converge.
The regularized generalized hypergeometric function
.. MATH::
{}_pF_q(a_1,\ldots,a_p;b_1,\ldots,b_q;z) = \sum_{k=0}^\infty \frac{(a_1)_k\dots(a_p)_k}{\Gamma(b_1+k)\dots\Gamma(b_q+k)} \frac {z^k} {k!}
is well-defined even when the lower parameters are nonpositive integers. Currently, this is only supported for some `p` and `q`.
EXAMPLES::
sage: CBF(1, pi/2).hypergeometric([], []) [+/- 7.72e-16] + [2.71828182845904 +/- 6.45e-15]*I
sage: CBF(1, pi).hypergeometric([1/4], [1/4]) [-2.7182818284590 +/- 7.11e-14] + [+/- 2.25e-14]*I
sage: CBF(1000, 1000).hypergeometric([10], [AA(sqrt(2))]) [9.79300951360e+454 +/- 5.01e+442] + [5.522579106816e+455 +/- 3.56e+442]*I sage: CBF(1000, 1000).hypergeometric([100], [AA(sqrt(2))]) [1.27967355557e+590 +/- 8.60e+578] + [-9.32333491987e+590 +/- 8.18e+578]*I
sage: CBF(0, 1).hypergeometric([], [1/2, 1/3, 1/4]) [-3.7991962344383 +/- 8.78e-14] + [23.878097177805 +/- 3.87e-13]*I
sage: CBF(0).hypergeometric([1], []) 1.000000000000000 sage: CBF(1, 1).hypergeometric([1], []) 1.000000000000000*I
sage: CBF(2+3*I).hypergeometric([1/4,1/3],[1/2]) [0.7871684267473 +/- 7.34e-14] + [0.2749254173721 +/- 9.23e-14]*I sage: CBF(2+3*I).hypergeometric([1/4,1/3],[1/2],regularized=True) [0.4441122268685 +/- 3.96e-14] + [0.1551100567338 +/- 5.75e-14]*I
sage: CBF(5).hypergeometric([2,3], [-5]) nan + nan*I sage: CBF(5).hypergeometric([2,3], [-5], regularized=True) [5106.925964355 +/- 5.41e-10]
sage: CBF(2016).hypergeometric([], [2/3]) [2.025642692328e+38 +/- 3.00e+25] sage: CBF(-2016).hypergeometric([], [2/3], regularized=True) [-0.0005428550847 +/- 5.00e-14]
sage: CBF(-7).hypergeometric([4], []) 0.0002441406250000000
sage: CBF(0, 3).hypergeometric([CBF(1,1)], [-4], regularized=True) [239.514000752841 +/- 8.03e-13] + [105.175157349015 +/- 6.28e-13]*I
TESTS::
sage: CBF(0, 1).hypergeometric([QQbar(sqrt(2)), RLF(pi)], [1r, 1/2]) [-8.7029449215408 +/- 6.17e-14] + [-0.8499070546106 +/- 5.21e-14]*I
""" cdef ComplexBall tmp, my_a, my_b, my_c regularized, prec(self)) regularized, prec(self)) self.value, regularized, prec(self)) raise NotImplementedError("regularized=True not yet supported in this case") cdef long s s = 1 self.value, -1, prec(self))
def hypergeometric_U(self, a, b): """ Return the Tricomi confluent hypergeometric function U(a, b, self) of this ball.
EXAMPLES::
sage: CBF(1000, 1000).hypergeometric_U(RLF(pi), -100) [-7.261605907166e-11 +/- 5.04e-24] + [-7.928136216391e-11 +/- 5.52e-24]*I sage: CBF(1000, 1000).hypergeometric_U(0, -100) 1.000000000000000 """
def erf(self): """ Return the error function with argument ``self``.
EXAMPLES::
sage: CBF(1, 1).erf() [1.316151281697947 +/- 7.26e-16] + [0.1904534692378347 +/- 6.03e-17]*I """
def erfc(self): """ Compute the complementary error function with argument ``self``.
EXAMPLES::
sage: CBF(20).erfc() [5.39586561160790e-176 +/- 6.73e-191] sage: CBF(100, 100).erfc() [0.00065234366376858 +/- 8.37e-18] + [-0.00393572636292141 +/- 7.21e-18]*I """
def airy(self): """ Return the Airy functions Ai, Ai', Bi, Bi' with argument ``self``, evaluated simultaneously.
EXAMPLES::
sage: CBF(10*pi).airy() ([1.2408955946101e-52 +/- 5.50e-66], [-6.965048886977e-52 +/- 5.23e-65], [2.288295683344e+50 +/- 5.10e+37], [1.2807602335816e+51 +/- 4.97e+37]) sage: ai, aip, bi, bip = CBF(1,2).airy() sage: (ai * bip - bi * aip) * CBF(pi) [1.0000000000000 +/- 1.25e-15] + [+/- 3.27e-16]*I
"""
def airy_ai(self): """ Return the Airy function Ai with argument ``self``.
EXAMPLES::
sage: CBF(1,2).airy_ai() [-0.2193862549814276 +/- 7.47e-17] + [-0.1753859114081094 +/- 4.06e-17]*I """
def airy_ai_prime(self): """ Return the Airy function derivative Ai' with argument ``self``.
EXAMPLES::
sage: CBF(1,2).airy_ai_prime() [0.1704449781789148 +/- 3.12e-17] + [0.387622439413295 +/- 1.06e-16]*I """
def airy_bi(self): """ Return the Airy function Bi with argument ``self``.
EXAMPLES::
sage: CBF(1,2).airy_bi() [0.0488220324530612 +/- 1.30e-17] + [0.1332740579917484 +/- 6.25e-17]*I """
def airy_bi_prime(self): """ Return the Airy function derivative Bi' with argument ``self``.
EXAMPLES::
sage: CBF(1,2).airy_bi_prime() [-0.857239258605362 +/- 3.47e-16] + [0.4955063363095674 +/- 9.22e-17]*I """
def bessel_J(self, nu): """ Return the Bessel function of the first kind with argument ``self`` and index ``nu``.
EXAMPLES::
sage: CBF(1, 1).bessel_J(1) [0.614160334922903 +/- 6.38e-16] + [0.365028028827088 +/- 3.96e-16]*I sage: CBF(100, -100).bessel_J(1/3) [1.108431870251e+41 +/- 5.53e+28] + [-8.952577603125e+41 +/- 2.93e+28]*I """
def bessel_J_Y(self, nu): """ Return the Bessel function of the first and second kind with argument ``self`` and index ``nu``, computed simultaneously.
EXAMPLES::
sage: J, Y = CBF(1, 1).bessel_J_Y(1) sage: J - CBF(1, 1).bessel_J(1) [+/- 3.75e-16] + [+/- 2.64e-16]*I sage: Y - CBF(1, 1).bessel_Y(1) [+/- 1.64e-14] + [+/- 1.62e-14]*I
""" my_nu.value, self.value, prec(self))
def bessel_Y(self, nu): """ Return the Bessel function of the second kind with argument ``self`` and index ``nu``.
EXAMPLES::
sage: CBF(1, 1).bessel_Y(1) [-0.6576945355913 +/- 5.29e-14] + [0.6298010039929 +/- 2.45e-14]*I sage: CBF(100, -100).bessel_Y(1/3) [-8.952577603125e+41 +/- 4.66e+28] + [-1.108431870251e+41 +/- 6.31e+28]*I """
def bessel_I(self, nu): """ Return the modified Bessel function of the first kind with argument ``self`` and index ``nu``.
EXAMPLES::
sage: CBF(1, 1).bessel_I(1) [0.365028028827088 +/- 3.96e-16] + [0.614160334922903 +/- 6.38e-16]*I sage: CBF(100, -100).bessel_I(1/3) [5.4362189595644e+41 +/- 6.40e+27] + [7.1989436985321e+41 +/- 2.92e+27]*I """
def bessel_K(self, nu): """ Return the modified Bessel function of the second kind with argument ``self`` and index ``nu``.
EXAMPLES::
sage: CBF(1, 1).bessel_K(0) [0.08019772694652 +/- 3.19e-15] + [-0.35727745928533 +/- 1.08e-15]*I sage: CBF(1, 1).bessel_K(1) [0.02456830552374 +/- 4.84e-15] + [-0.45971947380119 +/- 5.35e-15]*I sage: CBF(100, 100).bessel_K(QQbar(i)) [3.8693896656383e-45 +/- 2.76e-59] + [5.507100423418e-46 +/- 4.01e-59]*I """
def exp_integral_e(self, s): """ Return the image of this ball by the generalized exponential integral with index ``s``.
EXAMPLES::
sage: CBF(1+i).exp_integral_e(1) [0.00028162445198 +/- 2.79e-15] + [-0.17932453503936 +/- 2.12e-15]*I sage: CBF(1+i).exp_integral_e(QQbar(i)) [-0.10396361883964 +/- 3.78e-15] + [-0.16268401277783 +/- 3.69e-15]*I """
def ei(self): """ Return the exponential integral with argument ``self``.
EXAMPLES::
sage: CBF(1, 1).ei() [1.76462598556385 +/- 5.82e-15] + [2.38776985151052 +/- 4.29e-15]*I sage: CBF(0).ei() nan """
def si(self): """ Return the sine integral with argument ``self``.
EXAMPLES::
sage: CBF(1, 1).si() [1.10422265823558 +/- 2.48e-15] + [0.88245380500792 +/- 3.36e-15]*I sage: CBF(0).si() 0 """
def ci(self): """ Return the cosine integral with argument ``self``.
EXAMPLES::
sage: CBF(1, 1).ci() [0.882172180555936 +/- 4.85e-16] + [0.287249133519956 +/- 3.92e-16]*I sage: CBF(0).ci() nan + nan*I """
def shi(self): """ Return the hyperbolic sine integral with argument ``self``.
EXAMPLES::
sage: CBF(1, 1).shi() [0.88245380500792 +/- 3.36e-15] + [1.10422265823558 +/- 2.48e-15]*I sage: CBF(0).shi() 0 """
def chi(self): """ Return the hyperbolic cosine integral with argument ``self``.
EXAMPLES::
sage: CBF(1, 1).chi() [0.882172180555936 +/- 4.85e-16] + [1.28354719327494 +/- 1.07e-15]*I sage: CBF(0).chi() nan + nan*I """
def li(self, bint offset=False): """ Return the logarithmic integral with argument ``self``.
If ``offset`` is True, return the offset logarithmic integral.
EXAMPLES::
sage: CBF(1, 1).li() [0.61391166922120 +/- 6.40e-15] + [2.05958421419258 +/- 5.61e-15]*I sage: CBF(0).li() 0 sage: CBF(0).li(offset=True) [-1.045163780117493 +/- 5.54e-16] sage: li(0).n() 0.000000000000000 sage: Li(0).n() -1.04516378011749 """
def jacobi_theta(self, tau): r""" Return the four Jacobi theta functions evaluated at the argument ``self`` (representing `z`) and the parameter ``tau`` which should lie in the upper half plane.
The following definitions are used:
.. MATH::
\theta_1(z,\tau) = 2 q_{1/4} \sum_{n=0}^{\infty} (-1)^n q^{n(n+1)} \sin((2n+1) \pi z)
\theta_2(z,\tau) = 2 q_{1/4} \sum_{n=0}^{\infty} q^{n(n+1)} \cos((2n+1) \pi z)
\theta_3(z,\tau) = 1 + 2 \sum_{n=1}^{\infty} q^{n^2} \cos(2n \pi z)
\theta_4(z,\tau) = 1 + 2 \sum_{n=1}^{\infty} (-1)^n q^{n^2} \cos(2n \pi z)
where `q = \exp(\pi i \tau)` and `q_{1/4} = \exp(\pi i \tau / 4)`. Note that `z` is multiplied by `\pi`; some authors omit this factor.
EXAMPLES::
sage: CBF(3,-1/2).jacobi_theta(CBF(1/4,2)) ([-0.186580562274757 +/- 5.52e-16] + [0.93841744788594 +/- 2.48e-15]*I, [-1.02315311037951 +/- 4.10e-15] + [-0.203600094532010 +/- 7.33e-16]*I, [1.030613911309632 +/- 4.25e-16] + [0.030613917822067 +/- 1.89e-16]*I, [0.969386075665498 +/- 4.65e-16] + [-0.030613917822067 +/- 1.89e-16]*I)
sage: CBF(3,-1/2).jacobi_theta(CBF(1/4,-2)) (nan + nan*I, nan + nan*I, nan + nan*I, nan + nan*I)
sage: CBF(0).jacobi_theta(CBF(0,1)) (0, [0.91357913815612 +/- 3.96e-15], [1.086434811213308 +/- 8.16e-16], [0.913579138156117 +/- 8.89e-16])
"""
self.value, my_tau.value, prec(self))
def modular_j(self): """ Return the modular j-invariant with *tau* given by ``self``.
EXAMPLES::
sage: CBF(0,1).modular_j() [1728.0000000000 +/- 5.33e-11]
"""
def modular_eta(self): """ Return the Dedekind eta function with *tau* given by ``self``.
EXAMPLES::
sage: CBF(0,1).modular_eta() [0.768225422326057 +/- 9.18e-16] sage: CBF(12,1).modular_eta() [-0.768225422326057 +/- 9.18e-16]
"""
def modular_lambda(self): """ Return the modular lambda function with *tau* given by ``self``.
EXAMPLES::
sage: tau = CBF(sqrt(2),pi) sage: tau.modular_lambda() [-0.00022005123884157 +/- 6.41e-18] + [-0.0007978734645994 +/- 5.15e-17]*I sage: (tau + 2).modular_lambda() [-0.00022005123884157 +/- 6.41e-18] + [-0.0007978734645994 +/- 5.15e-17]*I sage: (tau / (1 - 2*tau)).modular_lambda() [-0.00022005123884 +/- 2.53e-15] + [-0.00079787346460 +/- 2.85e-15]*I
"""
def modular_delta(self): """ Return the modular discriminant with *tau* given by ``self``.
EXAMPLES::
sage: CBF(0,1).modular_delta() [0.0017853698506421 +/- 6.15e-17] sage: a, b, c, d = 2, 5, 1, 3 sage: tau = CBF(1,3) sage: ((a*tau+b)/(c*tau+d)).modular_delta() [0.20921376655 +/- 6.94e-12] + [1.5761192552 +/- 3.47e-11]*I sage: (c*tau+d)^12 * tau.modular_delta() [0.20921376654986 +/- 4.89e-15] + [1.5761192552253 +/- 4.45e-14]*I
"""
def eisenstein(self, long n): r""" Return the first ``n`` entries in the sequence of Eisenstein series `G_4(\tau), G_6(\tau), G_8(\tau), \ldots` where *tau* is given by ``self``. The output is a list.
EXAMPLES::
sage: a, b, c, d = 2, 5, 1, 3 sage: tau = CBF(1,3) sage: tau.eisenstein(4) [[2.1646498507193 +/- 6.30e-14], [2.0346794456073 +/- 2.44e-14], [2.0081609898081 +/- 3.67e-14], [2.0019857082706 +/- 4.60e-14]] sage: ((a*tau+b)/(c*tau+d)).eisenstein(3)[2] [331011.2004330 +/- 9.33e-8] + [-711178.1655746 +/- 7.51e-8]*I sage: (c*tau+d)^8 * tau.eisenstein(3)[2] [331011.20043304 +/- 7.62e-9] + [-711178.1655746 +/- 1.34e-8]*I
""" raise ValueError("n must be nonnegative")
def elliptic_p(self, tau, n=None): r""" Return the Weierstrass elliptic function with lattice parameter ``tau``, evaluated at ``self``. The function is doubly periodic in ``self`` with periods 1 and ``tau``, which should lie in the upper half plane.
If ``n`` is given, return a list containing the first ``n`` terms in the Taylor expansion at ``self``. In particular, with ``n`` = 2, compute the Weierstrass elliptic function together with its derivative, which generate the field of elliptic functions with periods 1 and ``tau``.
EXAMPLES::
sage: tau = CBF(1,4) sage: z = CBF(sqrt(2), sqrt(3)) sage: z.elliptic_p(tau) [-3.28920996772709 +/- 7.63e-15] + [-0.0003673767302933 +/- 6.04e-17]*I sage: (z + tau).elliptic_p(tau) [-3.28920996772709 +/- 7.97e-15] + [-0.000367376730293 +/- 6.51e-16]*I sage: (z + 1).elliptic_p(tau) [-3.28920996772709 +/- 7.63e-15] + [-0.0003673767302933 +/- 6.04e-17]*I
sage: z.elliptic_p(tau, 3) [[-3.28920996772709 +/- 7.62e-15] + [-0.0003673767302933 +/- 5.40e-17]*I, [0.002473055794309 +/- 5.01e-16] + [0.003859554040267 +/- 4.45e-16]*I, [-0.01299087561709 +/- 4.72e-15] + [0.00725027521915 +/- 4.32e-15]*I] sage: (z + 3 + 4*tau).elliptic_p(tau, 3) [[-3.28920996772709 +/- 8.4...e-15] + [-0.00036737673029 +/- 4.1...e-15]*I, [0.0024730557943 +/- 6.6...e-14] + [0.0038595540403 +/- 8.8...e-14]*I, [-0.01299087562 +/- 5.6...e-12] + [0.00725027522 +/- 3.5...e-12]*I]
""" cdef ComplexBall result cdef long nn cdef acb_ptr vec_r my_tau.value, prec(self)) else: raise ValueError("n must be nonnegative")
def elliptic_invariants(self): r""" Return the lattice invariants ``(g2, g3)``.
EXAMPLES::
sage: CBF(0,1).elliptic_invariants() ([189.07272012923 +/- 4.43e-12], [+/- 6.48e-12]) sage: CBF(sqrt(2)/2, sqrt(2)/2).elliptic_invariants() ([+/- 5.32e-12] + [-332.5338031465 +/- 5.03e-11]*I, [1254.4684215774 +/- 9.65e-11] + [1254.468421577 +/- 4.96e-10]*I) """
def elliptic_roots(self): r""" Return the lattice roots ``(e1, e2, e3)`` of `4 z^3 - g_2 z - g_3`.
EXAMPLES::
sage: e1, e2, e3 = CBF(0,1).elliptic_roots() sage: e1, e2, e3 ([6.8751858180204 +/- 6.18e-14], [+/- 1.20e-14], [-6.8751858180204 +/- 6.24e-14]) sage: g2, g3 = CBF(0,1).elliptic_invariants() sage: 4 * e1^3 - g2 * e1 - g3 [+/- 3.36e-11] """
def elliptic_k(self): """ Return the complete elliptic integral of the first kind evaluated at *m* given by ``self``.
EXAMPLES::
sage: CBF(2,3).elliptic_k() [1.04291329192852 +/- 5.9...e-15] + [0.62968247230864 +/- 3.4...e-15]*I
"""
def elliptic_e(self): """ Return the complete elliptic integral of the second kind evaluated at *m* given by ``self``.
EXAMPLES::
sage: CBF(2,3).elliptic_e() [1.472797144959 +/- 4.5...e-13] + [-1.231604783936 +/- 9.5...e-14]*I
"""
def elliptic_pi(self, m): """ Return the complete elliptic integral of the third kind evaluated at *m* given by ``self``.
EXAMPLES::
sage: CBF(2,3).elliptic_pi(CBF(1,1)) [0.27029997361983 +/- 1.31e-15] + [0.715676058329095 +/- 5.66e-16]*I
"""
def elliptic_f(self, m): r""" Return the incomplete elliptic integral of the first kind evaluated at *m*.
See :meth:`elliptic_k` for the corresponding complete integral
INPUT:
- ``m`` - complex ball
EXAMPLES::
sage: CBF(1,2).elliptic_f(CBF(0,1)) [0.6821522911854 +/- 2.96e-14] + [1.2482780628143 +/- 4.63e-14]*I
At parameter `\pi/2` it is a complete integral::
sage: phi = CBF(1,1) sage: (CBF.pi()/2).elliptic_f(phi) [1.5092369540513 +/- 6.62e-14] + [0.6251464152027 +/- 2.11e-14]*I sage: phi.elliptic_k() [1.50923695405127 +/- 5.07e-15] + [0.62514641520270 +/- 4.41e-15]*I
sage: phi = CBF(2, 3/7) sage: (CBF.pi()/2).elliptic_f(phi) [1.339358963909 +/- 5.02e-13] + [1.110436969072 +/- 1.37e-13]*I sage: phi.elliptic_k() [1.33935896390938 +/- 6.73e-15] + [1.11043696907194 +/- 6.41e-15]*I
"""
def elliptic_e_inc(self, m): r""" Return the incomplete elliptic integral of the second kind evaluated at *m*.
See :meth:`elliptic_e` for the corresponding complete integral
INPUT:
- ``m`` - complex ball
EXAMPLES::
sage: CBF(1,2).elliptic_e_inc(CBF(0,1)) [1.906576998914 +/- 5.01e-13] + [3.6896645289411 +/- 6.93e-14]*I
At parameter `\pi/2` it is a complete integral::
sage: phi = CBF(1,1) sage: (CBF.pi()/2).elliptic_e_inc(phi) [1.283840957898 +/- 3.23e-13] + [-0.5317843366915 +/- 7.79e-14]*I sage: phi.elliptic_e() [1.2838409578982 +/- 5.90e-14] + [-0.5317843366915 +/- 3.35e-14]*I
sage: phi = CBF(2, 3/7) sage: (CBF.pi()/2).elliptic_e_inc(phi) [0.787564350925 +/- 6.56e-13] + [-0.686896129145 +/- 4.60e-13]*I sage: phi.elliptic_e() [0.7875643509254 +/- 6.97e-14] + [-0.686896129145 +/- 3.11e-13]*I
"""
def elliptic_pi_inc(self, phi, m): r""" Return the Legendre incomplete elliptic integral of the third kind.
See: :meth:`elliptic_pi` for the complete integral.
INPUT:
- ``phi`` - complex ball
- ``m`` - complex ball
EXAMPLES::
sage: CBF(1,2).elliptic_pi_inc(CBF(0,1), CBF(2,-3)) [0.05738864021418 +/- 4.27e-15] + [0.55557494549951 +/- 5.71e-15]*I
At parameter `\pi/2` it is a complete integral::
sage: n = CBF(1,1) sage: m = CBF(-2/3, 3/5) sage: n.elliptic_pi_inc(CBF.pi()/2, m) [0.8934793755173 +/- 5.65e-14] + [0.9570786871075 +/- 1.98e-14]*I sage: n.elliptic_pi(m) [0.89347937551733 +/- 4.07e-15] + [0.95707868710750 +/- 1.23e-15]*I
sage: n = CBF(2, 3/7) sage: m = CBF(-1/3, 2/9) sage: n.elliptic_pi_inc(CBF.pi()/2, m) [0.296958874642 +/- 2.58e-13] + [1.318879533274 +/- 3.87e-13]*I sage: n.elliptic_pi(m) [0.29695887464189 +/- 4.98e-15] + [1.31887953327376 +/- 5.95e-15]*I """
def elliptic_rf(self, y, z): r""" Return the Carlson symmetric elliptic integral of the first kind evaluated at ``(self, y, z)``.
INPUT:
- ``y`` - complex ball
- ``z`` - complex ball
EXAMPLES::
sage: CBF(0,1).elliptic_rf(CBF(-1/2,1), CBF(-1,-1)) [1.469800396738515 +/- 3.70e-16] + [-0.2358791199824196 +/- 4.40e-17]*I
"""
def elliptic_rg(self, y, z): r""" Return the Carlson symmetric elliptic integral of the second kind evaluated at ``(self, y, z)``.
INPUT:
- ``y`` - complex ball
- ``z`` - complex ball
EXAMPLES::
sage: CBF(0,1).elliptic_rg(CBF(-1/2,1), CBF(-1,-1)) [0.1586786770922370 +/- 4.31e-17] + [0.2239733128130531 +/- 3.35e-17]*I
"""
def elliptic_rj(self, y, z, p): r""" Return the Carlson symmetric elliptic integral of the third kind evaluated at ``(self, y, z)``.
INPUT:
- ``y`` - complex ball
- ``z`` - complex ball
- ``p`` - complex bamm
EXAMPLES::
sage: CBF(0,1).elliptic_rj(CBF(-1/2,1), CBF(-1,-1), CBF(2)) [1.004386756285733 +/- 5.21e-16] + [-0.2451626834391645 +/- 3.50e-17]*I
"""
def elliptic_zeta(self, tau): r""" Return the value of the Weierstrass zeta function at ``(self, tau)``
EXAMPLES::
- ``tau`` - a complex ball with positive imaginary part
EXAMPLES::
sage: CBF(1,1).elliptic_zeta(CBF(1,3)) [3.2898676194970 +/- 5.93e-14] + [0.1365414361782 +/- 7.27e-14]*I """
def elliptic_sigma(self, tau): r""" Return the value of the Weierstrass sigma function at ``(self, tau)``
EXAMPLES::
- ``tau`` - a complex ball with positive imaginary part
EXAMPLES::
sage: CBF(1,1).elliptic_sigma(CBF(1,3)) [-0.543073363596 +/- 3.39e-13] + [3.635729118624 +/- 5.08e-13]*I
"""
def chebyshev_T(self, n): """ Return the Chebyshev function of the first kind of order ``n`` evaluated at ``self``.
EXAMPLES::
sage: CBF(1/3).chebyshev_T(20) [0.8710045668809 +/- 6.15e-14] sage: CBF(1/3).chebyshev_T(CBF(5,1)) [1.8429685451876 +/- 3.57e-14] + [0.20053614301799 +/- 7.05e-15]*I
"""
def chebyshev_U(self, n): """ Return the Chebyshev function of the second kind of order ``n`` evaluated at ``self``.
EXAMPLES::
sage: CBF(1/3).chebyshev_U(20) [0.6973126541184 +/- 2.83e-14] sage: CBF(1/3).chebyshev_U(CBF(5,1)) [1.7588496489342 +/- 5.99e-14] + [0.7497317165104 +/- 4.35e-14]*I
"""
def jacobi_P(self, n, a, b): r""" Return the Jacobi polynomial (or function) `P_n^{(a,b)}(z)` evaluated at ``self``.
EXAMPLES::
sage: CBF(5,-6).jacobi_P(8, CBF(1,2), CBF(2,3)) [-920983000.45982 +/- 2.22e-6] + [6069919969.92857 +/- 4.77e-6]*I
""" my_a.value, my_b.value, self.value, prec(self))
def gegenbauer_C(self, n, m): r""" Return the Gegenbauer polynomial (or function) `C_n^m(z)` evaluated at ``self``.
EXAMPLES::
sage: CBF(-10).gegenbauer_C(7, 1/2) [-263813415.6250000 +/- 9.57e-8]
""" my_m.value, self.value, prec(self))
def laguerre_L(self, n, m=0): r""" Return the Laguerre polynomial (or function) `L_n^m(z)` evaluated at ``self``.
EXAMPLES::
sage: CBF(10).laguerre_L(3) [-45.6666666666666 +/- 9.28e-14] sage: CBF(10).laguerre_L(3, 2) [-6.666666666667 +/- 4.15e-13] sage: CBF(5,7).laguerre_L(CBF(2,3), CBF(1,-2)) [5515.315030271 +/- 4.37e-10] + [-12386.942845271 +/- 5.47e-10]*I
""" my_m.value, self.value, prec(self))
def hermite_H(self, n): """ Return the Hermite function (or polynomial) of order ``n`` evaluated at ``self``.
EXAMPLES::
sage: CBF(10).hermite_H(1) 20.00000000000000 sage: CBF(10).hermite_H(30) [8.0574670961707e+37 +/- 3.28e+23]
"""
def legendre_P(self, n, m=0, type=2): r""" Return the Legendre function of the first kind `P_n^m(z)` evaluated at ``self``.
The ``type`` parameter can be either 2 or 3. This selects between different branch cut conventions. The definitions of the "type 2" and "type 3" functions are the same as those used by *Mathematica* and *mpmath*.
EXAMPLES::
sage: CBF(1/2).legendre_P(5) [0.08984375000000000 +/- 4.5...e-18] sage: CBF(1,2).legendre_P(CBF(2,3), CBF(0,1)) [0.10996180744364 +/- 7.45e-15] + [0.14312767804055 +/- 8.38e-15]*I sage: CBF(-10).legendre_P(5, 325/100) [-22104403.487377 +/- 6.81e-7] + [53364750.687392 +/- 7.25e-7]*I sage: CBF(-10).legendre_P(5, 325/100, type=3) [-57761589.914581 +/- 6.99e-7] + [+/- 5.14e-7]*I
""" raise ValueError("expected type 2 or 3") my_m.value, self.value, my_type - 2, prec(self))
def legendre_Q(self, n, m=0, type=2): r""" Return the Legendre function of the second kind `Q_n^m(z)` evaluated at ``self``.
The ``type`` parameter can be either 2 or 3. This selects between different branch cut conventions. The definitions of the "type 2" and "type 3" functions are the same as those used by *Mathematica* and *mpmath*.
EXAMPLES::
sage: CBF(1/2).legendre_Q(5) [0.55508089057168 +/- 2.79e-15] sage: CBF(1,2).legendre_Q(CBF(2,3), CBF(0,1)) [0.167678710 +/- 4.60e-10] + [-0.161558598 +/- 7.47e-10]*I sage: CBF(-10).legendre_Q(5, 325/100) [-83825154.36008 +/- 4.94e-6] + [-34721515.80396 +/- 5.40e-6]*I sage: CBF(-10).legendre_Q(5, 325/100, type=3) [-4.797306921692e-6 +/- 6.82e-19] + [-4.797306921692e-6 +/- 6.57e-19]*I
""" raise ValueError("expected type 2 or 3") my_m.value, self.value, my_type - 2, prec(self))
def spherical_harmonic(self, phi, n, m): r""" Return the spherical harmonic `Y_n^m(\theta,\phi)` evaluated at `\theta` given by ``self``. In the current implementation, ``n`` and ``m`` must be small integers.
EXAMPLES::
sage: CBF(1+I).spherical_harmonic(1/2, -3, -2) [0.80370071745224 +/- 4.02e-15] + [-0.07282031864711 +/- 4.69e-15]*I """ self.value, my_phi.value, prec(self))
CBF = ComplexBallField() |