Skip to content

Commit

Permalink
JuliaFormatter.jl
Browse files Browse the repository at this point in the history
  • Loading branch information
nathanaelbosch committed Feb 14, 2024
1 parent 895b666 commit efeaaad
Show file tree
Hide file tree
Showing 9 changed files with 71 additions and 43 deletions.
10 changes: 7 additions & 3 deletions src/alg_utils.jl
Original file line number Diff line number Diff line change
Expand Up @@ -13,9 +13,13 @@ OrdinaryDiffEq.isfsal(::AbstractEK) = false

for ALG in [:EK1, :DiagonalEK1]
@eval OrdinaryDiffEq._alg_autodiff(::$ALG{CS,AD}) where {CS,AD} = Val{AD}()
@eval OrdinaryDiffEq.alg_difftype(::$ALG{CS,AD,DiffType}) where {CS,AD,DiffType} = DiffType
@eval OrdinaryDiffEq.standardtag(::$ALG{CS,AD,DiffType,ST}) where {CS,AD,DiffType,ST} = ST
@eval OrdinaryDiffEq.concrete_jac(::$ALG{CS,AD,DiffType,ST,CJ}) where {CS,AD,DiffType,ST,CJ} = CJ
@eval OrdinaryDiffEq.alg_difftype(::$ALG{CS,AD,DiffType}) where {CS,AD,DiffType} =
DiffType
@eval OrdinaryDiffEq.standardtag(::$ALG{CS,AD,DiffType,ST}) where {CS,AD,DiffType,ST} =
ST
@eval OrdinaryDiffEq.concrete_jac(
::$ALG{CS,AD,DiffType,ST,CJ},
) where {CS,AD,DiffType,ST,CJ} = CJ
@eval OrdinaryDiffEq.get_chunksize(::$ALG{CS}) where {CS} = Val(CS)
@eval OrdinaryDiffEq.isimplicit(::$ALG) = true
end
Expand Down
46 changes: 31 additions & 15 deletions src/algorithms.jl
Original file line number Diff line number Diff line change
Expand Up @@ -3,7 +3,13 @@
########################################################################################
abstract type AbstractEK <: OrdinaryDiffEq.OrdinaryDiffEqAdaptiveAlgorithm end

