update constraint tutorial

Author: isaacsas

docs/src/model_creation/constraint_equations.md

@@ -1,4 +1,4 @@
-# [Constraint Equations and Events](@id constraint_equations)
+# [Coupled ODEs, Algebraic Equations, and Events](@id constraint_equations)
 In many applications one has additional algebraic or differential equations for
 non-chemical species that can be coupled to a chemical reaction network model.
 Catalyst supports coupled differential and algebraic equations, and currently
@@ -8,24 +8,81 @@ condition is reached, such as providing a drug treatment at a specified time, or
 turning off production of cells once the population reaches a given level.
 Catalyst supports the event representation provided by ModelingToolkit, see
 [here](https://docs.sciml.ai/ModelingToolkit/stable/basics/Events/), allowing
-for both continuous and discrete events (though only the latter are supported
-when converting to `JumpSystem`s currently).
+for both continuous and discrete events.
+
+In this tutorial we'll illustrate how to make use of coupled constraint (i.e.
+ODE/algebraic) equations and events. Let's consider a model of a cell with
+volume $V(t)$ that grows at a rate $\lambda$. For now we'll assume the cell can
+grow indefinitely. We'll also keep track of one protein $P(t)$, which is
+produced at a rate proportional $V$ and can be degraded.
+
+## [Coupling ODE constraints via the DSL](@id constraint_equations_dsl)
+The easiest way to include ODEs and algebraic equations is to just include them
+when using the DSL to specify a model. Here we include an ODE for $V(t)$ along
+with degradation and production reactions for $P(t)$:
+```@example ceq1
+rn = @reaction_network dynamic_cell begin
+    # the growth rate
+    @parameters λ = 1.0
+
+    # assume there is no protein initially
+    @species P(t) = 0.0
+
+    # set the initial volume to 1.0
+    @variables V(t) = 1.0
+
+    # the reactions
+    V, 0 --> P
+    1.0, P --> 0
+
+    # the coupled ODE for V(t)
+    @equations begin
+        D(V) ~ λ * V
+    end
+end
+```
+We can now create an `ODEProblem` from our model and solve it to see how $V(t)$
+and $P(t)$ evolve in time:
+```@example ceq1
+oprob = ODEProblem(growing_cell, [], (0.0, 1.0))
+sol = solve(oprob, Tsit5())
+plot(sol)
+```
+
+## Coupling ODE constraints via directly building a `ReactionSystem`
+As an alternative to the previous approach, we could have also constructed our
+`ReactionSystem` all at once using the symbolic interface:
+```@example ceq2
+using Catalyst, OrdinaryDiffEqTsit5, Plots
+t = default_t()
+D = default_time_deriv()
+
+@parameters λ = 1.0
+@variables V(t) = 1.0
+eq = D(V) ~ λ * V
+rx1 = @reaction $V, 0 --> P
+rx2 = @reaction 1.0, P --> 0
+@named growing_cell = ReactionSystem([rx1, rx2, eq], t)
+setdefaults!(growing_cell, [:P => 0.0])
+growing_cell = complete(growing_cell)
+
+oprob = ODEProblem(growing_cell, [], (0.0, 1.0))
+sol = solve(oprob, Tsit5())
+plot(sol)
+```
 
-In this tutorial we'll illustrate how to make use of constraint equations and
-events. Let's consider a model of a cell with volume $V(t)$ that grows at a rate
-$\lambda$. For now we'll assume the cell can grow indefinitely. We'll also keep
-track of one protein $P(t)$, which is produced at a rate proportional $V$ and
-can be degraded.
 
 ## [Coupling ODE constraints via extending a system](@id constraint_equations_coupling_constraints)
 
-There are several ways we can create our Catalyst model with the two reactions
-and ODE for $V(t)$. One approach is to use compositional modeling, create
-separate `ReactionSystem`s and `ODESystem`s with their respective components,
-and then extend the `ReactionSystem` with the `ODESystem`. Let's begin by
-creating these two systems.
+Finally, we could also construct our model by using compositional modeling. Here
+we create separate `ReactionSystem`s and `ODESystem`s with their respective
+components, and then extend the `ReactionSystem` with the `ODESystem`. Let's
+begin by creating these two systems.
 
-Here, to create differentials with respect to time (for our differential equations), we must import the time differential operator from Catalyst. We do this through `D = default_time_deriv()`. Here, `D(V)` denotes the differential of the variable `V` with respect to time.
+Here, to create differentials with respect to time (for our differential
+equations), we must import the time differential operator from Catalyst. We do
+this through `D = default_time_deriv()`. Here, `D(V)` denotes the differential
+of the variable `V` with respect to time.
 
 ```@example ceq1
 using Catalyst, OrdinaryDiffEqTsit5, Plots
@@ -72,27 +129,6 @@ sol = solve(oprob, Tsit5())
 plot(sol)
 ```
 
-## Coupling ODE constraints via directly building a `ReactionSystem`
-As an alternative to the previous approach, we could have constructed our
-`ReactionSystem` all at once by directly using the symbolic interface:
-```@example ceq2
-using Catalyst, OrdinaryDiffEqTsit5, Plots
-t = default_t()
-D = default_time_deriv()
-
-@parameters λ = 1.0
-@variables V(t) = 1.0
-eq = D(V) ~ λ * V
-rx1 = @reaction $V, 0 --> P
-rx2 = @reaction 1.0, P --> 0
-@named growing_cell = ReactionSystem([rx1, rx2, eq], t)
-setdefaults!(growing_cell, [:P => 0.0])
-growing_cell = complete(growing_cell)
-
-oprob = ODEProblem(growing_cell, [], (0.0, 1.0))
-sol = solve(oprob, Tsit5())
-plot(sol)
-```
 
 ## [Adding events](@id constraint_equations_events)
 Our current model is unrealistic in assuming the cell will grow exponentially