Merge pull request #159 from kalmarek/documentation

Documentation

.github/workflows/documentation.yml

@@ -0,0 +1,23 @@
+name: Documentation
+on:
+  push:
+    branches:
+      - master
+    tags: '*'
+  pull_request:
+
+jobs:
+  build:
+    runs-on: ubuntu-latest
+    steps:
+      - uses: actions/checkout@v3
+      - uses: julia-actions/setup-julia@latest
+        with:
+          version: '1.6'
+      - name: Install dependencies
+        run: julia --project=docs/ -e 'using Pkg; Pkg.develop(PackageSpec(path=pwd())); Pkg.instantiate()'
+      - name: Build and deploy
+        env:
+          GITHUB_TOKEN: ${{ secrets.GITHUB_TOKEN }} # If authenticating with GitHub Actions token
+          DOCUMENTER_KEY: ${{ secrets.DOCUMENTER_KEY }} # If authenticating with SSH deploy key
+        run: julia --project=docs/ docs/make.jl

.gitignore

@@ -2,4 +2,5 @@ deps/usr
 deps/build_*
 deps/build.log
 deps/deps.jl
+docs/build
 Manifest.toml

Project.toml

@@ -10,10 +10,16 @@ LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
 Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"
 Serialization = "9e88b42a-f829-5b0c-bbe9-9e923198166b"
 SpecialFunctions = "276daf66-3868-5448-9aa4-cd146d93841b"
+
+[extras]
+Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4"
 Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
 
 [compat]
 Arb_jll = "~200.2300"
 FLINT_jll = "~200.900.000"
 SpecialFunctions = "1.0, 2"
 julia = "1.6"
+
+[targets]
+test = ["Documenter", "Test"]

README.md

@@ -1,5 +1,8 @@
 # Arblib.jl
 
