Arbitrary precision complex balls

Arbitrary precision complex ball arithmetic is supplied by Arb which provides a ball representation which tracks error bounds rigorously. Complex numbers are represented in rectangular form $a+bi$ where $a,b$ are ArbFieldElem balls.

The corresponding field is constructed using the ComplexField constructor. This constructs the parent object for the Arb complex field.

The types of complex boxes in Nemo are given in the following table, along with the libraries that provide them and the associated types of the parent objects.

Library	Field	Element type	Parent type
Arb	$\mathbb{C}$ (boxes)	`ComplexFieldElem`	`ComplexField`

All the complex field types belong to the Field abstract type and the types of elements in this field, i.e. complex boxes in this case, belong to the FieldElem abstract type.

Complex ball functionality

The complex balls in Nemo provide all the field functionality defined by AbstractAlgebra:.

https://nemocas.github.io/AbstractAlgebra.jl/stable/field

Below, we document the additional functionality provided for complex balls.

Precision management

See Precision management.

Complex field constructors

In order to construct complex boxes in Nemo, one must first construct the Arb complex field itself. This is accomplished with the following constructor.

complex_field — Method

complex_field()

Return the field of complex numbers modelled via complex balls.

See precision and set_precision! on how to control the precision.

Here is an example of creating an Arb complex field and using the resulting parent object to coerce values into the resulting field.

Examples

julia> CC = complex_field()
Complex field

julia> a = CC("0.25")
0.25000000000000000000

julia> b = CC("0.1")
[0.100000000000000000 +/- 1.22e-20]

julia> c = CC(0.5)
0.50000000000000000000

julia> d = CC(12)
12.000000000000000000

Note that whilst one can coerce double precision floating point values into an Arb complex field, unless those values can be represented exactly in double precision the resulting ball can't be any more precise than the double precision supplied.

If instead, values can be represented precisely using decimal arithmetic then one can supply them to Arb using a string. In this case, Arb will store them to the precision specified when creating the Arb complex field.

If the values can be stored precisely as a binary floating point number, Arb will store the values exactly. See the function is_exact below for more information.

Constructors

onei — Method

onei(r::ComplexField)

Return exact one times $i$ in the given complex field.

Examples

julia> c = onei(CC)
1.0000000000000000000*im

Basic functionality

The following basic functionality is provided by the default Arb complex field implementation in Nemo, to support construction of generic rings over complex fields. Any custom complex field implementation in Nemo should provide analogues of these functions along with the usual arithmetic operations.

parent_type(::Type{ComplexFieldElem})

Gives the type of the parent object of an Arb complex field element.

elem_type(R::ComplexField)

Given the parent object for an Arb complex field, return the type of elements of the field.

mul!(c::ComplexFieldElem, a::ComplexFieldElem, b::ComplexFieldElem)

Multiply $a$ by $b$ and set the existing Arb complex field element $c$ to the result. This function is provided for performance reasons as it saves allocating a new object for the result and eliminates associated garbage collection.

deepcopy(a::ComplexFieldElem)

Return a copy of the Arb complex field element $a$, recursively copying the internal data. Arb complex field elements are mutable in Nemo so a shallow copy is not sufficient.

Given the parent object R for an Arb complex field, the following coercion functions are provided to coerce various elements into the Arb complex field. Developers provide these by overloading the call operator for the complex field parent objects.

R()

Coerce zero into the Arb complex field.

R(n::Integer)
R(f::ZZRingElem)
R(q::QQFieldElem)

Coerce an integer or rational value into the Arb complex field.

R(f::Float64)
R(f::BigFloat)

Coerce the given floating point number into the Arb complex field.

R(f::AbstractString)
R(f::AbstractString, g::AbstractString)

Coerce the decimal number, given as a string, into the Arb complex field. In each case $f$ is the real part and $g$ is the imaginary part.

R(f::ArbFieldElem)

Coerce the given Arb real ball into the Arb complex field.

R(f::ComplexFieldElem)

Take an Arb complex field element that is already in an Arb field and simply return it. A copy of the original is not made.

Here are some examples of coercing elements into the Arb complex field.

julia> RR = real_field()
Real field

julia> CC = complex_field()
Complex field

julia> a = CC(3)
3.0000000000000000000

julia> b = CC(QQ(2,3))
[0.6666666666666666666 +/- 8.48e-20]

julia> c = CC("3 +/- 0.0001")
[3.000 +/- 1.01e-4]

julia> d = CC("-1.24e+12345")
[-1.240000000000000000e+12345 +/- 1.16e+12326]

julia> f = CC("nan +/- inf")
nan

