Merge pull request #168 from ArnoStrouwen/formatter

reapply formatting

.JuliaFormatter.toml

@@ -1,2 +1,3 @@
 style = "sciml"
-format_markdown = true
+format_markdown = true
+format_docstrings = true

docs/pages.jl

@@ -4,5 +4,5 @@ pages = [
     "Home" => "index.md",
     "matrix_exponentials.md",
     "expv.md",
-    "arnoldi.md",
+    "arnoldi.md"
 ]

src/ExponentialUtilities.jl

@@ -38,11 +38,11 @@ include("krylov_phiv_error_estimate.jl")
 include("precompile.jl")
 
 export phi, phi!, KrylovSubspace, arnoldi, arnoldi!, lanczos!, ExpvCache, PhivCache,
-    expv, expv!, phiv, phiv!, kiops, expv_timestep, expv_timestep!, phiv_timestep,
-    phiv_timestep!,
-    StegrCache, get_subspace_cache, exponential!
+       expv, expv!, phiv, phiv!, kiops, expv_timestep, expv_timestep!, phiv_timestep,
+       phiv_timestep!,
+       StegrCache, get_subspace_cache, exponential!
 export ExpMethodHigham2005,
-    ExpMethodHigham2005Base, ExpMethodGeneric, ExpMethodNative,
-    ExpMethodDiagonalization
+       ExpMethodHigham2005Base, ExpMethodGeneric, ExpMethodNative,
+       ExpMethodDiagonalization
 
 end

src/StegrWork.jl

@@ -73,7 +73,7 @@ Allocate work arrays for diagonalization of real-symmetric tridiagonal
 matrices of sizes up to `n`×`n`.
 """
 function StegrWork(::Type{T}, n::Integer,
-    jobz::Char = 'V', range::Char = 'A') where {T}
+        jobz::Char = 'V', range::Char = 'A') where {T}
     n = convert(BlasInt, n)
     dv = Array{T}(undef, n)
     ev = Array{T}(undef, n)

src/arnoldi.jl

@@ -35,7 +35,7 @@ mutable struct KrylovSubspace{T, U, B, VType <: AbstractMatrix{T},
 end
 
 function KrylovSubspace{T, U}(n::Integer, maxiter::Integer = 30,
-    augmented::Integer = false) where {T, U}
+        augmented::Integer = false) where {T, U}
     V = Matrix{T}(undef, n + augmented, maxiter + 1)
     H = fill(zero(U), maxiter + 1, maxiter + !iszero(augmented))
     return KrylovSubspace{T, U, real(T), Matrix{T}, Matrix{U}}(maxiter, maxiter, augmented,
@@ -101,7 +101,7 @@ the dimension of `Ks` is smaller than `m`.
     Springer, Cham.
 """
 function arnoldi(A, b; m = min(30, size(A, 1)), ishermitian = LinearAlgebra.ishermitian(A),
-    kwargs...)
+        kwargs...)
     TA, Tb = eltype(A), eltype(b)
     T = promote_type(TA, Tb)
     n = length(b)
