GitHub

Arblib.jl

Arblib.jl is a Julia package for rigorous numerics based on ball arithmetic. It provides a thin and efficient wrapper around the parts of the C library FLINT that are concerned with real and complex numbers.

The purpose of the Arblib package extension is to make conversions between the types defined in IntervalArithmetic and Arblib straightforward. This allows for easy switching between the two packages, depending on which one suits a specific part of a computation best.

Some of the things that Arblib excels at are:

  1. Fast high precision computations, including linear algebra routines.
  2. A large library of special functions.
  3. Support for mutable arithmetic.
  4. Taylor series expansions.

Some of the things that IntervalArithmetic excels at are:

  1. Fast computations at Float32 and Float64 precision.
  2. Computations with wide intervals (FLINT is in general not optimized for this, though the situation is improving).
  3. Built-in safety features, such as decorations (see Decorations) and the "NG" flag (see Guarantee).

Conversion between Interval and Arblib types

The fundamental types in Arblib are Arb and Acb, corresponding to real and complex balls respectively. Conversion between Interval and Arb, as well as between Complex{<:Interval} and Acb, is done through the standard constructors. If no type is specified when calling interval, it defaults to Interval{BigFloat} for Arb and Complex{Interval{BigFloat}} for Acb.

julia> using IntervalArithmetic, ArblibPrecompiling IntervalArithmeticArblibExt...
    598.9 ms  ? IntervalArithmetic → IntervalArithmeticLinearAlgebraExt
Info Given IntervalArithmeticArblibExt was explicitly requested, output will be shown live 
┌ Warning: Module IntervalArithmetic with build ID fafbfcfd-70d9-b016-0000-00178308cb32 is missing from the cache.
│ This may mean IntervalArithmetic [d1acc4aa-44c8-5952-acd4-ba5d80a2a253] does not support precompilation but is imported by a module that does.
└ @ Base loading.jl:2541
    648.3 ms  ? IntervalArithmetic → IntervalArithmeticArblibExt