function ekargcheck(alg; diffusionmodel, pn_observation_noise, covariance_factorization, kwargs...)
function ekargcheck(
alg;
diffusionmodel,
pn_observation_noise,
covariance_factorization,
kwargs...,
)
if (isstatic(diffusionmodel) && diffusionmodel.calibrate) &&
(!isnothing(pn_observation_noise) && !iszero(pn_observation_noise))
throw(
Expand All @@ -14,20 +20,21 @@ function ekargcheck(alg; diffusionmodel, pn_observation_noise, covariance_factor
end
if alg == EK1
if diffusionmodel isa FixedMVDiffusion && diffusionmodel.calibrate
throw(
ArgumentError(
"The `EK1` algorithm does not support automatic global calibration of multivariate diffusion models. Either use a scalar diffusion model, or set `calibrate=false` and calibrate manually by optimizing `sol.pnstats.log_likelihood`. Or use a different solve, like `EK0` or `DiagonalEK1`.",
),
)
throw(
ArgumentError(
"The `EK1` algorithm does not support automatic global calibration of multivariate diffusion models. Either use a scalar diffusion model, or set `calibrate=false` and calibrate manually by optimizing `sol.pnstats.log_likelihood`. Or use a different solve, like `EK0` or `DiagonalEK1`.",
),
)
elseif diffusionmodel isa DynamicMVDiffusion
throw(
ArgumentError(
throw(
ArgumentError(
"The `EK1` algorithm does not support automatic calibration of local multivariate diffusion models. Either use a scalar diffusion model, or use a different solve, like `EK0` or `DiagonalEK1`.",
),
)
),
)
end
end
if diffusionmodel isa DynamicMVDiffusion && covariance_factorization == BlockDiagonalCovariance
if diffusionmodel isa DynamicMVDiffusion &&
covariance_factorization == BlockDiagonalCovariance
throw(
ArgumentError(
"Currenty the `DynamicMVDiffusion` does not work properly with the `BlockDiagonalCovariance`. Use `DenseCovariance` instead, or change the diffusionmodel to a scalar one and use `DynamicDiffusion`.",
Expand Down Expand Up @@ -111,7 +118,8 @@ struct EK0{PT,DT,IT,RT,CF} <: AbstractEK
) where {PT,DT,IT,RT,CF} = begin
ekargcheck(EK0; diffusionmodel, pn_observation_noise)
new{PT,DT,IT,RT,CF}(
prior, diffusionmodel, smooth, initialization, pn_observation_noise, covariance_factorization)
prior, diffusionmodel, smooth, initialization, pn_observation_noise,
covariance_factorization)
end
end

Expand Down Expand Up @@ -218,7 +226,11 @@ struct DiagonalEK1{CS,AD,DiffType,ST,CJ,PT,DT,IT,RT,CF} <: AbstractEK
standardtag=Val{true}(),
concrete_jac=nothing,
pn_observation_noise::RT=nothing,
covariance_factorization::CF=covariance_structure(DiagonalEK1, prior, diffusionmodel),
covariance_factorization::CF=covariance_structure(
DiagonalEK1,
prior,
diffusionmodel,
),
) where {PT,DT,IT,RT,CF} = begin
ekargcheck(DiagonalEK1; diffusionmodel, pn_observation_noise, covariance_factorization)
new{
Expand All @@ -243,7 +255,6 @@ struct DiagonalEK1{CS,AD,DiffType,ST,CJ,PT,DT,IT,RT,CF} <: AbstractEK
end
end


"""
ExpEK(; L, order=3, kwargs...)
Expand Down Expand Up @@ -323,7 +334,12 @@ function DiffEqBase.remake(thing::EK1{CS,AD,DT,ST,CJ}; kwargs...) where {CS,AD,D
)
end

function DiffEqBase.prepare_alg(alg::Union{EK1{0},DiagonalEK1{0}}, u0::AbstractArray{T}, p, prob) where {T}
function DiffEqBase.prepare_alg(
alg::Union{EK1{0},DiagonalEK1{0}},
u0::AbstractArray{T},
p,
prob,
) where {T}
# See OrdinaryDiffEq.jl: ./src/alg_utils.jl (where this is copied from).
# In the future we might want to make EK1 an OrdinaryDiffEqAdaptiveImplicitAlgorithm and
# use the prepare_alg from OrdinaryDiffEq; but right now, we do not use `linsolve` which
Expand Down
6 changes: 3 additions & 3 deletions src/blockdiagonals.jl
Original file line number Diff line number Diff line change
Expand Up @@ -30,7 +30,7 @@ end
_matmul!(
C::BlockDiagonal{T},
A::BlockDiagonal{T},
B::Adjoint{T, <:BlockDiagonal{T}},
B::Adjoint{T,<:BlockDiagonal{T}},
) where {T<:LinearAlgebra.BlasFloat} = begin
@assert length(C.blocks) == length(A.blocks) == length(B.parent.blocks)
@simd ivdep for i in eachindex(blocks(C))
Expand All @@ -41,7 +41,7 @@ end

_matmul!(
C::BlockDiagonal{T},
A::Adjoint{T, <:BlockDiagonal{T}},
A::Adjoint{T,<:BlockDiagonal{T}},
B::BlockDiagonal{T},
) where {T<:LinearAlgebra.BlasFloat} = begin
@assert length(C.blocks) == length(A.parent.blocks) == length(B.blocks)
Expand Down Expand Up @@ -82,7 +82,7 @@ LinearAlgebra.rmul!(B::BlockDiagonal, n::Number) = @simd ivdep for i in eachinde
rmul!(B.blocks[i], n)
end
LinearAlgebra.adjoint(B::BlockDiagonal) = Adjoint(B)
Base.:*(A::Adjoint{T, <:BlockDiagonal}, B::BlockDiagonal) where {T} = begin
Base.:*(A::Adjoint{T,<:BlockDiagonal}, B::BlockDiagonal) where {T} = begin
@assert length(A.parent.blocks) == length(B.blocks)
return BlockDiagonal([A.parent.blocks[i]' * B.blocks[i] for i in eachindex(B.blocks)])
end
27 changes: 16 additions & 11 deletions src/diffusions.jl
Original file line number Diff line number Diff line change
Expand Up @@ -6,24 +6,26 @@ isdynamic(diffusion::AbstractStaticDiffusion) = false
isstatic(diffusion::AbstractDynamicDiffusion) = false
isdynamic(diffusion::AbstractDynamicDiffusion) = true

apply_diffusion(Q::PSDMatrix{T, <:Matrix}, diffusion::Diagonal) where {T} = begin
apply_diffusion(Q::PSDMatrix{T,<:Matrix}, diffusion::Diagonal) where {T} = begin
d = size(diffusion, 1)
q = size(Q, 1) ÷ d - 1
return PSDMatrix(Q.R * sqrt.(kron(diffusion, I(q+1))))
return PSDMatrix(Q.R * sqrt.(kron(diffusion, I(q + 1))))
end
apply_diffusion(
Q::PSDMatrix{T, <:IsometricKroneckerProduct},
diffusion::Diagonal{T, <:FillArrays.Fill},
Q::PSDMatrix{T,<:IsometricKroneckerProduct},
diffusion::Diagonal{T,<:FillArrays.Fill},
) where {T} = begin
PSDMatrix(Q.R * sqrt.(diffusion.diag.value))
end
apply_diffusion(Q::PSDMatrix{T, <:BlockDiagonal}, diffusion::Diagonal) where {T} = begin
PSDMatrix(BlockDiagonal([
Q.R.blocks[i] * sqrt.(diffusion.diag[i]) for i in eachindex(Q.R.blocks)
]))
apply_diffusion(Q::PSDMatrix{T,<:BlockDiagonal}, diffusion::Diagonal) where {T} = begin
PSDMatrix(
BlockDiagonal([
Q.R.blocks[i] * sqrt.(diffusion.diag[i]) for i in eachindex(Q.R.blocks)
]),
)
end

apply_diffusion!(Q::PSDMatrix, diffusion::Diagonal{T, <:FillArrays.Fill}) where {T} =
apply_diffusion!(Q::PSDMatrix, diffusion::Diagonal{T,<:FillArrays.Fill}) where {T} =
rmul!(Q.R, sqrt.(diffusion.diag.value))
apply_diffusion!(
Q::PSDMatrix{T,<:BlockDiagonal},
Expand All @@ -33,10 +35,13 @@ apply_diffusion!(
rmul!(blocks(Q.R)[i], diffusion.diag[i])
end

apply_diffusion!(out::PSDMatrix, Q::PSDMatrix, diffusion::Diagonal{T,<:FillArrays.Fill}) where {T} =
apply_diffusion!(
out::PSDMatrix,
Q::PSDMatrix,
diffusion::Diagonal{T,<:FillArrays.Fill},
) where {T} =
rmul!(Q.R, sqrt.(diffusion.diag.value))


estimate_global_diffusion(diffusion::AbstractDynamicDiffusion, d, q, Eltype) = error()

"""
Expand Down
10 changes: 5 additions & 5 deletions src/filtering/markov_kernel.jl
Original file line number Diff line number Diff line change
Expand Up @@ -292,24 +292,24 @@ function compute_backward_kernel!(
_Kout = AffineNormalKernel(
Kout.A.blocks[i],
view(Kout.b, (i-1)*(q+1)+1:i*(q+1)),
PSDMatrix(Kout.C.R.blocks[i])
PSDMatrix(Kout.C.R.blocks[i]),
)
_xpred = Gaussian(
view(xpred.μ, (i-1)*(q+1)+1:i*(q+1)),
PSDMatrix(xpred.Σ.R.blocks[i])
PSDMatrix(xpred.Σ.R.blocks[i]),
)
_x = Gaussian(
view(x.μ, (i-1)*(q+1)+1:i*(q+1)),
PSDMatrix(x.Σ.R.blocks[i])
PSDMatrix(x.Σ.R.blocks[i]),
)
_K = AffineNormalKernel(
K.A.blocks[i],
ismissing(K.b) ? missing : view(K.b, (i-1)*(q+1)+1:i*(q+1)),
PSDMatrix(K.C.R.blocks[i])
PSDMatrix(K.C.R.blocks[i]),
)
_C_DxD = C_DxD.blocks[i]
compute_backward_kernel!(
_Kout, _xpred, _x, _K, C_DxD=_C_DxD, diffusion=diffusion
_Kout, _xpred, _x, _K, C_DxD=_C_DxD, diffusion=diffusion,
)
end
return Kout
Expand Down
8 changes: 6 additions & 2 deletions src/filtering/predict.jl
Original file line number Diff line number Diff line change
Expand Up @@ -66,7 +66,7 @@ function predict_cov!(
Qh::PSDMatrix,
C_DxD::AbstractMatrix,
C_2DxD::AbstractMatrix,
diffusion::Union{Number, Diagonal},
diffusion::Union{Number,Diagonal},
)
if iszero(diffusion)
fast_X_A_Xt!(Σ_out, Σ_curr, Ah)
Expand All @@ -85,7 +85,11 @@ function predict_cov!(
@warn "This is not yet implemented efficiently; TODO"
d = size(diffusion, 1)
q = D ÷ d - 1
_matmul!(view(R, D+1:2D, 1:D), Qh.R, sqrt.(kron(Eye(d)*diffusion, Eye(q + 1))))
_matmul!(
view(R, D+1:2D, 1:D),
Qh.R,
sqrt.(kron(Eye(d) * diffusion, Eye(q + 1))),
)
end
else
@.. R[D+1:2D, 1:D] = Qh.R
Expand Down
1 change: 0 additions & 1 deletion src/filtering/update.jl
Original file line number Diff line number Diff line change
Expand Up @@ -193,7 +193,6 @@ function update!(
return x_out, loglikelihood
end


function update!(
x_out::SRGaussian{T,<:BlockDiagonal},
x_pred::SRGaussian{T,<:BlockDiagonal},
Expand Down
4 changes: 2 additions & 2 deletions src/integrator_utils.jl
Original file line number Diff line number Diff line change
Expand Up @@ -70,11 +70,11 @@ diffusion estimates that are in there. Typically, `diffusion` is either a global
or the specified initial diffusion value if no calibration is desired.
"""
function set_diffusions!(solution::AbstractProbODESolution, diffusion)
if diffusion isa Diagonal{<:Number, <:FillArrays.Fill}
if diffusion isa Diagonal{<:Number,<:FillArrays.Fill}
@simd ivdep for i in eachindex(solution.diffusions)
solution.diffusions[i] = copy(diffusion)
end
elseif diffusion isa Diagonal{<:Number, <:Vector}
elseif diffusion isa Diagonal{<:Number,<:Vector}
@simd ivdep for d in solution.diffusions
copy!(d, diffusion)
end
Expand Down
2 changes: 1 addition & 1 deletion src/perform_step.jl
Original file line number Diff line number Diff line change
Expand Up @@ -222,7 +222,7 @@ function estimate_errors!(cache::AbstractODEFilterCache)

R = cache.C_Dxd

if local_diffusion isa Diagonal{<:Number, <:Vector}
if local_diffusion isa Diagonal{<:Number,<:Vector}
_Q = apply_diffusion(Qh, local_diffusion)
_matmul!(R, _Q.R, H')
error_estimate = view(cache.tmp, 1:d)
Expand Down

0 comments on commit efeaaad

Please sign in to comment.