You signed in with another tab or window. Reload to refresh your session.You signed out in another tab or window. Reload to refresh your session.You switched accounts on another tab or window. Reload to refresh your session.Dismiss alert
Copy file name to clipboardExpand all lines: docs/src/steady_state_functionality/examples/nullcline_plotting.md
+32-23Lines changed: 32 additions & 23 deletions
Display the source diff
Display the rich diff
Original file line number
Diff line number
Diff line change
@@ -1,28 +1,28 @@
1
1
# [Plotting Nullclines and Steady States in Phase Space](@id nullcline_plotting)
2
-
In this tutorial we will show how to extract a system's steady states and [nullclines](https://en.wikipedia.org/wiki/Nullcline), and how to plot these in [phase space](https://en.wikipedia.org/wiki/Phase_space). Generally, while nullclines are not directly needed for much analysis, plotting these can give some understanding of a systems steady state and stability properties. While nullclines can be "brute forced", we will here use [BifurcationKit.jl](https://github.com/bifurcationkit/BifurcationKit.jl) to rack their path exactly.
2
+
In this tutorial we will show how to extract a system's steady states and [nullclines](https://en.wikipedia.org/wiki/Nullcline), and how to plot these in [phase space](https://en.wikipedia.org/wiki/Phase_space). Generally, while nullclines are not directly needed for much analysis, plotting these can give some understanding of a system's steady state and stability properties. While nullclines can be "brute forced" (by solving their equations across grids of values), we will here use [BifurcationKit.jl](https://github.com/bifurcationkit/BifurcationKit.jl) to track their path exactly.
3
3
4
-
For a ordinary differential equation
4
+
For an ordinary differential equation
5
5
```math
6
6
\begin{aligned}
7
-
\frac{dx_1}{dt} &= f_1(x_1, x_2, ..., x_3) \\
8
-
\frac{dx_2}{dt} &= f_2(x_1, x_2, ..., x_3) \\
7
+
\frac{dx_1}{dt} &= f_1(x_1, x_2, ..., x_n) \\
8
+
\frac{dx_2}{dt} &= f_2(x_1, x_2, ..., x_n) \\
9
9
&\vdots\\
10
-
\frac{dx_n}{dt} &= f_3(x_1, x_2, ..., x_3) \\
10
+
\frac{dx_n}{dt} &= f_n(x_1, x_2, ..., x_n) \\
11
11
\end{aligned}
12
12
```
13
-
The nullclines are the curves
13
+
the nullclines are the curves
14
14
```math
15
15
\begin{aligned}
16
-
0 &= f_1(x_1, x_2, ..., x_3) \\
17
-
0 &= f_2(x_1, x_2, ..., x_3) \\
16
+
0 &= f_1(x_1, x_2, ..., x_n) \\
17
+
0 &= f_2(x_1, x_2, ..., x_n) \\
18
18
&\vdots\\
19
-
0 &= f_3(x_1, x_2, ..., x_3) \\
19
+
0 &= f_n(x_1, x_2, ..., x_n) \\
20
20
\end{aligned}
21
21
```
22
-
where the $i$'th nullclines is the curve along which $\frac{dx_i}{dt} = 0$. Generally, nullclines are primarily computed for models with 2 variable (as only here can they be easily plotted).
22
+
where the $i$'th nullclines is the curve along which $\frac{dx_i}{dt} = 0$. Generally, nullclines are primarily computed for models with 2 variables (as these can be easily plotted).
23
23
24
-
## [Computing nucclines and steady states for a bistable switch](@id nullcline_plotting_computation)
25
-
For our example we will use a simple bistable switch model, consisting of two species ($X$ and $Y$) which mutually inhibits each other through repressive hill functions. We will create our model [programmatically](@ref programmatic_CRN_construction).
24
+
## [Computing nullclines and steady states for a bistable switch](@id nullcline_plotting_computation)
25
+
For our example we will use a simple bistable switch model, consisting of two species ($X$ and $Y$) which mutually inhibit each other through repressive Hill functions. We will create our model [programmatically](@ref programmatic_CRN_construction).
26
26
```@example nullcline_plotting
27
27
using Catalyst
28
28
t = default_t()
@@ -50,11 +50,11 @@ Finally, we will compute the nullclines. We will first extract our model's stead
50
50
nlsys = convert(NonlinearSystem, bs_switch)
51
51
eqs = equations(nlsys)
52
52
```
53
-
Here, each equation is an expression of two variables ($X$ and $Y$). We wish to plot $Y$ as a function of $X$. To do this, we will substitute the species variable $X$ with the parameter $Xpar$. This will enable us to carry out bifurcation analysis of the equations' solutions as $Xpar$ is varied
53
+
Here, each equation is an expression of the two variables ($X$ and $Y$). We wish to plot $Y$ as a function of $X$. To do this, we will substitute the species variable $X$ with the parameter $Xpar$. This will enable us to carry out bifurcation analysis of the equations' solutions as $Xpar$ is varied
To input these into BifurcationKit we need to convert these into `NonlinearSystem`s.
60
60
```@example nullcline_plotting
@@ -64,7 +64,7 @@ nc_X_sys = complete(nc_X_sys)
64
64
nc_Y_sys = complete(nc_Y_sys)
65
65
nothing # hide
66
66
```
67
-
Finally, for out nullcline equations, we can use BifurcationKit's `continuation` function to track their solutions across all values of $Xpar$ (the workflow is similar to when we used `bifurcationdiagram`[here](@ref bifurcation_diagrams)).
67
+
Finally, for our nullcline equations, we can use BifurcationKit's `continuation` function to track their solutions across all values of $Xpar$ (the workflow is similar to when we used `bifurcationdiagram`[here](@ref bifurcation_diagrams)).
We are ready to create our plot. We will plot the steady states using `scatter`. We use the `steady_state_stability` function to [compute steady state stabilities](@ref steady_state_stability) (and use this to determine how to plot the steady state markers).
81
+
We are now ready to create our plot. We will plot the steady states using `scatter`. We use the `steady_state_stability` function to [compute steady state stabilities](@ref steady_state_stability) (and use this to determine how to plot the steady state markers). For the nullclines, we plot their `.x` values (the variable we solve for, $Y$) against their `.param` values (the parameter which value we vary, $Xpar$).
Here we can see how the steady states occur at the nullclines intersections.
99
104
100
105
!!! warn
101
-
BifurcationKit's `continuation` function will not detect disjoint branches, and above we can only use it because we know that each nullcline consists of a single branch. To handle disjoint nullclines, consider using [deflated continuation](@ref bifurcation_diagrams_disjoint_branches) (or possibly a package like [Contour.jl](https://github.com/JuliaGeometry/Contour.jl)).
106
+
BifurcationKit's `continuation` function will not detect disjoint branches, and here we can only use it because we know that each nullcline consists of a single branch. To handle nullclines with disjoint branches, consider using [deflated continuation](@ref bifurcation_diagrams_disjoint_branches) (or possibly a package like [Contour.jl](https://github.com/JuliaGeometry/Contour.jl)).
102
107
103
108
## [Plotting system directions in phase space](@id nullcline_plotting_directions)
104
-
One useful property of nullclines is that the sign of $dX/dt$ will only switch whenever the solution crosses the $dX/dt=0$ nullcline. This mean that, within each region defined by the nullclines, the direction of the solution remains constant. Below we use this to, for each such region, plot arrows showing the solution's direction.
105
-
109
+
One useful property of nullclines is that the sign of $dX/dt$ will only switch whenever the solution crosses the $dX/dt=0$ nullcline. This means that, within each region defined by the nullclines, the direction of the solution remains constant. Below we use this to, for each such region, plot arrows showing the solution's direction.
106
110
```@example nullcline_plotting
111
+
# Create a function which evaluates the ODE at a point in phase space.
107
112
nlprob = NonlinearProblem(complete(nlsys), [X => 0.0, Y => 0.0], ps)
0 commit comments