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/index.rst
+11-9Lines changed: 11 additions & 9 deletions
Display the source diff
Display the rich diff
Original file line number
Diff line number
Diff line change
@@ -7,7 +7,7 @@ The JiTCODE module
7
7
Overview
8
8
--------
9
9
10
-
JiTCODE (just-in-time compilation for ordinary differential equations) is an extension of SciPy’s `ODE`_ (`scipy.integrate.ode`) or `Solve IVP`_ (`scipy.integrate.solve_ivp`).
10
+
JiTCODE (just-in-time compilation for ordinary differential equations) is an extension of SciPy’s `ODE`_ (:class:`scipy.integrate.ode`) or `Solve IVP`_ (:func:`scipy.integrate.solve_ivp`).
11
11
Where the latter take a Python function as an argument, JiTCODE takes an iterable (or generator function or dictionary) of symbolic expressions, which it translates to C code, compiles on the fly, and uses as the function to feed into SciPy’s ODE or Solve IVP.
12
12
Symbolic expressions are mostly handled by `SymEngine`_, `SymPy`_’s compiled-backend-to-be (see `SymPy vs. SymEngine`_ for details).
13
13
@@ -40,7 +40,7 @@ As with SciPy’s ODE, this documentation assumes that the differential equation
40
40
41
41
.. math::
42
42
43
-
\dot{y} = f(t,y)
43
+
\dot{y} = f(t,y).
44
44
45
45
46
46
.. _example:
@@ -69,9 +69,11 @@ An explicit example
69
69
Details of the building process
70
70
-------------------------------
71
71
72
+
.. currentmodule:: jitcode
73
+
72
74
To generate the functions needed by SciPy’s ODE or Solve IVP, JiTCODE executes a series of distinct processing steps, each of which can be invoked through a command and tweaked as needed.
73
75
These commands execute previous steps as needed, i.e., if they have not been performed yet.
74
-
If you are happy with the default options, however, you do not need to bother with this and can just use the commands at the very end of the chain, namely `set_integrator`, `set_initial_value`, or `save_compiled`.
76
+
If you are happy with the default options, however, you do not need to bother with this and can just use the commands at the very end of the chain, namely :meth:`~_jitcode.jitcode.set_integrator`, :meth:`~_jitcode.jitcode.set_initial_value`, or :meth:`~_jitcode.jitcode.save_compiled`.
75
77
76
78
The following diagram details which command calls which other command when needed:
77
79
@@ -102,19 +104,19 @@ A more complicated example
102
104
Calculating Lyapunov exponents with `jitcode_lyap`
`jitcode_lyap` is a simple extension of `jitcode` that almost automatically handles calculating Lyapunov exponents by evolving tangent vectors [BGGS80]_.
106
-
It works just like `jitcode`, except that it generates and integrates additional differential equations for the tangent vectors.
107
-
After every call of `integrate`, the tangent vectors are orthonormalised, and the “local” Lyapunov exponents for this integration step are returned alongside with the system’s state.
107
+
:class:`~_jitcode.jitcode_lyap` is a simple extension of :class:`~_jitcode.jitcode` that almost automatically handles calculating Lyapunov exponents by evolving tangent vectors [BGGS80]_.
108
+
It works just like :class:`~_jitcode.jitcode`, except that it generates and integrates additional differential equations for the tangent vectors.
109
+
After every call of :meth:`~_jitcode.jitcode_lyap.integrate`, the tangent vectors are orthonormalised, and the “local” Lyapunov exponents for this integration step are returned alongside with the system’s state.
108
110
These can then be further processed to obtain the Lyapunov exponents.
109
111
The tangent vectors are initialised with random vectors, and you have to take care of the preiterations that the tangent vectors require to align themselves.
110
-
You also have to take care not to `integrate` for too long to avoid a numerical overflow in the tangent vectors.
112
+
You also have to take care not to :meth:`~_jitcode.jitcode_lyap.integrate` for too long to avoid a numerical overflow in the tangent vectors.
111
113
112
114
Estimates for the Lyapunov vectors are returned as well.
`jitcode_transversal_lyap` is a variant of `jitcode_lyap` that calculates Lyapunov exponents transversal to a user-defined synchronisation manifold.
129
+
:class:`~_jitcode.jitcode_transversal_lyap` is a variant of :class:`~_jitcode.jitcode_lyap` that calculates Lyapunov exponents transversal to a user-defined synchronisation manifold.
128
130
It automatically conflates the differential equations of a group of synchronised components into one equation on the synchronisation manifold.
129
131
Moreover, it transforms the equations for the tangent vectors such that the tangent vector is automatically orthogonal to the synchronisation manifold.
130
132
If you are interested in the mathematical details, please read `the accompanying paper`_.
Copy file name to clipboardExpand all lines: examples/SW_of_Roesslers.py
+6-6Lines changed: 6 additions & 6 deletions
Original file line number
Diff line number
Diff line change
@@ -12,25 +12,25 @@
12
12
\\dot{z}_i &= b + z_i (x_i -c) + k \\sum_{j=0}^N (z_j-z_i)
13
13
\\end{alignedat}
14
14
15
-
The control parameters shall be :math:`a = 0.165`, :math:`b = 0.2`, :math:`c = 10.0`, and :math:`k = 0.01`. The (frequency) parameter :math:`ω_i` shall be picked randomly from the uniform distribution on :math:`[0.8,1.0]` for each :math:`i`. :math:`A∈ℝ^{N×N}` shall be the adjacency matrix of a one-dimensional small-world network (which shall be provided by a function `small_world_network` in the following example code). So, the :math:`x` components are coupled diffusively with a small-world coupling topology, while the :math:`z` components are coupled diffusively to their mean field, with the coupling term being modulated with :math:`\\sin(t)`.
15
+
The control parameters shall be :math:`a = 0.165`, :math:`b = 0.2`, :math:`c = 10.0`, and :math:`k = 0.01`. The (frequency) parameter :math:`ω_i` shall be picked randomly from the uniform distribution on :math:`[0.8,1.0]` for each :math:`i`. :math:`A∈ℝ^{N×N}` shall be the adjacency matrix of a one-dimensional small-world network (which shall be provided by a function :any:`small_world_network` in the following example code). So, the :math:`x` components are coupled diffusively with a small-world coupling topology, while the :math:`z` components are coupled diffusively to their mean field, with the coupling term being modulated with :math:`\\sin(t)`.
16
16
17
17
Without further ado, here is the example code (`complete running example <https://raw.githubusercontent.com/neurophysik/jitcode/master/examples/SW_of_Roesslers.py>`_); highlighted lines will be commented upon below:
* The values of :math:`ω` are initialised globally (line 9). We cannot just define a function here, because the parameter is used twice for each oscillator. Moreover, if we were trying to calculate Lyapunov exponents or the Jacobian, the generator function would be called multiple times, and thus the value of the parameter would not be consistent (which would be desastrous).
27
+
* The values of :math:`ω` are initialised globally (line 10). We cannot just define a function here, because the parameter is used twice for each oscillator. Moreover, if we were trying to calculate Lyapunov exponents or the Jacobian, the generator function would be called multiple times, and thus the value of the parameter would not be consistent (which would be desastrous).
28
28
29
-
* Since we need :math:`\\sum_{j=0}^N x_j` to calculate the derivative of :math:`z` for every oscillator, it is prudent to only calculate this once. Therefore we define a helper symbol for this in lines 25 and 26, which we employ in line 34. (See the arguments of `jitcode` for details.)
29
+
* Since we need :math:`\\sum_{j=0}^N x_j` to calculate the derivative of :math:`z` for every oscillator, it is prudent to only calculate this once. Therefore we define a helper symbol for this in lines 26 and 27, which we employ in line 35. (See the arguments of `jitcode` for details.)
30
30
31
-
* In line 31, we implement :math:`\\sin(t)`. For this purpose we had to import `t` in line 1. Also, we need to use `symengine.sin` – in contrast to `math.sin` or `numpy.sin`.
31
+
* In line 32, we implement :math:`\\sin(t)`. For this purpose we had to import `t` in line 3. Also, we need to use `symengine.sin` – in contrast to :func:`math.sin` or :data:`numpy.sin`.
32
32
33
-
* As this is a large system, we use a generator function instead of a list to define :math:`f` (lines 28-35) and have the code automatically be split into chunks of 150 lines, corresponding to the equations of fifty oscillators, in line 43. (For this system, any reasonably sized multiple of 3 is a good choice for `chunk_size`; for other systems, the precise choice of the value may be more relevant.) See `large_systems` for the rationale.
33
+
* As this is a large system, we use a generator function instead of a list to define :math:`f` (lines 29-36) and have the code automatically be split into chunks of 150 lines, corresponding to the equations of fifty oscillators, in line 44. (For this system, any reasonably sized multiple of 3 is a good choice for `chunk_size`; for other systems, the precise choice of the value may be more relevant.) See `large_systems` for the rationale.
Copy file name to clipboardExpand all lines: examples/double_fhn_example.py
+6-6Lines changed: 6 additions & 6 deletions
Original file line number
Diff line number
Diff line change
@@ -5,13 +5,13 @@
5
5
6
6
.. math::
7
7
8
-
f(y) = \\left(
9
-
\\begin{matrix}
10
-
y_0 ( a-y_0 ) ( y_0-1) - y_1 + k (y_2 - y_0) \\\\
11
-
b_1 y_0 - c y_1 \\\\
12
-
y_2 ( a-y_2 ) ( y_2-1 ) - y_3 + k (y_0 - y_2)\\\\
8
+
f(y) =
9
+
\begin{bmatrix}
10
+
y_0 ( a-y_0 ) ( y_0-1) - y_1 + k (y_2 - y_0) \\
11
+
b_1 y_0 - c y_1 \\
12
+
y_2 ( a-y_2 ) ( y_2-1 ) - y_3 + k (y_0 - y_2)\\
13
13
b_2 y_2 - c y_3
14
-
\\end{matrix} \\right),
14
+
\end{bmatrix},
15
15
16
16
and :math:`a = -0.025794`, :math:`b_1 = 0.0065`, :math:`b_2 = 0.0135`, :math:`c = 0.02`, and :math:`k = 0.128`.
17
17
Then the following code integrates the above for 100000 time units (after discarding 2000 time units of transients), with :math:`y(t=0) = (1,2,3,4)`, and writes the results to :code:`timeseries.dat`:
Copy file name to clipboardExpand all lines: examples/double_fhn_transversal_lyap.py
+2-2Lines changed: 2 additions & 2 deletions
Original file line number
Diff line number
Diff line change
@@ -4,12 +4,12 @@
4
4
For instance, let’s interpret the system from `example` as two oscillators (which is what it is), one consisting of the first and second and one of the third and fourth component. Furthermore, let’s change the control parameters a bit to make the two oscillators identical. We can then calculate the transversal Lyapunov exponents to the synchronisation manifold as follows (important changes are highlighted):
Copy file name to clipboardExpand all lines: examples/lotka_volterra.py
+12-12Lines changed: 12 additions & 12 deletions
Original file line number
Diff line number
Diff line change
@@ -16,21 +16,21 @@
16
16
17
17
.. literalinclude:: ../examples/lotka_volterra.py
18
18
:start-after: example-st\u0061rt
19
-
:lines: 1-2
19
+
:lines: 1-3
20
20
:dedent: 1
21
21
22
22
… and defining the control parameters:
23
23
24
24
.. literalinclude:: ../examples/lotka_volterra.py
25
25
:start-after: example-st\u0061rt
26
-
:lines: 4-7
26
+
:lines: 5-8
27
27
:dedent: 1
28
28
29
29
The `y` that we imported from `jitcode` has to be used for the dynamical variables. However, to make our code use the same notation as the above equation, we can rename them:
30
30
31
31
.. literalinclude:: ../examples/lotka_volterra.py
32
32
:start-after: example-st\u0061rt
33
-
:lines: 9
33
+
:lines: 10
34
34
:dedent: 1
35
35
36
36
We here might as well have written `R,B = y(0),y(1)`, but the above is more handy for larger systems. Note that the following code is written such that the order of our dynamical variables does not matter.
@@ -39,7 +39,7 @@
39
39
40
40
.. literalinclude:: ../examples/lotka_volterra.py
41
41
:start-after: example-st\u0061rt
42
-
:lines: 11-14
42
+
:lines: 12-15
43
43
:dedent: 1
44
44
45
45
Note that there are other ways to define the derivative, e.g., used in the previous and following example.
@@ -48,35 +48,35 @@
48
48
49
49
.. literalinclude:: ../examples/lotka_volterra.py
50
50
:start-after: example-st\u0061rt
51
-
:lines: 16
51
+
:lines: 17
52
52
:dedent: 1
53
53
54
54
Before we start integrating, we have to choose an integrator and define initial conditions. We here choose the 5th-order Dormand–Prince method. Moreover the initial :math:`B` shall be 0.5, the initial :math:`R` shall be 0.2, and the starting time should be :math:`t=0`.
55
55
56
56
.. literalinclude:: ../examples/lotka_volterra.py
57
57
:start-after: example-st\u0061rt
58
-
:lines: 17-19
58
+
:lines: 18-20
59
59
:dedent: 1
60
60
61
61
We then define an array of time points where we want to observe the system, and a dictionary of empty lists that shall be filled with our results.
62
62
63
63
.. literalinclude:: ../examples/lotka_volterra.py
64
64
:start-after: example-st\u0061rt
65
-
:lines: 21-22
65
+
:lines: 22-23
66
66
:dedent: 1
67
67
68
68
Finally, we perform the actual integration by calling `ODE.integrate` for each of our `times`. After this, we use `ODE.y_dict` to access the current state as a dictionary and appending its values to the respective lists.
69
69
70
70
.. literalinclude:: ../examples/lotka_volterra.py
71
71
:start-after: example-st\u0061rt
72
-
:lines: 23-26
72
+
:lines: 24-27
73
73
:dedent: 1
74
74
75
75
We can now plot or analyse our data, but that’s beyond the scope of JiTCODE. So we just save it:
0 commit comments