Skip to content

Commit c1e8e45

Browse files
committed
update
1 parent cc47e15 commit c1e8e45

File tree

4 files changed

+49
-5
lines changed

4 files changed

+49
-5
lines changed

docs/src/catalyst_functionality/chemistry_related_functionality.md

Lines changed: 10 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -37,6 +37,16 @@ Finally, it is possible to check whether a species is a compound or not using th
3737
iscompound(CO2)
3838
```
3939

40+
!!! note
41+
Currently, compounds with components that depend on variables that are not `t` are not supported. E.g.
42+
```julia
43+
@variables x y
44+
@species O(x, y)
45+
@compound O2 ~ 2O
46+
```
47+
will currently throw an error.
48+
49+
4050
Compound components that are also compounds are allowed, e.g. we can create a carbonic acid compound (H₂CO₃) that consists of CO₂ and H₂O:
4151
```@example chem1
4252
@species H(t)

src/chemistry_functionality.jl

Lines changed: 8 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -47,7 +47,7 @@ function make_compound(expr)
4747

4848
# Extracts the composite species name, and a Vector{ReactantStruct} of its components.
4949
species_expr = expr.args[2] # E.g. :(CO2(t))
50-
species_expr = insert_independent_variable(species_expr, :t) # Add independent variable (e.g. goes from `CO2` to `CO2(t)`).
50+
species_expr = insert_independent_variable(species_expr, [:t]) # Add independent variable (e.g. goes from `CO2` to `CO2(t)`).
5151
species_name = find_varname_in_declaration(expr.args[2]) # E.g. :CO2
5252
composition = Catalyst.recursive_find_reactants!(expr.args[3], 1, Vector{ReactantStruct}(undef, 0))
5353

@@ -61,13 +61,20 @@ function make_compound(expr)
6161

6262
# Creates the found expressions that will create the compound species.
6363
# The `Expr(:escape, :(...))` is required so that teh expressions are evaluated in the scope the users use the macro in (to e.g. detect already exiting species).
64+
non_t_iv_error_check_expr = Expr(:escape, :(@species $species_expr)) # E.g. `@species CO2(t)`
6465
species_declaration_expr = Expr(:escape, :(@species $species_expr)) # E.g. `@species CO2(t)`
6566
compound_designation_expr = Expr(:escape, :($species_name = ModelingToolkit.setmetadata($species_name, Catalyst.CompoundSpecies, true))) # E.g. `CO2 = ModelingToolkit.setmetadata(CO2, Catalyst.CompoundSpecies, true)`
6667
components_designation_expr = Expr(:escape, :($species_name = ModelingToolkit.setmetadata($species_name, Catalyst.CompoundComponents, $components))) # E.g. `CO2 = ModelingToolkit.setmetadata(CO2, Catalyst.CompoundSpecies, [C, O])`
6768
coefficients_designation_expr = Expr(:escape, :($species_name = ModelingToolkit.setmetadata($species_name, Catalyst.CompoundCoefficients, $coefficients))) # E.g. `CO2 = ModelingToolkit.setmetadata(CO2, Catalyst.CompoundSpecies, [1, 2])`
6869

70+
# Currently, non-t independent variables are not supported for compounds. If there are any like these, we throw an error:
71+
non_t_iv_error_check_expr = Expr(:escape, :(issetequal(unique(reduce(vcat, arguments.(ModelingToolkit.unwrap.($components)))), [t]) || error("Currently, compounds depending on components that are not \"t\" are not supported.")))
72+
73+
println(non_t_iv_error_check_expr)
74+
6975
# Returns the rephrased expression.
7076
return quote
77+
$non_t_iv_error_check_expr
7178
$species_declaration_expr
7279
$compound_designation_expr
7380
$components_designation_expr

src/expression_utils.jl

Lines changed: 8 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -50,19 +50,23 @@ end
5050
# X(t), [metadata=true]
5151
# X(t) = 1.0, [metadata=true]
5252
# (In this example the independent variable t was inserted).
53+
# Extra dispatch is to ensure so that iv is a vector (so that we can handle both single or multiple iv insertion in one function).
54+
# - This way we don't have to make a check for how to create the final expression (which is different for symbol or vector of symbols).
5355
function insert_independent_variable(expr_in, iv)
56+
# If iv is a symbol (e.g. :t), we convert to vector (e.g. [:t]). This way we can handle single and multiple ivs using the same code, without needing to check at the end.
57+
(iv isa Symbol) && (iv = [iv])
5458
# If expr is a symbol, just attach the iv. If not we have to create a new expr and mutate it.
5559
# Because Symbols (a possible input) cannot be mutated, this function cannot mutate the input (would have been easier if Expr input was guaranteed).
56-
(expr_in isa Symbol) && (return :($(expr_in)($iv)))
60+
(expr_in isa Symbol) && (return Expr(:call, expr_in, iv...))
5761
expr = deepcopy(expr_in)
5862

5963
if expr.head == :(=) # Case: :(X = 1.0)
60-
expr.args[1] = :($(expr.args[1])($iv))
64+
expr.args[1] = Expr(:call, expr.args[1], iv...)
6165
elseif expr.head == :tuple
6266
if expr.args[1] isa Symbol # Case: :(X, [metadata=true])
63-
expr.args[1] = :($(expr.args[1])($iv))
67+
expr.args[1] = Expr(:call, expr.args[1], iv...)
6468
elseif (expr.args[1].head == :(=)) && (expr.args[1].args[1] isa Symbol) # Case: :(X = 1.0, [metadata=true])
65-
expr.args[1].args[1] = :($(expr.args[1].args[1])($iv))
69+
expr.args[1].args[1] = Expr(:call, expr.args[1].args[1], iv...)
6670
end
6771
end
6872
(expr == expr_in) && error("Failed to add independent variable $(iv) to expression: $expr_in")

test/miscellaneous_tests/compound_macro.jl

Lines changed: 23 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -67,6 +67,29 @@ let
6767
@test iscompound(rn.NH3_5)
6868
end
6969

70+
### Test Various Independent Variables ###
71+
let
72+
@variables t x y z
73+
@species C(t) H(x) N(x) O(t) P(t,x) S(x,y)
74+
75+
@compound CO2 ~ C + 2O
76+
@test issetequal(arguments(ModelingToolkit.unwrap(CO2)), [t])
77+
78+
# These 4 tests checks currently non-supported feature.
79+
@test_broken false
80+
@test_broken false
81+
@test_broken false
82+
@test_broken false
83+
# @compound NH4, [output=true] ~ N + 4H
84+
# @compound (H2O = 2.0) ~ 2H + N
85+
# @compound PH4 ~ P + 4H
86+
# @compound SO2 ~ S + 2O
87+
# @test issetequal(arguments(ModelingToolkit.unwrap(NH4)), [x])
88+
# @test issetequal(arguments(ModelingToolkit.unwrap(H2O)), [t, x])
89+
# @test issetequal(arguments(ModelingToolkit.unwrap(PH4)), [t, x])
90+
# @test issetequal(arguments(ModelingToolkit.unwrap(SO2)), [t, x, y])
91+
end
92+
7093
### Other Minor Tests ###
7194

7295
# Test base functionality in two cases.

0 commit comments

Comments
 (0)