1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
|
Equations
=========
An Equation is a set of single lines in a string:
(1) ``dx/dt = f : unit`` (differential equation)
(2) ``x = f : unit`` (equation)
(3) ``x = y`` (alias)
(4) ``x : unit`` (parameter)
The equations may be defined on multiple lines with the character \.
Comments using # may also be included.
Two special variables are defined: t (time) and xi (white noise).
Ultimately, it should be possible (using Sympy) to define equations implicitly,
e.g.: 'tau*dv/dt=-v : unit' (although it makes unit specification ambiguous).
An equation can be seen as a set of functions or code and a namespace to evaluate
them. A key part of object construction is the construction of the namespace
(after parsing).
Namespace construction
----------------------
The namespaces are stored in eq._namespace. Each equation (string) has a specific namespace.
Proposition for a simplification: there could be just one namespace per Equation object rather
than per string. Possible conflicts would be dealt with when equations are added (with prefix
as when inserting static variables, see below).
Variable substitution
~~~~~~~~~~~~~~~~~~~~~
These are simply string substitutions.
* Equation('dv/dt=-v/tau:volt',tau='taum')
The name of the variable (tau) is changed in the string to taum.
* Equation('dv/dt=-v/tau:volt',tau=None)
The name of the variable (tau) is changed in the string to a unique identifier.
Explicit namespace
~~~~~~~~~~~~~~~~~~
* Equation('dv/dt=-v/tau:volt',tau=2*ms)
The namespace is explicitly given: {'tau':2*ms}. In this case, Brian does not try
to build a namespace "magically", so the namespace must be exhaustive.
Units need not be passed.
Implicit namespace
~~~~~~~~~~~~~~~~~~
* Equation('dv/dt=-v/tau:volt')
The namespace is built from the globals and locals in the caller's frame.
For each identifier in the string, the name is looked up in:
1) locals of caller,
2) globals of caller,
3) globals of equations.py module (typically units).
Identifiers can be any Python object (for example functions).
Issues
~~~~~~
* Special variables (xi and t) are not taken into account
at this stage, i.e., they are integrated in the namespace if present.
This should probably be fixed and a warning should be raised.
A warning is raised for t at the preparation stage (see below).
* If identifiers are not found, then no error is raised. This is to allow
equations to be built in several pieces, which is useful in particular for
library objects.
* If an identifier is found whose name is the same as the name of a variable,
then no error is raised here and it is included in the namespace. This is difficult
to avoid in the case when equations are built in several pieces (e.g. the conflict
appears only when the pieces are put together). A warning is issued at the
preparation stage (see below).
Attributes after initialisation
-------------------------------
After initialisation, an Equation object contains:
* a namespace (_namespace)
* a dictionary of units for all variables (_units)
* a dictionary of strings corresponding to each variable (right hand side of each
equation), including parameters and aliases (_string). Parameters are defined as differential
equations with RHS 0*unit/second. All comments are removed and multiline strings are
concatenated.
* a list of variables of non-differential equations (_eq_names)
* a list of variables of differential equations, including parameters (_diffeq_names)
* a list of variables of differential equations, excluding parameters (_diffeq_names_nonzero)
* a dictionary of aliases (_alias), mapping a variable name to its alias
There is no explicit list of parameters, maybe it should be added.
Nothing more is done at initialisation time (no units checking, etc).
The reason is that the equation set might not be complete at this time, in the case when
equations are built in several pieces. Various checks are done using the prepare() method.
Finalisation (prepare())
------------------------
The Equation object is finalised by an explicit call to the prepare() method.
Finding Vm
~~~~~~~~~~
The first step is to find the name of the membrane potential variable (getVm()).
This is useful when the variable name for threshold or reset is not given (e.g. threshold=10*mV).
The method looks for one these names: 'v','V','vm','Vm'. If one is present, it is the
membrane potential variable. If none or more than one is present, no variable is found.
If it is found, the corresponding variable is swapped with the first variable in the
_diffeq_names list (note: not in the _diffeq_names_nonzero list). Otherwise, nothing happens.
This way, the first variable in the list is the membrane potential.
Possibly, a warning could be issued if it is not found. The problem it might issue
warnings too often. A better way would be to issue warnings only when threshold and reset
are used ambiguously (i.e., no Vm found and more than 1 variable).
Cleaning the namespace
~~~~~~~~~~~~~~~~~~~~~~
Then variables and t are removed from the namespace if present (N.B.: xi does not appear to be
removed), and warnings are issued using log_warn (method clean_namespace()).
Compiling functions
~~~~~~~~~~~~~~~~~~~
This is done by the compile_functions() method.
Python functions are created from the string definition of equations.
For each equation/differential equation, the list of identifiers is obtained from the string definition,
then only those referring to variables are kept. A Python lambda function of these remaining identifiers is then
compiled (using eval) and put in the _function dictionary.
Compiled functions are used for:
* checking units
* obtaining the list of arguments (this could be done independently)
* state updates
This step might be avoided and replaced by eval calls. It might actually be a little simpler because
arguments would be replaced by namespace. It seems to be faster with the current implementation,
but the string could be compiled with compile() (then evaluated in the relevant namespace).
Besides, with the way it is currently evaluated in the Euler update: ``f(*[S[var] for var in f.func_code.co_varnames])``,
it is not faster than direct evaluation in the namespace.
Checking units
~~~~~~~~~~~~~~
This is done by the check_units() method.
First, the static equations are ordered (see next section).
To check the units of a static equation, one calls the associated function (giving the RHS) where the
arguments are units (e.g., 1*volt for v, etc.) and adds the units of the LHS. A dimension error is raised
if it is not homogeneous. Currently, the message states "The differential equation is not homogeneous" but it
should be adapted to non-differential equations. One problem with this way of checking units is that the RHS function
may not be defined at the point it is checked.
Differential equations are checked in the same way, with two specificities: the units of RHS should be the units
of the variable divided by second (dx/dt), and noise (xi) has units of second**-.5 (this is put in the globals of
the function, which might not be a very clean way to do it).
Ordering static variables
~~~~~~~~~~~~~~~~~~~~~~~~~
It seems that this method (set_eq_order()) is already called by check_units() and therefore it is probably
not necessary to call it here.
This method computes the dependency graph of (static) equations on other static variables,
which must have no cycle (otherwise an error is raised). From that graph, an update list is built and
put in _eq_names. Then for each variable (static or differential), the list of dependent static variables is built
and sorted in update order. The result is put in the _dependencies dictionary.
This is a necessary step to calculate the RHS of any equation: it gives the ordered list of static variables
to calculate first before calculating the RHS.
Inserting static variables into differential equations
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
The value of static variables are then replaced by their string value (RHS) in all differential equations
(substitute_eq()). The previous step (ordering) ensures that the result if correct and does not depend on
static variables anymore.
To avoid namespace conflicts, all identifiers in the namespace of a static variable is augmented by a
prefix: name+'_' (e.g. 'x_y' for identifier y in equation 'x=2*y'). Then namespaces are merged.
It might not be optimal to do it in this way, because some of calculations will be done several times in
an update step. It might be better to keep the static variables separate.
Recompiling functions
~~~~~~~~~~~~~~~~~~~~~
Functions are then recompiled so that differential equations are now independent of static variables.
Checking free variables
~~~~~~~~~~~~~~~~~~~~~~~
Finally, the list of undefined identifiers is checked (free_variables()) and a warning is issued if
any is found.
Freezing
--------
Freezing is done by calling compile_functions(freeze=True). Each string expression is then frozen
with optimiser.freeze(), which replaces identifiers by their float value. This step does not necessarily
succeed, in which case a warning (not an error) is issued.
Adding Equation objects
-----------------------
Adding equations consists simply in merging the lists/dictionaries of variables, namespaces, strings, units
and functions. Conflicts raise an error.
This step must precede preparation of the object.
|