@@ -112,7 +112,8 @@ function arnoldi(A, b; m = min(30, size(A, 1)), ishermitian = LinearAlgebra.ishe
     # H is a small dense array in tridiagonal form
     H = zeros(U, (m + 1, m))
 
-    Ks = KrylovSubspace{T, U, real(T), typeof(V), typeof(H)}(m, m, false, zero(real(T)), false, V, H)
+    Ks = KrylovSubspace{T, U, real(T), typeof(V), typeof(H)}(
+        m, m, false, zero(real(T)), false, V, H)
 
     arnoldi!(Ks, A, b; m = m, ishermitian = ishermitian, kwargs...)
 end
@@ -211,8 +212,8 @@ end
 Take the `j` 'th step of the Lanczos iteration.
 """
 function arnoldi_step!(j::Integer, iop::Integer, A::AT,
-    V::AbstractMatrix{T}, H::AbstractMatrix{U},
-    n::Int = -1, p::Int = -1) where {AT, T, U}
+        V::AbstractMatrix{T}, H::AbstractMatrix{U},
+        n::Int = -1, p::Int = -1) where {AT, T, U}
     x, y = @view(V[:, j]), @view(V[:, j + 1])
     applyA!(y, A, x, V, j, n, p)
 
@@ -235,12 +236,11 @@ end
 Non-allocating version of `arnoldi`.
 """
 function arnoldi!(Ks::KrylovSubspace{T1, U}, A::AT, b;
-    tol::Real = 1e-7, m::Int = min(Ks.maxiter, size(A, 1)),
-    ishermitian::Bool = LinearAlgebra.ishermitian(A isa Tuple ? first(A) : A),
-    opnorm = nothing, iop::Int = 0,
-    init::Int = 0, t::Number = NaN, mu::Number = NaN,
-    l::Int = -1) where {T1 <: Number, U <: Number, AT}
-
+        tol::Real = 1e-7, m::Int = min(Ks.maxiter, size(A, 1)),
+        ishermitian::Bool = LinearAlgebra.ishermitian(A isa Tuple ? first(A) : A),
+        opnorm = nothing, iop::Int = 0,
+        init::Int = 0, t::Number = NaN, mu::Number = NaN,
+        l::Int = -1) where {T1 <: Number, U <: Number, AT}
     Ks.wasbreakdown = false
 
     ishermitian &&
@@ -277,10 +277,10 @@ end
 Take the `j`'th step of the Lanczos iteration.
 """
 function lanczos_step!(j::Integer, A,
-    V::AbstractMatrix{T},
-    u::AbstractVector{U},
-    v::AbstractVector{B},
-    n::Int = -1, p::Int = -1) where {B, T, U}
+        V::AbstractMatrix{T},
+        u::AbstractVector{U},
+        v::AbstractVector{B},
+        n::Int = -1, p::Int = -1) where {B, T, U}
     x, y = @view(V[:, j]), @view(V[:, j + 1])
     applyA!(y, A, x, V, j, n, p)
     α = u[j] = coeff(U, dot(x, y))
@@ -316,11 +316,10 @@ A variation of `arnoldi!` that uses the Lanczos algorithm for
 Hermitian matrices.
 """
 function lanczos!(Ks::KrylovSubspace{T1, U, B}, A::AT, b;
-    tol = 1e-7, m = min(Ks.maxiter, size(A, 1)),
-    opnorm = nothing,
-    init::Int = 0, t::Number = NaN, mu::Number = NaN,
-    l::Int = -1) where {T1 <: Number, U <: Number, B, AT}
-
+        tol = 1e-7, m = min(Ks.maxiter, size(A, 1)),
+        opnorm = nothing,
+        init::Int = 0, t::Number = NaN, mu::Number = NaN,
+        l::Int = -1) where {T1 <: Number, U <: Number, B, AT}
     Ks.wasbreakdown = false
 
     m > Ks.maxiter ? resize!(Ks, m) : Ks.m = m # might change if happy-breakdown occurs

src/exp_baseexp.jl

