Add Runic.jl formatter and allow default p for CoupledSDEs

.github/workflows/format.yml

@@ -0,0 +1,34 @@
+name: Format Check
+on:
+  push:
+    branches:
+      - 'main'
+    paths:
+      - '.github/workflows/FormatCheck.yml'
+      - '**.jl'
+  pull_request:
+    branches:
+      - 'main'
+    paths:
+      - '.github/workflows/FormatCheck.yml'
+      - '**.jl'
+    types:
+      - opened
+      - reopened
+      - synchronize
+      - ready_for_review
+
+jobs:
+  runic:
+    name: Runic
+    runs-on: ubuntu-latest
+    if: ${{ !github.event.pull_request.draft }}
+    steps:
+      - uses: actions/checkout@v6
+      - uses: julia-actions/setup-julia@v2
+        with:
+          version: '1'
+      - uses: julia-actions/cache@v2
+      - uses: fredrikekre/runic-action@v1
+        with:
+          version: '1'

Project.toml

@@ -1,7 +1,7 @@
 name = "DynamicalSystemsBase"
 uuid = "6e36e845-645a-534a-86f2-f5d4aa5a06b4"
 repo = "https://github.com/JuliaDynamics/DynamicalSystemsBase.jl.git"
-version = "3.15.4"
+version = "3.15.5"
 
 [deps]
 ForwardDiff = "f6369f11-7733-5829-9624-2563aa707210"

docs/make.jl

@@ -18,6 +18,7 @@ Downloads.download(
 )
 include("build_docs_with_style.jl")
 