┌ Warning: Module IntervalArithmetic with build ID fafbfcfd-70d9-b016-0000-00178308cb32 is missing from the cache.
│ This may mean IntervalArithmetic [d1acc4aa-44c8-5952-acd4-ba5d80a2a253] does not support precompilation but is imported by a module that does.
└ @ Base loading.jl:2541
julia> x = interval(π)Interval{Float64}(3.141592653589793, 3.1415926535897936, com)
julia> Arb(x)[3.141592653589793 +/- 5.61e-16]
julia> y = Arb(π)[3.1415926535897932384626433832795028841971693993751058209749445923078164062862 +/- 1.93e-77]
julia> interval(y) # Defaults to BigFloatInterval{BigFloat}(3.141592653589793238462643383279502884197169399375105820974944592307816406286163, 3.141592653589793238462643383279502884197169399375105820974944592307816406286233, com)
julia> interval(Float64, y) # Other type can be explicitly specifiedInterval{Float64}(3.141592653589793, 3.1415926535897936, com)
julia> z = complex(interval(2), interval(1 // 3))Interval{Float64}(2.0, 2.0, com) + im*Interval{Float64}(0.3333333333333333, 0.33333333333333337, com)
julia> Acb(z)2.0 + [0.3333333333333333 +/- 7.04e-17]im
julia> w = Acb(2, 1 // 3)2.0 + [0.33333333333333333333333333333333333333333333333333333333333333333333333333333 +/- 4.78e-78]im
julia> interval(w)Interval{BigFloat}(2.0, 2.0, com) + im*Interval{BigFloat}(0.3333333333333333333333333333333333333333333333333333333333333333333333333333261, 0.3333333333333333333333333333333333333333333333333333333333333333333333333333348, com)
julia> interval(Float64, w)Interval{Float64}(2.0, 2.0, com) + im*Interval{Float64}(0.3333333333333333, 0.33333333333333337, com)

Linear algebra

To use the optimized linear algebra routines from FLINT, the matrices should be converted to ArbMatrix or AcbMatrix (depending on whether they are real or complex). Basic methods can then be used directly:

julia> using IntervalArithmetic, Arblib, LinearAlgebra
julia> A = interval.(BigFloat, [1 2; 3 4])2×2 Matrix{Interval{BigFloat}}:
 Interval{BigFloat}(1.0, 1.0, com)  Interval{BigFloat}(2.0, 2.0, com)
 Interval{BigFloat}(3.0, 3.0, com)  Interval{BigFloat}(4.0, 4.0, com)
julia> B = ArbMatrix(A)2×2 Arblib.ArbMatrix:
 1.0  2.0
 3.0  4.0
julia> B * B2×2 Arblib.ArbMatrix:
  7.0  10.0
 15.0  22.0
julia> inv(B)2×2 Arblib.ArbMatrix:
 [-2.00{…71 digits…}00 ± 1.44e-76]   [1.00{…71 digits…}00 ± 7.78e-77]
  [1.50{…72 digits…}00 ± 9.07e-77]  [-0.50{…72 digits…}00 ± 3.68e-77]
julia> B \ B2×2 Arblib.ArbMatrix:
     [1.00{…71 digits…}00 ± 7.20e-77]  [+/- 6.91e-77]
 [+/- 3.24e-77]                            [1.00{…72 digits…}00 ± 5.19e-77]
julia> eigvals(AcbMatrix(B)) # Eigenvalues are only supported for AcbMatrix2-element Arblib.AcbVector:
 [-0.372281323269014329925305734109464659110132228991396183849938735282950360729 +/- 5.90e-76] + [+/- 2.51e-76]im
  [5.372281323269014329925305734109464659110132228991396183849938735282950360729 +/- 8.90e-76] + [+/- 2.51e-76]im

Full documentation about supported methods can be found in the FLINT documentation for arb_mat and acb_mat. Note that many of these methods are not wrapped in Arblib, but most of them can be used through a low level wrapper.

Special functions

Arblib wraps a large number of the special functions from SpecialFunctions.jl.

julia> using IntervalArithmetic, Arblib, SpecialFunctions
julia> x = interval(BigFloat, 2)Interval{BigFloat}(2.0, 2.0, com)
julia> interval(besselj0(Arb(x))) # Convert to Arb and then backInterval{BigFloat}(0.2238907791412356680518274546499486258251544822186076031283497060108539577680058, 0.2238907791412356680518274546499486258251544822186076031283497060108539577680145, com)
julia> interval(gamma(Arb(x)))Interval{BigFloat}(1.0, 1.0, com)
julia> interval(zeta(Arb(x)))Interval{BigFloat}(1.644934066848226436472415166646025189218949901206798437735558229370007470403185, 1.644934066848226436472415166646025189218949901206798437735558229370007470403219, com)

Note that FLINT implements several special functions that are not in SpecialFunctions.jl, such as the confluent hypergeometric functions and the Lerch transcendent; they can be used through the low-level wrapper:

julia> using IntervalArithmetic, Arblib
julia> z = interval(BigFloat, 2 + 3im)Interval{BigFloat}(2.0, 2.0, com) + im*Interval{BigFloat}(3.0, 3.0, com)
julia> s = interval(BigFloat, 2)Interval{BigFloat}(2.0, 2.0, com)
julia> a = interval(BigFloat, 4)Interval{BigFloat}(4.0, 4.0, com)
julia> interval(Arblib.dirichlet_lerch_phi!(Acb(), Acb(z), Acb(s), Acb(a)))Interval{BigFloat}(-0.0002662647805245549465036569063000041468826679484039490012317600272567739198984341, -0.0002662647805245549465036569063000041468826679484039490012317600272567739198722642, com) + im*Interval{BigFloat}(0.03373697955708857151314153116051267066240962948450506837346469744178310287937996, 0.03373697955708857151314153116051267066240962948450506837346469744178310287944149, com)