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
BifurcationKit computes bifurcation diagrams using the `bifurcationdiagram` function. From an initial point in the diagram, it tracks the solution (using a continuation algorithm) until the entire diagram is computed (BifurcationKit's continuation can be used for other purposes, however, this tutorial focuses on bifurcation diagram computation). The continuation settings are provided in a `ContinuationPar` structure. In this example, we will only specify three settings, `p_min` and `p_max` (which sets the minimum and maximum values over which the bifurcation parameter is varied) and `max_steps` (the maximum number of continuation steps to take as the bifurcation diagram is tracked). We wish to compute a bifurcation diagram over the interval *(2.0,20.0)*, and will use the following settings:
37
+
BifurcationKit computes bifurcation diagrams using the `bifurcationdiagram` function. From an initial point in the diagram, it tracks the solution (using a continuation algorithm) until the entire diagram is computed (BifurcationKit's continuation can be used for other purposes, however, this tutorial focuses on bifurcation diagram computation). The continuation settings are provided in a `ContinuationPar` structure. In this example, we will only specify three settings, `p_min` and `p_max` (which sets the minimum and maximum values over which the bifurcation parameter is varied) and `max_steps` (the maximum number of continuation steps to take as the bifurcation diagram is tracked). We wish to compute a bifurcation diagram over the interval $(2.0,20.0)$, and will use the following settings:
@@ -53,7 +53,7 @@ We can plot our bifurcation diagram using the Plots.jl package:
53
53
using Plots
54
54
plot(bif_dia; xguide="k1", yguide="X")
55
55
```
56
-
Here, the steady state concentration of *X* is shown as a function of *k1*'s value. Stable steady states are shown with thick lines, unstable ones with thin lines. The two fold bifurcation points are marked with "bp".
56
+
Here, the steady state concentration of $X$ is shown as a function of $k1$'s value. Stable steady states are shown with thick lines, unstable ones with thin lines. The two [fold bifurcation points](https://en.wikipedia.org/wiki/Saddle-node_bifurcation) are marked with "bp".
57
57
58
58
## Additional `ContinuationPar` options
59
59
Most of the options required by the `bifurcationdiagram` function are provided through the `ContinuationPar` structure. For full details, please read the [BifurcationKit documentation](https://bifurcationkit.github.io/BifurcationKitDocs.jl/dev/library/#BifurcationKit.ContinuationPar). However, a few common options, and how they affect the continuation computation, are described here:
@@ -75,18 +75,18 @@ nothing # hide
75
75
(however, in this case these additional settings have no significant effect on the result)
76
76
77
77
## Bifurcation diagrams with disjoint branches
78
-
Let's consider the previous case, but instead compute the bifurcation diagram over the interval *(2.0,15.0)*:
78
+
Let's consider the previous case, but instead compute the bifurcation diagram over the interval $(2.0,15.0)$:
Here, in the bistable region, we only see a single branch. The reason is that the continuation algorithm starts at our initial guess (here made at `k1 = 4.0` for `(X,Y) = (5.0,2.0)`) and tracks the diagram from there. However, with the upper bound set at `k1=15.0` the bifurcation diagram has a disjoint branch structure, preventing the full diagram from being computed by continuation alone. In this case it could be solved by increasing the bound from `k1=15.0`, however, this is not possible in all cases. In these cases, *deflation* can be used. This is described in the [BifurcationKit documentation](https://bifurcationkit.github.io/BifurcationKitDocs.jl/dev/tutorials/tutorials2/#Snaking-computed-with-deflation).
85
+
Here, in the bistable region, we only see a single branch. The reason is that the continuation algorithm starts at our initial guess (here made at $k1 = 4.0$ for $(X,Y) = (5.0,2.0)$) and tracks the diagram from there. However, with the upper bound set at $k1=15.0$ the bifurcation diagram has a disjoint branch structure, preventing the full diagram from being computed by continuation alone. In this case it could be solved by increasing the bound from $k1=15.0$, however, this is not possible in all cases. In these cases, *deflation* can be used. This is described in the [BifurcationKit documentation](https://bifurcationkit.github.io/BifurcationKitDocs.jl/dev/tutorials/tutorials2/#Snaking-computed-with-deflation).
86
86
87
87
88
88
## Systems with conservation laws
89
-
Some systems are under-determined determined at steady-state, so that for a given parameter set they have an infinite number of possible steady state solutions, preventing bifurcation diagrams from being computed. Similar to when we [compute steady states for fixed parameter values](@ref homotopy_continuation_conservation_laws), we can utilise Catalyst's ability to detect and eliminate conservation laws to resolve this issue. This requires us to provide information of the species concentrations at which we wish to compute the bifurcation diagram (to determine the values of conserved quantities). These are provided to the `BifurcationProblem` using the `u0` argument.
89
+
Some systems are under-determined at steady-state, so that for a given parameter set they have an infinite number of possible steady state solutions, preventing bifurcation diagrams from being computed. Similar to when we [compute steady states for fixed parameter values](@ref homotopy_continuation_conservation_laws), we can utilise Catalyst's ability to detect and eliminate conservation laws to resolve this issue. This requires us to provide information of the species concentrations at which we wish to compute the bifurcation diagram (to determine the values of conserved quantities). These are provided to the `BifurcationProblem` using the `u0` argument.
90
90
91
91
To illustrate this, we will create a simple model of a kinase that is produced and degraded (at rates $p$ and $d$). The kinase facilitates the phosphorylation of a protein ($X$), which is dephosphorylated at a constant rate. For this system, we will compute a bifurcation diagram, showing how the concentration of the phosphorylated protein ($Xp$) depends on the degradation rate of the kinase ($d$). We will set the total amount of protein ($X+Xp$) to $1.0$.
This bifurcation diagram does not contain any interesting features (such as bifurcation points), and only shows how the steady state concentration of $Xp$ is reduced as $d$ increases. For this example, we will note two facts:
110
-
- When providing the concentrations for computing the conserved quantities (in `u0`), we only have to designate the concentrations of species that are actually involved in conservation laws. For larger systems, determining which one are may, however, be difficult. In this case, it might still be wise to provide concentrations for all species.
111
-
- The steady state guess in `u_guess` does not actually have to fulfil the conserved concentrations provided in `u0`.
109
+
This bifurcation diagram does not contain any interesting features (such as bifurcation points), and only shows how the steady state concentration of $Xp$ is reduced as $d$ increases.
110
+
111
+
Finally, for additional clarity, we reiterate the purpose of the two `u` arguments used:
112
+
-`u_guess`: A guess of the initial steady states (which BifurcationKit uses to find its starting point). Typically, most trivial guesses work (e.g. setting all species concentrations to `1.0`). `u_guess`*does not* have to fulfil the conserved concentrations provided in `u0`.
113
+
-`u0`: Used to compute the concentrations of any conserved quantities (e.g. in our example $X + Xp = 1.0$). Technically, values are only required for species that are involved in conservation laws (in our case we do not need to provide a value for $K$). However, sometimes determining which species are actually involved in conservation laws can be difficult, and it might be easier to simply provide concentrations for all species.
112
114
113
115
!!! note
114
116
Computing bifurcation diagrams for [hierarchical models created by composing multiple systems](@ref compositional_modeling), that also contain conservations laws, is currently not supported. For these, please [flatten](@ref api_network_extension_and_modification) the system first.
0 commit comments