Merge pull request #786 from SciML/use_mtk_default_t_and_D

Use MTK's default `t` and `D` throughout repository

Author: TorkelE

docs/src/api/catalyst_api.md

@@ -36,8 +36,8 @@ models, and stochastic chemical kinetics jump process models.
 
 ```@example ex1
 using Catalyst, DifferentialEquations, Plots
+t = default_t()
 @parameters β γ
-@variables t
 @species S(t) I(t) R(t)
 
 rxs = [Reaction(β, [S,I], [I], [1,1], [2])

docs/src/catalyst_applications/simulation_structure_interfacing.md

@@ -104,8 +104,9 @@ Finally, we note that we cannot change the values of solution unknowns or parame
 Catalyst is built on an *intermediary representation* implemented by (ModelingToolkit.jl)[https://github.com/SciML/ModelingToolkit.jl]. ModelingToolkit is a modelling framework where one first declares a set of symbolic variables and parameters using e.g.
 ```@example ex2
 using ModelingToolkit
+t = default_t()
 @parameters σ ρ β
-@variables t x(t) y(t) z(t)
+@variables x(t) y(t) z(t)
 nothing # hide
 ```
 and then uses these to build systems of equations. Here, these symbolic variables (`x`, `y`, and `z`) and parameters (`σ`, `ρ`, and `β`) can be used to interface a `problem`, `integrator`, and `solution` object (like we did previously, but using Symbols, e.g. `:X`). Since Catalyst models are built on ModelingToolkit, these models also contain similar symbolic variables and parameters.

docs/src/catalyst_functionality/chemistry_related_functionality.md

@@ -10,7 +10,7 @@ While Catalyst has primarily been designed around the modelling of biological sy
 We will first show how to create compound species through [programmatic model construction](@ref programmatic_CRN_construction), and then demonstrate using the DSL. To create a compound species, use the `@compound` macro, first designating the compound, followed by its components (and their stoichiometries). In this example, we will create a CO₂ molecule, consisting of one C atom and two O atoms. First, we create species corresponding to the components:
 ```@example chem1
 using Catalyst
-@variables t
+t = default_t()
 @species C(t) O(t) 
 ```
 Next, we create the `CO2` compound species:
@@ -97,7 +97,8 @@ In all of these cases, the side to the left of the `~` must be enclosed within `
 ### Compounds with multiple independent variables
 While we generally do not need to specify independent variables for compound, if the components (together) have more than one independent variable, this have to be done:
 ```@example chem1
-@variables t s
+t = default_t()
+@variables s
 @species N(s) O(t) 
 @compound NO2(t,s) ~ N + 2O
 ```
@@ -117,7 +118,7 @@ which correctly finds the (rather trivial) solution `C + 2O --> CO2`. Here we no
 Let us consider a more elaborate example, the reaction between ammonia (NH₃) and oxygen (O₂) to form nitrogen monoxide (NO) and water (H₂O). Let us first create the components and the unbalanced reaction:
 ```@example chem2
 using Catalyst # hide
-@variables t
+t = default_t()
 @species N(t) H(t) O(t) 
 @compounds begin
     NH3 ~ N + 3H

docs/src/catalyst_functionality/compositional_modeling.md

@@ -51,8 +51,8 @@ plot(TreePlot(rn), method=:tree, fontsize=12, nodeshape=:ellipse)
 We could also have directly constructed `rn` using the same reaction as in
 `basern` as
 ```@example ex1
+t = default_t()
 @parameters k
-@variables t
 @species A(t), B(t), C(t)
 rxs = [Reaction(k, [A,B], [C])]
 @named rn = ReactionSystem(rxs, t; systems = [newrn, newestrn])
