Skip to content

Commit

Permalink
Merge pull request #853 from SciML/new_conserved_quantities_metadata
Browse files Browse the repository at this point in the history
Add `conservedquantity` metadata
  • Loading branch information
TorkelE committed May 21, 2024
2 parents 1808a3a + 334ccb2 commit a649acc
Show file tree
Hide file tree
Showing 2 changed files with 39 additions and 2 deletions.
20 changes: 19 additions & 1 deletion src/network_analysis.jl
Original file line number Diff line number Diff line change
Expand Up @@ -542,6 +542,23 @@ end

### Conservation Laws ###

# Implements the `conserved` parameter metadata.
struct ConservedParameter end
Symbolics.option_to_metadata_type(::Val{:conserved}) = ConservedParameter

"""
isconserved(p)
Checks if the input parameter (`p`) is a conserved quantity (i.e. have the `conserved`)
metadata.
"""
isconserved(x::Num, args...) = isconserved(Symbolics.unwrap(x), args...)
function isconserved(x, default = false)
p = Symbolics.getparent(x, nothing)
p === nothing || (x = p)
Symbolics.getmetadata(x, ConservedParameter, default)
end

"""
conservedequations(rn::ReactionSystem)
Expand Down Expand Up @@ -635,7 +652,8 @@ function cache_conservationlaw_eqs!(rn::ReactionSystem, N::AbstractMatrix, col_o
indepspecs = sts[indepidxs]
depidxs = col_order[(r + 1):end]
depspecs = sts[depidxs]
constants = MT.unwrap.(MT.scalarize((@parameters Γ[1:nullity])[1]))
constants = MT.unwrap.(MT.scalarize(only(
@parameters $(CONSERVED_CONSTANT_SYMBOL)[1:nullity] [conserved=true])))

conservedeqs = Equation[]
constantdefs = Equation[]
Expand Down
21 changes: 20 additions & 1 deletion test/network_analysis/conservation_laws.jl
Original file line number Diff line number Diff line change
Expand Up @@ -142,4 +142,23 @@ let
sprob3.g(g3, u3, sprob3.p, 1.0)
@test isapprox(g1, g2[istsidxs, :])
@test isapprox(g2[istsidxs, :], g3)
end
end

### ConservedParameter Metadata Tests ###

# Checks that `conserved` metadata is added correctly to parameters.
# Checks that the `isconserved` getter function works correctly.
let
# Creates ODESystem with conserved quantities.
rs = @reaction_network begin
(k1,k2), X1 <--> X2
(k1,k2), Y1 <--> Y2
end
osys = convert(ODESystem, rs; remove_conserved = true)

# Checks that the correct parameters have the `conserved` metadata.
@test Catalyst.isconserved(osys.Γ[1])
@test Catalyst.isconserved(osys.Γ[2])
@test !Catalyst.isconserved(osys.k1)
@test !Catalyst.isconserved(osys.k2)
end

0 comments on commit a649acc

Please sign in to comment.