+[![][docs-stable-img]][docs-stable-url]
+[![][docs-dev-img]][docs-dev-url]
+
 This package is a thin, efficient wrapper around [Arb](http://arblib.org) - a C library for arbitrary-precision ball arithmetic.
 
 The package is currently in early development. More features and
@@ -254,3 +257,8 @@ Enabling a threaded version of flint can be done by setting the
 environment variable `NEMO_THREADED=1`. Note that this should be
 set before `Arblib.jl` is loaded. To set the actual number of threads,
 use `Arblib.flint_set_num_threads($numberofthreads)`.
+
+[docs-dev-img]: https://img.shields.io/badge/docs-dev-blue.svg
+[docs-dev-url]: https://Kalmarek.github.io/Arblib.jl/dev/
+[docs-stable-img]: https://img.shields.io/badge/docs-stable-blue.svg
+[docs-stable-url]: https://Kalmarek.github.io/Arblib.jl/stable

docs/Project.toml

@@ -0,0 +1,7 @@
+[deps]
+BenchmarkTools = "6e4b80f9-dd63-53aa-95a3-0cdb28fa8baf"
+Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4"
+
+[compat]
+BenchmarkTools = "1"
+Documenter = "1"

docs/make.jl

@@ -0,0 +1,27 @@
+using Documenter, Arblib
+
+DocMeta.setdocmeta!(Arblib, :DocTestSetup, :(using Arblib); recursive = true)
+
+makedocs(
+    sitename = "Arblib.jl",
+    modules = [Arblib],
+    pages = [
+        "index.md",
+        "Low level wrapper" => [
+            "Types" => "wrapper-types.md",
+            "Methods" => "wrapper-methods.md",
+            "Floating point wrapper" => "wrapper-fpwrap.md",
+        ],
+        "High level interface" => [
+            "Types" => "interface-types.md",
+            "Ball methods" => "interface-ball.md",
+            "Integration" => "interface-integration.md",
+            "Series" => "interface-series.md",
+            "Mutable arithmetic" => "interface-mutable.md",
+        ],
+        "Rigorous numerics" => "rigorous.md",
+    ],
+    warnonly = [:missing_docs],
+)
+
+deploydocs(repo = "github.com/kalmarek/Arblib.jl")

docs/src/index.md

@@ -0,0 +1,20 @@
+# Arblib.jl Documentation
+
+This package is a wrapper around [Arb](http://arblib.org) - a C
+library for arbitrary-precision ball arithmetic. Other wrappers of Arb
+for Julia include [Nemo](https://github.com/Nemocas/Nemo.jl) and
+[ArbNumerics.jl](https://github.com/JeffreySarnoff/ArbNumerics.jl).
+
+The **goal** of Arblib.jl is to supply a **low lever wrapper** of the
+methods in Arb as well as a **high level interface**. The low level
+wrapper should allow for writing methods using mutability and with
+performance very close to that of those written in C. The high level
+interface should make it easy to use in generic Julia code, similarly
+to how `BigFloat` is a wrapper around the MPFR library. In addition it
+should be possible to seamlessly switch between the high level
+interface and the low level wrapper when needed.
+
+The above goals can be put into contrast with Nemo, whose high level
+interface is made for use in the
+[AbstractAlgebra.jl](https://github.com/Nemocas/AbstractAlgebra.jl)
+universe and not general Julia code.

docs/src/interface-ball.md

@@ -0,0 +1,37 @@
+# Ball methods
+The following methods are useful for explicitly dealing with the ball
+representation of `Arb` and related values.
+
+## Construction
+For constructing balls the methods below are useful. Note that there
+is no `setinterval` method, this is instead accomplished with `Arb((a,
+b))` for constructing a ball containing the interval ``[a, b]``.
+
+``` @docs
+setball
+add_error
+```
+
+## Destruction
+For extracting information about the ball representation the following
+methods are useful.
+
+``` @docs
+radius
+midpoint
+lbound
+ubound
+abs_lbound
+abs_ubound
+getinterval
+getball
+```
+
+## Union and intersection
+The `Base.union` and `Base.intersect` methods are overloaded to
+compute the union and intersection of balls.
+
+``` @docs
+union(::Arb, ::Arb)
+intersect(::Arb, ::Arb)
+```

docs/src/interface-integration.md

@@ -0,0 +1,4 @@
+``` @docs
+Arblib.integrate
+Arblib.integrate!
+```

docs/src/interface-mutable.md

@@ -0,0 +1,85 @@
+# Mutable arithmetic
+The high level interface can be combined with the low level wrapper to
+allow for efficient computations using mutable arithmetic.
+
+In the future it would be nice to have an interface to
+[MutableArithmetics.jl](https://github.com/jump-dev/MutableArithmetics.jl),
+see [#118](https://github.com/kalmarek/Arblib.jl/issues/118).
+
+The following methods are useful for mutating part of a value
+
+``` @docs
+Arblib.radref
+Arblib.midref
+Arblib.realref
+Arblib.imagref
+Arblib.ref
+```
+
+## Examples
+
+Compare computing ``\sqrt{x^2 + y^2}`` using mutable arithmetic with
+the default.
+
+``` @repl
+using Arblib, BenchmarkTools
+
+x = Arb(1 // 3)
+y = Arb(1 // 5)
+res = zero(x)
+
+f(x, y) = sqrt(x^2 + y^2)
+f!(res, x, y) = begin
+    Arblib.sqr!(res, x)
+    Arblib.fma!(res, res, y, y)
+    return Arblib.sqrt!(res, res)
+end
+
+@benchmark f($x, $y) samples=10000 evals=500
+@benchmark f!($res, $x, $y) samples=10000 evals=500
+```
+
+Set the radius of the real part of an `Acb`.
+
+``` @repl
+using Arblib
+
+z = Acb(1, 2)
+Arblib.set!(Arblib.radref(Arblib.realref(z)), 1e-10)
+z
+```
+
+Compare a naive implementation of polynomial evaluation using
+mutable arithmetic with one not using using it.
+
+``` @repl
+using Arblib, BenchmarkTools
+
+p = ArbPoly(1:10)
+x = Arb(1 // 3)
+res = zero(x)
+
+function eval(p, x)
+    res = zero(x)
+    xi = one(x)
+    for i in 0:Arblib.degree(p)
+        res += p[i] * xi
+        xi *= x
+    end
+    return res
+end
+
+function eval!(res, p, x)
+    Arblib.zero!(res)
+    xi = one(x)
+    for i in 0:Arblib.degree(p)
+        Arblib.addmul!(res, Arblib.ref(p, i), xi)
+        Arblib.mul!(xi, xi, x)
+    end
+    return res
+end
+
+@benchmark eval($p, $x) samples = 10000 evals = 30
+@benchmark eval!($res, $p, $x) samples = 10000 evals = 30
+@benchmark $p($x) samples = 10000 evals = 30 # Arb implementation for reference
+```

docs/src/interface-series.md

@@ -0,0 +1,25 @@
+# Series
+Taylor series arithmetic allows for the computation of truncated
+Taylor series of functions and is a form of higher order automatic
+differentiation. See e.g.
+[`TaylorSeries.jl`](https://github.com/JuliaDiff/TaylorSeries.jl) and
+[`TaylorDiff.jl`](https://github.com/JuliaDiff/TaylorDiff.jl) for
+implementations of Taylor series in Julia.
+
+The Arb library has good support for computing with polynomials as
+Taylor expansions. The types [`ArbSeries`](@ref) and
+[`AcbSeries`](@ref) are intended to make this easy to use from Julia.
+They are given by an [`ArbPoly`](@ref)/[`AcbPoly`](@ref) together with
+the length of the expansion.
+
+## Example
+```@repl 1
+using Arblib
+
+x0 = Arb(1 // 3, prec = 64)
+x = ArbSeries((x0, 1), degree = 5)
+
+sin(x)
+
+sin(x)^2 + cos(x)^2
+```

docs/src/interface-types.md

@@ -0,0 +1,61 @@
+# Types
+
+The package defines a number of types for the high level interface.
+
+## Basic
+These types directly map to corresponding Arb types.
+
+``` @docs
+Mag
+Arf
+Arb
+Acb
+ArbVector
+AcbVector
+ArbPoly
+AcbPoly
+ArbMatrix
+AcbMatrix
+```
+
+## Series
+The package defines two series types, which are wrapper for the
+polynomial types with a specified degree.
+
+``` @docs
+ArbSeries
+AcbSeries
+```
+
+## Ref
+In addition to these there are a number of `Ref` types, which allow
+for non-allocating access in a number of cases.
+
+``` @docs
+MagRef
+ArfRef
+ArbRef
+AcbRef
+ArbRefVector
+AcbRefVector
+ArbRefMatrix
+AcbRefMatrix
+```
+
+## Correspondence between types
+We have the following table for the correspondence with between the
+[Low level wrapper types](wrapper-types.md) and the high level
+interface types.
+
+| Arb      | Wrapper          | High level  | Ref            |
+|----------|------------------|-------------|----------------|
+| `mag_t`  | `mag_struct`     | `Mag`       | `MagRef`       |
+| `arf_t`  | `arf_struct`     | `Arf`       | `ArfRef`       |
+| `arb_t`  | `arb_struct`     | `Arb`       | `ArbRef`       |
+| `acb_t`  | `acb_struct`     | `Acb`       | `AcbRef`       |
+| `arb_t*` | `arb_vec_struct` | `ArbVector` | `ArbRefVector` |
+| `acb_t*` | `acb_vec_struct` | `AcbVector` | `AcbRefVector` |
+| `arb_poly_t` | `arb_poly_struct` | `ArbPoly` or `ArbSeries` | |
+| `acb_poly_t` | `acb_poly_struct` | `AcbPoly` or `AcbSeries` | |
+| `arb_mat_t` | `arb_mat_struct` | `ArbMatrix` | `ArbRefMatrix` |
+| `acb_mat_t` | `acb_mat_struct` | `AcbMatrix` | `AcbRefMatrix` |

docs/src/rigorous.md

@@ -0,0 +1,48 @@
+# Rigorous numerics
+Arb is made for rigorous numerics and any functions which do not
+produce rigorous results are clearly marked as such. This is not the
+case with Julia in general and you therefore have to be careful when
+interacting with the ecosystem if you want your results to be
+completely rigorous. Below we discuss some things to be extra careful
+with.
+
+## Implicit promotion
+Julia automatically promotes types in many cases and in particular you
+have to watch out for temporary non-rigorous values. For example
+`2(π * Arb(1 // 3))` is okay, but not `2π * Arb(1 // 3)`
+
+``` repl
+x = 2(π * Arb(1 // 3))
+y = 2π * Arb(1 // 3)
+Arblib.overlaps(x, y)
+```
+
+## Non-rigorous algorithms
+Standard numerical algorithms typically return (hopefully good)
+approximations. These algorithms can then not directly be used in
+rigorous numerical computations unless the error can be bounded.
+
+For example Julias built in methods for solving linear systems doesn't
+produce rigorous results. Instead you would have to use the solves
+provided by Arb, such as `Arblib.solve!`.
+
+Other examples would include integration and solving of differential
+equations.
+
+## Implementation details
+In some cases the implementation in Julia implicitly makes certain
+assumptions to improve performance and this can lead to issues.
+
+For example, prior to Julia version 1.8 the `minimum` and `maximum`
+methods in Julia checked for `NaN` results (on which is short fuses)
+using `x == x`, which works for most numerical types but not for `Arb`
+(`x == x` is only true if the radius is zero). See
+<https://github.com/JuliaLang/julia/issues/36287> and in particular
+<https://github.com/JuliaLang/julia/issues/45932> for more details.
+Since Julia version 1.8 the `minimum` and `maximum` methods work
+correctly for `Arb`, for earlier versions of Julia it only works
+correctly in some cases.
+
+These types of problems are the hardest to find since they are not
+clear from the documentation but you have to read the implementation,
+`@which` and `@less` are your friends in these cases.