@@ -114,7 +114,7 @@ ability to substitute the value of these variables into the DSL (see
 [Interpolation of Julia Variables](@ref dsl_description_interpolation_of_variables)). To make the repressilator we now make
 three genes, and then compose them together
 ```@example ex1
-@variables t
+t = default_t()
 @species G3₊P(t)
 @named G1 = repressed_gene(; R=ParentScope(G3₊P))
 @named G2 = repressed_gene(; R=ParentScope(G1.P))
@@ -130,7 +130,7 @@ plot(TreePlot(repressilator), method=:tree, fontsize=12, nodeshape=:ellipse)
 In building the repressilator we needed to use two new features. First, we
 needed to create a symbolic variable that corresponds to the protein produced by
 the third gene before we created the corresponding system. We did this via
-`@variables t, G3₊P(t)`. We also needed to set the scope where each repressor
+`@variables G3₊P(t)`. We also needed to set the scope where each repressor
 lived. Here `ParentScope(G3₊P)`, `ParentScope(G1.P)`, and `ParentScope(G2.P)`
 signal Catalyst that these variables will come from parallel systems in the tree
 that have the same parent as the system being constructed (in this case the

docs/src/catalyst_functionality/constraint_equations.md

@@ -23,19 +23,22 @@ 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:
+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.
 
 ```@example ceq1
 using Catalyst, DifferentialEquations, Plots
+t = default_t()
+D = default_time_deriv()
 
 # set the growth rate to 1.0
 @parameters λ = 1.0
 
 # set the initial volume to 1.0
-@variables t V(t) = 1.0
+@variables V(t) = 1.0
 
 # build the ODESystem for dV/dt
-D = Differential(t)
 eq = [D(V) ~ λ * V]
 @named osys = ODESystem(eq, t)
 
@@ -70,10 +73,11 @@ 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, DifferentialEquations, Plots
+t = default_t()
+D = default_time_deriv()
 
 @parameters λ = 1.0
-@variables t V(t) = 1.0
-D = Differential(t)
+@variables V(t) = 1.0
 eq = D(V) ~ λ * V
 rx1 = @reaction $V, 0 --> P
 rx2 = @reaction 1.0, P --> 0
@@ -102,11 +106,12 @@ advanced_simulations) tutorial.
 Let's first create our equations and unknowns/species again
 ```@example ceq3
 using Catalyst, DifferentialEquations, Plots
+t = default_t()
+D = default_time_deriv()
 
 @parameters λ = 1.0
-@variables t V(t) = 1.0
+@variables V(t) = 1.0
 @species P(t) = 0.0
-D = Differential(t)
 eq = D(V) ~ λ * V
 rx1 = @reaction $V, 0 --> $P
 rx2 = @reaction 1.0, $P --> 0
@@ -126,8 +131,8 @@ plot(sol; plotdensity = 1000)
 ```
 We can also model discrete events. Similar to our example with continuous events, we start by creating reaction equations, parameters, variables, and unknowns. 
 ```@example ceq3
+t = default_t()
 @parameters k_on switch_time k_off
-@variables t
 @species A(t) B(t)
 
 rxs = [(@reaction k_on, A --> B), (@reaction k_off, B --> A)]

docs/src/catalyst_functionality/dsl_description.md

@@ -54,7 +54,7 @@ the model. We do this by creating a mapping from each symbolic variable
 representing a chemical species to its initial value
 ```@example tut2
 # define the symbolic variables
-@variables t
+t = default_t()
 @species X(t) Y(t) Z(t) XY(t) Z1(t) Z2(t)
 
 # create the mapping
@@ -610,8 +610,8 @@ parameters outside of the macro, which can then be used within expressions in
 the DSL (see the [Programmatic Construction of Symbolic Reaction Systems](@ref programmatic_CRN_construction)
 tutorial for details on the lower-level symbolic interface). For example,
 ```@example tut2
+t = default_t()
 @parameters k α
-@variables t
 @species A(t)
 spec = A
 par = α

docs/src/catalyst_functionality/example_networks/hodgkin_huxley_equation.md

@@ -15,6 +15,8 @@ We begin by importing some necessary packages.
 using ModelingToolkit, Catalyst, NonlinearSolve
 using DifferentialEquations, Symbolics
 using Plots
+t = default_t()
+D = default_time_deriv()
 ```
 
 We'll build a simple Hodgkin-Huxley model for a single neuron, with the voltage,
@@ -51,7 +53,7 @@ We now declare the symbolic variable, `V(t)`, that will represent the
 transmembrane potential
 
 ```@example hh1
-@variables t V(t)
+@variables V(t)
 nothing # hide
 ```
 
