equations are now supported, rather than throwing an error.

Author: TorkelE

ext/CatalystHomotopyContinuationExtension/homotopy_continuation_extension.jl

@@ -43,9 +43,8 @@ function Catalyst.hc_steady_states(rs::ReactionSystem, ps; filter_negative=true,
 end
 
 # For a given reaction system, parameter values, and initial conditions, find the polynomial that HC solves to find steady states.
-function steady_state_polynomial(rs_in::ReactionSystem, ps, u0)
-    any(!(eq isa Reaction) for eq in rs_in.eqs) && error("This feature is currently not supported for reaction systems containing (non-reaction) equations.")
-    rs = @set rs_in.rxs = [Catalyst.expand_registered_functions(rx) for rx in rs_in.rxs]
+function steady_state_polynomial(rs::ReactionSystem, ps, u0)
+    rs = expand_registered_functions(rs)
     ns = convert(NonlinearSystem, rs; remove_conserved = true)
     pre_varmap = [symmap_to_varmap(rs,u0)..., symmap_to_varmap(rs,ps)...]
     conservationlaw_errorcheck(rs, pre_varmap)

src/registered_functions.jl

@@ -132,7 +132,15 @@ function expand_registered_functions(expr)
     end
     return expr
 end
-# If applied to a reaction, return a reaction with its rate modified.
+# If applied to a Reaction, return a reaction with its rate modified.
 function expand_registered_functions(rx::Reaction)
     Reaction(expand_registered_functions(rx.rate), rx.substrates, rx.products, rx.substoich, rx.prodstoich, rx.netstoich, rx.only_use_rate)
 end
+# If applied to a Equation, returns it with it applied to lhs and rhs
+function expand_registered_functions(eq::Equation)
+    return expand_registered_functions(eq.lhs) ~ expand_registered_functions(eq.rhs)
+end
+# If applied to a ReactionSystem, applied function to all Reactions and other Equations, and return updated system.
+function expand_registered_functions(rs::ReactionSystem)
+    @set rs.eqs = [Catalyst.expand_registered_functions(eq) for eq in rs.eqs]
+end