Single cell simulation¶
-
class
myokit.
Simulation
(model, protocol=None, apd_var=None)¶ Runs single cell simulations using the CVODE solver (see [1]); CVODE uses an implicit multi-step method to achieve high accuracy and stability with adaptive step sizes.
The model passed to the simulation is cloned and stored internally, so changes to the original model object will not affect the simulation. A protocol can be passed in as
protocol
or set later usingset_protocol()
.Simulations maintain an internal state consisting of
- the current simulation time
- the current state
- the default state
When a simulation is created, the simulation time is set to 0 and both the current and the default state are copied from the model. After each call to
Simulation.run()
the time variable and current state are updated, so that each successive call to run continues where the previous simulation left off. Areset()
method is provided that will set the time back to 0 and revert the current state to the default state. To change the time or state manually, useset_time()
andset_state()
.A pre-pacing method
pre()
is provided that doesn’t affect the simulation time but will update the current and the default state. This allows you to pre-pace, run a simulation, reset to the pre-paced state, run another simulation etc.To get action potential duration (APD) measurements, the simulation can be run with threshold crossing detection. To enable this, the membrane potential variable must be specified when the simulation is created using the
apd_var
argument. This can be either a variable object or a string containing the variable’s fully qualified name. When running a simulation a threshold value can be passed in. In addition to the usual simulation log the run method will then return a list of all times at whichapd_var
crossed the threshold. Please note this is an APD calculated as the time between the crossing of a fixed threshold, it does not calculate dynamic thresholds like “90% of max(V) - min(V)”.The simulation provides four inputs a model variable can be bound to:
time
- This input provides the simulation time.
pace
- This input provides the current value of the pacing variable. This is determined using the protocol passed into the Simulation.
evaluations
- This input provides the number of rhs evaluations used at each point in time and can be used to gain some insight into the solver’s behaviour.
realtime
- This input provides the elapsed system time at each logged point.
No variable labels are required for this simulation type.
[1] SUNDIALS: Suite of nonlinear and differential/algebraic equation solvers. Hindmarsh, Brown, Woodward, et al. (2005) ACM Transactions on Mathematical Software.
-
default_state
()¶ Returns the default state.
-
eval_derivatives
(y=None)¶ Evaluates and returns the state derivatives.
The state to evaluate for can be given as
y
. If no state is given the current simulation state is used.
-
last_number_of_evaluations
()¶ Returns the number of rhs evaluations performed by the solver during the last simulation.
-
last_number_of_steps
()¶ Returns the number of steps taken by the solver during the last simulation.
-
last_state
()¶ If the last simulation resulted in an error, this will return the last state reached during that simulation. In all other cases, this method will return
None
.
-
pre
(duration, progress=None, msg='Pre-pacing Simulation')¶ This method can be used to perform an unlogged simulation, typically to pre-pace to a (semi-)stable orbit.
After running this method
- The simulation time is not affected
- The current state and the default state are updated to the final state reached in the simulation.
Calls to
reset()
after usingpre()
will set the current state to this new default state.To obtain feedback on the simulation progress, an object implementing the
myokit.ProgressReporter
interface can be passed in. passed in asprogress
. An optional description of the current simulation to use in the ProgressReporter can be passed in as msg.
-
reset
()¶ Resets the simulation:
- The time variable is set to 0
- The state is set to the default state
-
run
(duration, log=None, log_interval=None, log_times=None, apd_threshold=None, progress=None, msg='Running simulation')¶ Runs a simulation and returns the logged results. Running a simulation has the following effects:
- The internal state is updated to the last state in the simulation.
- The simulation’s time variable is updated to reflect the time
- elapsed during the simulation.
The number of time units to simulate can be set with
duration
.The method returns a
myokit.DataLog
dictionary that maps variable names to lists of logged values. The variables to log can be indicated using thelog
argument. There are several options for its value:None
(default), to log all states.- An integer flag or a combination of flags. Options:
myokit.LOG_NONE
,myokit.LOG_STATE
,myokit.LOG_BOUND
,myokit.LOG_INTER
,myokit.LOG_DERIV
ormyokit.LOG_ALL
. - A sequence of variable names. To log derivatives, use “dot(membrane.V)”.
- A
myokit.DataLog
object. In this case, the new data will be appended to the existing log.
For detailed information about the
log
argument, see the functionmyokit.prepare_log()
.By default, every step the solver takes is logged. This is usually advantageous, since more points are added exactly at the times the system gets more interesting. However, if equidistant points are required a
log_interval
can be set. Alternatively, thelog_times
argument can be used to specify logging times directly.To obtain accurate measurements of the action potential (AP) duration, the argument
apd_threshold
can be set to a fixed threshold level used to define the AP. This functionality is only available for simulations created with a validapd_var
argument. If apd measurements are enabled, the value returned by this method has the form(log, apds)
.To obtain feedback on the simulation progress, an object implementing the
myokit.ProgressReporter
interface can be passed in. passed in asprogress
. An optional description of the current simulation to use in the ProgressReporter can be passed in asmsg
.The
duration
argument cannot be negative, and special care needs to be taken when very small (positive) values are used. Ifduration
is zero or so small thatsimulation.time() + duration == simulation.time()
, then the method returns without updating the internal states or time. However, for some extremely short durations (approx2 epsilon * time
), the simulation will try to run but the underlying CVODE engine may return aCV_TOO_CLOSE
error, causing amyokit.SimulationError
to be raised.
-
set_constant
(var, value)¶ Changes a model constant. Only literal constants (constants not dependent on any other variable) can be changed.
The constant
var
can be given as aVariable
or a string containing a variable qname. Thevalue
should be given as a float.
-
set_default_state
(state)¶ Allows you to manually set the default state.
-
set_fixed_form_protocol
(times=None, values=None)¶ Configures this simulation to run with a predetermined protocol instead of the usual event-based mechanism.
A 1D time-series should be given as input. During the simulation, the value of the pacing variable will be determined by linearly interpolating between the two nearest points in the series. If the simulation time is outside the bounds of the time-series, the first or last value in the series will be used.
Setting a predetermined protocol clears any previously set (event-based or pre-determined) protocol. To clear all protocols, call this method with times=None. When a simulation is run without any protocol, the value of any variables bound to pace will be set to 0.
Arguments:
times
- A non-decreasing array of times. If any times appear more than once, only the value at the highest index will be used.
values
- An array of values for the pacing variable. Must have the same size
as
times
.
-
set_max_step_size
(dtmax=None)¶ Sets a maximum step size. To let the solver pick any step size it likes use
dtmax = None
.
-
set_min_step_size
(dtmin=None)¶ Sets a minimum step size. To let the solver pick any step size it likes use
dtmin = None
.
-
set_protocol
(protocol=None)¶ Sets the pacing
Protocol
used by this simulation.To run without pacing call this method with
protocol = None
. In this case, the value of any variables bound to pace will be set to 0.
-
set_state
(state)¶ Sets the current state.
-
set_time
(time=0)¶ Sets the current simulation time.
-
set_tolerance
(abs_tol=1e-06, rel_tol=0.0001)¶ Sets the solver tolerances. Absolute tolerance is set using
abs_tol
, relative tolerance usingrel_tol
. For more information on these values, see the Sundials CVODE documentation.
-
state
()¶ Returns the current state.
-
time
()¶ Returns the current simulation time.