@@ -75,7 +77,6 @@ I = I₀ * sin(2*pi*t/30)^2
 # get the gating variables to use in the equation for dV/dt
 @unpack m,n,h = hhrn
 
-Dₜ = Differential(t)
 eqs = [Dₜ(V) ~ -1/C * (ḡK*n^4*(V-EK) + ḡNa*m^3*h*(V-ENa) + ḡL*(V-EL)) + I/C]
 @named voltageode = ODESystem(eqs, t)
 nothing # hide

docs/src/catalyst_functionality/example_networks/smoluchowski_coagulation_equation.md

@@ -52,7 +52,7 @@ end
 We'll store the reaction rates in `pars` as `Pair`s, and set the initial condition that only monomers are present at ``t=0`` in `u₀map`.
 ```julia
 # unknown variables are X, pars stores rate parameters for each rx
-@variables t
+t = default_t()
 @species k[1:nr] (X(t))[1:N]
 pars = Pair.(collect(k), kv)
 

docs/src/catalyst_functionality/parametric_stoichiometry.md

@@ -35,8 +35,8 @@ rx.substrates[1],rx.substoich[1]
 We could have equivalently specified our systems directly via the Catalyst
 API. For example, for `revsys` we would could use
 ```@example s1
+t = default_t()
 @parameters k₊, k₋, m, n
-@variables t
 @species A(t), B(t)
 rxs = [Reaction(k₊, [A], [B], [m], [m*n]),
        Reaction(k₋, [B], [A])]
@@ -177,9 +177,10 @@ calculate and plot the average amount of protein (which is also plotted in the
 MomentClosure.jl
 [tutorial](https://augustinas1.github.io/MomentClosure.jl/dev/tutorials/geometric_reactions+conditional_closures/)).
 ```@example s1
+t = default_t()
 function getmean(jprob, Nsims, tv)
     Pmean = zeros(length(tv))
-    @variables t, P(t)
+    @variables P(t)
     for n in 1:Nsims
         sol = solve(jprob, SSAStepper())
         Pmean .+= sol(tv, idxs=P)

docs/src/catalyst_functionality/programmatic_CRN_construction.md

@@ -15,12 +15,12 @@ and then define symbolic variables for each parameter and species in the system
 (the latter corresponding to a `variable` or `unknown` in ModelingToolkit
 terminology)
 ```@example ex
+t = default_t()
 @parameters α K n δ γ β μ
-@variables t
 @species m₁(t) m₂(t) m₃(t) P₁(t) P₂(t) P₃(t)
 nothing    # hide
 ```
-*Note, each species is declared as a function of time!*
+Note: each species is declared as a function of time. Here, we first import the *time independent variable*, and stores it in `t`, using `t = default_t()`, and then use it to declare out species.
 
 !!! note
        For users familiar with ModelingToolkit, chemical species must be declared
@@ -117,8 +117,8 @@ Reaction(rate, nothing, [P₁,...,Pₙ], nothing, [β₁,...,βₙ])
 Finally, we note that the rate constant, `rate` above, does not need to be a
 constant or fixed function, but can be a general symbolic expression:
 ```julia
+t = default_t()
 @parameters α, β
-@variables t
 @species A(t), B(t)
 rx = Reaction(α + β*t*A, [A], [B])
 ```
@@ -133,7 +133,7 @@ reactions using the [`@reaction`](@ref) macro.
 
 For example, the repressilator reactions could also have been constructed like
 ```julia
-@variables t
+t = default_t()
 @species P₁(t) P₂(t) P₃(t)
 rxs = [(@reaction hillr($P₃,α,K,n), ∅ --> m₁),
        (@reaction hillr($P₁,α,K,n), ∅ --> m₂),
@@ -162,8 +162,9 @@ rx = @reaction hillr(P,α,K,n), A --> B
 ```
 is equivalent to
 ```julia
+t = default_t()
 @parameters P α K n
-@variables t A(t) B(t)
+@variables A(t) B(t)
 rx = Reaction(hillr(P,α,K,n), [A], [B])
 ```
 Here `(P,α,K,n)` are parameters and `(A,B)` are species.