Skip to content

Commit

Permalink
Finish everything for tagging 0.1 (#1)
Browse files Browse the repository at this point in the history
* fix missing tanh test

* finish docs of Process

* rename notimevar to no timeder

* discuss default lhs

* correct doc reference page

* add correct stuff in docstrings

* finish docs

* GOOD EXPORTS

* note similarity with ODESystem construction

* ignore all errors to build docs
  • Loading branch information
Datseris authored Feb 20, 2024
1 parent 2b8a844 commit e5a0bc7
Show file tree
Hide file tree
Showing 7 changed files with 97 additions and 48 deletions.
8 changes: 2 additions & 6 deletions docs/make.jl
Original file line number Diff line number Diff line change
Expand Up @@ -10,16 +10,12 @@ Downloads.download(
include("build_docs_with_style.jl")

pages = [
"Introduction" => "index.md",
"Overarching tutorial" => "tutorial.md",
"Contents" => "contents.md",
"Animations, GUIs, Visuals" => "visualizations.md",
"Contributor Guide" => "contributors_guide.md",
"Documentation" => "index.md",
]

build_docs_with_style(pages, ProcessBasedModelling;
authors = "George Datseris <datseris.george@gmail.com>",
# We need to remove the cross references because we don't list here
# the whole `DynamicalSystem` API...
warnonly = [:doctest, :missing_docs, :cross_references],
warnonly = true,
)
32 changes: 16 additions & 16 deletions docs/src/index.md
Original file line number Diff line number Diff line change
Expand Up @@ -12,7 +12,7 @@ In ProcessBasedModelling.jl, each variable is governed by a "process".
Conceptually this is just an equation that _defines_ the given variable.
To couple the variable with the process it is governed by, a user either defines simple equations of the form "variable = expression", or creates an instance of [`Process`](@ref) if the left-hand-side of the equation needs to be anything more complex. In either case, the variable and the expression are both _symbolic expressions_ created via ModellingToolkit.jl (more specifically, via Symbolics.jl).

Once all the processes about the physical system are collected, they are given as a `Vector` to the [`processes_to_mtkmodel`](@ref) central function. This function also defines what quantifies as a "process" in more specificity.
Once all the processes about the physical system are collected, they are given as a `Vector` to the [`processes_to_mtkmodel`](@ref) central function, similarly to how one gives a `Vector` of `Equation`s to e.g., `ModelingToolkit.ODESystem`. This function also defines what quantifies as a "process" in more specificity.

## Example

Expand All @@ -30,6 +30,7 @@ symbolically using ModelingToolkit.jl (**MTK**). We define
using ModelingToolkit
using OrdinaryDiffEq: Tsit5
@variables t # independent variable
@variables z(t) = 0.0
@variables x(t) # no default value
@variables y(t) = 0.0
Expand All @@ -38,8 +39,6 @@ ProcessBasedModelling.jl (**PBM**) strongly recommends that all defined variable

To make the equations we want, we can use MTK directly, and call
```@example MAIN
@variables t
eqs = [
Differential(t)(z) ~ x^2 - z
Differential(x) ~ 0.1y
Expand Down Expand Up @@ -71,28 +70,26 @@ end

Interestingly, the error is wrong. ``x`` is defined and has an equation, at least on the basis of our scientific reasoning. However ``y`` that ``x`` introduced does not have an equation. Moreover, in our experience these errors messages become increasingly less useful when a model has many equations and/or variables, as many variables get cited as "missing" from the variable map even when only one should be.

**PBM** resolves these problems and always gives accurate error messages. This is because on top of the variable map that MTK constructs automatically, **PBM** requires the user to implicitly create a map of variables to processes that govern said variables. **PBM** creates the map automatically, the only thing the user has to do is to define the equations in terms of what [`processes_to_mtkmodel`](@ref) wants (which are either [`Process`](@ref)es or `Equation`s as above).
**PBM** resolves these problems and always gives accurate error messages. This is because on top of the variable map that MTK constructs automatically, **PBM** requires the user to implicitly provide a map of variables to processes that govern said variables. **PBM** creates the map automatically, the only thing the user has to do is to define the equations in terms of what [`processes_to_mtkmodel`](@ref) wants (which are either [`Process`](@ref)es or `Equation`s as above).
Here is what the user defines to make the same system of equations:

```@example MAIN
processes = [
ExpRelaxation(z, x^2), # introduces x variable
ExpRelaxation(z, x^2), # introduces x variable
TimeDerivative(x, 0.1*y), # introduces y variable
y ~ z-x, # can be an equation because LHS is single variable
y ~ z - x, # can be an equation because LHS is single variable
]
```

Internally, all of these

which is the given to
which is then given to
```@example MAIN
model = processes_to_mtkmodel(processes)
model = processes_to_mtkmodel(processes; name = :example)
equations(model)
```

Notice that the resulting **MTK** model is not `structural_simplify`-ed, to allow composing it with other models.
Notice that the resulting **MTK** model is not `structural_simplify`-ed, to allow composing it with other models. By default `t` is taken as the independent variable.

Now, in contrast to before, if we "forgot" a process, **PBM** will react accordingly. For example, we forgot the 2nd process, then the construction will error informatively, telling us exactly which variable is missing, and because of which processes it is missing:
Now, in contrast to before, if we "forgot" a process, **PBM** will react accordingly. For example, if we forgot the 2nd process, then the construction will error informatively, telling us exactly which variable is missing, and because of which processes it is missing:
```@example MAIN
try
model = processes_to_mtkmodel(processes[[1, 3]])
Expand All @@ -101,7 +98,7 @@ catch e
end
```

If instead we "forgot" the ``y`` process, **PBM** will not error, but instead warn, and make ``y`` a named parameter:
If instead we "forgot" the ``y`` process, **PBM** will not error, but instead warn, and make ``y`` equal to a named parameter:
```@example MAIN
model = processes_to_mtkmodel(processes[1:2])
equations(model)
Expand All @@ -112,6 +109,7 @@ parameters(model)
```

Lastly, [`processes_to_mtkmodel`](@ref) also allows the concept of "default" processes, that can be used for introduced "process-less" variables.
Default processes like `processes` given as a 2nd argument to [`process_to_mtkmodel`](@ref).
For example,

```@example MAIN
Expand Down Expand Up @@ -163,14 +161,16 @@ This API describes how you can implement your own `Process` subtype, if the [exi

```@docs
Process
rhs
timescale
NoTimeVariability
ProcessBasedModelling.lhs_variable
ProcessBasedModelling.rhs
ProcessBasedModelling.NoTimeDerivative
ProcessBasedModelling.lhs
```

## Utility functions

```@docs
default_value
has_variable
new_derived_named_parameter
```
52 changes: 41 additions & 11 deletions src/API.jl
Original file line number Diff line number Diff line change
@@ -1,26 +1,36 @@
"""
A process subtype `p::Process` extends the functions:
- `lhs_variable(p)` which returns the variable the process describes (left-hand-side variable)
A process subtype `p::Process` extends the following unexported functions:
- `lhs_variable(p)` which returns the variable the process describes
(left-hand-side variable). There is a default implementation
`lhs_variable(p) = p.variable` if the field exists.
- `rhs(p)` which is the right-hand-side expression, i.e., the "actual" process.
- (optional) `timescale`, which defaults to [`NoTimeVariability`](@ref).
- (optional) `timescale`, which defaults to [`NoTimeDerivative`](@ref).
- (optional) `lhs(p)` which returns the left-hand-side. Let `τ = timescale(p)`.
Then `lhs(p)` behavior depends on `τ` as follows:
- Just `lhs_variable(p)` if `τ == NoTimeVariability()`.
Then default `lhs(p)` behaviour depends on `τ` as follows:
- Just `lhs_variable(p)` if `τ == NoTimeDerivative()`.
- `Differential(t)(p)` if `τ == nothing`.
- `τ_var*Differential(t)(p)` if `τ isa Union{Real, Num}`. If real,
a new named parameter `τ_var` is created that has the prefix `:τ_` and then the
lhs-variable name and has default value `τ`. Else if `Num`, `τ_var = τ` as given.
- Explicitly extend `lhs_variable` if the above do not suit you.
"""
abstract type Process end

"""
NoTimeVariability()
ProcessBasedModelling.NoTimeDerivative()
Singleton value that is the default output of the [`timescale`](@ref) function
for variables that do not vary in time autonomously (i.e., no d/dt derivative).
for variables that do not vary in time autonomously, i.e., they have no d/dt derivative
and hence the concept of a "timescale" does not apply to them.
"""
struct NoTimeDerivative end

"""
struct NoTimeVariability end
ProcessBasedModelling.lhs_variable(p::Process)
Return the variable (a single symbolic variable) corresponding to `p`.
"""
function lhs_variable(p::Process)
if !hasfield(typeof(p), :variable)
error("`lhs_variable` not defined for process $(nameof(typeof(p))).")
Expand All @@ -29,32 +39,52 @@ function lhs_variable(p::Process)
end
end

timescale(::Process) = NoTimeVariability()
"""
ProcessBasedModelling.timescale(p::Process)
Return the timescale associated with `p`. See [`Process`](@ref) for more.
"""
timescale(::Process) = NoTimeDerivative()

"""
ProcessBasedModelling.lhs(p::Process)
Return the left-hand-side of the equation that `p` represents as an `Expression`.
If [`timescale`](@ref) is implemented for `p`, typically `lhs` does not need to be as well.
See [`Process`](@ref) for more.
"""
function lhs(p::Process)
τ = timescale(p)
v = lhs_variable(p)
if isnothing(τ) # time variability exists but timescale is nonexistent (unity)
return Differential(t)(v)
elseif τ isa NoTimeVariability || iszero(τ) # no time variability
elseif τ isa NoTimeDerivative || iszero(τ) # no time variability
return v
else # τ is either Num or Real
τvar = new_derived_named_parameter(v, τ, "τ", false)
return τvar*Differential(t)(v)
end
end

"""
ProcessBasedModelling.rhs(p::Process)
Return the right-hand-side of the equation that `p` represents as an `Expression`.
See [`Process`](@ref) for more.
"""
function rhs(p::Process)
error("Right-hand side (`rhs`) is not defined for process $(nameof(typeof(p))).")
end

# Extensions for `Equation`:
rhs(e::Equation) = e.rhs
lhs(e::Equation) = lhs_variable(e)
function lhs_variable(e::Equation)
x = e.lhs
# we first check whether `x` is a variable
if !is_variable(x)
throw(ArgumentError("In given equation $(e), the left-hand-side does "*
"not represent a variable."))
"not represent a single variable."))
end
return x
end
12 changes: 2 additions & 10 deletions src/ProcessBasedModelling.jl
Original file line number Diff line number Diff line change
Expand Up @@ -17,22 +17,14 @@ include("utils.jl")
include("make.jl")
include("processes_basic.jl")


# TODO: In MAKE, make it so that if a variable does not have a process
# a constant process is created for it if it has a default value.
# Add a keyword `use_default` which would warn if no process but there is default
# and otherwise would error.
# TODO: Make an "addition process" that adds to processes
# It checks whether they target the same variable
# TODO: Package should compose with ODESystem
# so that component-based modelling can be utilized as well.


# TODO: Perhaps not don't export `t, rhs`?
# TODO: Perhaps not don't export `t`?
export t
export lhs_variable, rhs, timescale, NoTimeVariability
export Process, ParameterProcess, TimeDerivative, ExpRelaxation
export processes_to_mtkmodel
export new_derived_named_parameter
export hs_variable, default_value

end
13 changes: 11 additions & 2 deletions src/make.jl
Original file line number Diff line number Diff line change
Expand Up @@ -18,8 +18,17 @@ but they don't themselves have an assigned process.
It is expected that downstream packages that use ProcessBasedModelling.jl to make a
field-specific library implement a 1-argument version of `processes_to_mtkmodel`,
or provide a wrapper function for it, and add a default value for `default`.
## Keyword arguments
- `type = ODESystem`: the model type to make
- `name = nameof(type)`: the name of the model
- `independent = t`: the independent variable (default: `@variables t`).
`t` is also exported by ProcessBasedModelling.jl for convenience.
"""
function processes_to_mtkmodel(_processes, _default = []; type = ODESystem, name = nameof(type))
function processes_to_mtkmodel(_processes, _default = [];
type = ODESystem, name = nameof(type), independent = t
)
processes = expand_multi_processes(_processes)
default = default_dict(_default)
# Setup: obtain lhs-variables so we can track new variables that are not
Expand Down Expand Up @@ -72,7 +81,7 @@ function processes_to_mtkmodel(_processes, _default = []; type = ODESystem, name
end
end
end
sys = type(eqs, t; name)
sys = type(eqs, independent; name)
return sys
end
# version without given processes
Expand Down
9 changes: 6 additions & 3 deletions src/utils.jl
Original file line number Diff line number Diff line change
@@ -1,6 +1,3 @@
export print_system_info, has_variable, default_value
export new_derived_named_parameter, @named_parameters

"""
has_variable(eq, var)
Expand All @@ -13,6 +10,12 @@ function has_variable(eq::Equation, var)
end
has_variable(eqs, var) = any(eq -> has_variable(eq, var), eqs)

"""
default_value(x)
Return the default value of a symbolic variable `x` or `nothing`
if it doesn't have any. Return `x` if `x` is not a symbolic variable.
"""
default_value(x) = x
default_value(x::Num) = default_value(x.val)
function default_value(x::ModelingToolkit.SymbolicUtils.Symbolic)
Expand Down
19 changes: 19 additions & 0 deletions test/runtests.jl
Original file line number Diff line number Diff line change
Expand Up @@ -25,6 +25,25 @@ using OrdinaryDiffEq
return absorbed_shortwave - emitted_longwave
end

# make a new type of process
struct TanhProcess <: Process
variable
driver_variable
left
right
scale
reference
end
function ProcessBasedModelling.rhs(p::TanhProcess)
x = p.driver_variable
(; left, right, scale, reference) = p
return tanh_expression(x, left, right, scale, reference)
end
function tanh_expression(T, left, right, scale, reference)
return left + (right - left)*(1 + tanh(2(T - reference)/(scale)))*0.5
end


processes = [
TanhProcess(α, T, 0.7, 0.289, 10.0, 274.5),
TanhProcess(ε, T, 0.5, 0.41, 2.0, 288.0),
Expand Down

0 comments on commit e5a0bc7

Please sign in to comment.