-build_docs_with_style(pages,
-    DynamicalSystemsBase, SciMLBase, StateSpaceSets, StochasticSystemsBase;
+build_docs_with_style(
+    pages,
+    DynamicalSystemsBase, SciMLBase, StateSpaceSets, StochasticSystemsBase
 )

ext/StochasticSystemsBase.jl

@@ -1,11 +1,11 @@
 module StochasticSystemsBase
 
 using DynamicalSystemsBase: DynamicalSystemsBase, SciMLBase, correct_state, CoupledODEs,
-                            CoupledSDEs, StateSpaceSets, isinplace, _delete, set_parameter!,
-                            set_state!, dynamic_rule, isdeterministic, current_state,
-                            DynamicalSystemsBase, _set_parameter!, u_modified!,
-                            additional_details, referrenced_sciml_prob, DEFAULT_DIFFEQ,
-                            SVector, SMatrix, current_parameters
+    CoupledSDEs, StateSpaceSets, isinplace, _delete, set_parameter!,
+    set_state!, dynamic_rule, isdeterministic, current_state,
+    DynamicalSystemsBase, _set_parameter!, u_modified!,
+    additional_details, referrenced_sciml_prob, DEFAULT_DIFFEQ,
+    SVector, SMatrix, current_parameters
 using SciMLBase: SDEProblem, AbstractSDEIntegrator, __init, SDEFunction, step!
 using StochasticDiffEq: SOSRA
 using LinearAlgebra

ext/src/CoupledSDEs.jl

@@ -7,8 +7,8 @@ import DynamicalSystemsBase: CoupledSDEs
 # SOSRA is only applicable for additive  noise and must be adaptive
 # default sciml tolerance is 1e-2
 const DEFAULT_SDE_SOLVER = SOSRA()
-const DEFAULT_STOCH_DIFFEQ_KWARGS = (abstol=1e-2, reltol=1e-2, dt=0.1)
-const DEFAULT_STOCH_DIFFEQ = (alg=DEFAULT_SDE_SOLVER, DEFAULT_STOCH_DIFFEQ_KWARGS...)
+const DEFAULT_STOCH_DIFFEQ_KWARGS = (abstol = 1.0e-2, reltol = 1.0e-2, dt = 0.1)
+const DEFAULT_STOCH_DIFFEQ = (alg = DEFAULT_SDE_SOLVER, DEFAULT_STOCH_DIFFEQ_KWARGS...)
 
 # Function from user `@xlxs4`, see
 # https://github.com/JuliaDynamics/jl/pull/153
@@ -25,45 +25,46 @@ end
 ###########################################################################################
 
 function DynamicalSystemsBase.CoupledSDEs(
-    f,
-    u0,
-    p=SciMLBase.NullParameters();
-    g=nothing,
-    noise_strength=1.0,
-    covariance=nothing,
-    t0=0.0,
-    diffeq=DEFAULT_STOCH_DIFFEQ,
-    noise_prototype=nothing,
-    noise_process=nothing,
-    seed=UInt64(0)
-)
+        f,
+        u0,
+        p = SciMLBase.NullParameters();
+        g = nothing,
+        noise_strength = 1.0,
+        covariance = nothing,
+        t0 = 0.0,
+        diffeq = DEFAULT_STOCH_DIFFEQ,
+        noise_prototype = nothing,
+        noise_process = nothing,
+        seed = UInt64(0)
+    )
     IIP = isinplace(f, 4) # from SciMLBase
     if !isnothing(g)
         @assert IIP == isinplace(g, 4) "f and g must both be in-place or out-of-place"
     end
 
     noise_type, cov = find_noise_type(g, u0, p, t0, noise_process, covariance, noise_prototype, IIP)
     g, noise_prototype = construct_diffusion_function(
-        g, cov, noise_prototype, noise_strength, length(u0), IIP)
+        g, cov, noise_prototype, noise_strength, length(u0), IIP
+    )
 
     s = correct_state(Val{IIP}(), u0)
     T = eltype(s)
     prob = SDEProblem{IIP}(
         f,
         g,
         s,
-        (T(t0), T(1e11)),
+        (T(t0), T(1.0e11)),
         p;
-        noise_rate_prototype=noise_prototype,
-        noise=noise_process,
-        seed=seed
+        noise_rate_prototype = noise_prototype,
+        noise = noise_process,
+        seed = seed
     )
     return CoupledSDEs(prob, diffeq, noise_type)
 end
 
 function DynamicalSystemsBase.CoupledSDEs(
-    prob::SDEProblem, diffeq=DEFAULT_STOCH_DIFFEQ, noise_type=nothing; special_kwargs...
-)
+        prob::SDEProblem, diffeq = DEFAULT_STOCH_DIFFEQ, noise_type = nothing; special_kwargs...
+    )
     if haskey(special_kwargs, :diffeq)
         throw(
             ArgumentError(
@@ -77,7 +78,7 @@ function DynamicalSystemsBase.CoupledSDEs(
     if prob.tspan === (nothing, nothing)
         # If the problem was made via MTK, it is possible to not have a default timespan.
         U = eltype(prob.u0)
-        prob = SciMLBase.remake(prob; tspan=(U(0), U(Inf)))
+        prob = SciMLBase.remake(prob; tspan = (U(0), U(Inf)))
     end
     if isnothing(noise_type)
         noise_type, _ = find_noise_type(prob, IIP)
@@ -89,43 +90,45 @@ function DynamicalSystemsBase.CoupledSDEs(
         solver;
         remaining...,
         # Integrators are used exclusively iteratively. There is no reason to save anything.
-        save_start=false,
-        save_end=false,
-        save_everystep=false,
+        save_start = false,
+        save_end = false,
+        save_everystep = false,
         # DynamicalSystems.jl operates on integrators and `step!` exclusively,
         # so there is no reason to limit the maximum time evolution
-        maxiters=Inf
+        maxiters = Inf
     )
-    return CoupledSDEs{IIP,D,typeof(integ),P}(
+    return CoupledSDEs{IIP, D, typeof(integ), P}(
         integ, deepcopy(prob.p), diffeq, noise_type
     )
 end
 # This preserves the referrenced MTK system and the originally passed diffeq kwargs
 CoupledSDEs(ds::CoupledSDEs, diffeq) = CoupledSDEs(SDEProblem(ds), merge(ds.diffeq, diffeq))
 
 """
-    CoupledSDEs(ds::CoupledODEs, p; kwargs...)
+    CoupledSDEs(ds::CoupledODEs, p = current_parameters(ds); kwargs...)
 
-Converts a [`CoupledODEs`
-](https://juliadynamics.github.io/DynamicalSystems.jl/stable/tutorial/#CoupledODEs)
-system into a [`CoupledSDEs`](@ref).
+Converts a [`CoupledODEs`](@ref) into a [`CoupledSDEs`](@ref).
+While `p` is optional, it is likely to change as the
+diffusion (noise) function `g` is added.
 """
 function DynamicalSystemsBase.CoupledSDEs(
-    ds::CoupledODEs,
-    p; # the parameter is likely changed as the diffusion function g is added.
-    g=nothing,
-    noise_strength=1.0,
-    covariance=nothing,
-    diffeq=DEFAULT_STOCH_DIFFEQ,
-    noise_prototype=nothing,
-    noise_process=nothing,
-    seed=UInt64(0)
-)
+        ds::CoupledODEs,
+        p = current_parameters(ds); # the parameter is likely changed as the diffusion function g is added.
+        g = nothing,
+        noise_strength = 1.0,
+        covariance = nothing,
+        diffeq = DEFAULT_STOCH_DIFFEQ,
+        noise_prototype = nothing,
+        noise_process = nothing,
+        seed = UInt64(0)
+    )
     prob = referrenced_sciml_prob(ds)
     # we want the symbolic jacobian to be transfered over
     # dynamic_rule(ds) takes the deepest nested f wich does not have a jac field
-    return CoupledSDEs(prob.f, current_state(ds), p;
-        g, noise_strength, covariance, diffeq, noise_prototype, noise_process, seed)
+    return CoupledSDEs(
+        prob.f, current_state(ds), p;
+        g, noise_strength, covariance, diffeq, noise_prototype, noise_process, seed
+    )
 end
 
 """
@@ -135,12 +138,13 @@ Converts a [`CoupledSDEs`](@ref) into a [`CoupledODEs`](@ref) by extracting the
 deterministic part of `ds`.
 """
 function DynamicalSystemsBase.CoupledODEs(
-    sys::CoupledSDEs; diffeq=DEFAULT_DIFFEQ, t0=0.0)
+        sys::CoupledSDEs; diffeq = DEFAULT_DIFFEQ, t0 = 0.0
+    )
     prob = referrenced_sciml_prob(sys)
     # we want the symbolic jacobian to be transfered over
     # dynamic_rule(ds) takes the deepest nested f wich does not have a jac field
     return CoupledODEs(
-        prob.f, SVector{length(sys.integ.u)}(sys.integ.u), sys.p0; diffeq=diffeq, t0=t0
+        prob.f, SVector{length(sys.integ.u)}(sys.integ.u), sys.p0; diffeq = diffeq, t0 = t0
     )
 end
 
@@ -150,7 +154,7 @@ function DynamicalSystemsBase.additional_details(ds::CoupledSDEs)
     return [
         "SDE solver" => string(nameof(typeof(solver))),
         "SDE kwargs" => remaining,
-        "Noise type" => ds.noise_type
+        "Noise type" => ds.noise_type,
     ]
 end
 
@@ -159,7 +163,7 @@ end
 ###########################################################################################
 
 SciMLBase.isinplace(::CoupledSDEs{IIP}) where {IIP} = IIP
-StateSpaceSets.dimension(::CoupledSDEs{IIP,D}) where {IIP,D} = D
+StateSpaceSets.dimension(::CoupledSDEs{IIP, D}) where {IIP, D} = D
 DynamicalSystemsBase.current_state(ds::CoupledSDEs) = current_state(ds.integ)
 DynamicalSystemsBase.isdeterministic(ds::CoupledSDEs) = false
 
@@ -193,7 +197,7 @@ If this is not the case, returns `nothing`.
 Note: The diffusion matrix ``Σ`` is the square root of the noise covariance matrix ``Q`` (see
 [`covariance_matrix`](@ref)), defined via the Cholesky decomposition ``Q = Σ Σ^\\top``.
 """
-function diffusion_matrix(ds::CoupledSDEs{IIP,D})::AbstractMatrix where {IIP,D}
+function diffusion_matrix(ds::CoupledSDEs{IIP, D})::AbstractMatrix where {IIP, D}
     if ds.noise_type[:invertible]
         diffusion = diffusion_function(ds)
         A = diffusion(zeros(D), current_parameters(ds), 0.0)
@@ -219,7 +223,7 @@ See also [`diffusion_matrix`](@ref).
 """
 function covariance_matrix(ds::CoupledSDEs)::AbstractMatrix
     A = diffusion_matrix(ds)
-    (A == nothing) ? nothing : A * A'
+    return (A == nothing) ? nothing : A * A'
 end
 
 jacobian(sde::CoupledSDEs) = DynamicalSystemsBase.jacobian(CoupledODEs(sde))
@@ -241,8 +245,8 @@ Here `D` is the system dimension and `IIP` indicated whether the function `g` is
 Returns `g, noise_prototype`.
 """
 function construct_diffusion_function(
-    g, covariance, noise_prototype, noise_strength, D, IIP
-)
+        g, covariance, noise_prototype, noise_strength, D, IIP
+    )
     if isnothing(g) # diagonal additive noise
         cov = isnothing(covariance) ? LinearAlgebra.I(D) : covariance
         size(cov) != (D, D) &&
@@ -258,11 +262,11 @@ function construct_diffusion_function(
             end
         else
             if isdiag(cov)
-                g = (u, p, t) -> SVector{length(diag(A)),eltype(A)}(
+                g = (u, p, t) -> SVector{length(diag(A)), eltype(A)}(
                     diag(noise_strength .* A)
                 )
             else
-                g = (u, p, t) -> SMatrix{size(A)...,eltype(A)}(
+                g = (u, p, t) -> SMatrix{size(A)..., eltype(A)}(
                     noise_strength .* A
                 )
                 noise_prototype = zeros(size(A))

ext/src/classification.jl

@@ -6,7 +6,7 @@ Checks whether the function g depends on u based on 10 random points around the
 function is_state_independent(g, u, p, t)
     rdm_states = [u .+ rand(eltype(u), length(u)) .- 0.5 for _ in 1:10]
     val = map(u -> g(u, p, t), rdm_states)
-    length(unique(val)) == 1
+    return length(unique(val)) == 1
 end
 
 """
@@ -15,14 +15,14 @@ Checks whether g depends explicitly on time for select points on the interval [t
 function is_time_independent(g, u, p, t0)
     trange = t0 .+ [0.0, 0.101, 1.01, 10.1, 101.0]
     val = map(t -> g(u, p, t), trange)
-    length(unique(val)) == 1
+    return length(unique(val)) == 1
 end
 
 """
 Checks whether a matrix x is invertible by verifying that it has nonzero determinant.
 """
-function is_invertible(x; tol=1e-10)
-    F = lu(x, check=false)
+function is_invertible(x; tol = 1.0e-10)
+    F = lu(x, check = false)
     det = abs(prod(diag(F.U)))
     return det > tol
 end
@@ -38,7 +38,7 @@ function is_linear(f, x, y, c)
 end
 
 function diffusion_function(g, IIP, noise_prototype)
-    function diffusion(u, p, t)
+    return function diffusion(u, p, t)
         if IIP
             du = deepcopy(isnothing(noise_prototype) ? u : noise_prototype)
             g(du, u, p, t)
@@ -50,7 +50,7 @@ function diffusion_function(g, IIP, noise_prototype)
 end
 
 function diffusion_function(ds::CoupledSDEs{IIP}) where {IIP}
-    diffusion_function(ds.integ.g, IIP, referrenced_sciml_prob(ds).noise_rate_prototype)
+    return diffusion_function(ds.integ.g, IIP, referrenced_sciml_prob(ds).noise_rate_prototype)
 end
 
 """
@@ -100,11 +100,13 @@ function find_noise_type(g, u0, p, t0, noise, covariance, noise_prototype, IIP)
         islinear = true
         if !state_independent
             for i in 1:10
-                check = is_linear(u -> diffusion(u, p, t0),
-                    u0 + i .* rand(D), u0 + i .* rand(D), 2.0)
+                check = is_linear(
+                    u -> diffusion(u, p, t0),
+                    u0 + i .* rand(D), u0 + i .* rand(D), 2.0
+                )
                 check ? nothing : islinear = false
             end
-        end        
+        end
 
         # Previous formulation:
         #islinear = !state_independent ?
@@ -127,13 +129,16 @@ function find_noise_type(g, u0, p, t0, noise, covariance, noise_prototype, IIP)
         end
     end
 
-    noise_type = (additive=isadditive, autonomous=isautonomous,
-        linear=islinear, invertible=isinvertible)
+    noise_type = (
+        additive = isadditive, autonomous = isautonomous,
+        linear = islinear, invertible = isinvertible,
+    )
     return noise_type, covariance
 end
 
 function find_noise_type(prob::SDEProblem, IIP)
-    find_noise_type(
+    return find_noise_type(
         prob.g, prob.u0, prob.p, prob.tspan[1], prob.noise,
-        nothing, prob.noise_rate_prototype, IIP)
-end
+        nothing, prob.noise_rate_prototype, IIP
+    )
+end