julia> g = CC(RR(3))
3.0000000000000000000

In addition to the above, developers of custom complex field types must ensure that they provide the equivalent of the function base_ring(R::ComplexField) which should return Union{}. In addition to this they should ensure that each complex field element contains a field parent specifying the parent object of the complex field element, or at least supply the equivalent of the function parent(a::ComplexFieldElem) to return the parent object of a complex field element.

Basic manipulation

isfinite — Method

isfinite(x::ComplexFieldElem)

Return true if $x$ is finite, i.e. its real and imaginary parts have finite midpoint and radius, otherwise return false.

is_exact — Method

is_exact(x::ComplexFieldElem)

Return true if $x$ is exact, i.e. has its real and imaginary parts have zero radius, otherwise return false.

isinteger — Method

isinteger(x::ComplexFieldElem)

Return true if $x$ is an exact integer, otherwise return false.

accuracy_bits — Method

accuracy_bits(x::ComplexFieldElem)

Return the relative accuracy of $x$ measured in bits, capped between typemax(Int) and -typemax(Int).

Examples

julia> a = CC("1.2 +/- 0.001")
[1.20 +/- 1.01e-3]

julia> b = CC(3)
3.0000000000000000000

julia> isreal(a)
true

julia> isfinite(b)
true

julia> isinteger(b)
true

julia> c = real(a)
[1.20 +/- 1.01e-3]

julia> d = imag(b)
0

julia> f = accuracy_bits(a)
9

Containment

It is often necessary to determine whether a given exact value or box is contained in a given complex box or whether two boxes overlap. The following functions are provided for this purpose.

overlaps — Method

overlaps(x::ComplexFieldElem, y::ComplexFieldElem)

Returns true if any part of the box $x$ overlaps any part of the box $y$, otherwise return false.

contains — Method

contains(x::ComplexFieldElem, y::ComplexFieldElem)

Returns true if the box $x$ contains the box $y$, otherwise return false.

contains — Method

contains(x::ComplexFieldElem, y::Integer)

Returns true if the box $x$ contains the given integer value, otherwise return false.

contains — Method

contains(x::ComplexFieldElem, y::ZZRingElem)

Returns true if the box $x$ contains the given integer value, otherwise return false.

contains — Method

contains(x::ComplexFieldElem, y::QQFieldElem)

Returns true if the box $x$ contains the given rational value, otherwise return false.

The following functions are also provided for determining if a box intersects a certain part of the complex number plane.

contains_zero — Method

contains_zero(x::ComplexFieldElem)

Returns true if the box $x$ contains zero, otherwise return false.

Examples

julia> x = CC("1 +/- 0.001")
[1.00 +/- 1.01e-3]

julia> y = CC("3")
3.0000000000000000000

julia> overlaps(x, y)
false

julia> contains(x, y)
false

julia> contains(y, 3)
true

