Merge pull request #3 from SciML/master

updating all commits to local repo

.codecov.yml

@@ -1,3 +1 @@
 comment: false
-ignore:
-  - "src/graphs.jl"         # ignore graphs for testing until GraphViz jll works robustly

.github/workflows/Documentation.yml

@@ -14,7 +14,7 @@ jobs:
       - uses: actions/checkout@v2
       - uses: julia-actions/setup-julia@latest
         with:
-          version: '1.4'
+          version: '1'
       - name: Install dependencies
         run: julia --project=docs/ -e 'using Pkg; Pkg.develop(PackageSpec(path=pwd())); Pkg.instantiate()'
       - name: Build and deploy

Project.toml

@@ -1,27 +1,33 @@
 name = "Catalyst"
 uuid = "479239e8-5488-4da2-87a7-35f2df7eef83"
-version = "6.8.0"
+version = "6.12.1"
 
 [deps]
 DataStructures = "864edb3b-99cc-5e75-8d2d-829cb0a9cfe8"
+DocStringExtensions = "ffbed154-4ef7-542d-bbb7-c09d3a79fcae"
 Latexify = "23fbe1c1-3f47-55db-b15f-69d7ec21a316"
 MacroTools = "1914dd2f-81c6-5fcd-8719-6d5c9610ff09"
 ModelingToolkit = "961ee093-0014-501f-94e3-6117800e7a78"
 Parameters = "d96e819e-fc66-5662-9728-84c9c7592b0a"
 Reexport = "189a3867-3050-52da-a836-e630ba90ab69"
+Requires = "ae029012-a4dd-5104-9daa-d747884805df"
 
 [compat]
 DataStructures = "0.18"
-Latexify = "0.14"
+DocStringExtensions = "0.8"
+Latexify = "0.14, 0.15"
 MacroTools = "0.5.5"
-ModelingToolkit = "5.6.1"
+ModelingToolkit = "5.14.0"
 Parameters = "0.12"
 Reexport = "0.2, 1.0"
-julia = "1.3"
+Requires = "1.0"
+julia = "1.5"
 
 [extras]
 DiffEqBase = "2b5f629d-d688-5b77-993f-72d75c75574e"
 DiffEqJump = "c894b116-72e5-5b58-be3c-e6d8d4ac2b12"
+Graphviz_jll = "3c863552-8265-54e4-a6dc-903eb78fde85"
+LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
 OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed"
 Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"
 SafeTestsets = "1bc83da4-3b8d-516f-aca4-4fe02f6d838f"
@@ -30,7 +36,6 @@ Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2"
 SteadyStateDiffEq = "9672c7b4-1e72-59bd-8a11-6ac3964bc41f"
 StochasticDiffEq = "789caeaf-c7a9-5a7d-9973-96adeb23e2a0"
 Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
-UnPack = "3a884ed6-31ef-47d7-9d2a-63182c4928ed"
 
 [targets]
-test = ["DiffEqBase", "DiffEqJump", "OrdinaryDiffEq", "Random", "SafeTestsets", "StableRNGs", "Statistics", "SteadyStateDiffEq", "StochasticDiffEq", "Test", "UnPack"]
+test = ["DiffEqBase", "DiffEqJump", "Graphviz_jll", "LinearAlgebra", "OrdinaryDiffEq", "Random", "SafeTestsets", "StableRNGs", "Statistics", "SteadyStateDiffEq", "StochasticDiffEq", "Test"]

README.md