@@ -5,7 +5,7 @@ The same as `ExpMethodHigham2005` but follows `Base.exp` closer.
 """
 struct ExpMethodHigham2005Base end
 function alloc_mem(A::StridedMatrix{T},
-    method::ExpMethodHigham2005Base) where {T <: LinearAlgebra.BlasFloat}
+        method::ExpMethodHigham2005Base) where {T <: LinearAlgebra.BlasFloat}
     n = LinearAlgebra.checksquare(A)
     A2 = Matrix{T}(undef, n, n)
     P = Matrix{T}(undef, n, n)
@@ -21,7 +21,7 @@ end
 ## Non-allocating version of `LinearAlgebra.exp!`. Modifies `A` to
 ## become (approximately) `exp(A)`.
 function exponential!(A::StridedMatrix{T}, method::ExpMethodHigham2005Base,
-    cache = alloc_mem(A, method)) where {T <: LinearAlgebra.BlasFloat}
+        cache = alloc_mem(A, method)) where {T <: LinearAlgebra.BlasFloat}
     X = A
     n = LinearAlgebra.checksquare(A)
     # if ishermitian(A)

src/exp_generic.jl

@@ -3,7 +3,7 @@ intlog2(x::T) where {T <: Integer} = T(8 * sizeof(T) - leading_zeros(x - one(T))
 intlog2(x) = x > typemax(UInt) ? ceil(Int, log2(x)) : intlog2(ceil(UInt, x)) % Int
 
 function naivemul!(C::StridedMatrix{T}, A::StridedMatrix{T},
-    B::StridedMatrix{T}) where {T <: LinearAlgebra.BlasFloat}
+        B::StridedMatrix{T}) where {T <: LinearAlgebra.BlasFloat}
     mul!(C, A, B)
 end
 function naivemul!(C, A, B)
@@ -37,7 +37,7 @@ _const(A) = A
 _const(A::Array) = Base.Experimental.Const(A)
 # Separated to make it easier to test.
 @generated function naivemul!(C::AbstractMatrix{T}, A, B, Maxis, Naxis, ::Val{MU},
-    ::Val{NU}) where {T, MU, NU}
+        ::Val{NU}) where {T, MU, NU}
     nrem_body = quote
         m = first(Maxis) - 1
         while m < M - $(MU - 1)
@@ -124,7 +124,7 @@ struct ExpMethodGeneric{T} end
 ExpMethodGeneric() = ExpMethodGeneric{Val(13)}();
 
 function exponential!(x, method::ExpMethodGeneric{Vk},
-    cache = alloc_mem(x, method)) where {Vk}
+        cache = alloc_mem(x, method)) where {Vk}
     nx = opnorm(x, 1)
     if isnan(nx) || nx > 4611686018427387904 # i.e. 2^62 since it would cause overflow in 2^s
         # This should (hopefully) ensure that the result is Inf or NaN depending on
@@ -197,7 +197,7 @@ function exp_generic_core!(y1, y2, y3, x, ::Val{13})
 end
 
 function exp_pade_p!(y1::AbstractMatrix{T}, y2::AbstractMatrix{T}, x::AbstractMatrix,
-    ::Val{13}, ::Val{13}) where {T}
+        ::Val{13}, ::Val{13}) where {T}
     N = size(x, 1) # check square is in `exp_generic`
     y1 .= x .* 1.5440497506703088e-17
     for c in (2.8101705462199623e-15, 2.529153491597966e-13, 1.48377004840414e-11,

src/kiops.jl

@@ -48,8 +48,8 @@ References:
   - Niesen, J. and Wright, W.M., 2012. Algorithm 919: A Krylov subspace algorithm for evaluating the ``φ``-functions appearing in exponential integrators. ACM Transactions on Mathematical Software (TOMS), 38(3), p.22
 """
 function kiops(tau_out, A, u; mmin::Int = 10, mmax::Int = 128, m::Int = min(mmin, mmax),
-    tol::Real = 1e-7, opnorm = LinearAlgebra.opnorm(A, Inf), iop::Int = 2,
-    ishermitian::Bool = LinearAlgebra.ishermitian(A), task1::Bool = false)
+        tol::Real = 1e-7, opnorm = LinearAlgebra.opnorm(A, Inf), iop::Int = 2,
+        ishermitian::Bool = LinearAlgebra.ishermitian(A), task1::Bool = false)
     n, ppo = size(u, 1), size(u, 2)
     p = ppo - 1
 