julia> contains(x, ZZ(1)//2)
false

julia> contains_zero(x)
false

Comparison

Nemo provides a full range of comparison operations for Arb complex boxes.

In addition to the standard comparisons, we introduce an exact equality. This is distinct from arithmetic equality implemented by ==, which merely compares up to the minimum of the precisions of its operands.

isequal — Method

isequal(x::ComplexFieldElem, y::ComplexFieldElem)

Return true if the boxes $x$ and $y$ are precisely equal, i.e. their real and imaginary parts have the same midpoints and radii.

A full range of ad hoc comparison operators is provided. These are implemented directly in Julia, but we document them as though only == were provided.

Function
`==(x::ComplexFieldElem, y::Integer)`
`==(x::Integer, y::ComplexFieldElem)`
`==(x::ComplexFieldElem, y::ZZRingElem)`
`==(x::ZZRingElem, y::ComplexFieldElem)`
`==(x::ArbFieldElem, y::ZZRingElem)`
`==(x::ZZRingElem, y::ArbFieldElem)`
`==(x::ComplexFieldElem, y::Float64)`
`==(x::Float64, y::ComplexFieldElem)`

Examples

julia> x = CC("1 +/- 0.001")
[1.00 +/- 1.01e-3]

julia> y = CC("3")
3.0000000000000000000

julia> z = CC("4")
4.0000000000000000000

julia> isequal(x, deepcopy(x))
true

julia> x == 3
false

julia> ZZ(3) == y
true

julia> z != 1.23
true

Absolute value

Examples

julia> x = CC("-1 +/- 0.001")
[-1.00 +/- 1.01e-3]

julia> a = abs(x)
[1.00 +/- 1.01e-3]

Shifting

Examples

julia> x = CC("-3 +/- 0.001")
[-3.00 +/- 1.01e-3]

julia> a = ldexp(x, 23)
[-2.52e+7 +/- 4.26e+4]

julia> b = ldexp(x, -ZZ(15))
[-9.16e-5 +/- 7.78e-8]

Miscellaneous operations

trim — Method

trim(x::ComplexFieldElem)

Return an ComplexFieldElem box containing $x$ but which may be more economical, by rounding off insignificant bits from midpoints.

unique_integer — Method

unique_integer(x::ComplexFieldElem)

Return a pair where the first value is a boolean and the second is an ZZRingElem integer. The boolean indicates whether the box $x$ contains a unique integer. If this is the case, the second return value is set to this unique integer.

Examples

julia> x = CC("-3 +/- 0.001", "0.1")
[-3.00 +/- 1.01e-3] + [0.100000000000000000 +/- 1.22e-20]*im

julia> a = trim(x)
[-3.00 +/- 1.01e-3] + [0.100000000000000000 +/- 1.22e-20]*im

julia> b, c = unique_integer(x)
(false, 0)

julia> d = conj(x)
[-3.00 +/- 1.01e-3] + [-0.100000000000000000 +/- 1.22e-20]*im

julia> f = angle(x)
[3.1083 +/- 3.95e-5]

Constants

const_pi — Method

const_pi(r::ComplexField)

Return $\pi = 3.14159\ldots$ as an element of $r$.

Examples

CC = complex_field()
set_precision!(ComplexField, 200) do
  a = const_pi(CC)
end

Mathematical and special functions

rsqrt — Method

rsqrt(x::ComplexFieldElem)

Return the reciprocal of the square root of $x$, i.e. $1/\sqrt{x}$.

cispi — Method

cispi(x::ComplexFieldElem)

Return the exponential of $\pi i x$.

root_of_unity — Method

root_of_unity(C::ComplexField, k::Int)

Return $\exp(2\pi i/k)$.

log_sinpi — Method

log_sinpi(x::ComplexFieldElem)

Return $\log\sin(\pi x)$, constructed without branch cuts off the real line.

gamma — Method

gamma(x::ComplexFieldElem)

Return the Gamma function evaluated at $x$.

lgamma — Method

lgamma(x::ComplexFieldElem)

Return the logarithm of the Gamma function evaluated at $x$.

rgamma — Method

rgamma(x::ComplexFieldElem)

Return the reciprocal of the Gamma function evaluated at $x$.

digamma — Method

digamma(x::ComplexFieldElem)

Return the logarithmic derivative of the gamma function evaluated at $x$, i.e. $\psi(x)$.

zeta — Method

zeta(x::ComplexFieldElem)

Return the Riemann zeta function evaluated at $x$.

barnes_g — Method

barnes_g(x::ComplexFieldElem)

Return the Barnes $G$-function, evaluated at $x$.

log_barnes_g — Method

log_barnes_g(x::ComplexFieldElem)

Return the logarithm of the Barnes $G$-function, evaluated at $x$.

erf — Method

erf(x::ComplexFieldElem)

Return the error function evaluated at $x$.

erfi — Method

erfi(x::ComplexFieldElem)

Return the imaginary error function evaluated at $x$.

exp_integral_ei — Method

exp_integral_ei(x::ComplexFieldElem)

Return the exponential integral evaluated at $x$.

sin_integral — Method

sin_integral(x::ComplexFieldElem)

Return the sine integral evaluated at $x$.

cos_integral — Method

cos_integral(x::ComplexFieldElem)

Return the exponential cosine integral evaluated at $x$.

sinh_integral — Method

sinh_integral(x::ComplexFieldElem)

Return the hyperbolic sine integral evaluated at $x$.

cosh_integral — Method

cosh_integral(x::ComplexFieldElem)

Return the hyperbolic cosine integral evaluated at $x$.

dedekind_eta — Method

dedekind_eta(x::ComplexFieldElem)

Return the Dedekind eta function $\eta(\tau)$ at $\tau = x$.

modular_weber_f — Method

modular_weber_f(x::ComplexFieldElem)

Return the modular Weber function $\mathfrak{f}(\tau) = \frac{\eta^2(\tau)}{\eta(\tau/2)\eta(2\tau)},$ at $x$ in the complex upper half plane.

modular_weber_f1 — Method

modular_weber_f1(x::ComplexFieldElem)

Return the modular Weber function $\mathfrak{f}_1(\tau) = \frac{\eta(\tau/2)}{\eta(\tau)},$ at $x$ in the complex upper half plane.

modular_weber_f2 — Method

modular_weber_f2(x::ComplexFieldElem)

Return the modular Weber function $\mathfrak{f}_2(\tau) = \frac{\sqrt{2}\eta(2\tau)}{\eta(\tau)}$ at $x$ in the complex upper half plane.

j_invariant — Method

j_invariant(x::ComplexFieldElem)

Return the $j$-invariant $j(\tau)$ at $\tau = x$.

j_invariant(E::EllipticCurve) -> FieldElem

Compute the j-invariant of $E$.

modular_lambda — Method

modular_lambda(x::ComplexFieldElem)

Return the modular lambda function $\lambda(\tau)$ at $\tau = x$.

modular_delta — Method

modular_delta(x::ComplexFieldElem)

Return the modular delta function $\Delta(\tau)$ at $\tau = x$.

eisenstein_g — Method

eisenstein_g(k::Int, x::ComplexFieldElem)

Return the non-normalized Eisenstein series $G_k(\tau)$ of $\mathrm{SL}_2(\mathbb{Z})$. Also defined for $\tau = i \infty$.

hilbert_class_polynomial — Method

hilbert_class_polynomial(D::Int, R::ZZPolyRing)

Return in the ring $R$ the Hilbert class polynomial of discriminant $D$, which is only defined for $D < 0$ and $D \equiv 0, 1 \pmod 4$.

elliptic_k — Method

elliptic_k(x::ComplexFieldElem)

Return the complete elliptic integral $K(x)$.

elliptic_e — Method

elliptic_e(x::ComplexFieldElem)

Return the complete elliptic integral $E(x)$.

agm — Method

agm(x::ComplexFieldElem)

Return the arithmetic-geometric mean of $1$ and $x$.

agm — Method

agm(x::ComplexFieldElem, y::ComplexFieldElem)

Return the arithmetic-geometric mean of $x$ and $y$.

polygamma — Method

polygamma(s::ComplexFieldElem, a::ComplexFieldElem)

Return the generalised polygamma function $\psi(s,z)$.

zeta — Method

zeta(s::ComplexFieldElem, a::ComplexFieldElem)

Return the Hurwitz zeta function $\zeta(s,a)$.

rising_factorial — Method

rising_factorial(x::ComplexFieldElem, n::Int)

Return the rising factorial $x(x + 1)\ldots (x + n - 1)$.

rising_factorial2 — Method

rising_factorial2(x::ComplexFieldElem, n::Int)

Return a tuple containing the rising factorial $x(x + 1)\ldots (x + n - 1)$ and its derivative.

polylog — Method

polylog(s::Union{ComplexFieldElem,Int}, a::ComplexFieldElem)

Return the polylogarithm Li$_s(a)$.

log_integral — Method

log_integral(x::ComplexFieldElem)

Return the logarithmic integral, evaluated at $x$.

log_integral_offset — Method

log_integral_offset(x::ComplexFieldElem)

Return the offset logarithmic integral, evaluated at $x$.

exp_integral_e — Method

exp_integral_e(s::ComplexFieldElem, x::ComplexFieldElem)

Return the generalised exponential integral $E_s(x)$.

gamma — Method

gamma(s::ComplexFieldElem, x::ComplexFieldElem)

Return the upper incomplete gamma function $\Gamma(s,x)$.

gamma_regularized — Method

gamma_regularized(s::ComplexFieldElem, x::ComplexFieldElem)

Return the regularized upper incomplete gamma function $\Gamma(s,x) / \Gamma(s)$.

gamma_lower — Method

gamma_lower(s::ComplexFieldElem, x::ComplexFieldElem)

Return the lower incomplete gamma function $\gamma(s,x) / \Gamma(s)$.

gamma_lower_regularized — Method

gamma_lower_regularized(s::ComplexFieldElem, x::ComplexFieldElem)

Return the regularized lower incomplete gamma function $\gamma(s,x) / \Gamma(s)$.

airy_ai — Method

airy_ai(x::ComplexFieldElem)

Return the Airy function $\operatorname{Ai}(x)$.

airy_ai_prime — Method

airy_ai_prime(x::ComplexFieldElem)

Return the derivative of the Airy function $\operatorname{Ai}^\prime(x)$.

airy_bi — Method

airy_bi(x::ComplexFieldElem)

Return the Airy function $\operatorname{Bi}(x)$.

airy_bi_prime — Method

airy_bi_prime(x::ComplexFieldElem)

Return the derivative of the Airy function $\operatorname{Bi}^\prime(x)$.

bessel_j — Method

bessel_j(nu::ComplexFieldElem, x::ComplexFieldElem)

Return the Bessel function $J_{\nu}(x)$.

bessel_y — Method

bessel_y(nu::ComplexFieldElem, x::ComplexFieldElem)

Return the Bessel function $Y_{\nu}(x)$.

bessel_i — Method

bessel_i(nu::ComplexFieldElem, x::ComplexFieldElem)

Return the Bessel function $I_{\nu}(x)$.

bessel_k — Method

bessel_k(nu::ComplexFieldElem, x::ComplexFieldElem)

Return the Bessel function $K_{\nu}(x)$.

hypergeometric_1f1 — Method

hypergeometric_1f1(a::ComplexFieldElem, b::ComplexFieldElem, x::ComplexFieldElem)

Return the confluent hypergeometric function ${}_1F_1(a,b,x)$.

hypergeometric_1f1_regularized — Method

hypergeometric_1f1_regularized(a::ComplexFieldElem, b::ComplexFieldElem, x::ComplexFieldElem)

Return the regularized confluent hypergeometric function ${}_1F_1(a,b,x) / \Gamma(b)$.

hypergeometric_u — Method

hypergeometric_u(a::ComplexFieldElem, b::ComplexFieldElem, x::ComplexFieldElem)

Return the confluent hypergeometric function $U(a,b,x)$.

hypergeometric_2f1 — Method

hypergeometric_2f1(a::ComplexFieldElem, b::ComplexFieldElem, c::ComplexFieldElem, x::ComplexFieldElem; flags=0)

Return the Gauss hypergeometric function ${}_2F_1(a,b,c,x)$.

jacobi_theta — Method

jacobi_theta(z::ComplexFieldElem, tau::ComplexFieldElem)

Return a tuple of four elements containing the Jacobi theta function values $\theta_1, \theta_2, \theta_3, \theta_4$ evaluated at $z, \tau$.

weierstrass_p — Method

weierstrass_p(z::ComplexFieldElem, tau::ComplexFieldElem)

Return the Weierstrass elliptic function $\wp(z,\tau)$.

Examples

julia> s = CC(1, 2)
1.0000000000000000000 + 2.0000000000000000000*im

julia> z = CC("1.23", "3.45")
[1.230000000000000000 +/- 2.00e-19] + [3.450000000000000000 +/- 3.91e-19]*im

julia> a = sin(z)^2 + cos(z)^2
[1.000000000000000 +/- 4.92e-16] + [+/- 4.12e-16]*im

julia> b = zeta(z)
[0.685803329024164062 +/- 6.30e-19] + [-0.038574782404586856 +/- 7.54e-19]*im

julia> c = bessel_j(s, z)
[0.63189634741402481 +/- 4.85e-18] + [0.00970090757446076 +/- 4.66e-18]*im

julia> d = hypergeometric_1f1(s, s+1, z)
[-1.3355297330012291 +/- 5.83e-17] + [-0.1715020340928697 +/- 4.97e-17]*im

Linear dependence

lindep — Method

lindep(A::Vector{ComplexFieldElem}, bits::Int)

Find a small linear combination of the entries of the array $A$ that is small (using LLL). The entries are first scaled by the given number of bits before truncating the real and imaginary parts to integers for use in LLL. This function can be used to find linear dependence between a list of complex numbers. The algorithm is heuristic only and returns an array of Nemo integers representing the linear combination.

lindep — Method

lindep(A::Matrix{ComplexFieldElem}, bits::Int)

Find a (common) small linear combination of the entries in each row of the array $A$, that is small (using LLL). It is assumed that the complex numbers in each row of the array share the same linear combination. The entries are first scaled by the given number of bits before truncating the real and imaginary parts to integers for use in LLL. This function can be used to find a common linear dependence shared across a number of lists of complex numbers. The algorithm is heuristic only and returns an array of Nemo integers representing the common linear combination.

Examples

julia> # These are two of the roots of x^5 + 3x + 1

julia> a = CC(1.0050669478588622428791051888364775253, -0.93725915669289182697903585868761513585)
[1.0050669478588623029 +/- 2.25e-20] - [0.93725915669289183718 +/- 1.50e-21]*im

julia> b = CC(-0.33198902958450931620250069492231652319)
-[0.33198902958450932088 +/- 4.15e-22]

julia> V1 = [CC(1), a, a^2, a^3, a^4, a^5]; # We recover the polynomial from one root....

julia> W = lindep(V1, 20)
6-element Vector{ZZRingElem}:
 1
 3
 0
 0
 0
 1

julia> V2 = [CC(1), b, b^2, b^3, b^4, b^5]; # ...or from two

julia> Vs = [transpose(V1); transpose(V2)];

julia> X = lindep(Vs, 20)
6-element Vector{ZZRingElem}:
 1
 3
 0
 0
 0
 1