Add documentation

docs/Project.toml

@@ -1,6 +1,4 @@
 [deps]
-LinearOperators = "5c8ed15e-5a4c-59e4-a42b-c7e8811fb125"
 Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4"
-
-[compat]
-Documenter = "~0.21"
+Krylov = "ba0b0d4f-ebba-5204-a429-3ac8c609bfb7"
+LinearOperators = "5c8ed15e-5a4c-59e4-a42b-c7e8811fb125"

docs/make.jl

@@ -4,12 +4,12 @@ makedocs(
   modules = [Krylov],
   doctest = true,
   linkcheck = true,
-  # strict = true,
-  assets = ["assets/style.css"],
-  format = Documenter.HTML(prettyurls = get(ENV, "CI", nothing) == "true"),
+  strict = true,
+  format = Documenter.HTML(assets = ["assets/style.css"], prettyurls = get(ENV, "CI", nothing) == "true"),
   sitename = "Krylov.jl",
   pages = ["Home" => "index.md",
            "API" => "api.md",
+           "Solvers" => "solvers.md",
            "Reference" => "reference.md",
           ]
 )

docs/src/api.md

@@ -1,23 +1,17 @@
-# API
-
-```@contents
-Pages = ["api.md"]
-```
-
 ## Types
 
 ```@docs
-KrylovStats
-SimpleStats
-LanczosStats
-SymmlqStats
+Krylov.KrylovStats
+Krylov.SimpleStats
+Krylov.LanczosStats
+Krylov.SymmlqStats
 ```
 
 ## Utilities
 
 ```@docs
-roots_quadratic
-sym_givens
-to_boundary
-vec2str
+Krylov.roots_quadratic
+Krylov.sym_givens
+Krylov.to_boundary
+Krylov.vec2str
 ```

docs/src/solvers.md

@@ -0,0 +1,28 @@
+All methods require that A is an `AbstractLinearOperator`.
+But a variant allows you to give A as an `AbstractMatrix`. Thereafter A is automatically wrapped in a `LinearOperator`.
+
+## Krylov methods
+
+```@docs
+cg
+cr
+symmlq
+cg_lanczos
+cg_lanczos_shift_seq
+minres
+minres_qlp
+diom
+dqgmres
+usymqr
+bilq
+cgs
+cgls
+crls
+cgne
+crmr
+lslq
+lsqr
+lsmr
+craig
+craigmr
+```

src/Krylov.jl

@@ -14,6 +14,7 @@ mutable struct SimpleStats{T} <: KrylovStats{T}
   status :: String
 end
 
+"Type for statistics returned by Lanczos solvers"
 mutable struct LanczosStats{T} <: KrylovStats{T}
   solved :: Bool
   residuals :: Array{T}
@@ -23,6 +24,7 @@ mutable struct LanczosStats{T} <: KrylovStats{T}
   status :: String
 end
 
+"Type for statistics returned by SYMMLQ"
 mutable struct SymmlqStats{T} <: KrylovStats{T}
   solved :: Bool
   residuals :: Array{T}

src/cg_lanczos.jl

@@ -15,7 +15,7 @@ export cg_lanczos, cg_lanczos_shift_seq
 """The Lanczos version of the conjugate gradient method to solve the
 symmetric linear system
 
-  Ax = b
+    Ax = b
 
 The method does _not_ abort if A is not definite.
 
@@ -118,14 +118,13 @@ end
 """The Lanczos version of the conjugate gradient method to solve a family
 of shifted systems
 
-  (A + αI) x = b  (α = α₁, α₂, ...)
+    (A + αI) x = b  (α = α₁, ..., αₙ)
 
 The method does _not_ abort if A + αI is not definite.
 
 A preconditioner M may be provided in the form of a linear operator and is
 assumed to be symmetric and positive definite.
 """
-
 function cg_lanczos_shift_seq(A :: AbstractLinearOperator, b :: AbstractVector{Tb},
                               shifts :: AbstractVector{Ts}; M :: AbstractLinearOperator=opEye(),
                               atol :: Tb=√eps(Tb), rtol :: Tb=√eps(Tb), itmax :: Int=0,

src/cgls.jl

@@ -5,7 +5,7 @@
 #
 # equivalently, of the normal equations
 #
-#  A'Ax = A'b.
+#  AᵀAx = Aᵀb.
 #
 # CGLS is formally equivalent to applying the conjugate gradient method
 # to the normal equations but should be more stable. It is also formally
@@ -24,16 +24,16 @@ export cgls
 
 """Solve the regularized linear least-squares problem
 