@@ -22,6 +22,34 @@ be used with solvers throughout the broader [SciML](https://sciml.ai) ecosystem,
 including higher-level SciML packages (e.g., for sensitivity analysis, parameter
 estimation, machine learning applications, etc.).
 
+## New in 6.11
+*1.* Plain text arrows "<--" and "<-->" for backward and reversible reactions are
+   available if using Julia 1.6:
+```julia
+rn = @reaction_network begin 
+  (k1,k2), A + B <--> C
+  k3, 0 <-- C
+end k1 k2 k3
+```
+*2.* Reaction networks can be named
+```julia
+rn = @reaction_network Reversible_Reaction begin
+  k1, A --> B
+  k2, B --> A
+  end k1 k2 
+nameof(rn) == :Reversible_Reaction
+```
+Note, empty networks can no longer be created with parameters, i.e. only
+```julia
+rn = @reaction_network          # uses a randomly generated name
+rn = @reaction_network MyName   # is named MyName
+```
+are allowed.
+
+*3.* Compositional modeling with generated `ODESystem`s, see
+[here](https://github.com/SciML/ModelingToolkit.jl/blob/master/test/reactionsystem_components.jl)
+for an example that composes three gene modules to make the repressilator.
+
 ## Tutorials and Documentation
 
 For information on using the package, [see the stable

docs/src/api/catalyst_api.md

@@ -11,8 +11,9 @@ be converted to other `ModelingToolkit.AbstractSystem`s, including a
 `ModelingToolkit.ODESystem`, `ModelingToolkit.SDESystem`, or
 `ModelingToolkit.JumpSystem`.
 
-An empty network can be generated using [`@reaction_network`](@ref) with no arguments or
-the [`make_empty_network`](@ref) function. These can then be extended
+An empty network can be generated using [`@reaction_network`](@ref) with no
+arguments (or one argument to name the system), or the
+[`make_empty_network`](@ref) function. These can then be extended
 programmatically using [`addspecies!`](@ref), [`addparam!`](@ref), and
 [`addreaction!`](@ref).
 
@@ -52,6 +53,7 @@ dependents
 dependants
 substoichmat
 prodstoichmat
+netstoichmat
 ```
 
 ## Functions to extend a Network

docs/src/tutorials/advanced.md

@@ -12,8 +12,8 @@ rn = @reaction_network begin
 end
 ```
 In some cases, it may be necessary or desirable to register functions with
-ModelingToolkit before their use in Catalyst, see the discussion
-[here](https://mtk.sciml.ai/stable/tutorials/symbolic_functions/#Registering-Functions-1).
+Symbolics.jl before their use in Catalyst, see the discussion
+[here](https://symbolics.juliasymbolics.org/dev/manual/functions/).
 
 #### Ignoring mass action kinetics
 While generally one wants the reaction rate to use the law of mass action, so

docs/src/tutorials/basics.md

@@ -65,12 +65,13 @@ rn = @reaction_network begin
   1.0, XY ← X + Y      
 end
 ```
-*Note, currently Julia's parser does not support `<--`, `<->` or `<-->`, so that
-`-->` is the only supported plain text arrow.*
+*On Julia 1.6 and up the plain text arrows `<--` (for backward reactions) and
+`<-->` (for reversible reactions) also work. Note, these are not available on
+earlier Julia versions.*
 
 #### Using bi-directional arrows
 Bi-directional unicode arrows can be used to designate a reaction that goes two
-ways. These two models are equivalent:
+ways. These three models are equivalent:
 ```julia
 rn = @reaction_network begin
   2.0, X + Y → XY             

docs/src/tutorials/generated_systems.md

@@ -4,13 +4,13 @@ has a number of fields that can be accessed directly or via the [Catalyst.jl
 API](@ref) (the recommended route). Below we list these components, with the recommended
 API method listed first:
 
-* [`species(rn)`](@ref), `states(rn)` and `rn.states` is a vector of all the chemical
+* [`species(rn)`](@ref) and `states(rn)` is a vector of all the chemical
   species within the system, each represented as a `ModelingToolkit.Term`.
-* [`params(rn)`](@ref), `parameters(rn)` and `rn.ps` is a vector of all the parameters
+* [`params(rn)`](@ref) and `parameters(rn)` is a vector of all the parameters
   within the system, each represented as a `ModelingToolkit.Sym`.
-* [`reactions(rn)`](@ref), `equations(rn)` and `rn.eqs` is a vector of all the
+* [`reactions(rn)`](@ref) and `equations(rn)` is a vector of all the
   `Reaction`s within the system.
-* `independent_variable(rn)` and `rn.iv` are the independent variable of the
+* `independent_variable(rn)` and `ModelingToolkit.get_iv(rn)` are the independent variable of the
   system, usually `t` for time, represented as a `ModelingToolkit.Sym`.
 
 Each `Reaction` within `reactions(rn)` has a number of subfields. For `rx` a
@@ -33,7 +33,7 @@ Each `Reaction` within `reactions(rn)` has a number of subfields. For `rx` a
   non-filled arrows and should ignore mass action kinetics. `false` by default.
 
 Empty `ReactionSystem`s can be generated via [`make_empty_network`](@ref) or
-[`@reaction_network`](@ref) with no arguments. `ReactionSystem`s can be
-programmatically extended using [`addspecies!`](@ref), [`addparam!`](@ref),
-[`addreaction!`](@ref), [`@add_reactions`](@ref), or composed using `merge` and
-`merge!`.
+[`@reaction_network`](@ref) with no arguments (giving one argument to the latter
+will specify a system name). `ReactionSystem`s can be programmatically extended
+using [`addspecies!`](@ref), [`addparam!`](@ref), [`addreaction!`](@ref),
+[`@add_reactions`](@ref), or composed using `merge` and `merge!`.

src/Catalyst.jl

@@ -1,11 +1,29 @@
+"""
+$(DocStringExtensions.README)
+"""
 module Catalyst
 
+using DocStringExtensions
 using Reexport, ModelingToolkit
-using ModelingToolkit: Symbolic, value, istree
+using ModelingToolkit: Symbolic, value, istree, get_states, get_ps, get_iv, get_systems, get_eqs
+const DEFAULT_IV = (@parameters t)[1]
 @reexport using ModelingToolkit
 import MacroTools
 import Base: (==), merge!, merge
-using Latexify
+using Latexify, Requires
+
+# as used in Catlab
+const USE_GV_JLL = Ref(false)
+function __init__()
+    @require Graphviz_jll="3c863552-8265-54e4-a6dc-903eb78fde85" begin
+      USE_GV_JLL[] = true
+      let cfg = joinpath(Graphviz_jll.artifact_dir, "lib", "graphviz", "config6")
+        if !isfile(cfg)
+          Graphviz_jll.dot(path -> run(`$path -c`))
+        end
+      end
+    end
+  end
 
 const ExprValues = Union{Expr,Symbol,Float64,Int}
 
@@ -19,7 +37,7 @@ export @reaction_network, @add_reactions
 include("networkapi.jl")
 export species, params, reactions, speciesmap, paramsmap, numspecies, numreactions, numparams
 export make_empty_network, addspecies!, addparam!, addreaction!
-export dependants, dependents, substoichmat, prodstoichmat
+export dependants, dependents, substoichmat, prodstoichmat, netstoichmat
 
 # for Latex printing of ReactionSystems
 include("latexify_recipes.jl")

src/graphs.jl

@@ -36,11 +36,11 @@ as_attributes(d::AbstractDict) =
 @with_kw_noshow struct Graph <: Expression
   name::String
   directed::Bool
-  prog::String="dot"
-  stmts::Vector{Statement}=Statement[]
-  graph_attrs::Attributes=Attributes()
-  node_attrs::Attributes=Attributes()
-  edge_attrs::Attributes=Attributes()
+  prog::String = "dot"
+  stmts::Vector{Statement} = Statement[]
+  graph_attrs::Attributes = Attributes()
+  node_attrs::Attributes = Attributes()
+  edge_attrs::Attributes = Attributes()
 end
 
 Graph(name::String, stmts::Vector{Statement}; kw...) =
@@ -53,11 +53,11 @@ Digraph(name::String, stmts::Vararg{Statement}; kw...) =
   Graph(; name=name, directed=true, stmts=collect(stmts), kw...)
 
 @with_kw_noshow struct Subgraph <: Statement
-  name::String="" # Subgraphs can be anonymous
-  stmts::Vector{Statement}=Statement[]
-  graph_attrs::Attributes=Attributes()
-  node_attrs::Attributes=Attributes()
-  edge_attrs::Attributes=Attributes()
+  name::String = "" # Subgraphs can be anonymous
+  stmts::Vector{Statement} = Statement[]
+  graph_attrs::Attributes = Attributes()
+  node_attrs::Attributes = Attributes()
+  edge_attrs::Attributes = Attributes()
 end
 
 Subgraph(stmts::Vector{Statement}; kw...) = Subgraph(; stmts=stmts, kw...)
@@ -106,11 +106,16 @@ For bindings to the Graphviz C API, see the the package
 GraphViz.jl is unmaintained.
 """
 function run_graphviz(io::IO, graph::Graph; prog::Union{String,Nothing}=nothing,
-                      format::String="json0")
+  format::String="json0")
   if isnothing(prog)
     prog = graph.prog
   end
-  @assert prog in ("dot","neato","fdp","sfdp","twopi","circo")
+  @assert prog in ("dot", "neato", "fdp", "sfdp", "twopi", "circo")
+  if USE_GV_JLL[]
+    print("here!!!")
+    fun = getfield(Graphviz_jll, Symbol(prog))
+    prog = fun(identity)
+  end
   open(`$prog -T$format`, io, write=true) do gv
     pprint(gv, graph)
   end
@@ -138,11 +143,11 @@ function pprint(io::IO, graph::Graph, n::Int)
   print(io, graph.directed ? "digraph " : "graph ")
   print(io, graph.name)
   println(io, " {")
-  pprint_attrs(io, graph.graph_attrs, n+2; pre="graph", post=";\n")
-  pprint_attrs(io, graph.node_attrs, n+2; pre="node", post=";\n")
-  pprint_attrs(io, graph.edge_attrs, n+2; pre="edge", post=";\n")
+  pprint_attrs(io, graph.graph_attrs, n + 2; pre="graph", post=";\n")
+  pprint_attrs(io, graph.node_attrs, n + 2; pre="node", post=";\n")
+  pprint_attrs(io, graph.edge_attrs, n + 2; pre="edge", post=";\n")
   for stmt in graph.stmts
-    pprint(io, stmt, n+2, directed=graph.directed)
+    pprint(io, stmt, n + 2, directed=graph.directed)
     println(io)
   end
   indent(io, n)
@@ -153,16 +158,16 @@ function pprint(io::IO, subgraph::Subgraph, n::Int; directed::Bool=false)
   indent(io, n)
   if isempty(subgraph.name)
     println(io, "{")
-  else
+        else
     print(io, "subgraph ")
     print(io, subgraph.name)
     println(io, " {")
   end
-  pprint_attrs(io, subgraph.graph_attrs, n+2; pre="graph", post=";\n")
-  pprint_attrs(io, subgraph.node_attrs, n+2; pre="node", post=";\n")
-  pprint_attrs(io, subgraph.edge_attrs, n+2; pre="edge", post=";\n")
+  pprint_attrs(io, subgraph.graph_attrs, n + 2; pre="graph", post=";\n")
+  pprint_attrs(io, subgraph.node_attrs, n + 2; pre="node", post=";\n")
+  pprint_attrs(io, subgraph.edge_attrs, n + 2; pre="edge", post=";\n")
   for stmt in subgraph.stmts
-    pprint(io, stmt, n+2, directed=directed)
+    pprint(io, stmt, n + 2, directed=directed)
     println(io)
   end
   indent(io, n)
@@ -185,7 +190,7 @@ function pprint(io::IO, node::NodeID, n::Int)
   if !isempty(node.anchor)
     print(io, ":")
     print(io, node.anchor)
-  end
+end
 end
 
 function pprint(io::IO, edge::Edge, n::Int; directed::Bool=false)
@@ -204,7 +209,7 @@ function pprint_attrs(io::IO, attrs::Attributes, n::Int=0;
                       pre::String="", post::String="")
   if !isempty(attrs)
     indent(io, n)
-    print(io, pre)
+        print(io, pre)
     print(io, " [")
     for (i, (key, value)) in enumerate(attrs)
       if (i > 1) print(io, ",") end
@@ -227,16 +232,16 @@ indent(io::IO, n::Int) = print(io, " "^n)
 # https://github.com/mehalter/Petri.jl
 #######################################################################
 
-graph_attrs = Attributes(:rankdir=>"LR")
-node_attrs  = Attributes(:shape=>"plain", :style=>"filled", :color=>"white")
-edge_attrs  = Attributes(:splines=>"splines")
+graph_attrs = Attributes(:rankdir => "LR")
+node_attrs  = Attributes(:shape => "plain", :style => "filled", :color => "white")
+edge_attrs  = Attributes(:splines => "splines")
 
 function edgify(δ, i, reverse::Bool)
     attr = Attributes()
     return map(δ) do p        
         val = String(p[1].f.name)
       weight = "$(p[2])"
-      attr = Attributes(:label=>weight, :labelfontsize=>"6")
+      attr = Attributes(:label => weight, :labelfontsize => "6")
       return Edge(reverse ? ["rx_$i", "$val"] :
                             ["$val", "rx_$i"], attr)
     end
@@ -245,7 +250,7 @@ end
 # make distinguished edge based on rate constant
 function edgifyrates(rxs, specs)
     es = Edge[]
-    for (i,rx) in enumerate(rxs)
+    for (i, rx) in enumerate(rxs)
         deps = rx.rate isa Number ? Any[] : get_variables(rx.rate, specs) 
         for dep in deps
             val = String(dep.f.name)
@@ -277,19 +282,21 @@ Notes:
 function Graph(rn::ReactionSystem)
     rxs   = reactions(rn)
     specs = species(rn)
-    statenodes = [Node(string(s.f.name), Attributes(:shape=>"circle", :color=>"#6C9AC3")) for s in specs]
-    transnodes = [Node(string("rx_$i"), Attributes(:shape=>"point", :color=>"#E28F41", :width=>".1")) for (i,r) in enumerate(rxs)]
+    statenodes = [Node(string(s.f.name), Attributes(:shape => "circle", :color => "#6C9AC3")) for s in specs]
+    transnodes = [Node(string("rx_$i"), Attributes(:shape => "point", :color => "#E28F41", :width => ".1")) for (i, r) in enumerate(rxs)]
 
     stmts = vcat(statenodes, transnodes)
-    edges = map(enumerate(rxs)) do (i,r)
-      vcat(edgify(zip(r.substrates,r.substoich), i, false),
-           edgify(zip(r.products,r.prodstoich), i, true))
+    edges = map(enumerate(rxs)) do (i, r)
+      vcat(edgify(zip(r.substrates, r.substoich), i, false),
+           edgify(zip(r.products, r.prodstoich), i, true))
     end
     es = edgifyrates(rxs, specs)
     (!isempty(es)) && push!(edges, es)
 
-    stmts = vcat(stmts, collect(flatten(edges)))
-    g = Digraph("G", stmts; graph_attrs=graph_attrs, node_attrs=node_attrs, edge_attrs=edge_attrs)
+    stmts2 = Vector{Statement}()
+    append!(stmts2, stmts)
+    append!(stmts2, collect(flatten(edges)))
+    g = Digraph("G", stmts2; graph_attrs=graph_attrs, node_attrs=node_attrs, edge_attrs=edge_attrs)
     return g
 end