File: clocks.py

package info (click to toggle)
brian 2.9.0-2
  • links: PTS, VCS
  • area: main
  • in suites: forky, sid
  • size: 6,872 kB
  • sloc: python: 51,820; cpp: 2,033; makefile: 108; sh: 72
file content (233 lines) | stat: -rw-r--r-- 7,607 bytes parent folder | download | duplicates (2)
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
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
"""
Clocks for the simulator.
"""

__docformat__ = "restructuredtext en"

import numpy as np

from brian2.core.names import Nameable
from brian2.core.variables import Variables
from brian2.groups.group import VariableOwner
from brian2.units.allunits import second
from brian2.units.fundamentalunits import Quantity, check_units
from brian2.utils.logger import get_logger

__all__ = ["Clock", "defaultclock"]

logger = get_logger(__name__)


def check_dt(new_dt, old_dt, target_t):
    """
    Check that the target time can be represented equally well with the new
    dt.

    Parameters
    ----------
    new_dt : float
        The new dt value
    old_dt : float
        The old dt value
    target_t : float
        The target time

    Raises
    ------
    ValueError
        If using the new dt value would lead to a difference in the target
        time of more than `Clock.epsilon_dt` times ``new_dt`` (by default,
        0.01% of the new dt).

    Examples
    --------
    >>> from brian2 import *
    >>> check_dt(float(17*ms), float(0.1*ms), float(0*ms))  # For t=0s, every dt is fine
    >>> check_dt(float(0.05*ms), float(0.1*ms), float(10*ms))  # t=10*ms can be represented with the new dt
    >>> check_dt(float(0.2*ms), float(0.1*ms), float(10.1*ms))  # t=10.1ms cannot be represented with dt=0.2ms # doctest: +ELLIPSIS
    Traceback (most recent call last):
    ...
    ValueError: Cannot set dt from 100. us to 200. us, the time 10.1 ms is not a multiple of 200. us.
    """
    old_t = np.int64(np.round(target_t / old_dt)) * old_dt
    new_t = np.int64(np.round(target_t / new_dt)) * new_dt
    error_t = target_t
    if abs(new_t - old_t) / new_dt > Clock.epsilon_dt:
        old = str(old_dt * second)
        new = str(new_dt * second)
        t = str(error_t * second)
        raise ValueError(
            f"Cannot set dt from {old} to {new}, the "
            f"time {t} is not a multiple of {new}."
        )


class Clock(VariableOwner):
    """
    An object that holds the simulation time and the time step.

    Parameters
    ----------
    dt : float
        The time step of the simulation as a float
    name : str, optional
        An explicit name, if not specified gives an automatically generated name

    Notes
    -----
    Clocks are run in the same `Network.run` iteration if `~Clock.t` is the
    same. The condition for two
    clocks to be considered as having the same time is
    ``abs(t1-t2)<epsilon*abs(t1)``, a standard test for equality of floating
    point values. The value of ``epsilon`` is ``1e-14``.
    """

    def __init__(self, dt, name="clock*"):
        # We need a name right away because some devices (e.g. cpp_standalone)
        # need a name for the object when creating the variables
        Nameable.__init__(self, name=name)
        self._old_dt = None
        self.variables = Variables(self)
        self.variables.add_array(
            "timestep", size=1, dtype=np.int64, read_only=True, scalar=True
        )
        self.variables.add_array(
            "t",
            dimensions=second.dim,
            size=1,
            dtype=np.float64,
            read_only=True,
            scalar=True,
        )
        self.variables.add_array(
            "dt",
            dimensions=second.dim,
            size=1,
            values=float(dt),
            dtype=np.float64,
            read_only=True,
            constant=True,
            scalar=True,
        )
        self.variables.add_constant("N", value=1)
        self._enable_group_attributes()
        self.dt = dt
        logger.diagnostic(f"Created clock {self.name} with dt={self.dt}")

    @check_units(t=second)
    def _set_t_update_dt(self, target_t=0 * second):
        new_dt = self.dt_
        old_dt = self._old_dt
        target_t = float(target_t)
        if old_dt is not None and new_dt != old_dt:
            self._old_dt = None
            # Only allow a new dt which allows to correctly set the new time step
            check_dt(new_dt, old_dt, target_t)

        new_timestep = self._calc_timestep(target_t)
        # Since these attributes are read-only for normal users, we have to
        # update them via the variables object directly
        self.variables["timestep"].set_value(new_timestep)
        self.variables["t"].set_value(new_timestep * new_dt)
        logger.diagnostic(f"Setting Clock {self.name} to t={self.t}, dt={self.dt}")

    def _calc_timestep(self, target_t):
        """
        Calculate the integer time step for the target time. If it cannot be
        exactly represented (up to 0.01% of dt), round up.

        Parameters
        ----------
        target_t : float
            The target time in seconds

        Returns
        -------
        timestep : int
            The target time in integers (based on dt)
        """
        new_i = np.int64(np.round(target_t / self.dt_))
        new_t = new_i * self.dt_
        if new_t == target_t or np.abs(new_t - target_t) / self.dt_ <= Clock.epsilon_dt:
            new_timestep = new_i
        else:
            new_timestep = np.int64(np.ceil(target_t / self.dt_))
        return new_timestep

    def __repr__(self):
        return f"Clock(dt={self.dt!r}, name={self.name!r})"

    def _get_dt_(self):
        return self.variables["dt"].get_value().item()

    @check_units(dt_=1)
    def _set_dt_(self, dt_):
        self._old_dt = self._get_dt_()
        self.variables["dt"].set_value(dt_)

    @check_units(dt=second)
    def _set_dt(self, dt):
        self._set_dt_(float(dt))

    dt = property(
        fget=lambda self: Quantity(self.dt_, dim=second.dim),
        fset=_set_dt,
        doc="""The time step of the simulation in seconds.""",
    )
    dt_ = property(
        fget=_get_dt_,
        fset=_set_dt_,
        doc="""The time step of the simulation as a float (in seconds)""",
    )

    @check_units(start=second, end=second)
    def set_interval(self, start, end):
        """
        set_interval(self, start, end)

        Set the start and end time of the simulation.

        Sets the start and end value of the clock precisely if
        possible (using epsilon) or rounding up if not. This assures that
        multiple calls to `Network.run` will not re-run the same time step.
        """
        self._set_t_update_dt(target_t=start)
        end = float(end)
        self._i_end = self._calc_timestep(end)
        if self._i_end > 2**40:
            logger.warn(
                "The end time of the simulation has been set to "
                f"{str(end*second)}, which based on the dt value of "
                f"{str(self.dt)} means that {self._i_end} "
                "time steps will be simulated. This can lead to "
                "numerical problems, e.g. the times t will not "
                "correspond to exact multiples of "
                "dt.",
                "many_timesteps",
            )

    #: The relative difference for times (in terms of dt) so that they are
    #: considered identical.
    epsilon_dt = 1e-4


class DefaultClockProxy:
    """
    Method proxy to access the defaultclock of the currently active device
    """

    def __getattr__(self, name):
        if name == "_is_proxy":
            return True
        from brian2.devices.device import active_device

        return getattr(active_device.defaultclock, name)

    def __setattr__(self, key, value):
        from brian2.devices.device import active_device

        setattr(active_device.defaultclock, key, value)


#: The standard clock, used for objects that do not specify any clock or dt
defaultclock = DefaultClockProxy()