@@ -259,8 +259,8 @@ function kiops(tau_out, A, u; mmin::Int = 10, mmax::Int = 128, m::Int = min(mmin
 end
 
 Base.@propagate_inbounds function kiops_update_solution!(tau_now, tau, tau_out, w, l, V, F,
-    H, beta, j, n, step, numSteps,
-    reject, ireject)
+        H, beta, j, n, step, numSteps,
+        reject, ireject)
     # Yep, got the required tolerance update
     reject = reject + ireject
     step = step + 1

src/krylov_phiv_adaptive.jl

@@ -44,12 +44,12 @@ end
 Non-allocating version of `expv_timestep`.
 """
 function expv_timestep!(u::AbstractVector{T}, t::tType, A, b::AbstractVector{T};
-    kwargs...) where {T <: Number, tType <: Real}
+        kwargs...) where {T <: Number, tType <: Real}
     expv_timestep!(reshape(u, length(u), 1), [t], A, b; kwargs...)
     return u
 end
 function expv_timestep!(U::AbstractMatrix{T}, ts::Vector{tType}, A, b::AbstractVector{T};
-    kwargs...) where {T <: Number, tType <: Real}
+        kwargs...) where {T <: Number, tType <: Real}
     B = reshape(b, length(b), 1)
     phiv_timestep!(U, ts, A, B; kwargs...)
 end
@@ -95,19 +95,19 @@ end
 Non-allocating version of `phiv_timestep`.
 """
 function phiv_timestep!(u::AbstractVector{T}, t::tType, A, B::AbstractMatrix{T};
-    kwargs...) where {T <: Number, tType <: Real}
+        kwargs...) where {T <: Number, tType <: Real}
     phiv_timestep!(reshape(u, length(u), 1), [t], A, B; kwargs...)
     return u
 end
 function phiv_timestep!(U::AbstractMatrix{T}, ts::Vector{tType}, A, B::AbstractMatrix{T};
-    tau::Real = 0.0,
-    m::Int = min(10, size(A, 1)), tol::Real = 1e-7,
-    opnorm = LinearAlgebra.opnorm(A, Inf), iop::Int = 0,
-    correct::Bool = false, caches = nothing, adaptive = false,
-    delta::Real = 1.2,
-    ishermitian::Bool = LinearAlgebra.ishermitian(A),
-    gamma::Real = 0.8, NA::Int = 0,
-    verbose = false) where {T <: Number, tType <: Real}
+        tau::Real = 0.0,
+        m::Int = min(10, size(A, 1)), tol::Real = 1e-7,
+        opnorm = LinearAlgebra.opnorm(A, Inf), iop::Int = 0,
+        correct::Bool = false, caches = nothing, adaptive = false,
+        delta::Real = 1.2,
+        ishermitian::Bool = LinearAlgebra.ishermitian(A),
+        gamma::Real = 0.8, NA::Int = 0,
+        verbose = false) where {T <: Number, tType <: Real}
     # Choose initial timestep
     opnorm = opnorm isa Number ? opnorm : opnorm(A, Inf) # backward compatibility
     abstol = tol * opnorm
@@ -246,7 +246,7 @@ function phiv_timestep!(U::AbstractMatrix{T}, ts::Vector{tType}, A, B::AbstractM
 end
 # Helper functions for phiv_timestep!
 function _phiv_timestep_adapt(m, tau, epsilon, m_old, tau_old, epsilon_old, q, kappa,
-    gamma, omega, maxtau, n, p, NA, iop, Hnorm, verbose)
+        gamma, omega, maxtau, n, p, NA, iop, Hnorm, verbose)
     # Compute new m and tau (Algorithm 4)
     if tau_old > tau
         q = log(tau / tau_old) / log(epsilon / epsilon_old) - 1

src/krylov_phiv_error_estimate.jl

