Documentation

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"]

docs/make.jl

@@ -2,6 +2,7 @@ using Documenter, Arblib
 
 makedocs(
     sitename = "Arblib",
+    modules = [Arblib],
     pages = [
         "index.md",
         "Low level wrapper" => [
@@ -18,4 +19,5 @@ makedocs(
         ],
         "Rigorous numerics" => "rigorous.md",
     ],
+    warnonly = [:missing_docs, :cross_references],
 )

docs/src/interface-series.md

@@ -16,9 +16,7 @@ the length of the expansion.
 ```@repl 1
 using Arblib
 
-setprecision(Arb, 64)
-
-x0 = Arb(1 // 3)
+x0 = Arb(1 // 3, prec = 64)
 x = ArbSeries((x0, 1), degree = 5)
 
 sin(x)

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.

src/ArbCall/ArbFunction.jl

@@ -16,6 +16,8 @@ declaration in C. Note that note all types of function declarations
 are supported, only those which are relevant for the Arb library.
 
 ```jldoctest
+julia> using Arblib
+
 julia> Arblib.ArbCall.ArbFunction("void arb_zero(arb_t x)")
 Arblib.ArbCall.ArbFunction{Nothing}("arb_zero", Arblib.ArbCall.Carg[Arblib.ArbCall.Carg{Arb}(:x, false)])
 ```

src/ArbCall/Carg.jl

@@ -7,6 +7,8 @@ corresponding Julia type is `T`, `name` is the name of the argument,
 `isconst` is true if the argument is declared as a const.
 
 ```jldoctest
+julia> using Arblib
+
 julia> Arblib.ArbCall.Carg("const arb_t x")
 Arblib.ArbCall.Carg{Arb}(:x, true)
 ```
@@ -48,7 +50,7 @@ end
     jltype(ca::Carg{T})
 
 The most general Julia type for which we allow automatic conversion to
-the [`ctype`](@ref) of `ca`.
+the [`Arblib.ArbCall.ctype`](@ref) of `ca`.
 
 These conversations should be done without any loss of information,
 for example for floating point numbers we only allow conversion from

src/ArbCall/parse.jl

@@ -73,9 +73,9 @@ end
     generate_file(filename, title, sections; verbose, manual_overrides)
 
 Given a title and a list of sections, as returned by
-[`parse_arbdoc`](@ref), return a string with the Julia code to load
-all of those functions. If a file name is given then also write the
-string to that file.
+[`Arblib.ArbCall.parse_arbdoc`](@ref), return a string with the Julia
+code to load all of those functions. If a file name is given then also
+write the string to that file.
 
 For each function it checks that it is correctly parsed. Those that do
 not parse due to types which are not supported and likely never will

src/calc_integrate.jl

@@ -158,14 +158,18 @@ particular, there are three differences from `integrate`:
 
 # Examples
 ```jldoctest
+julia> using Arblib
+
 julia> Arblib.integrate!(Arblib.sin!, Acb(0), Acb(0), Acb(10)) # Integrate sin from 0 to 10
 [1.83907152907645245225886394782406483451993016513316854683595373104879258687 +/- 5.15e-75]
 
 julia> Arblib.integrate!(Arblib.inv!, Acb(0), Acb(1, -5), Acb(1, 5)) # Integrate 1/z from 1 - 5i to 1 + 5i
 [+/- 2.02e-75] + [2.74680153389003172172254385288992229730199919179940161793956671182574846633 +/- 2.83e-75]im
 
 julia> # Integrate √z from 1 - 5im to 1 + 5im, taking into account the branch cut at (-∞, 0]
+
 julia> f! = (res, z; analytic = false) -> Arblib.sqrt_analytic!(res, z, analytic);
+
 julia> Arblib.integrate!(f!, Acb(0), Acb(1, -5), Acb(1, 10), check_analytic = true, prec = 64)
 [-9.0064084416559764 +/- 6.53e-17] + [23.8636067095598007 +/- 6.98e-17]im
 ```
@@ -302,20 +306,24 @@ See also: [`integrate!`](@ref).
 
 # Examples
 ```jldoctest
+julia> using Arblib
+
 julia> Arblib.integrate(sin, 0, 10) # Integrate sin from 0 to 10
 [1.83907152907645245225886394782406483451993016513316854683595373104879258687 +/- 5.15e-75]
 
 julia> Arblib.integrate(z -> 1/z, Acb(1, -5), Acb(1, 5)) # Integrate 1/z from 1 - 5i to 1 + 5i
 [+/- 2.02e-75] + [2.74680153389003172172254385288992229730199919179940161793956671182574846633 +/- 2.83e-75]im
 
 julia> # Integrate √z from 1 - 5im to 1 + 5im, taking into account the branch cut at (-∞, 0]
+
 julia> f = (z; analytic = false) -> begin
            if analytic && Arblib.contains_nonpositive(real(z))
                return Acb(NaN, prec = precision(z))
            else
                return sqrt(z)
            end
        end;
+
 julia> Arblib.integrate(f, Acb(1, -5), Acb(1, 10), check_analytic = true, prec = 64)
 [-9.0064084416559764 +/- 7.40e-17] + [23.8636067095598007 +/- 9.03e-17]im
 ```

test/runtests.jl

@@ -1,8 +1,11 @@
 ENV["NEMO_THREADED"] = 1
 
 using Arblib, Test, LinearAlgebra, Random, Serialization, SpecialFunctions
+using Documenter
 
 @testset "Arblib" begin
+    doctest(Arblib)
+
     include("ArbCall/runtests.jl")
 
     include("arb_types.jl")