-  minimize ‖b - Ax‖₂² + λ ‖x‖₂²
+    minimize ‖b - Ax‖₂² + λ ‖x‖₂²
 
 using the Conjugate Gradient (CG) method, where λ ≥ 0 is a regularization
 parameter. This method is equivalent to applying CG to the normal equations
 
-  (A'A + λI) x = A'b
+    (AᵀA + λI) x = Aᵀb
 
 but is more stable.
 
-CGLS produces monotonic residuals ‖r‖₂ but not optimality residuals ‖A'r‖₂.
+CGLS produces monotonic residuals ‖r‖₂ but not optimality residuals ‖Aᵀr‖₂.
 It is formally equivalent to LSQR, though can be slightly less accurate,
 but simpler to implement.
 """

src/cgne.jl

@@ -10,7 +10,7 @@
 # and is equivalent to applying the conjugate gradient method
 # to the linear system
 #
-#  AA'y = b.
+#  AAᵀy = b.
 #
 # This method is also known as Craig's method, CGME, and other
 # names, and is described in
@@ -31,21 +31,21 @@ export cgne
 
 """Solve the consistent linear system
 
-  Ax + √λs = b
+    Ax + √λs = b
 
 using the Conjugate Gradient (CG) method, where λ ≥ 0 is a regularization
 parameter. This method is equivalent to applying CG to the normal equations
 of the second kind
 
-  (AA' + λI) y = b
+    (AAᵀ + λI) y = b
 
 but is more stable. When λ = 0, this method solves the minimum-norm problem
 
   min ‖x‖₂  s.t. Ax = b.
 
 When λ > 0, it solves the problem
 
-  min ‖(x,s)‖₂  s.t. Ax + √λs = b.
+    min ‖(x,s)‖₂  s.t. Ax + √λs = b.
 
 CGNE produces monotonic errors ‖x-x*‖₂ but not residuals ‖r‖₂.
 It is formally equivalent to CRAIG, though can be slightly less accurate,

src/craig.jl

@@ -35,7 +35,7 @@ export craig
 
 """Find the least-norm solution of the consistent linear system
 
-  Ax + √λs = b
+    Ax + √λs = b
 
 using the Golub-Kahan implementation of Craig's method, where λ ≥ 0 is a
 regularization parameter. This method is equivalent to CGNE but is more
@@ -45,14 +45,14 @@ For a system in the form Ax = b, Craig's method is equivalent to applying
 CG to AAᵀy = b and recovering x = Aᵀy. Note that y are the Lagrange
 multipliers of the least-norm problem
 
-  minimize ‖x‖  subject to Ax = b.
+    minimize ‖x‖  subject to Ax = b.
 
 Preconditioners M⁻¹ and N⁻¹ may be provided in the form of linear operators and are
 assumed to be symmetric and positive definite.
 Afterward CRAIG solves the symmetric and quasi-definite system
 
-  [ -N   Aᵀ ] [ x ]   [ 0 ]
-  [  A   M  ] [ y ] = [ b ],
+    [ -N   Aᵀ ] [ x ]   [ 0 ]
+    [  A   M  ] [ y ] = [ b ],
 
 which is equivalent to applying CG to (M + AN⁻¹Aᵀ)y = b.
 

src/craigmr.jl

@@ -30,28 +30,28 @@ export craigmr
 
 """Solve the consistent linear system
 
-  Ax + √λs = b
+    Ax + √λs = b
 
 using the CRAIG-MR method, where λ ≥ 0 is a regularization parameter.
 This method is equivalent to applying the Conjugate Residuals method
 to the normal equations of the second kind
 
-  (AAᵀ + λI) y = b
+    (AAᵀ + λI) y = b
 
 but is more stable. When λ = 0, this method solves the minimum-norm problem
 
-  min ‖x‖₂  s.t.  x ∈ argmin ‖Ax - b‖₂.
+    min ‖x‖₂  s.t.  x ∈ argmin ‖Ax - b‖₂.
 
 When λ > 0, this method solves the problem
 
-  min ‖(x,s)‖₂  s.t. Ax + √λs = b.
+    min ‖(x,s)‖₂  s.t. Ax + √λs = b.
 
 Preconditioners M⁻¹ and N⁻¹ may be provided in the form of linear operators and are
 assumed to be symmetric and positive definite.
 Afterward CRAIGMR solves the symmetric and quasi-definite system
 
-  [ -N   Aᵀ ] [ x ]   [ 0 ]
-  [  A   M  ] [ y ] = [ b ],
+    [ -N   Aᵀ ] [ x ]   [ 0 ]
+    [  A   M  ] [ y ] = [ b ],
 
 which is equivalent to applying MINRES to (M + AN⁻¹Aᵀ)y = b.
 

src/crls.jl

@@ -23,12 +23,12 @@ export crls
 
 """Solve the linear least-squares problem
 
-  minimize ‖b - Ax‖₂² + λ ‖x‖₂²
+    minimize ‖b - Ax‖₂² + λ ‖x‖₂²
 
 using the Conjugate Residuals (CR) method. This method is equivalent to
 applying MINRES to the normal equations
 
-  (AᵀA + λI) x = Aᵀb.
+    (AᵀA + λI) x = Aᵀb.
 
 This implementation recurs the residual r := b - Ax.
 

src/crmr.jl

@@ -30,21 +30,21 @@ export crmr
 
 """Solve the consistent linear system
 
-  Ax + √λs = b
+    Ax + √λs = b
 
 using the Conjugate Residual (CR) method, where λ ≥ 0 is a regularization
 parameter. This method is equivalent to applying CR to the normal equations
 of the second kind
 
-  (AAᵀ + λI) y = b
+    (AAᵀ + λI) y = b
 
 but is more stable. When λ = 0, this method solves the minimum-norm problem
 
-  min ‖x‖₂  s.t.  x ∈ argmin ‖Ax - b‖₂.
+    min ‖x‖₂  s.t.  x ∈ argmin ‖Ax - b‖₂.
 
 When λ > 0, this method solves the problem
 
-  min ‖(x,s)‖₂  s.t. Ax + √λs = b.
+    min ‖(x,s)‖₂  s.t. Ax + √λs = b.
 
 CGMR produces monotonic residuals ‖r‖₂.
 It is formally equivalent to CRAIG-MR, though can be slightly less accurate,

src/krylov_utils.jl

@@ -1,6 +1,9 @@
-# Display an array in the form
-# [ -3.0e-01 -5.1e-01  1.9e-01 ... -2.3e-01 -4.4e-01  2.4e-01 ]
-# with ndisp/2 elements on each side.
+"""Display an array in the form
+
+    [ -3.0e-01 -5.1e-01  1.9e-01 ... -2.3e-01 -4.4e-01  2.4e-01 ]
+
+with (ndisp - 1)/2 elements on each side.
+"""
 function vec2str(x :: Array{T,1}; ndisp :: Int=7) where T <: Number
   n = length(x);
   if n <= ndisp

src/lsmr.jl

@@ -27,12 +27,12 @@ export lsmr
 
 """Solve the regularized linear least-squares problem
 
-  minimize ‖b - Ax‖₂² + λ² ‖x‖₂²
+    minimize ‖b - Ax‖₂² + λ² ‖x‖₂²
 
 using the LSMR method, where λ ≥ 0 is a regularization parameter.
 LSQR is formally equivalent to applying MINRES to the normal equations
 
-  (AᵀA + λ² I) x = Aᵀb
+    (AᵀA + λ² I) x = Aᵀb
 
 (and therefore to CRLS) but is more stable.
 
@@ -43,16 +43,16 @@ Preconditioners M and N may be provided in the form of linear operators and are
 assumed to be symmetric and positive definite. If `sqd` is set to `true`,
 we solve the symmetric and quasi-definite system
 
-  [ E   Aᵀ ] [ r ]   [ b ]
-  [ A  -F  ] [ x ] = [ 0 ],
+    [ E   Aᵀ ] [ r ]   [ b ]
+    [ A  -F  ] [ x ] = [ 0 ],
 
 where E = M⁻¹  and F = N⁻¹.
 
 If `sqd` is set to `false` (the default), we solve the symmetric and
 indefinite system
 
-  [ E   Aᵀ ] [ r ]   [ b ]
-  [ A   0  ] [ x ] = [ 0 ].
+    [ E   Aᵀ ] [ r ]   [ b ]
+    [ A   0  ] [ x ] = [ 0 ].
 
 In this case, `N` can still be specified and indicates the norm
 in which `x` should be measured.