Skip to content

Commit

Permalink
Merge pull request #3085 from ArnoStrouwen/docs3
Browse files Browse the repository at this point in the history
Start rewriting programmatically generating systems docs
  • Loading branch information
ChrisRackauckas authored Oct 3, 2024
2 parents c0d1112 + 39a63f9 commit 85d8d10
Showing 1 changed file with 26 additions and 18 deletions.
44 changes: 26 additions & 18 deletions docs/src/tutorials/programmatically_generating.md
Original file line number Diff line number Diff line change
@@ -1,42 +1,44 @@
# [Programmatically Generating and Scripting ODESystems](@id programmatically)

In the following tutorial we will discuss how to programmatically generate `ODESystem`s.
This is for cases where one is writing functions that generating `ODESystem`s, for example
if implementing a reader which parses some file format to generate an `ODESystem` (for example,
SBML), or for writing functions that transform an `ODESystem` (for example, if you write a
function that log-transforms a variable in an `ODESystem`).
In the following tutorial, we will discuss how to programmatically generate `ODESystem`s.
This is useful for functions that generate `ODESystem`s, for example
when you implement a reader that parses some file format, such as SBML, to generate an `ODESystem`.
It is also useful for functions that transform an `ODESystem`, for example
when you write a function that log-transforms a variable in an `ODESystem`.

## The Representation of a ModelingToolkit System

ModelingToolkit is built on [Symbolics.jl](https://symbolics.juliasymbolics.org/dev/),
a symbolic Computer Algebra System (CAS) developed in Julia. As such, all CAS functionality
is available on ModelingToolkit systems, such as symbolic differentiation, Groebner basis
is also available to be used on ModelingToolkit systems, such as symbolic differentiation, Groebner basis
calculations, and whatever else you can think of. Under the hood, all ModelingToolkit
variables and expressions are Symbolics.jl variables and expressions. Thus when scripting
a ModelingToolkit system, one simply needs to generate Symbolics.jl variables and equations
as demonstrated in the Symbolics.jl documentation. This looks like:

```@example scripting
using Symbolics
using ModelingToolkit: t_nounits as t, D_nounits as D
@variables x(t) y(t) # Define variables
using ModelingToolkit # reexports Symbolics
@variables t x(t) y(t) # Define variables
D = Differential(t)
eqs = [D(x) ~ y
D(y) ~ x] # Define an array of equations
```

However, ModelingToolkit has many higher-level features which will make scripting ModelingToolkit systems more convenient.
For example, as shown in the next section, defining your own independent variables and differentials is rarely needed.

## The Non-DSL (non-`@mtkmodel`) Way of Defining an ODESystem

Using `@mtkmodel` is the preferred way of defining ODEs with MTK. However, let us
look at how we can define the same system without `@mtkmodel`. This is useful for
defining PDESystem etc.
Using `@mtkmodel`, like in the [getting started tutorial](@ref getting_started),
is the preferred way of defining ODEs with MTK.
However generating the contents of a `@mtkmodel` programmatically can be tedious.
Let us look at how we can define the same system without `@mtkmodel`.

```@example scripting
using ModelingToolkit
using ModelingToolkit: t_nounits as t, D_nounits as D
@variables x(t) # independent and dependent variables
@parameters τ # parameters
@variables x(t) = 0.0 # independent and dependent variables
@parameters τ = 3.0 # parameters
@constants h = 1 # constants
eqs = [D(x) ~ (h - x) / τ] # create an array of equations
Expand All @@ -45,10 +47,16 @@ eqs = [D(x) ~ (h - x) / τ] # create an array of equations
# Perform the standard transformations and mark the model complete
# Note: Complete models cannot be subsystems of other models!
fol_model = structural_simplify(model)
fol = structural_simplify(model)
prob = ODEProblem(fol, [], (0.0, 10.0), [])
using DifferentialEquations: solve
sol = solve(prob)
using Plots
plot(sol)
```

As you can see, generating an ODESystem is as simple as creating the array of equations
As you can see, generating an ODESystem is as simple as creating an array of equations
and passing it to the `ODESystem` constructor.

## Understanding the Difference Between the Julia Variable and the Symbolic Variable
Expand Down

0 comments on commit 85d8d10

Please sign in to comment.