@@ -84,7 +84,7 @@ function expv!(w::AbstractVector{T}, t::Number, A, b::AbstractVector{T},
     @. V[:, 1] = b / Ks.beta
 
     ε = atol + rtol * Ks.beta
-    verbose && @printf("Initial norm: β₀ %e, stopping threshold: %e\n", Ks.beta,ε)
+    verbose && @printf("Initial norm: β₀ %e, stopping threshold: %e\n", Ks.beta, ε)
 
     α = @diagview(H)
     β = @diagview(H, -1)
@@ -99,7 +99,7 @@ function expv!(w::AbstractVector{T}, t::Number, A, b::AbstractVector{T},
         #   Saad, Y. (1992). Analysis of some Krylov subspace
         #   approximations. SIAM Journal on Numerical Analysis.
         σ = β[j] * Ks.beta * abs(cache.v[j])
-        verbose && @printf("iter %d, α[%d] %e, β[%d] %e, σ %e\n", j, j, α[j], j, β[j],σ)
+        verbose && @printf("iter %d, α[%d] %e, β[%d] %e, σ %e\n", j, j, α[j], j, β[j], σ)
         if σ < ε
             Ks.m = j
             break

src/phi.jl

@@ -18,7 +18,7 @@ package for computing matrix exponentials. ACM Transactions on Mathematical
 Software (TOMS), 24(1), 130-156. Theorem 1).
 """
 function phi(z::T, k::Integer; cache = nothing,
-    expmethod = ExpMethodHigham2005Base()) where {T <: Number}
+        expmethod = ExpMethodHigham2005Base()) where {T <: Number}
     # Construct the matrix
     if cache == nothing
         cache = fill(zero(T), k + 1, k + 1)
@@ -59,9 +59,9 @@ end
 Non-allocating version of `phiv_dense`.
 """
 function phiv_dense!(w::AbstractMatrix{T}, A::AbstractMatrix{T},
-    v::AbstractVector{T}, k::Integer;
-    cache = nothing,
-    expmethod = ExpMethodHigham2005Base()) where {T <: Number}
+        v::AbstractVector{T}, k::Integer;
+        cache = nothing,
+        expmethod = ExpMethodHigham2005Base()) where {T <: Number}
     @assert size(w, 1)==size(A, 1)==size(A, 2)==length(v) "Dimension mismatch"
     @assert size(w, 2)==k + 1 "Dimension mismatch"
     m = length(v)
@@ -118,7 +118,7 @@ end
 Non-allocating version of `phi` for non-diagonal matrix inputs.
 """
 function phi!(out::Vector{Matrix{T}}, A::AbstractMatrix{T}, k::Integer; caches = nothing,
-    expmethod = ExpMethodHigham2005Base()) where {T <: Number}
+        expmethod = ExpMethodHigham2005Base()) where {T <: Number}
     m = size(A, 1)
     @assert length(out) == k + 1&&all(P -> size(P) == (m, m), out) "Dimension mismatch"
     if caches == nothing
@@ -142,7 +142,7 @@ function phi!(out::Vector{Matrix{T}}, A::AbstractMatrix{T}, k::Integer; caches =
     return out
 end
 function phi!(out::Vector{Diagonal{T, V}}, A::Diagonal{T, V}, k::Integer;
-    caches = nothing) where {T <: Number, V <: AbstractVector{T}}
+        caches = nothing) where {T <: Number, V <: AbstractVector{T}}
     for i in 1:size(A, 1)
         phiz = phi(A[i, i], k; cache = caches)
         for j in 1:(k + 1)

src/precompile.jl

@@ -32,11 +32,11 @@
     end
 
     precomp_ms = [
-        #ExpMethodHigham2005(),
-        ExpMethodHigham2005Base(),
-        #ExpMethodGeneric(),
-        #ExpMethodNative(),
-        #ExpMethodDiagonalization(),
+    #ExpMethodHigham2005(),
+        ExpMethodHigham2005Base()
+    #ExpMethodGeneric(),
+    #ExpMethodNative(),
+    #ExpMethodDiagonalization(),
     ]
 
     Txs = [Float64]