####################################################################################
#
# STEPS - STochastic Engine for Pathway Simulation
# Copyright (C) 2007-2023 Okinawa Institute of Science and Technology, Japan.
# Copyright (C) 2003-2006 University of Antwerp, Belgium.
#
# See the file AUTHORS for details.
# This file is part of STEPS.
#
# STEPS is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License version 3,
# as published by the Free Software Foundation.
#
# STEPS is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program. If not, see <http://www.gnu.org/licenses/>.
#
#################################################################################
###
import collections
import copy
from enum import Enum
import inspect
import itertools
import numbers
import warnings
from steps import stepslib
from . import geom as ngeom
from . import utils as nutils
try:
COMPLEX_FILTER_MAX_VALUE = stepslib._py_COMPLEX_FILTER_MAX_VALUE
except AttributeError as ex:
# If the installed STEPS version does not support complexes
COMPLEX_FILTER_MAX_VALUE = 0
stepslib._py_Complex = None
stepslib._py_ComplexReac = None
__all__ = [
'Model',
'VolumeSystem',
'SurfaceSystem',
'VesicleSurfaceSystem',
'RaftSurfaceSystem',
'Vesicle',
'Raft',
'Species',
'LinkSpecies',
'Complex',
'Channel',
'Reaction',
'Diffusion',
'Endocytosis',
'Exocytosis',
'RaftEndocytosis',
'RaftGen',
'RaftDis',
'VesicleBind',
'VesicleUnbind',
'Current',
'OhmicCurr',
'GHKCurr',
'Location',
'In',
'Out',
'Surf',
'VesSurf',
'RaftSurf',
'ReactionManager',
'ReactionElement',
'ReactingElement',
'ReactionSide',
'SubUnit',
'SubUnitState',
'SubUnitSelector',
'ComplexState',
'ComplexSelector',
'NoOrdering',
'StrongOrdering',
'RotationalSymmetryOrdering',
'XDepFunc',
'VDepFunc',
'CompDepFunc',
'VDepRate',
'CompDepRate',
'CompDepCond',
'CompDepP',
'CompDepDcst',
'IMMOBILIZING',
'MOBILIZING',
'NO_EFFECT',
]
###################################################################################################
# Enums
try:
Immobilization = stepslib._py_Immobilization
except AttributeError:
# Older versions of cython do not export this enum
class Immobilization(Enum):
IMMOBILIZING = 0
MOBILIZING = 1
NO_EFFECT = 2
# Aliases
IMMOBILIZING = Immobilization.IMMOBILIZING
MOBILIZING = Immobilization.MOBILIZING
NO_EFFECT = Immobilization.NO_EFFECT
###################################################################################################
# Model
[docs]@nutils.FreezeAfterInit
class Model(nutils.UsableObject, nutils.StepsWrapperObject):
"""Top-level container for the objects in a kinetic model
Should be used as a context manager for the declaration of Volume systems, reactions,
diffusion rules, etc::
mdl = Model()
with mdl:
vsys = VolumeSystem.Create()
... # Declare other objects in mdl
After having declared children objects in the ``with mdl:`` block, they can be accessed
as attributes of ``mdl`` with their name (see :py:class:`steps.API_2.utils.NamedObject` and
:py:func:`steps.API_2.utils.NamedObject.Create`)::
mdl.vsys
"""
def __init__(self, *args, _createObj=True, **kwargs):
super().__init__(*args, **kwargs)
self.stepsModel = self._createStepsObj() if _createObj else None
self.volSysConstraints = []
def _SetUpMdlDeps(self, geom):
"""Set up structures that depend on objects declared in the model."""
# Start with vesicles because they can update rafts
for pl in self._getChildrenOfType(Vesicle):
pl._SetUpMdlDeps(self, geom)
for pl in self._getChildrenOfType(Raft):
pl._SetUpMdlDeps(self, geom)
for pl in self._getChildrenOfType(VesicleSurfaceSystem):
pl._SetUpMdlDeps(geom)
for pl in self._getChildrenOfType(RaftSurfaceSystem):
pl._SetUpMdlDeps(geom)
def _getStepsObjects(self):
"""Return a list of the steps objects that this named object holds."""
return [self.stepsModel]
@classmethod
def _FromStepsObject(cls, obj):
"""Create the interface object from a STEPS object."""
mdl = cls(_createObj=False)
mdl.stepsModel = obj
with mdl:
for spec in obj.getAllSpecs():
Species._FromStepsObject(spec, mdl)
for chan in obj.getAllChans():
Channel._FromStepsObject(chan, mdl)
for lspec in obj.getAllLinkSpecs():
LinkSpecies._FromStepsObject(lspec, mdl)
for ves in obj.getAllVesicles():
Vesicle._FromStepsObject(ves, mdl)
for raft in obj.getAllRafts():
Raft._FromStepsObject(raft, mdl)
for vssys in obj.getAllVesSurfsyss():
VesicleSurfaceSystem._FromStepsObject(vssys, mdl)
for rssys in obj.getAllRaftsyss():
RaftSurfaceSystem._FromStepsObject(rssys, mdl)
for vsys in obj.getAllVolsyss():
VolumeSystem._FromStepsObject(vsys, mdl)
for ssys in obj.getAllSurfsyss():
SurfaceSystem._FromStepsObject(ssys, mdl)
return mdl
def _addVolSysConstraint(self, reac, volSys, surfSys, loc):
self.volSysConstraints.append((reac, volSys, surfSys, loc))
def _createStepsObj(self):
return stepslib._py_Model()
###################################################################################################
# Space systems
@nutils.FreezeAfterInit
class SpaceSystem(nutils.UsableObject, nutils.UsingObjects(Model), nutils.StepsWrapperObject):
"""Base class for volume and surface system classes
:meta private:
"""
def __init__(self, *args, _createObj=True, **kwargs):
super().__init__(*args, **kwargs)
(mdl,) = self._getUsedObjects()
self.stepsSys = self._createStepsObj(mdl) if _createObj else None
self.locations = []
def _addLocation(self, l):
self.locations.append(l)
def _getStepsObjects(self):
"""Return a list of the steps objects that this named object holds."""
return [self.stepsSys]
[docs]@nutils.FreezeAfterInit
class VolumeSystem(SpaceSystem):
"""A container that groups reactions and diffusion rules located in a volume
Should be used as a context manager for the declaration of reactions, diffusion rules, etc::
...
with mdl:
...
vsys = VolumeSystem.Create()
with vsys:
diff1 = Diffusion.Create(...)
... # Declare other objects in vsys
After having declared children objects in the ``with vsys:`` block, they can be accessed
as attributes of ``vsys`` with their name (see :py:class:`steps.API_2.utils.NamedObject` and
:py:func:`steps.API_2.utils.NamedObject.Create`)::
vsys.diff1
"""
def __init__(self, *args, **kwargs):
super().__init__(*args, **kwargs)
def _createStepsObj(self, mdl):
return stepslib._py_Volsys(self.name, mdl.stepsModel)
@classmethod
def _FromStepsObject(cls, obj, mdl):
"""Create the interface object from a STEPS object."""
vsys = cls(_createObj=False, name=obj.getID())
vsys.stepsSys = obj
with vsys:
for diff in obj.getAllDiffs():
Diffusion._FromStepsObject(diff, mdl)
for reac in obj.getAllReacs():
# need to convert to actual steps classes
Reaction._FromStepsObject(reac, mdl)
for bind in obj.getAllVesBinds():
VesicleBind._FromStepsObject(bind, mdl)
for unbind in obj.getAllVesUnbinds():
VesicleUnbind._FromStepsObject(unbind, mdl)
return vsys
[docs]@nutils.FreezeAfterInit
class SurfaceSystem(SpaceSystem):
"""A container that groups reactions and diffusion rules located in a surface
Should be used as a context manager for the declaration of surface reactions, surface
diffusion rules, etc::
...
with mdl:
...
ssys = SurfaceSystem.Create()
with ssys:
diff2 = Diffusion.Create(...)
... # Declare other objects in ssys
After having declared children objects in the ``with ssys:`` block, they can be accessed
as attributes of ``ssys`` with their name (see :py:class:`steps.API_2.utils.NamedObject` and
:py:func:`steps.API_2.utils.NamedObject.Create`)::
ssys.diff2
"""
def __init__(self, *args, **kwargs):
super().__init__(*args, **kwargs)
def _createStepsObj(self, mdl):
return stepslib._py_Surfsys(self.name, mdl.stepsModel)
def _addLocation(self, l):
super()._addLocation(l)
if isinstance(l, ngeom.Patch):
if l.innerComp is not None:
l.innerComp.addSystem(self, Location.IN)
if l.outerComp is not None:
l.outerComp.addSystem(self, Location.OUT)
@classmethod
def _FromStepsObject(cls, obj, mdl):
"""Create the interface object from a STEPS object."""
ssys = cls(_createObj=False, name=obj.getID())
ssys.stepsSys = obj
with ssys:
for diff in obj.getAllDiffs():
Diffusion._FromStepsObject(diff, mdl)
for sreac in obj.getAllSReacs():
Reaction._FromStepsObject(sreac, mdl)
for vdsreac in obj.getAllVDepSReacs():
Reaction._FromStepsObject(vdsreac, mdl)
# TODO not urgent: implement _FromStepsObject for GHK and Ohmic currents
# for ghkc in obj.getAllGHKcurrs():
# GHKCurr._FromStepsObject(ghkc, mdl)
# for ohmc in obj.getAllOhmicCurrs():
# OhmicCurr._FromStepsObject(ohmc, mdl)
return ssys
[docs]@nutils.FreezeAfterInit
class VesicleSurfaceSystem(SurfaceSystem):
"""A container that groups reactions and diffusion rules located in a vesicle surface
Should be used as a context manager for the declaration of vesicle surface reactions, vesicle surface
diffusion rules, etc::
...
with mdl:
...
vssys = VesicleSurfaceSystem.Create()
with vssys:
diff3 = Diffusion.Create(...)
... # Declare other objects in vssys
After having declared children objects in the ``with vssys:`` block, they can be accessed
as attributes of ``vssys`` with their name (see :py:class:`steps.API_2.utils.NamedObject` and
:py:func:`steps.API_2.utils.NamedObject.Create`)::
vssys.diff3
"""
def __init__(self, *args, **kwargs):
super().__init__(*args, **kwargs)
def _createStepsObj(self, mdl):
return stepslib._py_VesSurfsys(self.name, mdl.stepsModel)
def _addLocation(self, l):
super()._addLocation(l)
def _SetUpMdlDeps(self, geom):
"""Set up structures that depend on objects declared in the model."""
# Potentially all compartments can contain any vesicle
for comp in geom.ALL(ngeom.Compartment):
comp.addSystem(self, Location.OUT)
@classmethod
def _FromStepsObject(cls, obj, mdl):
"""Create the interface object from a STEPS object."""
vssys = cls(_createObj=False, name=obj.getID())
vssys.stepsSys = obj
with vssys:
for diff in obj.getAllVesSDiffs():
Diffusion._FromStepsObject(diff, mdl)
for vsreac in obj.getAllVesSReacs():
Reaction._FromStepsObject(vsreac, mdl)
for exc in obj.getAllExocytosis():
Exocytosis._FromStepsObject(exc, mdl)
return vssys
[docs]@nutils.FreezeAfterInit
class RaftSurfaceSystem(SurfaceSystem):
"""A container that groups reactions and diffusion rules located in a raft surface
Should be used as a context manager for the declaration of raft surface reactions, raft endocytosis
rules, etc::
...
with mdl:
...
rssys = RaftSurfaceSystem.Create()
with rssys:
rend1 = RaftEndocytosis.Create(...)
... # Declare other objects in rssys
After having declared children objects in the ``with rssys:`` block, they can be accessed
as attributes of ``rssys`` with their name (see :py:class:`steps.API_2.utils.NamedObject` and
:py:func:`steps.API_2.utils.NamedObject.Create`)::
rssys.rend1
"""
def __init__(self, *args, **kwargs):
super().__init__(*args, **kwargs)
def _createStepsObj(self, mdl):
return stepslib._py_Raftsys(self.name, mdl.stepsModel)
def _addLocation(self, l):
super()._addLocation(l)
def _SetUpMdlDeps(self, geom):
"""Set up structures that depend on objects declared in the model."""
# Potentially all patches can contain any raft
for patch in geom.ALL(ngeom.Patch):
patch.addSystem(self)
if patch.innerComp is not None:
patch.innerComp.addSystem(self, Location.IN)
if patch.outerComp is not None:
patch.outerComp.addSystem(self, Location.OUT)
@classmethod
def _FromStepsObject(cls, obj, mdl):
"""Create the interface object from a STEPS object."""
rssys = cls(_createObj=False, name=obj.getID())
rssys.stepsSys = obj
with rssys:
for rsreac in obj.getAllRaftSReacs():
Reaction._FromStepsObject(rsreac, mdl)
for end in obj.getAllRaftEndocytosiss():
RaftEndocytosis._FromStepsObject(end, mdl)
for dis in obj.getAllRaftDiss():
RaftDis._FromStepsObject(dis, mdl)
return rssys
###################################################################################################
# Complexes and channels
####################################
# State Ordering Functions
[docs]def StrongOrdering(state):
"""Strong ordering function for complex states
:param state: A complex state, i.e. a sequence of :py:class:`SubUnitState`
:type state: Tuple[:py:class:`SubUnitState`]
:returns: The equivalent unique ordering
:rtype: Tuple[:py:class:`SubUnitState`]
Take a tuple of :py:class:`SubUnitState` as input and return the unique equivalent based on
the type of ordering chosen.
StrongOrdering corresponds to no symmetry in the complex.
"""
return tuple(state)
[docs]def NoOrdering(state):
"""No ordering function for complex states
:param state: A complex state, i.e. a sequence of :py:class:`SubUnitState`
:type state: Tuple[:py:class:`SubUnitState`]
:returns: The equivalent unique ordering
:rtype: Tuple[:py:class:`SubUnitState`]
Take a tuple of :py:class:`SubUnitState` as input and return the unique equivalent based on
the type of ordering chosen.
NoOrdering considers that the subunits have no physical position, a state is just a set of
SubUnitStates.
"""
return tuple(sorted(state, key=lambda x: x.name))
[docs]def RotationalSymmetryOrdering(state):
"""Rotational symmetry ordering function for complex states
:param state: A complex state, i.e. a sequence of :py:class:`SubUnitState`
:type state: Tuple[:py:class:`SubUnitState`]
:returns: The equivalent unique ordering
:rtype: Tuple[:py:class:`SubUnitState`]
Take a tuple of :py:class:`SubUnitState` as input and return the unique equivalent based on
the type of ordering chosen.
RotationalSymmetryOrdering represents central symmetry: state [a,b,c] is equivalent to
[b,c,a] and [c,a,b].
"""
# Number the states based on sorting their names
srt = sorted(enumerate(state), key=lambda x: x[1].name)
num = 0
s2n = {state[srt[0][0]]: num}
for i, s in enumerate(srt[1:]):
if s[1] != srt[i][1]:
num += 1
s2n[state[s[0]]] = num
# Compute all rotations
n = len(state)
rots = [tuple(state[(i + shift) % n] for i in range(n)) for shift in range(n)]
# select the one that maximizes the number formed by the sequence of digits
return max(rots, key=lambda state: sum(s2n[s] * (num + 1) ** p for p, s in enumerate(state)))
####################################
# Complexes
[docs]@nutils.FreezeAfterInit
class Complex(nutils.UsingObjects(Model)):
"""Multi-state complex class
A complex is composed of one or more subunits, each subunit has a finite (preferably low)
number of states. Subunit states can be used inside reactions so as to simplify the declaration
of reactions for multi-subunits complexes.
This class can also be used as a way to automatically declare complex states as species.
:param subUnits: The list of :py:class:`SubUnit` of the complex. Can also be a list of
:py:class:`SubUnitState` in which case a :py:class:`SubUnit` will automatically be created
and the complex will be composed of only this :py:class:`SubUnit`.
:type subUnits: List[:py:class:`SubUnit`] or List[:py:class:`SubUnitState`]
:param statesAsSpecies: If True, all the possible states of the complex are declared as
:py:class:`Species`. If not, and if your STEPS implementation supports it, the complex will be
declared as a native STEPS multi-state complex.
:type statesAsSpecies: bool
:param order: A callable that defines a custom state ordering. It should map a tuple of
:py:class:`SubUnitState` into a unique ordering. It defaults to and :py:func:`NoOrdering`,
which considers that subunits have no spatial position. See :py:func:`StrongOrdering`,
:py:func:`NoOrdering`, and :py:func:`RotationalSymmetryOrdering` for details.
:type order: Callable[[Tuple[:py:class:`SubUnitState`]], Tuple[:py:class:`SubUnitState`]]
Declaration::
S1A, S1B, S1C, S2A, S2B = SubUnitState.Create()
SU1, SU2 = SubUnit.Create(
[S1A, S1B, S1C],
[S2A, S2B]
)
Comp = Complex.Create([SU1, SU1, SU2])
Usage in reactions:
- Centered on complexes::
# Get a single 'Comp' complex, when several complexes of the same type are used in
# the same reaction, it allows to determine which complex undergoes which changes.
C = Comp.get()
# Ca binding to the first subunit of complex Comp in state S1A and changing its
# state to S1B
C[S1A, :, :] + Ca <r[1]> C[S1B, :, :]
- Centered on subunit states::
with Comp[...]:
# Ca Binding to any subunit in state S1A and changing it to state S1B
S1A + Ca <r[1]> S1B
# Ca Binding to first subunit in state S1A and changing it to state S1B
S1A[0] + Ca <r[1]> S1B[0]
with Comp[...] as C1, Comp[..., S2B] as C2:
S1A[C1] + S1B[C2] <r[1]> S1B[C1] + S1C[C2]
"""
_elemStr = 'Complex'
_SIMPATH_ONLY_CHILDREN = True
def __init__(self, subUnits, *args, statesAsSpecies=False, order=NoOrdering, _createObj=True, **kwargs):
"""
Declare a complex from a list of subunits.
If statesAsSpecies is True, automatically declare each complex state as a steps Species.
'order' allows the user to define a custom ordering, it should be callable and take
a state (tuple of SubUnitStates) as argument and return a unique version of this state
according to the custom ordering. It defaults to NoOrdering, which considers that subunits
have no spatial position.
"""
# Creating from a list of states instead of a list of subunits
if all(isinstance(su, SubUnitState) for su in subUnits):
subUnits = [SubUnit(subUnits)]
self._subUnits = subUnits
super().__init__(*args, **kwargs)
for i, sub in enumerate(self._subUnits):
sub._addComplex(self)
# Add children
self._addChildren(sub)
for sus in sub._states:
self._addChildren(sus)
if order is not NoOrdering and not statesAsSpecies:
raise Exception(
f'Steps complexes need to be declared with "order=NoOrdering", no '
f'other orders are accepted.'
)
if not statesAsSpecies and stepslib._py_Complex is None:
warnings.warn(
f'Your version of STEPS does not support multi-state complexes, '
f'Complex {self.name} will instead be declared with '
f'statesAsSpecies = True.'
)
statesAsSpecies = True
self._order = order
(mdl,) = self._getUsedObjects()
self._statesAsSpecies = statesAsSpecies
self.stepsComplex = None
self._compStates = {}
if _createObj:
if statesAsSpecies:
# Create all complex states as steps species
self.stepsComplex = self._createStepsStates(mdl)
else:
self._initializePoolStructs()
self.stepsComplex = self._createStepsObj(mdl)
def _areStatesAsSpecies(self):
return self._statesAsSpecies
def _getAdditionalChildren(self):
"""
Return a list of NamedObjects that should be added as children of the used object
along with the current object.
"""
res = set()
for su in self._subUnits:
res |= su._states
return res
def _initializePoolStructs(self):
"""Initialize structures related to real complexes."""
allSus = set()
for su in set(self._subUnits):
allSus |= su._states
self._poolElems = sorted(allSus, key=lambda x: x.name)
self._sus2PoolInd = {sus: self._poolElems.index(sus) for sus in allSus}
def _getSUSInd(self, sus):
"""Return the index os the sus SubUnitState."""
if sus not in self._sus2PoolInd:
raise Exception(f'SubUnitState {sus} is not associated with Complex {self}.')
return self._sus2PoolInd[sus]
def _createStepsObj(self, mdl):
# Not implemented in steps yet
return stepslib._py_Complex(self.name, mdl.stepsModel, len(self._subUnits), len(self._poolElems))
def _getStepsObjects(self):
"""Return a list of the steps objects that this named object holds."""
if self._statesAsSpecies:
return [obj for state, obj in self._compStates.items()]
else:
return [self.stepsComplex]
def _getAllElems(self, loc=None):
"""
Return all elements that can directly be accessed by name from the parent physical
location.
"""
return [state for state, obj in self._compStates.items()]
def _createStepsStates(self, mdl):
"""Create steps species for each complex state."""
for comb in itertools.product(*[sorted(su._states, key=lambda x: x.name) for su in self._subUnits]):
cs = ComplexState(self, comb)
if cs not in self._compStates:
spec = stepslib._py_Spec(cs.name, mdl.stepsModel)
cs._setStepsObjects(spec)
self._compStates[cs] = spec
return None
[docs] def __iter__(self):
"""Iterator over all the possible states of the complex.
:returns: An iterator that iterates over all unique (as defined by the custom ordering)
:py:class:`ComplexState` of the complex.
:rtype: Iterator[:py:class:`ComplexState`]
Usage::
Comp1 = Complex.Create(...)
for state in Comp1:
# Do something with the state
...
:meta public:
"""
return iter(self.__getitem__(...))
def _simPathWalkExpand(self):
"""Return an iterable of the elements that should be part of ResultSelector paths."""
if self._statesAsSpecies:
return self.__iter__()
else:
return [self.__getitem__(...)]
def _simPathCombinerClass(self):
"""Return the class that needs to be used to combine expanded elements."""
if self._statesAsSpecies:
return nutils.SumSimPathComb
else:
return None
[docs] def __getitem__(self, subStates, _selName=None):
"""Square-bracket operator for defining states or sets of states.
:param subStates: A sequence of one or more :py:class:`SubUnitState` or
:py:class:`SubUnitSelector` that corresponds to one or several :py:class:`ComplexState`.
:type subStates: Tuple[Union[:py:class:`SubUnitState`, :py:class:`SubUnitSelector`], ...]
:returns: A single complex state or a set of complex states matching the arguments.
:rtype: Union[:py:class:`ComplexState`, :py:class:`ComplexSelector`]
Return the :py:class:`ComplexState` corresponding to the combination of
:py:class:`SubUnitState` given as argument. If the state is fully determined (the number and
type of SubUnitStates must match the :py:class:`SubUnit` that have been declared when the
complex was created), return the corresponding :py:class:`ComplexState`. Otherwise, return
a :py:class:`ComplexSelector` representing all the possible states.
The slice operator ('``:``') as well as ellipsis ('``...``') can be used in the same way as
when slicing numpy arrays. See https://numpy.org/doc/1.18/reference/arrays.indexing.html
for details.
Valid examples::
Comp[S1A, S1B, S2A] # Valid, returns the corresponding ComplexState
Comp[S1A, :, :] # Valid, returns a ComplexSelector in which the first SubUnit is
# in state S1A and the other subunits are in any state
Comp[S1A, ...] # Valid, the ellipsis operator ('...') is equivalent to filling
# the remaining subUnits with ':'
Comp[..., S2A] # Valid, returns a ComplexSelector in which the last SubUnit is in
# state S2A and the other subunits are in any state
Comp[S1A|S1B, :, :] # Valid, first subunit is either in state S1A or S1B and the other
# subunits are in any state
Comp[~S1C, :, :] # Valid, first subunit is not in state S1C and the other subunits
# are in any state
See :py:class:`SubUnitState` and :py:class:`SubUnitSelector` for the syntax relative to
the last two examples.
Invalid examples::
Comp[S2A, S1B, S1A] # Invalid, the order of states does not match the one used for
# complex creation
Comp[S1A] # Invalid, only one subunit is specified while the complex has 3
Comp[..., S1A, ...] # Invalid, cannot use more than one ellipsis operator
:meta public:
"""
if isinstance(subStates, ComplexState):
subStates = tuple(subStates)
if not isinstance(subStates, tuple):
subStates = (subStates,)
# If ellipsis is used, expand it first
if len(subStates) > len(self._subUnits):
raise KeyError(
'The number of specified states is higher than the number of subunits ' 'in the complex.'
)
elif len(subStates) <= len(self._subUnits):
if any(s is Ellipsis for s in subStates):
if subStates.count(Ellipsis) > 1:
raise KeyError('Cannot use the ellipsis operator ("...") more than once.')
# If the ellipsis operator was used correctly, fill the missing slots with ':'
i = subStates.index(Ellipsis)
m = len(self._subUnits) - len(subStates) + 1
subStates = subStates[:i] + (slice(None),) * m + subStates[i + 1 :]
elif len(subStates) < len(self._subUnits):
raise KeyError(
'The number of states does not match the number of subunits. '
'A complex selector needs to either specify all subunits states '
'or use the Ellipsis operator ("...").'
)
# If the state is fully declared
if all(isinstance(s, SubUnitState) for s in subStates):
cs = ComplexState(self, subStates, _parentSelName=_selName)
cs._stepsObj = self._compStates[cs] if self._statesAsSpecies else None
return cs
else:
suSelectors = []
for i, su in enumerate(subStates):
if isinstance(su, SubUnitState):
if su._subUnit is not self._subUnits[i]:
raise KeyError(
f'SubUnitState {su} is not associated to SubUnit ' f'{self._subUnits[i]}.'
)
suSelectors.append(SubUnitSelector(self._subUnits[i], [su], _pos=i))
elif isinstance(su, SubUnitSelector):
if su._isEmpty():
raise KeyError(
f'SubUnitSelector {su} is empty. Cannot build a complex '
f'selector from an empty SubUnitSelector.'
)
su2 = copy.copy(su)
su2._pos = i
suSelectors.append(su2)
elif su == slice(None):
suSelectors.append(SubUnitSelector(self._subUnits[i], self._subUnits[i]._states, _pos=i))
else:
raise KeyError(f'Unexpected entry in complex selector expression: {su}.')
return ComplexSelector(self, suSelectors, name=_selName)
[docs] def get(self):
"""Get a specific instance of the complex (i.e. a single molecule)
:returns: An object that represents a specific instance of a complex, to be used in
reactions. This object behaves like a :py:class:`Complex`.
:rtype: :py:class:`_ComplexSelectorSpawner`
"""
return _ComplexSelectorSpawner(ComplexSelector(self, []))
def __eq__(self, other):
return isinstance(other, self.__class__) and self.name == other.name
def __hash__(self):
return hash(self.name)
def _solverStr(self):
"""Return the string that is used as part of method names for this specific object."""
return Complex._elemStr
[docs]@nutils.FreezeAfterInit
class SubUnit(nutils.NamedObject):
"""Class for representing complex subunits
:param subStates: List of states in which the subunit can be.
:type subStates: List[:py:class:`SubUnitState`]
Does not necessarily represent a biological subunit, it can be used to group substates of
a complex. For example, if a complex can bind to species, it could have a binding
:py:class:`SubUnit` with :py:class:`SubUnitState` representing the species being bound.
A :py:class:`SubUnit` can only be asscoiated to a single :py:class:`Complex`.
"""
_elemStr = 'SU'
def __init__(self, subStates, *args, **kwargs):
"""Initialize the SubUnit from a list of SubUnitState objects."""
super().__init__(*args, **kwargs)
self._orderedStates = subStates
self._states = set(subStates)
for state in self._states:
state._setSubUnit(self)
self._complex = None
[docs] def __iter__(self):
"""Iterator over all the possible states of the subunit.
:returns: An iterator that iterates over all :py:class:`SubUnitState` in the subunit
state. The order is the same as the one defined during the creation of the subunit.
:rtype: Iterator[:py:class:`SubUnitState`]
Usage::
>>> su = SubUnit.Create([h0, h1])
>>> for subState in su:
... print(subState.name)
...
h0
h1
:meta public:
"""
return iter(self._orderedStates)
def _addComplex(self, comp):
if self._complex is not None and self._complex != comp:
raise Exception(f'A SubUnit can only be associated to one complex.')
self._complex = comp
def _areStatesAsSpecies(self):
return self._complex._areStatesAsSpecies()
def __hash__(self):
return id(self)
def _solverStr(self):
"""Return the string that is used as part of method names for this specific object."""
return SubUnit._elemStr
class _ComplexSelectorSpawner:
"""
Utility class for spawning complex selectors with identical names.
Having identical names means that they will represent the same instance of a complex in
reactions.
Usage:
C = Comp.get() # Get the _ComplexSelectorSpawner
C[S1A, ...] <r[1]> C[S2A, ...] # Signifies that the complexSelectors on both sides are
# related to the same complex.
"""
def __init__(self, cs):
self._cs = cs
def __getitem__(self, key):
return self._cs._complex.__getitem__(key, _selName=self._cs.name)
def __lshift__(self, other):
return self._cs._complex.__getitem__(..., _selName=self._cs.name) << other
####################################
# Channels
[docs]@nutils.FreezeAfterInit
class Channel(Complex):
"""Ion channel class
:parameters: See :py:class:`Complex` for a description of the parameters.
Channels behave like complexes. This simplifies the declaration of reactions between channel
states in cases like the Hodgkin-Huxley model of Na+ channel, in which gating variable m and h
are independent and can thus be declared as separate SubUnits. All combinations are then
automatically computed for creating the channel states.
"""
def __init__(self, subUnits, *args, statesAsSpecies=True, **kwargs):
if not statesAsSpecies:
raise NotImplementedError()
super().__init__(subUnits, *args, statesAsSpecies=True, **kwargs)
@classmethod
def _FromStepsObject(cls, obj, mdl):
"""Create the interface object from a STEPS object."""
chanStates = obj.getAllChanStates()
allSus = [SubUnitState() for state in chanStates]
channel = cls(allSus, statesAsSpecies=True, _createObj=False, name=obj.getID())
channel.stepsComplex = obj
for state, sus in zip(chanStates, allSus):
cs = ComplexState(channel, (sus,))
cs._setStepsObjects(state)
channel._compStates[cs] = state
return channel
def _createStepsStates(self, mdl):
"""Create steps ChanState objects for each channel state."""
chan = stepslib._py_Chan(self.name, mdl.stepsModel)
# Compute the channel state from the subunit states
for comb in itertools.product(*[sorted(su._states, key=lambda x: x.name) for su in self._subUnits]):
cs = ComplexState(self, comb)
if cs not in self._compStates:
chanState = stepslib._py_ChanState(cs.name, mdl.stepsModel, chan)
cs._setStepsObjects(chanState)
self._compStates[cs] = chanState
return chan
###################################################################################################
# Reaction elements
[docs]@nutils.FreezeAfterInit
class ReactionElement(nutils.NamedObject):
"""Base class for objects that can be included in a reaction
Can be combined to other :py:class:`ReactionElement` with the ``+`` operator to form a
:py:class:`ReactionSide`.
"""
def __init__(self, *args, **kwargs):
super().__init__(*args, **kwargs)
def _toReactionSide(self):
return ReactionSide([ReactingElement(self, 1)])
[docs] def __add__(self, other):
"""Add two elements and form a reaction side
:param other: The other element
:type other: :py:class:`ReactionElement`
:returns: A reaction side containing itself and the other element it was added to.
:rtype: :py:class:`ReactionSide`
Usage::
...
S1 + S2 <r[1]> ...
S1 + 2*S2 <r[2]> ...
S1 + (S2 + S3) <r[3]> ...
...
:meta public:
"""
if not isinstance(other, ReactionElement):
raise TypeError(f'Cannot add {other} to a reaction.')
reacSelf = self if isinstance(self, ReactingElement) else ReactingElement(self, 1)
reacOther = other if isinstance(other, ReactingElement) else ReactingElement(other, 1)
return ReactionSide([reacSelf, reacOther])
[docs] def __rmul__(self, stoich):
"""Specify a stoichiometric coefficient
:param stoich: The stoichimetric coefficient
:type stoich: int
:returns: A reacting element object that uses stoich as a stoichiometric coefficient.
:rtype: :py:class:`ReactingElement`
Usage::
...
2*S1 <r[1]> ...
# Equivalent to:
S1 + S1 <r[1]> ...
...
:meta public:
"""
self._checkStoech(stoich)
return ReactingElement(self, stoich)
def _toReactionElem(self):
return self
def _checkStoech(self, stoich):
if not isinstance(stoich, numbers.Integral) or stoich <= 0:
raise TypeError(
f'{stoich} is not an acceptable stoichiometric coefficient, it needs '
f'to be a strictly positive integer.'
)
def _makeCopy(self):
"""
Make a copy of the reaction element.
For most elements, like species, self should be returned because there should not be more
than one instance of this element.
"""
return self
# Location specifiers
@property
def i(self):
"""Get a version of the element located in the inside compartment
Shorthand for :py:func:`In`.
:type: :py:class:`ReactingElement`, read-only
"""
return In(self)
@property
def o(self):
"""Get a version of the element located in the outside compartment
Shorthand for :py:func:`Out`.
:type: :py:class:`ReactingElement`, read-only
"""
return Out(self)
@property
def s(self):
"""Get a version of the element located on the surface patch
Shorthand for :py:func:`Surf`.
:type: :py:class:`ReactingElement`, read-only
"""
return Surf(self)
def _simPathAutoMetaData(self):
"""Return a dictionary with string keys and string or numbers values."""
return {'obj_type': self._solverStr(), 'obj_id': self.name}
@property
def v(self):
"""Get a version of the element located on the surface of a vesicle
Shorthand for :py:func:`VesSurf`.
:type: :py:class:`ReactingElement`, read-only
"""
return VesSurf(self)
@property
def r(self):
"""Get a version of the element located on the surface of a raft
Shorthand for :py:func:`RaftSurf`.
:type: :py:class:`ReactingElement`, read-only
"""
return RaftSurf(self)
[docs]@nutils.FreezeAfterInit
class ReactingElement(ReactionElement):
"""A reaction element with a stoichiometric coefficient and a location
All elements in a :py:class:`ReactionSide` are :py:class:`ReactingElement`. While
:py:class:`ReactionElement` is the base class for all objects that could be included in a
reaction, :py:class:`ReactingElement` is a wrapper that is actually included in reactions.
It contains additional information about stoichiometry and location of the element.
.. note::
This class should not be instantiated by the user, it is only documented for clarity.
"""
def __init__(self, elem, stoich, loc=None, *args, **kwargs):
super().__init__(*args, **kwargs)
self._elem = elem._toReactionElem()
self.name = self._elem.name
self.stoich = stoich
self.loc = loc
def _toReactionSide(self):
return ReactionSide([self])
def _SetLoc(self, loc):
self.loc = loc
def __rmul__(self, stoich):
"""Return a ReactingElement object that uses stoich as a stoichiometric coefficient."""
self._checkStoech(stoich)
return ReactingElement(self._elem, stoich * self.stoich, self.loc)
def __repr__(self):
if self.loc is None:
return '{}{}'.format(self.stoich if self.stoich > 1 else '', str(self._elem))
else:
s = f'{self.stoich} ' if self.stoich > 1 else ''
pos = {
Location.IN: 'i', Location.OUT: 'o', Location.SURF: 's', Location.VESSURF: 'v',
Location.RAFTSURF: 'r'
}[self.loc]
return f'{s}{self._elem}.{pos}'
def _getStepsObjects(self):
"""Return a list of the steps objects that this named object holds."""
return self._elem._getStepsObjects()
def __eq__(self, other):
return isinstance(other, ReactingElement) and (self._elem, self.stoich, self.loc) == (
other._elem,
other.stoich,
other.loc,
)
def __hash__(self):
return hash((self._elem, self.stoich, self.loc))
####################################
# Location specifiers
[docs]class Location(Enum):
"""The different locations that reaction elements can be in"""
IN = 1
SURF = 2
OUT = 3
VESSURF = 4
RAFTSURF = 5
ALL_LOCATIONS = [None, Location.IN, Location.SURF, Location.OUT, Location.VESSURF, Location.RAFTSURF]
ALL_SURF_LOCATIONS = [Location.SURF, Location.VESSURF, Location.RAFTSURF]
def _locSpecifier(elem, loc):
"""Return a localized ReactingElement with location 'loc'."""
if isinstance(elem, ReactionElement) and not isinstance(elem, ReactingElement):
return ReactingElement(elem, 1, loc)
elif isinstance(elem, ReactionSide):
elem._SetLoc(loc)
return elem
else:
raise TypeError(f'Expected a ReactionElement or a ReactionSide, got {elem} instead.')
[docs]def In(elem):
"""Get a version of an element located in the inside compartment
:param elem: The element
:type elem: Union[:py:class:`ReactionElement`, :py:class:`ReactionSide`]
:returns: A version of the element located in the inside compartment
:rtype: Union[:py:class:`ReactingElement`, :py:class:`ReactionSide`]
"""
return _locSpecifier(elem, Location.IN)
[docs]def Out(elem):
"""Get a version of an element located in the outside compartment
:param elem: The element
:type elem: Union[:py:class:`ReactionElement`, :py:class:`ReactionSide`]
:returns: A version of the element located in the outside compartment
:rtype: Union[:py:class:`ReactingElement`, :py:class:`ReactionSide`]
"""
return _locSpecifier(elem, Location.OUT)
[docs]def Surf(elem):
"""Get a version of an element located on the surface patch
:param elem: The element
:type elem: Union[:py:class:`ReactionElement`, :py:class:`ReactionSide`]
:returns: A version of the element located on the surface patch
:rtype: Union[:py:class:`ReactingElement`, :py:class:`ReactionSide`]
"""
return _locSpecifier(elem, Location.SURF)
[docs]def VesSurf(elem):
"""Get a version of an element located on the surface of a vesicle
:param elem: The element
:type elem: Union[:py:class:`ReactionElement`, :py:class:`ReactionSide`]
:returns: A version of the element located on the surface of a vesicle
:rtype: Union[:py:class:`ReactingElement`, :py:class:`ReactionSide`]
"""
return _locSpecifier(elem, Location.VESSURF)
[docs]def RaftSurf(elem):
"""Get a version of an element located on the surface of a raft
:param elem: The element
:type elem: Union[:py:class:`ReactionElement`, :py:class:`ReactionSide`]
:returns: A version of the element located on the surface of a raft
:rtype: Union[:py:class:`ReactingElement`, :py:class:`ReactionSide`]
"""
return _locSpecifier(elem, Location.RAFTSURF)
####################################
# Vesicles
class _RaftVesLocation(nutils.UsingObjects(Model),
nutils.StepsWrapperObject,
nutils.ParameterizedObject
):
def __init__(self, *args, **kwargs):
super().__init__(*args, **kwargs)
self.systems = []
self.sysNames = []
def _SetUpMdlDeps(self, mdl, geom):
"""Set up the structures that will allow species or reaction access from locations."""
for sname, loc in self.sysNames:
s = getattr(mdl, sname)
self.systems.append((s, loc))
# Add the reactions and the corresponding reactants to children
for name, c in s.children.items():
if _RaftVesLocation._canBeChild(self, c):
if self._canBeChild(c):
self._addChildren(c)
for re in c._getAllElems(loc):
if re.name not in self.children:
self._addChildren(re)
s._addLocation(self)
# Add all species as children since users can set any species on vesicles and rafts
for spec in mdl.ALL(Species):
if spec.name not in self.children:
self._addChildren(spec)
def _canBeChild(self, c):
"""Return whether c can be a child of self."""
return isinstance(c, (Reaction, Diffusion, Exocytosis, RaftEndocytosis, RaftDis))
def addSystem(self, sys, _loc=None):
"""Add a vesicle or raft surface system to the location
:param sys: The volume or surface system to be added, or its name.
:type sys: :py:class:`VesicleSurfaceSystem`, :py:class:`RaftSurfaceSystem` or `str`
:returns: None
"""
if isinstance(sys, SpaceSystem):
self.sysNames.append((sys.name, _loc))
elif isinstance(sys, str):
self.sysNames.append((sys, _loc))
else:
raise TypeError(
f'Expected a VesicleSurfaceSystem or a RaftSurfaceSystem, got {sys} instead.'
)
class _VesicleSelection(nutils.SolverPathObject):
"""Class that represents a subselection of vesicles with additional information"""
def __init__(self, ves, inside):
self._ves = ves
self._inside = inside
def _solverStr(self):
"""Return the string that is used as part of method names for this specific object."""
return self._ves._solverStr() + ('Inner' if self._inside else 'Surface')
def _solverId(self):
"""Return the id that is used to identify an object in the solver."""
return self._ves._solverId()
def _simPathAutoMetaData(self):
"""Return a dictionary with string keys and string or numbers values."""
return {'ves_loc': self.loc}
def __getattr__(self, name):
return getattr(self._ves, name)
@property
def children(self):
return self._ves.children
@property
def ves(self):
return self._ves
@property
def loc(self):
return 'in' if self._inside else 'surf'
def __repr__(self):
return f"{self.ves}('{self.loc}')"
[docs]@nutils.FreezeAfterInit
class Vesicle(_RaftVesLocation):
"""Represents a vesicle object
Vesicles are described by a spherical diameter and diffusion rate.
It models the many behaviours and interactions of these complex biological
entities. A vesicle may be formed by endocytosis or by direct user input,
may exit by exocytosis or user input, may transport species on its surface or
luminally, and may undergo many interactions with its environment. These
interactions are individually described in each corresponding class.
:param diameter: Diameter of the vesicle (in m)
:type diameter: float
:param dcst: Default diffusion coefficient of the vesicle (in m^2 s^-1)
:param vesSurfSys: Optional vesicle surface system to be associated with the vesicle
:type vesSurfSys: :py:class:`VesicleSurfaceSystem`
"""
_locStr = 'Vesicle'
def __init__(self, diameter, dcst=0, vesSurfSys=None, _createObj=True, **kwargs):
super().__init__(**kwargs)
(mdl,) = self._getUsedObjects()
self._setParameter('Diameter', diameter, nutils.Units('m'))
self._setParameter('Dcst', dcst, nutils.Units('m^2 s^-1'))
if _createObj:
self.stepsVesicle = stepslib._py_Vesicle(self.name, mdl.stepsModel, self.Diameter, self.Dcst)
else:
self.stepsVesicle = None
if vesSurfSys is not None:
self.addSystem(vesSurfSys)
def _SetUpMdlDeps(self, mdl, geom):
"""Set up the structures that will allow species or reaction access from locations."""
super()._SetUpMdlDeps(mdl, geom)
# Add all relevant link species as children
for vsys in mdl.ALL(VolumeSystem):
for vbind in vsys.ALL(VesicleBind, VesicleUnbind):
for ves, ls in vbind._getVesLinkSpecPairs():
if ves == ves and ls.name not in self.children:
self._addChildren(ls)
def _canBeChild(self, c):
"""Return whether c can be a child of self."""
return super()._canBeChild(c)
def _getStepsObjects(self):
"""Return a list of the steps objects that this named object holds."""
return [self.stepsVesicle]
def _solverStr(self):
"""Return the string that is used as part of method names for this specific object."""
return Vesicle._locStr
[docs] def addSystem(self, sys, _loc=None):
"""Add a vesicle surface system to the vesicle
:param sys: The vesicle surface system to be added, or its name.
:type sys: Union[:py:class:`steps.API_2.model.VesicleSurfaceSystem`, str]
:returns: None
"""
super().addSystem(sys, Location.VESSURF)
if _loc is None:
if isinstance(sys, SpaceSystem):
sys = sys.name
self.stepsVesicle.addVesSurfsys(sys)
[docs] def __call__(self, loc):
"""Get a version of the vesicle with added information
Parentheses notation (function call) can be used on a vesicle to specify additional
information. This should never be needed during model declaration but is frequently needed
for simulation control and data saving (see :py:class:`steps.API_2.sim.SimPath`).
:param loc: If `'surf'`, the simulation path represents the surface of the selected vesicles,
if `'in'`, it represents the inside of the selected vesicles
:type inside: str
:returns: An object that represents the vesicle with the added information
:meta public:
"""
if loc not in {'surf', 'in'}:
raise ValueError(f"Location string can only be 'surf' or 'in', got {loc} instead.")
return _VesicleSelection(self, loc == 'in')
@classmethod
def _FromStepsObject(cls, obj, mdl):
"""Create the interface object from a STEPS object."""
ves = cls(obj.getDiameter(), obj.getDcst(), _createObj=False, name=obj.getID())
ves.stepsVesicle = obj
for sysname in obj.getVesSurfsys():
ves.addSystem(sysname)
return ves
@property
@nutils.ParameterizedObject.RegisterGetter(units=nutils.Units('m'))
def Diameter(self):
"""Vesicle diameter (in m)
:type: Union[float, :py:class:`steps.API_2.nutils.Parameter`], read-only
"""
pass
@property
@nutils.ParameterizedObject.RegisterGetter(units=nutils.Units('m^2 s^-1'))
def Dcst(self):
"""Diffusion constant (in m^2 s^-1)
:type: Union[float, :py:class:`steps.API_2.nutils.Parameter`]
"""
pass
@Dcst.setter
@nutils.ParameterizedObject.RegisterSetter(units=nutils.Units('m^2 s^-1'))
def Dcst(self, val):
if isinstance(val, numbers.Number):
self.stepsVesicle.setDcst(val)
else:
raise TypeError(f'{val} is not a valid diffusion constant.')
[docs]@nutils.FreezeAfterInit
class Raft(_RaftVesLocation):
"""Represents a membrane raft
The membrane analogy of vesicles, rafts effectively group surface species
within a defined radius and exist in and may diffuse within patches. Rafts
and the species they contain undergo special interactions such as raft
endocytosis, raft generation and raft dissolution.
Rafts can represent lipid rafts or other clusters. They occupy space
defined by a diameter projected onto the surface, and can be mobile.
:param diameter: Diameter of the raft (in m)
:type diameter: float
:param dcst: Default diffusion coefficient of the raft (in m^2 s^-1)
:param raftSurfSys: Optional raft surface system to be associated with the raft
:type raftSurfSys: :py:class:`RaftSurfaceSystem`
"""
_locStr = 'Raft'
def __init__(self, diameter, dcst=0, raftSurfSys=None, _createObj=True, **kwargs):
super().__init__(**kwargs)
(mdl,) = self._getUsedObjects()
self._setParameter('Diameter', diameter, nutils.Units('m'))
self._setParameter('Dcst', dcst, nutils.Units('m^2 s^-1'))
if _createObj:
self.stepsRaft = stepslib._py_Raft(self.name, mdl.stepsModel, self.Diameter, self.Dcst)
else:
self.stepsRaft = None
if raftSurfSys is not None:
self.addSystem(raftSurfSys)
def _canBeChild(self, c):
"""Return whether c can be a child of self."""
if not super()._canBeChild(c):
return False
if ((isinstance(c, Reaction) and not c._isRaftSurfReac()) or
isinstance(c, (Diffusion, Exocytosis))
):
return False
return True
def _getStepsObjects(self):
"""Return a list of the steps objects that this named object holds."""
return [self.stepsRaft]
def _solverStr(self):
"""Return the string that is used as part of method names for this specific object."""
return Raft._locStr
[docs] def addSystem(self, sys, _loc=None):
"""Add a raft surface system to the raft
:param sys: The raft surface system to be added, or its name.
:type sys: Union[:py:class:`steps.API_2.model.RaftSurfaceSystem`, str]
:returns: None
"""
super().addSystem(sys, Location.RAFTSURF)
if _loc is None:
if isinstance(sys, SpaceSystem):
sys = sys.name
self.stepsRaft.addRaftsys(sys)
@classmethod
def _FromStepsObject(cls, obj, mdl):
"""Create the interface object from a STEPS object."""
raft = cls(obj.getDiameter(), obj.getDcst(), _createObj=False, name=obj.getID())
raft.stepsRaft = obj
for sysname in obj.getRaftsys():
raft.addSystem(sysname)
return raft
@property
@nutils.ParameterizedObject.RegisterGetter(units=nutils.Units('m'))
def Diameter(self):
"""Raft diameter (in m)
:type: Union[float, :py:class:`steps.API_2.nutils.Parameter`], read-only
"""
pass
@property
@nutils.ParameterizedObject.RegisterGetter(units=nutils.Units('m^2 s^-1'))
def Dcst(self):
"""Raft diffusion constant (in m^2 s^-1)
:type: Union[float, :py:class:`steps.API_2.nutils.Parameter`]
"""
pass
@Dcst.setter
@nutils.ParameterizedObject.RegisterSetter(units=nutils.Units('m^2 s^-1'))
def Dcst(self, val):
if isinstance(val, numbers.Number):
self.stepsRaft.setDcst(val)
else:
raise TypeError(f'{val} is not a valid diffusion constant.')
####################################
[docs]@nutils.FreezeAfterInit
class Species(
nutils.UsingObjects(Model),
ReactionElement,
nutils.StepsWrapperObject,
nutils.ParameterizedObject
):
"""A chemical species
Species can be involved in :py:class:`Reaction`, :py:class:`Diffusion`, or other transport
mechanisms.
:param valence: Optional, the valence of the species (defaults to 0)
:type valence: int
"""
_elemStr = 'Spec'
def __init__(self, *args, valence=0, _createObj=True, **kwargs):
super().__init__(*args, **kwargs)
(mdl,) = self._getUsedObjects()
self.stepsSpecies = stepslib._py_Spec(self.name, mdl.stepsModel) if _createObj else None
if valence != 0:
self.valence = valence
def _getStepsObjects(self):
"""Return a list of the steps objects that this named object holds."""
return [self.stepsSpecies]
@classmethod
def _FromStepsObject(cls, obj, mdl):
"""Create the interface object from a STEPS object."""
spec = cls(_createObj=False, name=obj.getID())
spec.stepsSpecies = obj
val = obj.getValence()
if val != 0:
spec.valence = val
return spec
@property
@nutils.ParameterizedObject.RegisterGetter()
def valence(self):
"""The valence of the species (defaults to 0)
:type: Union[int, :py:class:`steps.API_2.utils.Parameter`]
"""
return self.stepsSpecies.getValence()
@valence.setter
@nutils.ParameterizedObject.RegisterSetter()
def valence(self, val):
self.stepsSpecies.setValence(val)
def __repr__(self):
return self.name
def _solverStr(self):
"""Return the string that is used as part of method names for this specific object."""
return Species._elemStr
[docs]@nutils.FreezeAfterInit
class LinkSpecies(
nutils.UsingObjects(Model),
ReactionElement,
nutils.StepsWrapperObject,
nutils.ParameterizedObject
):
"""A chemical species that can form links between vesicles
A link species is a special kind of chemical species that links two
vesicles (by forming a bond with another link species). A link species is
formed by a vesicle binding event and exists within a specified upper and
lower bound of length for the duration of its existence. A link species
may diffuse on a vesicle surface, but only within this length bound. Link
species are destroyed by a vesicle unbinding event.
:param dcst: Optional, diffusion coefficient of the vesicle (in m^2 s^-1, defaults to 0)
:type dcst: float
"""
_elemStr = 'LinkSpec'
def __init__(self, dcst=0, _createObj=True, **kwargs):
super().__init__(**kwargs)
(mdl,) = self._getUsedObjects()
self._setParameter('Dcst', dcst, nutils.Units('m^2 s^-1'))
self.stepsSpecies = stepslib._py_LinkSpec(self.name, mdl.stepsModel, self.Dcst) if _createObj else None
def _getStepsObjects(self):
"""Return a list of the steps objects that this named object holds."""
return [self.stepsSpecies]
@classmethod
def _FromStepsObject(cls, obj, mdl):
"""Create the interface object from a STEPS object."""
spec = cls(obj.getDcst(), _createObj=False, name=obj.getID())
spec.stepsSpecies = obj
return spec
@property
@nutils.ParameterizedObject.RegisterGetter(units=nutils.Units('m^2 s^-1'))
def Dcst(self):
"""Diffusion constant (in m^2 s^-1)
:type: float, read-only
"""
pass
def _solverStr(self):
"""Return the string that is used as part of method names for this specific object."""
return LinkSpecies._elemStr
[docs]@nutils.FreezeAfterInit
class SubUnitState(ReactionElement):
"""State of a :py:class:`SubUnit`
Can take part in reactions (see :py:class:`Complex`) and be combined with ``|`` and ``~`` to
form a :py:class:`SubUnitSelector`.
In reactions, further information can be provided with the square bracket operator, see
:py:func:`SubUnitState.__getitem__`.
A :py:class:`SubUnitState` can only be associated to a single :py:class:`SubUnit`.
.. warning::
When used in reactions in conjunction with the `|` operator, the expression must be wrapped
in parentheses to avoid incorrect operator precedence.
"""
_elemStr = 'SUS'
def __init__(self, _orig=None, *args, **kwargs):
super().__init__(*args, **kwargs)
self._subUnit = None
# Each subunit state creates its own modified version, for use in reaction computation
if _orig is None:
self._modVers = SubUnitState(name=f'special%mod', _orig=self)
self._orig = None
else:
self._modVers = None
self._orig = _orig
def _setSubUnit(self, su):
if self._subUnit is None:
self._subUnit = su
if self._modVers is not None:
self._modVers._setSubUnit(su)
else:
raise Exception('A SubUnitState can only be associated to one subunit.')
def _getModif(self):
return self._modVers if self._modVers is not None else self
def _getOrig(self):
return self if self._orig is None else self._orig
def _toReactionElem(self):
return self._toSelector()
def _areStatesAsSpecies(self):
return self._subUnit._areStatesAsSpecies()
[docs] def __or__(self, other):
"""Get the union between two :py:class:`SubUnitState`
:param other: The other state(s)
:type other: Union[:py:class:`SubUnitState`, :py:class:`SubUnitSelector`]
:returns: An object representing the union between the current state and the state(s)
from 'other'.
:rtype: :py:class:`SubUnitSelector`
Usage::
(S1A | S1B)
# It can also be chained
(S1A | S1B | S1C)
If the associated :py:class:`SubUnit` was defined with ``S1A``, ``S1B``, and ``S1C`` as
states, ``(S1A | S1B | S1C)`` is equivalent to ``:`` (all SubUnitStates).
.. warning::
Do not forget to wrap the overall expression in parentheses in order to avoid issues
with operator precedence when used in reaction declaration.
:meta public:
"""
return self._toSelector() | other
[docs] def __invert__(self):
"""Get all :py:class:`SubUnitState` except itself
:returns: An object representing all states of the associated :py:class:`SubUnit`
except itself
:rtype: :py:class:`SubUnitSelector`
Usage::
~S1A # All subunits except S1A
If the associated :py:class:`SubUnit` was defined with ``S1A``, ``S1B``, and ``S1C`` as
states, ``~S1A`` is equivalent to ``(S1B | S1C)``
:meta public:
"""
return SubUnitSelector(self._subUnit, [s for s in self._subUnit._states if s is not self])
def _toSelector(self):
return SubUnitSelector(self._subUnit, [self])
[docs] def __getitem__(self, key):
"""Specify additional information
Transform the subunit state into a :py:class:`SubUnitSelector` that only selects the
current state, and call ``__getitem__`` on it.
See :py:func:`SubUnitSelector.__getitem__` for details.
:meta public:
"""
return SubUnitSelector(self._subUnit, [self]).__getitem__(key)
def __iter__(self):
raise NotImplementedError()
def __repr__(self):
return self.name
def __eq__(self, other):
return isinstance(other, SubUnitState) and self.name == other.name
def __hash__(self):
return hash(self.name)
def _solverStr(self):
"""Return the string that is used as part of method names for this specific object."""
if self._areStatesAsSpecies():
return ''
else:
return SubUnitState._elemStr
def _solverId(self):
"""Return the id that is used to identify an object in the solver."""
if self._areStatesAsSpecies():
return tuple()
else:
return (self._subUnit._complex._getSUSInd(self),)
def _solverModifier(self):
"""Return None or a function that will modify the output from the solver."""
if self._areStatesAsSpecies():
return lambda x: (x, self)
else:
return None
[docs]@nutils.FreezeAfterInit
class SubUnitSelector(ReactionElement):
"""Set of :py:class:`SubUnitState` associated to a :py:class:`SubUnit`
Represents a selection among the states of a :py:class:`SubUnit`.
.. note::
Should not be instantiated by a user, the class is only documented for clarity.
"""
_OR = 0
_AND = 1
def __init__(self, su, states, *args, _pos=None, _compSel=None, _ident=None, **kwargs):
super().__init__(*args, **kwargs)
self._subUnit = su
if any(sus._subUnit != su for sus in states):
raise Exception(f'Some subunit states in {states} are not associated with subunit {su}.')
self._states = frozenset(states)
self._pos = _pos
self._compSel = _compSel
self._id = _ident
def _setOp(self, op, other):
if isinstance(other, SubUnitState):
other = other._toSelector()
if isinstance(other, SubUnitSelector):
if other._subUnit is not self._subUnit:
raise TypeError('Cannot combine subunit states from different subunits.')
if other._compSel != self._compSel:
raise Exception('Cannot combine subunit states from different complex selectors.')
if other._pos != self._pos:
raise Exception('Cannot combine subunit states from different positions.')
tmpId = self._id if self._id is not None else other._id
if op == SubUnitSelector._OR:
return SubUnitSelector(
self._subUnit,
self._states | other._states,
_pos=self._pos,
_compSel=self._compSel,
_ident=tmpId,
)
elif op == SubUnitSelector._AND:
return SubUnitSelector(
self._subUnit,
self._states & other._states,
_pos=self._pos,
_compSel=self._compSel,
_ident=tmpId,
)
else:
raise NotImplementedError()
elif isinstance(other, (ReactionSide, ReactionElement)):
# Common syntax error
raise SyntaxError(
'SubUnitSelectors with "|" or "&" operator need to be wrapped by '
'parentheses when used in reactions.'
)
else:
raise TypeError(f'Expected a SubUnitState or a SubUnitSelector, got {other} instead.')
[docs] def __getitem__(self, key):
"""Specify additional information
:param key: The information(s) (it is possible to provide several within the same square
brackets)
:type key: Union[int, str, :py:class:`ComplexSelector`]
:returns: The subunit selector, annotated with the information
:rtype: :py:class:`SubUnitSelector`
This method allows the usage of square brackets notation for providing the following
additional informations:
- Position of the associated :py:class:`SubUnit` in the complex::
with Comp[...]:
# Ca binds to the first subunit (index 0) when in state S1A
S1A[0] + Ca <r[1]> S1B[0]
- Association with a specific complex within a reaction::
with Comp[...] as C1, Comp[...] as C2:
# Two instances of Comp react, S1A in the first one is modified into S1B
# while S1B in the second is modified to S1C.
S1A[C1] + S1B[C2] <r[1]> S1B[C1] + S1C[C2]
- Custom identifier (useful when a reaction is ambiguous)::
with Comp[...]:
# A subunit in state S1A reacts with another subunit is state S1B. Since
# both states are associated to the same subunit (and this subunit is
# present twice in the complex), the following reaction could be read in
# two different ways:
# S1A + S1B >r[1]> S1B + S1C
# We could have (S1A, S1B) -> (S1B, S1C) or (S1A, S1B) -> (S1C, S1B). In
# the first case both subunit changed while in the second, only the first
# one changed. Custom identifiers allow us to specify this reaction
# without ambiguities:
S1A['a'] + S1B['b'] >r[1]> S1B['a'] + S1C['b']
Note that it is possible to provide several information at the same time. For example::
with Comp[...] as C1, Comp[...] as C2:
# Two instances of Comp react, in the first one (S1A, S1B) -> (S1B, S1C) while
# in the second one S1B -> S1A
S1A[C1, 'a'] + S1B[C1, 'b'] + S1B[C2] >r[1]> S1B[C1, 'a'] + S1C[C1, 'b'] + S1A[C2]
:meta public:
"""
if not isinstance(key, tuple):
key = (key,)
if len(key) == 0 or len(key) > 3:
raise KeyError(
'SubUnitStates can only be indexed with one complex selector, one '
'integer indicating their position, or one identifier string.'
)
res = copy.copy(self)
for v in key:
if isinstance(v, _ComplexSelectorSpawner):
v = v[...]
if isinstance(v, ComplexSelector):
if res._compSel is not None:
raise Exception('The SubUnitSelector is already associated with a complex selector.')
res._compSel = v
elif isinstance(v, numbers.Integral):
if res._pos is not None:
raise Exception('The SubUnitSelector is already associated with a position.')
res._pos = v
elif isinstance(v, str):
if res._id is not None:
raise Exception('The SubUnitSelector is already associated with an identifier.')
res._id = v
else:
raise TypeError(
f'Only complex selectors, strings, and integers can be used as a '
f'key, got {v} instead.'
)
return res
def _isEmpty(self):
return len(self._states) == 0
def _isFull(self):
"""Return True if the SubUnitSelector covers all states of its corresponding subunit."""
return len(self._states) == len(self._subUnit._states)
def _checkPos(self, comp):
"""Raise an exception is the position is not compatible with complex 'comp'."""
n = len(comp._subUnits)
if self._pos is not None and not (0 <= self._pos < n):
raise Exception(f'{self._pos} is not a valid position, it needs to be >= 0 and < {n}')
[docs] def __or__(self, other):
"""Get the union between the states of two :py:class:`SubUnitSelector`
:param other: The other subunit selector, has to be associated to the same
:py:class:`SubUnit`. Can also be a single :py:class:`SubUnitState`.
:type other: Union[:py:class:`SubUnitState`, :py:class:`SubUnitSelector`]
:returns: A subunit selector representing the union between the current selector states and
the state(s) from 'other'.
:rtype: :py:class:`SubUnitSelector`
See :py:class:`SubUnitState.__or__` for usage examples.
.. warning::
Do not forget to wrap the overall expression in parentheses in order to avoid issues
with operator precedence when used in reaction declaration.
:meta public:
"""
return self._setOp(SubUnitSelector._OR, other)
[docs] def __and__(self, other):
"""Get the intersection between the states of two :py:class:`SubUnitSelector`
:param other: The other subunit selector, has to be associated to the same
:py:class:`SubUnit`. Can also be a single :py:class:`SubUnitState`.
:type other: Union[:py:class:`SubUnitState`, :py:class:`SubUnitSelector`]
:returns: A subunit selector representing the intersection between the current selector
states and the state(s) from 'other'.
:rtype: :py:class:`SubUnitSelector`
Usage::
AorB = S1A | S1B
BorC = S1B | S1C
(AorB & BorC) # Equivalent to S1B
(AorB & (~S1A) & (~S1C)) # Equivalent to S1B
If the intersection is empty, it returns an empty :py:class:`SubUnitSelector` which
might raise exceptions if later used in a reaction declaration.
.. warning::
Do not forget to wrap the overall expression in parentheses in order to avoid issues
with operator precedence when used in reaction declaration.
:meta public:
"""
return self._setOp(SubUnitSelector._AND, other)
[docs] def __invert__(self):
"""Get all :py:class:`SubUnitState` except the ones in self
:returns: An object representing all states of the associated :py:class:`SubUnit`
minus the ones present in self
:rtype: :py:class:`SubUnitSelector`
Usage::
~(S1A | S1C) # All subunits except S1A and S1C
If the associated :py:class:`SubUnit` was defined with ``S1A``, ``S1B``, and ``S1C`` as
states, ``~(S1A | S1C)`` is equivalent to ``S1B``
:meta public:
"""
return SubUnitSelector(self._subUnit, self._subUnit._states - self._states)
def __iter__(self):
return iter(sorted(self._states, key=lambda x: x.name))
[docs] def __contains__(self, sus):
"""Test whether a subunit state is present in the subunit selector
:param sus: The subunit state
:type sus: :py:class:`SubUnitState`
:returns: True if the given subunit state is in the current selector, False otherwise.
:rtype: bool
Usage::
>>> AorB = S1A | S1B
>>> S1A in AorB
True
>>> S1C in AorB
False
:meta public:
"""
return isinstance(sus, SubUnitState) and sus._subUnit is self._subUnit and sus in self._states
def __repr__(self):
if self._states == self._subUnit._states:
return ':'
elif len(self._states) <= len(self._subUnit._states) / 2:
return '|'.join(map(str, self._states))
else:
# Express it as a negation if it is more concise
if len(self._states) == len(self._subUnit._states) - 1:
return f'~{~self}'
else:
return f'~({~self})'
def __eq__(self, other):
return isinstance(other, SubUnitSelector) and (
self._states,
self._compSel._complex if self._compSel is not None else None,
) == (other._states, other._compSel._complex if other._compSel is not None else None)
def __hash__(self):
return hash((self._states, self._compSel._complex if self._compSel is not None else None))
def _makeCopy(self):
"""Make a copy of the reaction element."""
return copy.copy(self)
[docs]@nutils.FreezeAfterInit
class ComplexState(ReactionElement):
"""Fully specified state of a complex
Can be used in the same way as a :py:class:`ReactionElement`.
.. note::
A :py:class:`ComplexState` object should not be directly instantiated by the user, instead
it is obtained by using the square bracket notation on a :py:class:`Complex` object::
Comp[S1B, S1B, S2A] # Returns the complex state for which the two S1 subunits are
# in state S1B and the S2 subunit is in state S2A.
See :py:func:`Complex.__getitem__` for details.
"""
_SIMPATH_ONLY_CHILDREN = True
def __init__(
self, comp, state, *args, _stepsObj=None, _parentCompSel=None, _parentSelName=None, **kwargs
):
self._comp = comp
# check that the subunit states are compatible with the subunits
for su, st in zip(self._comp._subUnits, state):
if st not in su._states:
raise TypeError(f'{st} is not a state of subunit {su}.')
self._state = self._reOrderState(self._comp._order(state))
self._stepsObj = _stepsObj
self._parentCompSel = _parentCompSel
self._parentSelName = _parentSelName
super().__init__(*args, name=self._getStateName(), **kwargs)
def _getReferenceObject(self):
"""
Return the object this object was derived from. Useful for getting the complex associated
with a complex selector, etc.
"""
return self._comp
def _getComplexName(self):
return self._comp.name
def _getStateName(self):
"""Return the name of the state."""
return self._comp.name + '_' + '_'.join(map(str, self._state))
def _setStepsObjects(self, obj):
self._stepsObj = obj
def _areStatesAsSpecies(self):
return self._comp._areStatesAsSpecies()
def _getStepsObjects(self):
if self._stepsObj is None and self._comp._areStatesAsSpecies():
return [self._comp._compStates[self]]
else:
return [self._stepsObj]
def _reOrderState(self, state):
"""
Transforms the output of the ordering function in such a way that subunit states are
associated to the correct subunits, i.e. reorder subunit states to match subunit order.
Conserve relative ordering of subunit states linked to identical subunits.
"""
subunits = copy.copy(self._comp._subUnits)
res = [None] * len(subunits)
for sus in state:
i = subunits.index(sus._subUnit)
res[i] = sus
subunits[i] = None
if any(sus is None for sus in res) or any(su is not None for su in subunits):
raise Exception(
f'An invalid complex state was used. There is probably an issue with the ordering function.'
)
return tuple(res)
def _toComplexSelector(self):
"""
Return a complex selector describing this complex state, keep the link to the original
complexSelector that was used to instantiante the state, if it exists.
"""
stateVals = zip(range(len(self._state)), self._comp._subUnits, self._state)
subSels = [SubUnitSelector(su, [s], _pos=i) for i, su, s in stateVals]
if self._parentCompSel is not None:
res = copy.copy(self._parentCompSel)
res._initSubSels(subSels)
return res
elif self._parentSelName is not None:
return ComplexSelector(self._comp, subSels, name=self._parentSelName)
else:
return ComplexSelector(self._comp, subSels)
def _toUnorderedVect(self):
"""
Return a subunit state pool vector representing the complex state
"""
vec = [0] * len(self._comp._poolElems)
for sus in self:
vec[self._comp._sus2PoolInd[sus]] += 1
return tuple(vec)
def _toUnorderedFilter(self):
"""
Return a filter representing the complex state
"""
vec = [(0, 0) for i in range(len(self._comp._poolElems))]
for sus in self:
_min, _max = vec[self._comp._sus2PoolInd[sus]]
vec[self._comp._sus2PoolInd[sus]] = (_min + 1, _max + 1)
return (tuple(vec),)
[docs] def Count(self, sus):
"""Count the number of subunits that match a given subunit selector
:param sus: A subunit state or a subunit selector
:type sus: Union[:py:class:`SubUnitSelector`, :py:class:`SubUnitState`]
:returns: The number of subunits in the complex state that match ``sus``.
:rtype: int
Usage::
>>> state = Comp[S1A, S1B, S2A]
>>> state.Count(S1B)
1
>>> state.Count(S1C)
0
>>> state.Count(S1A | S1B)
2
"""
if isinstance(sus, SubUnitState):
sus = sus._toSelector()
if not isinstance(sus, SubUnitSelector):
raise TypeError(f'Expected a SubUnitState or a SubUnitSelector, got {sus} instead.')
return sum(1 for s in self._state if s in sus)
[docs] def __iter__(self):
"""Iterator over all the subunit states of the complex state.
:returns: An iterator that iterates over all :py:class:`SubUnitState` in the complex
state. The order is the same as the one defined during the creation of the associated
:py:class:`Complex`.
:rtype: Iterator[:py:class:`SubUnitState`]
Usage::
>>> state = Comp[S1A, S1B, S2A]
>>> for subState in state:
... print(subState.name)
...
S1A
S1B
S2A
:meta public:
"""
return iter(self._state)
def __eq__(self, other):
return isinstance(other, ComplexState) and other._comp == self._comp and other._state == self._state
def __hash__(self):
return hash((self._comp, self._state))
[docs] def __or__(self, other):
"""Get a complex selector representing the union between the states of self and other
:param other: The other complex state(s), has to be associated to the same
:py:class:`Complex`.
:type other: Union[:py:class:`ComplexState`, :py:class:`ComplexSelector`]
:returns: An object representing the union between the current state and the complex
state(s) from 'other'.
:rtype: :py:class:`ComplexSelector`
Usage::
Comp[S1A, S1A, S2A] | Comp[S1B, S1B, S2A] # Represents only these two states
# Note that Comp[S1A | S1B, S1A | S1B, S2A] is equivalent to more than two states
Comp[S1A, S1A, S2A] | Comp[S1A, S1B, S2A] | Comp[S1B, S1B, S2A]
# If the custom order of the complex (See :py:class:`Complex`) is set to
# StrongOrdering, it is equivalent to 4 states:
Comp[S1A, S1A, S2A] | Comp[S1A, S1B, S2A] | Comp[S1B, S1A, S2A] | Comp[S1B, S1B, S2A]
.. warning::
Do not forget to wrap the overall expression in parentheses in order to avoid issues
with operator precedence when used in reaction declaration.
:meta public:
"""
if isinstance(other, ComplexState):
return self._toComplexSelector().__or__(other._toComplexSelector())
elif isinstance(other, ComplexSelector):
return self._toComplexSelector().__or__(other)
else:
raise TypeError('ComplexStates can only be combined with ComplexSelectors or ' 'ComplexStates.')
def __repr__(self):
return f'{self._comp.name}[' + ', '.join(map(str, self)) + ']'
def _solverStr(self):
"""Return the string that is used as part of method names for this specific object."""
if self._areStatesAsSpecies():
return Species._elemStr
else:
return Complex._elemStr
def _solverId(self):
"""Return the id that is used to identify an object in the solver."""
if self._areStatesAsSpecies():
return (self.name,)
else:
return (
self._comp.name,
self._toUnorderedFilter(),
)
def _solverModifier(self):
"""Return None or a function that will modify the output from the solver."""
if self._areStatesAsSpecies():
def modif(x):
if isinstance(x, tuple):
return x[0] * self.Count(x[1]._toSelector())
else:
return x
return modif
else:
return None
def _simPathCheckParent(self):
"""
Determines whether the object needs to be in the parent's children in order to be added to
a simulation path
"""
return self._areStatesAsSpecies()
class _ComplexReactionElement(ReactionElement):
def _areStatesAsSpecies(self):
pass
[docs]@nutils.FreezeAfterInit
class ComplexSelector(_ComplexReactionElement, nutils.MultiUsable):
"""Partially specified state of a complex
Can be used in the same way as a :py:class:`ReactionElement`.
.. note::
A :py:class:`ComplexSelector` object should not be directly instantiated by the user,
instead it is obtained by using the square bracket notation on a :py:class:`Complex`
object in a way that does not fully specify the whole state::
Comp[S1B, S1B, :] # Returns the complex state for which the two S1 subunits are
# in state S1B and the S2 subunit is in any state.
See :py:func:`Complex.__getitem__` for details.
"""
_OR = 0
_AND = 1
def __init__(self, comp, subSels, *args, **kwargs):
super().__init__(*args, **kwargs)
self._complex = comp
self._initSubSels(subSels)
def _getReferenceObject(self):
"""
Return the object this object was derived from. Useful for getting the complex associated
with a complex selector, etc.
"""
return self._complex
def _areStatesAsSpecies(self):
return self._complex._areStatesAsSpecies()
def _getComplexName(self):
return self._complex.name
def _initSubSels(self, subSels):
self._subSels = set(
[tuple(SubUnitSelector(su, su._states, _pos=i) for i, su in enumerate(self._complex._subUnits))]
)
for ss in subSels:
self._addSubSel(ss)
def _addSubSel(self, ss):
newSubSels = set()
for row in self._subSels:
newRow = row[: ss._pos] + (row[ss._pos] & ss,) + row[ss._pos + 1 :]
if not newRow[ss._pos]._isEmpty():
newSubSels.add(newRow)
self._subSels = newSubSels
def _isEmpty(self):
"""Return True if no state can match the ComplexSelector."""
return len(self._subSels) == 0
def _toUnorderedFilter(self):
"""
Return a filter vector representing the complex selector
If necessary, the returned vector contains several chained filter vectors.
"""
totvec = []
for row in self._subSels:
filt = []
for sus in self._complex._poolElems:
minv, maxv, realMax = 0, 0, 0
for ss in row:
if sus in ss:
maxv += 1
if len(ss._states) == 1:
minv += 1
if sus in ss._subUnit._states:
realMax += 1
filt.append((minv, maxv if maxv < realMax else COMPLEX_FILTER_MAX_VALUE))
totvec.append(tuple(filt))
return tuple(totvec)
def _setOp(self, op, other):
"""
Return the complex selector resulting from applying the operation op to self and other.
"""
if isinstance(other, ComplexState):
other = other._toComplexSelector()
if not isinstance(other, ComplexSelector):
raise TypeError('ComplexSelector objects can only be combined with other ' '_Complexselectors.')
if other._complex is not self._complex:
raise Exception('Cannot combine _Complexselectors associated to different complexes.')
cs = copy.copy(self)
if op == ComplexSelector._OR:
cs._subSels = cs._subSels | other._subSels
elif op == ComplexSelector._AND:
cs._subSels = set()
for ss1 in self._subSels:
for ss2 in other._subSels:
inter = tuple(su1 & su2 for su1, su2 in zip(ss1, ss2))
if not any(su._isEmpty() for su in inter):
cs._subSels.add(inter)
else:
raise NotImplementedError()
return cs
def _getAllStates(self):
"""Return the set of all states that match the complex selector."""
states = set()
for ss in self._subSels:
for comb in itertools.product(*[sus._states for sus in ss]):
states.add(ComplexState(self._complex, comb, _parentCompSel=self))
return sorted(states, key=lambda x: x.name)
def _getMatchingStates(self, rcs):
"""
Return the pairs of states needed to declare a reaction in which self is on the lhs and
rcs on the rhs.
"""
l2rstates = {}
for lss in self._subSels:
for comb in itertools.product(*[sorted(sus._states, key=lambda x: x.name) for sus in lss]):
# For each state in the lhs complex, try to find the corresponding state
# in the rhs complex
destComb = None
moreGeneralError = 0
for rss in rcs._subSels:
tmpDestComb = []
ok = True
# For each subunit
for i, s in enumerate(comb):
# If there are identifiers and they are not matching, just skip
if rss[i]._id != lss[i]._id:
ok = False
break
rstates = rss[i]._states
lstates = lss[i]._states
if len(lstates) < len(rstates):
# If the destination state is more general than the source state
moreGeneralError += 1
ok = False
break
if len(rstates) == 1:
# If the destination state modifies the subunit
tmpDestComb.append(next(iter(rstates)))
elif s in rstates:
# If the destination state can contain the source state
tmpDestComb.append(s)
else:
# Otherwise, stop considering
ok = False
break
if ok:
if destComb is None:
destComb = tmpDestComb
elif ComplexState(rcs._complex, destComb) != ComplexState(rcs._complex, tmpDestComb):
# Only raise an exception if the states corresponding to the
# combinations are different
raise Exception(
f'ComplexSelector {self} can be matched with ComplexSelector '
f'{rcs} in more than one way. The reaction is ambiguous. State '
f'{comb} could at least be matched to {destComb} or to '
f'{tmpDestComb}. Try adding position or identification to '
f'SubUnitStates.'
)
if destComb is None:
# Try to give a specific error message
if moreGeneralError == len(rcs._subSels):
raise Exception(
f'ComplexSelector {self} cannot be matched with ComplexSelector {rcs} '
f'because the latter is more general than the former. This can happen '
f'when the sides of a reaction using SubUnitStates are not balanced.'
)
else:
# If not possible, give a generic one
raise Exception(
f'ComplexSelector {self} cannot be matched with ComplexSelector {rcs} '
f'for state {comb}. The reaction is undefined.'
)
lstate = ComplexState(self._complex, comb, _parentCompSel=self)
rstate = ComplexState(rcs._complex, destComb, _parentCompSel=rcs)
l2rstates.setdefault(lstate, set()).add(rstate)
allPairs = []
for ls, rss in l2rstates.items():
for rs in sorted(rss, key=lambda s: s.name):
allPairs.append((ls, rs, 1 / len(rss)))
return allPairs
[docs] def __or__(self, other):
"""Get a complex selector representing the union between the states of self and other
:param other: The other complex state(s), has to be associated to the same
:py:class:`Complex`.
:type other: Union[:py:class:`ComplexSelector`, :py:class:`ComplexState`]
:returns: A complex selector representing the union between the states in the current
selector and the state(s) from 'other'.
:rtype: :py:class:`ComplexSelector`
Usage:
States that have the first two subunits as ``S1A`` *or* the first as ``S1A`` and third
as ``S2A``::
Comp[S1A, S1A, :] | Comp[S1A, :, S2A]
Beware, the above expression is **NOT** equivalent to ``Comp[S1A, :, :]``::
>>> for state in (Comp[S1A, S1A, :] | Comp[S1A, :, S2A]):
... print(state.name)
...
Comp_S1A_S1A_S2A
Comp_S1A_S1A_S2B
Comp_S1A_S1B_S2A
Comp_S1A_S1C_S2A
>>> for state in Comp[S1A, :, :]:
... print(state.name)
...
Comp_S1A_S1A_S2A
Comp_S1A_S1A_S2B
Comp_S1A_S1B_S2A
Comp_S1A_S1B_S2B
Comp_S1A_S1C_S2A
Comp_S1A_S1C_S2B
In general, an OR operation between two :py:class:`ComplexSelector` is **NOT**
equivalent to doing OR operations between their :py:class:`SubUnitSelector`.
.. warning::
Do not forget to wrap the overall expression in parentheses in order to avoid issues
with operator precedence when used in reaction declaration.
:meta public:
"""
return self._setOp(ComplexSelector._OR, other)
[docs] def __and__(self, other):
"""Get a complex selector representing the intersection between the states of self and other
:param other: The other complex state(s), has to be associated to the same
:py:class:`Complex`.
:type other: Union[:py:class:`ComplexSelector`, :py:class:`ComplexState`]
:returns: A complex selector representing the intersection between the states in the
current selector and the state(s) from 'other'.
:rtype: :py:class:`ComplexSelector`
Usage:
States that have the first two subunits as ``S1A`` *and* the first as ``S1A`` and third
as ``S2A``::
Comp[S1A, S1A, :] & Comp[S1A, :, S2A]
The above expression corresponds to only one state: ``Comp[S1A, S1A, S2A]``. Note that,
in constrast to :py:func:`__or__`, this does correspond to a series of AND operations
between the :py:class:`SubUnitSelector`.
If the intersection is empty, it returns an empty :py:class:`ComplexSelector` which
might raise exceptions if later used in a reaction declaration.
.. warning::
Do not forget to wrap the overall expression in parentheses in order to avoid issues
with operator precedence when used in reaction declaration.
:meta public:
"""
return self._setOp(ComplexSelector._AND, other)
[docs] def __lshift__(self, subSels):
"""Get a complex selector representing the insertion of one or several subunit selector(s).
This operator is mainly used with complexes that have a :py:func:`StrongOrdering`, meaning
that subunits have specific positions. With this ordering, ``Comp[S1A, S1B, S2A]`` is
not the same as ``Comp[S1B, S1A, S2A]`` (the order of the first two subunits was changed).
Since, with this order function, the specific positions of the subunits matter, it would be
cumbersome to write ``Comp[S1A, :, :] | Comp[:, S1A, :]`` to express "a SU1 subunit is in
state S1A, in any position". And it would become harder and harder as the number of
identical subunits grows. The ``<<`` operator helps expressing this type of selectors.
:param subSels: One or several subunit states or selectors. If several, they should be
provided in the form of a :py:class:`ReactionSide` (``S1A + S1B`` for example).
:type subSels: Union[:py:class:`SubUnitState`, :py:class:`SubUnitSelector`,
:py:class:`ReactionSide`]
:returns: A complex selector to which subunit selectors present in subSels were
added in all possible positions.
:rtype: :py:class:`ComplexSelector`
Assuming we have a complex that has 3 identical subunits, declared like so::
with mdl:
S3A, S3B = SubUnitState.Create()
S3U = SubUnit.Create([S3A, S3B])
Comp2 = Complex([S3U, S3U, S3U], order=StrongOrdering, statesAsSpecies=True)
Usage of the ``<<`` operator::
(Comp2[...] << S3A) # Equivalent to:
# (Comp2[S3A, : , : ] |
# Comp2[ : , S3A, : ] |
# Comp2[ : , : , S3A])
(Comp2[...] << S3A + S3B) # Equivalent to:
# (Comp2[S3A, S3B, : ] |
# Comp2[S3A, : , S3B] |
# Comp2[S3B, S3A, : ] |
# Comp2[ : , S3A, S3B] |
# Comp2[S3B, : , S3A] |
# Comp2[ : , S3B, S3A])
.. warning::
Do not forget to wrap the overall expression in parentheses in order to avoid issues
with operator precedence when used in reaction declaration.
:meta public:
"""
if isinstance(subSels, (ReactionElement, ReactingElement)):
subSels = subSels._toReactionSide()
if not isinstance(subSels, ReactionSide):
raise TypeError(f'Expected reaction elements, got {subSels} instead.')
if not all(isinstance(el._elem, SubUnitSelector) for el in subSels.reacElems):
raise TypeError(f'Some elements in "{subSels}" are not subunit states or selectors.')
return self._mergeSubSels(subSels, enforceNoLocation=True)
[docs] def __iter__(self):
"""Iterator over all the possible states of the complex selector.
:returns: An iterator that iterates over all unique (as defined by the custom ordering)
:py:class:`ComplexState` of the complex selector.
:rtype: Iterator[:py:class:`ComplexState`]
Usage::
>>> selector = Comp[S1A, :, S2A]
>>> for state in selector:
... print(state.name)
...
Comp_S1A_S1A_S2A
Comp_S1A_S1B_S2A
Comp_S1A_S1C_S2A
:meta public:
"""
return iter(self._getAllStates())
[docs] def __len__(self):
"""Return the number of states covered by the complex selector.
:returns: The number of unique states covered by the complex selector.
:rtype: int
Usage::
>>> print(len(Comp[S1A, :, S2A]))
3
:meta public:
"""
return len(self._getAllStates())
def _simPathWalkExpand(self):
"""Return an iterable of the elements that should be part of ResultSelector paths."""
if self._areStatesAsSpecies():
return self.__iter__()
else:
return [self]
def _simPathCombinerClass(self):
"""Return the class that needs to be used to combine expanded elements."""
if self._areStatesAsSpecies():
return nutils.SumSimPathComb
else:
return None
def _simPathCheckParent(self):
"""
Determines whether the object needs to be in the parent's children in order to be added to
a simulation path
"""
return False
def _mergeSubSels(self, subSels, enforceNoLocation=False):
"""
Return a complex selector to which SubUnitSelectors present in subSels were added in all
possible positions if they were unordered or in their position if specified.
"""
if enforceNoLocation and subSels._hasLocation():
raise Exception('Location cues are not accepted when building a ComplexSelector.')
cs = self
for re in subSels.reacElems:
for s in range(re.stoich):
newcs = None
ss = re._elem
if ss._pos is None:
# If the subunitselector has no specified positions, generate the different
# positions available
for p in cs._getFreeSlots(ss):
ss2 = copy.copy(ss)
ss2._pos = p
if newcs is None:
newcs = cs._applysSubSel(ss2)
else:
newcs |= cs._applysSubSel(ss2)
else:
newcs = cs._applysSubSel(ss)
if newcs is None or newcs._isEmpty():
raise Exception(f'Complex selector {cs} is incompatible with SubUnitSelector ' f'{ss}.')
cs = newcs
return cs
def _applysSubSel(self, ss):
res = copy.copy(self)
res._subSels = set()
for row in self._subSels:
if len(row[ss._pos]._states) > len(ss._states) and not (row[ss._pos] & ss)._isEmpty():
newSus = row[ss._pos] & ss
res._subSels.add(row[: ss._pos] + (newSus,) + row[ss._pos + 1 :])
return res
def _getFreeSlots(self, ss):
"""
Return the subunit positions where ss could be applied. Only considers the positions that
are completely undetermined.
"""
for i in range(len(self._complex._subUnits)):
for st in self._subSels:
if st[i]._subUnit is ss._subUnit and st[i]._isFull():
yield i
break
def _getFullyDeterminedState(self):
"""
If the ComplexSelector only represents a single fully determined state, return this state,
otherwise, return None.
"""
fullState = None
for row in self._subSels:
if any(len(ss._states) != 1 for ss in row):
return None
state = tuple(next(iter(ss._states)) for ss in row)
if fullState is not None and state != fullState:
return None
fullState = state
return ComplexState(self._complex, fullState, _parentCompSel=self)
def _getUpdateVects(self, rcs):
"""
Return a list of pool vectors tuples corresponding to the complex reaction.
"""
if len(rcs._subSels) > 1:
raise Exception(
f'Cannot declare a reaction with composite complex selectors on the right hand '
f'side: {rcs}'
)
rrow = next(iter(rcs._subSels))
comp = self._complex
n = len(comp._poolElems)
filt = self._toUnorderedFilter()
updates = set()
for lrow in self._subSels:
# For each pair of SubUnitSelector rows
lds, rds = [], []
for lss, rss in zip(lrow, rrow):
if lss != rss and not lss._isFull():
if len(rss._states) == 1:
rsus = next(iter(rss))
lds.append(list(lss._states))
rds.append(rsus)
elif any(lsus not in rss for lsus in lss):
# Should never happen
raise Exception()
# For each possible state of the considered SubUnits
for comb in itertools.product(*lds):
req = [0] * n
upd = [0] * n
for lsus, rsus in zip(comb, rds):
req[comp._sus2PoolInd[lsus]] += 1
upd[comp._sus2PoolInd[lsus]] += -1
upd[comp._sus2PoolInd[rsus]] += 1
updates.add((tuple(req), tuple(upd)))
react = [0] * n
return [(filt, react, sorted(updates))]
@staticmethod
def _removeMoreSpecific(lst, func):
"""
Takes a list of identical size iterables and returns a simplified list by removing the
elements that are strictly more specific than other elements. The function 'func' should
take two elements and return True if the first one is more specific than the other.
"""
finallst = []
while len(lst) > 0:
delInd = 0
for i, vect1 in enumerate(lst):
totLst = finallst + lst[0:i] + lst[i + 1 :]
if any(func(vect1, vect2) for vect2 in totLst):
# If vect1 is always more specific than at least one other vector, remove vect1
break
delInd += 1
finallst += lst[0:delInd]
lst = lst[delInd + 1 :]
return finallst
def __repr__(self):
return ' | '.join(f'{self._complex.name}[' + ', '.join(map(str, row)) + ']' for row in self._subSels)
def __eq__(self, other):
return isinstance(other, ComplexSelector) and self.name == other.name
def __hash__(self):
return hash(self.name)
def _solverStr(self):
"""Return the string that is used as part of method names for this specific object."""
if self._areStatesAsSpecies():
return Species._elemStr
else:
return Complex._elemStr
def _solverId(self):
"""Return the id that is used to identify an object in the solver."""
if self._areStatesAsSpecies():
raise Exception('_solverId should never be called when statesAsSpecied is True.')
return (
self._complex.name,
self._toUnorderedFilter(),
)
@nutils.FreezeAfterInit
class _ComplexReactants(_ComplexReactionElement):
"""
Represents the combination of a complex selectors and a list of reactants involved in a
reaction. The reactants are SubUnitSelectors and can be associated with any compatible position
in the complex selector.
"""
def __init__(self, filt, reactants, *args, **kwargs):
super().__init__(*args, **kwargs)
self._filter = filt
self._reactants = []
for re in reactants.reacElems:
self._reactants += [re._elem] * re.stoich
def _areStatesAsSpecies(self):
return self._filter._areStatesAsSpecies()
def _hasSameReactants(self, cr):
"""Return True if cr has the same reactants in the same stoichiometry."""
return set((r, self._reactants.count(r)) for r in self._reactants) == set(
(r, self._reactants.count(r)) for r in cr._reactants
)
def _getMatchingStates(self, rcs):
"""
Return the pairs of states needed to declare a reaction in which self is on the lhs and
rcs on the rhs.
"""
# Compute the pairs of reactants
reactantPairs = self._getMatchingReactants(rcs)
# Order them from the smallest to the biggest to ensure correct positions computation
reactantPairs = sorted(reactantPairs, key=lambda x: (len(x[0]._states), hash(x[0])))
groupedStatePairs = {}
for filtSubSels in self._filter._subSels:
# Find all possible positions for the pairs of reactants
for posReacPairs in _ComplexReactants._getAllReactantsPositions(filtSubSels, reactantPairs):
# Get resulting src and dst rows
src = filtSubSels
dst = filtSubSels
reactantPos = set()
for lrp, rrp, p in posReacPairs:
lrp = copy.copy(lrp)
lrp._pos = p
rrp = copy.copy(rrp)
rrp._pos = p
src = src[0:p] + (lrp & src[p],) + src[p + 1 :]
if (rrp & dst[p])._isEmpty():
dst = dst[0:p] + (rrp,) + dst[p + 1 :]
else:
dst = dst[0:p] + (rrp & dst[p],) + dst[p + 1 :]
reactantPos.add(p)
# Generate states from src
for comb in itertools.product(*[sorted(sus._states, key=lambda x: x.name) for sus in src]):
# Compute dst state
dstState, dstModState = [], []
for srcState, srcSubSel, dstSubSel, p in zip(comb, src, dst, range(len(comb))):
if len(dstSubSel._states) == 1:
sus = next(iter(dstSubSel._states))
dstState.append(sus)
# Tag the destination state as modified, useful if destination is the same as src
dstModState.append(sus._getModif() if p in reactantPos else sus)
elif srcSubSel == dstSubSel:
dstState.append(srcState)
# Tag the destination state as modified, useful if destination is the same as src
dstModState.append(srcState._getModif() if p in reactantPos else srcState)
else:
raise Exception(
f'Cannot tranform state {srcState} associated with '
f'SubUnitSelector {srcSubSel} to SubUnitSelector '
f'{dstSubSel}, reaction is ambiguous.'
)
# Get the corresponding generic state (identical if using ordered Complexes)
srcos = tuple(comb)
dstos = tuple(dstState)
dstoms = tuple(dstModState)
srcgs = ComplexState(self._filter._complex, srcos, _parentCompSel=self)
dstgs = ComplexState(rcs._filter._complex, dstos, _parentCompSel=rcs)
# Add the state pairs to the pair counting structure
key = (srcgs, dstgs)
groupedStatePairs.setdefault(key, {}).setdefault(srcos, set()).add(dstoms)
if len(groupedStatePairs) == 0:
raise Exception(
f'It was not possible to find a combination of subunit states position compatible with the '
f'complex selector.'
)
allPairs = []
# Count the pairs to determine the required rate multiplier
for groupedPair, pairs in groupedStatePairs.items():
# Count total number of subpairs
totSubPairs = sum(len(dstStates) for srcos, dstStates in pairs.items())
rateMult = totSubPairs / len(pairs.keys())
allPairs.append(groupedPair + (rateMult,))
return allPairs
@staticmethod
def _getAllReactantsPositions(filtSubSels, reactantPairs):
"""
Compute all the possible positions for reactant pairs given a filter.
"""
if len(reactantPairs) == 0:
yield []
else:
lrp, rrp = reactantPairs[0]
# Get the possible positions for the first pair
if lrp._pos is None:
positions = list(_ComplexReactants._getFreeSlotsForReactantsPair(filtSubSels, lrp, rrp))
else:
if (filtSubSels[lrp._pos] & lrp)._isEmpty():
positions = []
else:
positions = [lrp._pos]
# Recursively compute the other positions
for p in positions:
empty = lrp & ~lrp
newFilt = filtSubSels[0:p] + (empty,) + filtSubSels[p + 1 :]
for subRPpos in _ComplexReactants._getAllReactantsPositions(newFilt, reactantPairs[1:]):
if len(subRPpos) == len(reactantPairs[1:]) :
yield [(lrp, rrp, p)] + subRPpos
@staticmethod
def _getFreeSlotsForReactantsPair(filtSubSels, lrp, rrp):
"""
Return the positions in filtSubSels for which both lrp and rrp could fit. Only considers
the positions that are completely undetermined.
"""
for i, ss in enumerate(filtSubSels):
unpos = copy.copy(ss)
unpos._pos = None
if lrp._subUnit is ss._subUnit and not (lrp & unpos)._isEmpty():
yield i
def _getMatchingReactants(self, rcs):
"""
Check that the reaction is balanced and return the pairs of reactants needed to
Generate the lhs and rhs states.
"""
# Group by subunit
lhsSU2React = {}
for r in self._reactants:
lhsSU2React.setdefault(r._subUnit, []).append(r)
rhsSU2React = {}
for r in rcs._reactants:
rhsSU2React.setdefault(r._subUnit, []).append(r)
# Check that the same subUnits are covered by both sides
if set(lhsSU2React.keys()) != set(rhsSU2React.keys()):
raise Exception(
f'The Subunits involved in the LHS of the reaction ({lhsSU2React.keys()}) '
f'are different from those involved in the RHS ({rhsSU2React.keys()})'
)
reactantPairs = []
# For each subunit
for su, lhsReacts in lhsSU2React.items():
# First perform basic checks
if su not in rhsSU2React:
raise Exception(f'SubUnitSelectors {lhsReacts} have no match in the RHS of the reaction.')
rhsReacts = rhsSU2React[su]
if len(lhsReacts) != len(rhsReacts):
raise Exception(
f'The number of SubUnitSelectors associated to SubUnit {su} should '
f'be the same on both sides of the reaction.'
)
# Try to find an unambiguous matching
for lhsr in lhsReacts:
tmpDest = []
for rhsr in rhsReacts:
# If there are identifiers or positions and they are not matching, just skip
if (
lhsr._id == rhsr._id
and lhsr._pos == rhsr._pos
and len(lhsr._states) >= len(rhsr._states)
):
tmpDest.append(rhsr)
# Check that the matching is possible
if len(tmpDest) == 0:
raise Exception(
f'SubunitSelector {lhsr} from the LHS of the reaction could '
f'not be matched with any SubUnitSelector from the RHS.'
)
if len(set(tmpDest)) > 1:
raise Exception(
f'The reaction is ambiguous: SubUnitSelector {lhsr} from the LHS '
f'could be matched with any of {tmpDest} from the RHS.'
)
# Add the pair
reactantPairs.append((lhsr, tmpDest[0]))
# Remove therhs one that just got added
rhsReacts.remove(tmpDest[0])
# If there are SubUnitSelector from the rhs that were not added
if len(rhsReacts) > 0:
raise Exception(
f'SubUnitSelectors from the RHS {rhsReacts} could not be matched '
f'with any SubUnitSelectors from the LHS.'
)
# TODO Not urgent: Throw exception if reactantPairs is empty?
return reactantPairs
def _getUpdateVects(self, rcs):
"""
Return a list of pool vectors tuples corresponding to the complex reaction.
"""
comp = self._filter._complex
n = len(comp._poolElems)
allFilts = self._filter._toUnorderedFilter()
allVects = []
for filt in allFilts:
allReactants = [[0] * n]
allUpdates = [([0] * n, [0] * n)]
# Get matching reactants
for lhsr, rhsr in self._getMatchingReactants(rcs):
newReactants = []
newUpdates = []
for lsus in sorted(lhsr._states, key=lambda x: x.name):
for reac, (req, upd) in zip(allReactants, allUpdates):
nr = copy.copy(reac)
nreq, nupd = copy.copy(req), copy.copy(upd)
ind = comp._sus2PoolInd[lsus]
nr[ind] += 1
nreq[ind] += 1
nupd[ind] += -1
if len(rhsr._states) == 1:
rsus = next(iter(rhsr._states))
elif lsus in rhsr._states:
rsus = lsus
else:
# This should never happen
raise Exception()
nupd[comp._sus2PoolInd[rsus]] += 1
newReactants.append(nr)
newUpdates.append((nreq, nupd))
allReactants = newReactants
allUpdates = newUpdates
# Update the filters with the reactants
for reac, (req, upd) in zip(allReactants, allUpdates):
newFilt = []
for i, r in enumerate(reac):
# SubUnitStates from the filter can also be reactants
newFilt.append((max(filt[i][0], r), filt[i][1]))
allVects.append(([newFilt], reac, [(req, upd)]))
# Remove duplicates or more specific
def func(a, b):
return all(
ra == rb and a[0][0][i][0] >= b[0][0][i][0] for i, ra, rb in zip(range(len(a[2])), a[2], b[2])
)
return ComplexSelector._removeMoreSpecific(allVects, func)
def _getReferenceObject(self):
"""
Return the object this object was derived from. Useful for getting the complex associated
with a complex selector, etc.
"""
return self._filter._getReferenceObject()
@property
def _complex(self):
return self._filter._complex
def __eq__(self, other):
return isinstance(other, _ComplexReactants) and self._filter == other._filter
def __hash__(self):
return hash(self._filter)
###################################################################################################
# Reactions
[docs]@nutils.FreezeAfterInit
class ReactionSide:
"""Unordered collection of reactants
Represents one of the sides of a reaction. Reactants can be added with `+`.
Reactants are unordered.
.. note::
This class should not be instantiated by the user, it is only documented for clarity.
"""
def __init__(self, elems, *args, **kwargs):
super().__init__(*args, **kwargs)
self.reacElems = list(elems)
def _add(self, other):
res = copy.copy(self)
if other is None or self.reacElems is None:
raise Exception('Cannot add reactant to None.')
if isinstance(other, ReactingElement):
res.reacElems = self.reacElems + [other]
elif isinstance(other, ReactionElement):
res.reacElems = self.reacElems + [ReactingElement(other, 1)]
elif isinstance(other, ReactionSide):
res.reacElems = self.reacElems + other.reacElems
else:
raise TypeError(f'Cannot add {other} to a reaction side.')
return res
[docs] def __add__(self, other):
"""Get the reaction side resulting from the addition of other to self
:param other: The object to be added to the reaction side
:type other: Union[:py:class:`ReactionElement`, :py:class:`ReactionSide`]
:returns: The reaction side resulting from the addition of other to self
:rtype: :py:class:`ReactionSide`
:meta public:
"""
return self._add(other)
def __radd__(self, other):
return self._add(other)
def __or__(self, other):
"""
Raise an exception for the common case in which a SubUnitSelector made from '|' operator
was not wrapped in parentheses.
"""
if isinstance(other, SubUnitState):
raise SyntaxError(
'SubUnitSelectors with "|" operator need to be wrapped by '
'parentheses when used in reactions.'
)
else:
raise SyntaxError('Cannot use "|" operator in reactions except for combining SubUnitStates.')
def _SetLoc(self, loc):
for s in self.reacElems:
s.loc = loc
def _GetStepsElems(self, loc=None):
res = []
for re in self.reacElems:
if loc is None or re.loc == loc:
res += re._getStepsObjects() * re.stoich
return res
def _getElemsOfType(self, *cls):
for e in self.reacElems:
for c in cls:
if isinstance(e._elem, c):
yield e
break
def _hasElemsOfType(self, *cls):
try:
next(self._getElemsOfType(*cls))
return True
except StopIteration:
return False
def _isEmpty(self):
return len(self.reacElems) == 0
def _getAllElems(self, loc=None):
res = set()
for re in self.reacElems:
if loc is None or re.loc == loc:
res.add(re._elem)
if isinstance(re._elem, ComplexState):
res.add(re._elem._comp)
if isinstance(re._elem, ComplexSelector):
res.add(re._elem._complex)
return res
def _isSurfaceSide(self):
if len(self.reacElems) == 0:
return None
elif all(s.loc is None for s in self.reacElems):
return False
elif all(s.loc is not None for s in self.reacElems):
return True
else:
raise Exception(
'Cannot mix localized and non localized reacElems in the same reaction side.'
)
def _hasSingleLocation(self):
return len(set(s.loc for s in self.reacElems)) < 2
def _hasVolumeLocation(self):
return all(s.loc not in ALL_SURF_LOCATIONS for s in self.reacElems)
def _hasLocation(self):
return any(s.loc is not None for s in self.reacElems)
def _getLocation(self):
if len(self.reacElems) > 0 and self._hasSingleLocation():
return self.reacElems[0].loc
else:
return None
def __repr__(self):
if len(self.reacElems) > 0:
return ' + '.join(map(str, self.reacElems))
else:
return 'Void'
def __iter__(self):
return iter(self.reacElems)
def __len__(self):
return len(self.reacElems)
def _toReactionSide(self):
return self
def _getFlatList(self):
"""Return an equivalent list of reacting elements that all have stoich = 1"""
res = []
for re in self:
res += [ReactingElement(re._elem, 1, re.loc)] * re.stoich
return res
def _makeCopy(self):
"""Make a copy of the reaction side."""
res = []
for e in self.reacElems:
e2 = copy.copy(e)
e2._elem = e._elem._makeCopy()
res.append(e2)
return ReactionSide(res)
def _getFinal(self, globCompSels):
"""
Return a ReactionSide that only contains ComplexSelectors and Species. Transform all
SubUnitSelectors into ComplexSelectors.
"""
# Split the ReactionSide into ComplexSelectors, SubUnitSelectors, and Species
compSelElems = list(self._getElemsOfType(ComplexSelector))
specElems = list(self._getElemsOfType(Species, LinkSpecies))
compStateElems = list(self._getElemsOfType(ComplexState))
subSelElems = list(self._getElemsOfType(SubUnitSelector))
# Transform complex states in complex selectors if they are real steps complexes
tmpCS = []
for cs in compStateElems:
if not cs._elem._comp._areStatesAsSpecies():
cs._elem = cs._elem._toComplexSelector()
compSelElems.append(cs)
else:
tmpCS.append(cs)
compStateElems = tmpCS
# Automatically attribute ComplexSelectors to SubUnitSelectors
compSels2subSels = {}
for re in subSelElems:
ss = re._elem
if ss._compSel is None:
if len(globCompSels) == 0:
raise Exception(
f'The reaction was declared outside of a ComplexSelector and '
f'SubUnitSelector {ss} was not associated with any '
f'ComplexSelector.'
)
potCompSels = [c for c in globCompSels if ss._subUnit in c._complex._subUnits]
if len(potCompSels) == 0:
raise Exception(
f'SubUnitSelector {ss} does not match any of the complex '
f'selectors: {globCompSels}.'
)
if len(potCompSels) > 1:
raise Exception(
f'SubUnitSelector {ss} does not have a clear complex selector.'
f' It could be matched with any of these: {potCompSels}.'
)
compSel = potCompSels[0]
else:
compSel = ss._compSel
compSels2subSels.setdefault(compSel, ReactionSide([]))
# Check that the positions cues are valid for this complex
ss._checkPos(compSel._complex)
# Clear compSel, it should not be used from this point on
ss._compSel = None
compSels2subSels[compSel] += re
# Generate ComplexReactants objects
resCs = ReactionSide([])
for cs, subSels in compSels2subSels.items():
if not subSels._hasSingleLocation():
raise Exception(f'SubUnits from the same complex have different locations: {subSels}.')
if cs in [e._elem for e in compSelElems]:
raise Exception(
f'Cannot mix ComplexSelector and SubunitSelectors for the same '
f'ComplexSelector in a reaction side. Either use one complex '
f'selector or only SubUnitSelectors.'
)
compReacts = _ComplexReactants(cs, subSels)
resCs += ReactionSide([ReactingElement(compReacts, 1, subSels._getLocation())])
# Reassemble the ReactionSide
finalSide = (
ReactionSide(compSelElems) + resCs + ReactionSide(specElems) + ReactionSide(compStateElems)
)
# Check that each ComplexSelector is only present once
allCompSels = list(finalSide._getElemsOfType(ComplexSelector))
if len(allCompSels) > 0:
compSels, stoichs = zip(*[(re._elem, re.stoich) for re in allCompSels])
for cs, st in zip(compSels, stoichs):
if st > 1 or compSels.count(cs) > 1:
raise Exception(
f'Complex selector {cs} was used more than once in the reaction side.'
)
return finalSide
@nutils.FreezeAfterInit
class _SubReaction(nutils.NamedObject):
"""
Wrapper for a single steps reaction object.
"""
_elemStrDict = {
stepslib._py_Reac: 'Reac',
stepslib._py_SReac: 'SReac',
stepslib._py_VDepSReac: 'VDepSReac',
stepslib._py_VesSReac: 'VesSReac',
stepslib._py_RaftSReac: 'RaftSReac',
stepslib._py_ComplexReac: 'ComplexReac',
stepslib._py_ComplexSReac: 'ComplexSReac',
}
def __init__(self, parent, stepsReac, lhs=None, rateMult=1, *args, **kwargs):
super().__init__(*args, **kwargs)
self._parent = parent
self._stepsReac = stepsReac
cls = stepsReac.__class__
if cls not in _SubReaction._elemStrDict:
raise TypeError(f'Unknown steps reaction type {cls}.')
self._elemStr = _SubReaction._elemStrDict[cls]
self._lhs = lhs
self._rateMult = rateMult
def _getRateUnits(self):
"""Return rate units of the reaction"""
return self._parent._getRateUnits()
def _getStepsObjects(self):
"""Return a list of the steps objects that this named object holds."""
return [self._stepsReac]
def _solverStr(self):
"""Return the string that is used as part of method names for this specific object."""
return self._elemStr
def _simPathAutoMetaData(self):
"""Return a dictionary with string keys and string or numbers values."""
return {'obj_type': 'Reac', 'obj_id': self.name}
def _solverSetValue(self, valName, v):
"""
Return the value that should actually be set in the solver when value 'v' is given by
the user.
"""
if isinstance(v, CompDepFunc):
if self._lhs is None:
raise NotImplementedError()
v = v(self._lhs)
if valName == 'K':
return self._rateMult * v
else:
return v
@nutils.FreezeAfterInit
class _SubReactionList(nutils.SolverPathObject, nutils.ParameterizedObject, list):
"""
Represents a unidirectional reaction consisting in a list of single steps reaction objects.
"""
def __init__(self, parent, reacName, LRP, compEvs, isFwd=True, *args, **kwargs):
super().__init__(*args, **kwargs)
self._parent = parent
self._reacName = reacName
self._LRP = LRP
self._compEvs = compEvs
self._isFwd = isFwd
# Add the reaction name as parameter if the reaction is not anonymous
if not isinstance(self._parent, AnonymousReaction):
self._setParameter('Reaction ID', self._parent.name + ("['fwd']" if self._isFwd else "['bkw']"))
# Add the complex selector as parameter, if any
_, compSels, *_ = self._parent._usedObjects
if len(compSels) > 0:
self._setParameter('Complex', nutils._ParamList(compSels))
self.K = 0
self._reacElems = {side: {loc: set() for loc in ALL_LOCATIONS} for side in 'LR'}
def append(self, item):
if not isinstance(item, _SubReaction):
raise TypeError('Only _SubReaction objects can be added to _SubReactionList.')
super().append(item)
def _getRateUnits(self):
"""Return rate units of the reaction"""
units = set()
# some reactions are volume-volume regardless of their content
isVolVolReac = self._parent._isVesicleSurfReac()
for l, r, rm in self._LRP:
l = l._getFlatList()
molExp = 1 - len(l)
unit = ''
if molExp != 0:
if isVolVolReac or len(l) == 0 or any(s.loc != Location.SURF for s in l):
unit = f'M'
else:
unit = f'(mol m^-2)'
unit += f'^{molExp} ' if molExp != 1 else ' '
units.add(unit + 's^-1')
if len(units) > 1:
raise Exception(f'Inconsistent rate units in reaction {self}: {units}')
return nutils.Units(list(units)[0])
@property
@nutils.ParameterizedObject.RegisterGetter(units=_getRateUnits)
def K(self):
"""Reaction rate(s)
If the reaction is reversible, this value corresponds to a tuple composed of the forward
and backward rate. The rate can be a special function that depends e.g. on membrane potential
(:py:class:`VDepRate`) or on complex states (:py:class:`CompDepRate`).
Usage examples::
# Bidirectional reactions
SA + SB <r[1]> 2 * SC
r[1].K = forward_rate, backward_rate
# Unidirectional reactions
SC >r[2]> SD
r[1].K = single_rate
The reaction rate(s) default(s) to 0.
:type: Union[float, :py:class:`steps.API_2.utils.Parameter`, :py:class:`XDepFunc`]
"""
pass
@K.setter
@nutils.ParameterizedObject.RegisterSetter(units=_getRateUnits)
def K(self, rate):
if not isinstance(rate, (numbers.Number, XDepFunc)):
raise TypeError(f'{rate} cannot be used as a rate.')
@property
@nutils.ParameterizedObject.RegisterGetter(units=None)
def Dependencies(self):
"""Reaction dependencies
A set of species that need to be present in order for the reaction to happen but their amount
will not have any effect on the reaction rate.
If the reaction is reversible, this value corresponds to a tuple composed of the forward
and backward dependencies. The dependencies can be a special function that depends e.g. on
complex states (:py:class:`CompDepFunc`).
Note that for now, only species on the surfaces of rafts and vesicles can be used as dependencies.
The dependencies should be given in the shape of a reaction side, for example: `SA.v + SC.v`
Usage examples::
# Bidirectional reactions
SA.v + SB.o <r[1]> SA.v + SB.v
r[1].K = fwd_rate, bkw_rate
# SC is a dependency for the forward reaction, SD for the backward:
r[1].Dependencies = SC.v, SD.v
# Unidirectional reactions
SC.v >r[2]> SD.v
r[1].K = single_rate
r[1].Dependencies = SA.v
The dependencies are empty by default.
:type: Union[:py:class:`steps.API_2.ReactionElement`, :py:class:`steps.API_2.ReactionSide`,
:py:class:`CompDepFunc`]
"""
# Return a default value if the property was not set
return None
@Dependencies.setter
@nutils.ParameterizedObject.RegisterSetter(units=None)
def Dependencies(self, deps):
if not isinstance(deps, (ReactionElement, ReactionSide, CompDepFunc, type(None))):
raise TypeError(f'{deps} cannot be used as a dependency.')
@property
@nutils.ParameterizedObject.RegisterGetter(units=None)
def AntiDependencies(self):
"""Reaction anti-dependencies
A set of species that, if present, will prevent the reaction from happenning.
If the reaction is reversible, this value corresponds to a tuple composed of the forward
and backward anti-dependencies. The anti-dependencies can be a special function that depends e.g. on
complex states (:py:class:`CompDepFunc`).
Note that for now, only species on the surfaces of rafts can be used as anti-dependencies.
The anti-dependencies should be given in the shape of a reaction side, for example: `SA.r + SC.r`
In that example, the reaction will be prevented if either `SA` or `SC` are present in the raft.
If e.g. `2*SA.R + SC.r` is used, the reaction will be prevented is at least two `SA` are present
or at least one `SC`.
Usage examples::
# Bidirectional reactions
SA.r + SB.s <r[1]> SA.r + SB.r
r[1].K = fwd_rate, bkw_rate
# SC will prevent the forward reaction, SD the backward:
r[1].AntiDependencies = SC.r, SD.r
# Unidirectional reactions
SC.r >r[2]> SD.r
r[1].K = single_rate
r[1].AntiDependencies = SA.r
The anti-dependencies are empty by default.
:type: Union[:py:class:`steps.API_2.ReactionElement`, :py:class:`steps.API_2.ReactionSide`,
:py:class:`CompDepFunc`]
"""
# Return a default value if the property was not set
return None
@AntiDependencies.setter
@nutils.ParameterizedObject.RegisterSetter(units=None)
def AntiDependencies(self, deps):
if not isinstance(deps, (ReactionElement, ReactionSide, CompDepFunc, type(None))):
raise TypeError(f'{deps} cannot be used as an anti-dependency.')
@property
@nutils.ParameterizedObject.RegisterGetter(units=None)
def Immobilization(self):
"""Reaction immobilization flag
The immobilization flag is only applicable to reactions involving rafts or vesicles, it can
either be :py:data:`IMMOBILIZING`, :py:data:`NO_EFFECT` or :py:data:`MOBILIZING`.
If the reaction is reversible, this value corresponds to a tuple composed of the forward
and backward immobilization flags. Immobilization flags can be a special function that depends e.g. on
complex states (:py:class:`CompDepFunc`).
Usage examples::
# Bidirectional reactions
SA.v + SB.o <r[1]> SA.v + SB.v
r[1].K = fwd_rate, bkw_rate
# Forward reaction immobilizes, backward reaction does not:
r[1].Immobilization = IMMOBILIZING, NO_EFFECT
# Unidirectional reactions
SC.v >r[2]> SD.v
r[1].K = single_rate
r[1].Immobilization = IMMOBILIZING
The immobilization flag defaults to :py:data:`NO_EFFECT`.
:type: Union[:py:data:`steps.API_2.model.IMMOBILIZING`, :py:data:`steps.API_2.model.NO_EFFECT`,
:py:data:`steps.API_2.model.MOBILIZING`, :py:class:`CompDepFunc`]
"""
# Return a default value if the property was not set
return NO_EFFECT
@Immobilization.setter
@nutils.ParameterizedObject.RegisterSetter(units=None)
def Immobilization(self, val):
if not isinstance(val, Immobilization):
raise TypeError(f'Expected an immobilization flag, got {val} instead.')
@property
@nutils.ParameterizedObject.RegisterGetter(units=nutils.Units('m'))
def MaxDistance(self):
"""Maximum distance (in m) between vesicle and patch for the reaction to occur
Note that only reactions involving vesicles can have a maximum distance parameter.
If the reaction is reversible, this value corresponds to a tuple composed of the forward
and backward maximum distances. Maximum distances can be a special function that depends e.g. on
complex states (:py:class:`CompDepFunc`).
Usage examples::
# Bidirectional reactions
SA.v + SB.s <r[1]> SA.v + SB.v
r[1].K = fwd_rate, bkw_rate
# Forward reaction has a maximum distance, backward reaction does not:
r[1].MaxDistance = 0.1e-6, None
# Unidirectional reactions
SC.v + SD.s >r[2]> 2 * SD.v
r[1].K = single_rate
r[1].MaxDistance = 5e-7
The maximum distance defaults to None, i.e. no maximum distance.
:type: Union[float, :py:class:`CompDepFunc`]
"""
# Return a default value if the property was not set
return None
@MaxDistance.setter
@nutils.ParameterizedObject.RegisterSetter(units=nutils.Units('m'))
def MaxDistance(self, md):
if md is not None and (isinstance(md, bool) or not isinstance(md, (numbers.Number, CompDepFunc))):
raise TypeError(f'{md} cannot be used as maximum distance.')
def __repr__(self):
if self._isFwd:
return f'{self._parent.lhs} → {self._parent.rhs}'
else:
return f'{self._parent.rhs} → {self._parent.lhs}'
@classmethod
def _getDisplayName(cls):
"""Name that will be used as column title during parameter export
Defaults to class name.
"""
return Reaction._getDisplayName()
def _getStepsObjects(self):
"""Return a list of the steps objects that this named object holds."""
res = []
for sr in self:
res += sr._getStepsObjects()
return res
def _simPathWalkExpand(self):
"""Return an iterable of the elements that should be part of ResultSelector paths."""
return self
def _getAllElems(self, loc, sides='LR'):
elems = set()
for side in sides:
elems |= self._reacElems[side][loc]
return elems
def _decl(self):
"""
Return a string that describes the declaration of the object. It returns its name by default
but it can contain e.g. the file and line at which it was declared.
"""
return self._parent._decl()
def _addToSteps(self):
"""
Declare and add all steps reactions corresponding to the (lhs, rhs) pairs in LRP.
"""
mdl, globCompSels, volSys, surfSys = self._parent._usedObjects
for i, reac in enumerate(self._LRP):
l, r, rm = reac
# Add reac elems
for loc in ALL_LOCATIONS:
self._reacElems['L'][loc] |= l._getAllElems(loc)
self._reacElems['R'][loc] |= r._getAllElems(loc)
# Try to declare complex reac, if applicable
for j, compEvents in enumerate(self._compEvs):
# Add complex reac elems
for ce in compEvents:
for loc in ce.getLocations():
for side in ce.affectedSides():
self._reacElems[side][loc].add(ce.comp)
# Declare reaction
name = self._reacName + (f'_{i}' if len(self._LRP) > 1 else '')
if len(compEvents) > 0:
name_c = name + (f'_{j}' if len(self._compEvs) > 1 else '')
res = self._declareComplexReac(
name_c, l, r, compEvents, mdl, volSys, surfSys, rm * self.K, rm,
self.Dependencies, self.AntiDependencies, self.Immobilization, self.MaxDistance
)
else:
res = self._declareSimpleReac(
name, l, r, mdl, volSys, surfSys, rm * self.K, rm, self.Dependencies,
self.AntiDependencies, self.Immobilization, self.MaxDistance
)
# Add the subreactions to self
self.append(res)
def _declareComplexReac(self, name, simpLHS, simpRHS, compEvents, mdl, volSys, surfSys, rate, rateMult,
deps, antideps, immobilization, max_dist):
"""Declare complex reaction involving real complexes."""
if deps is not None:
raise Exception(f'Cannot use dependencies with complex reactions.')
nutils._print(f'\tAdding STEPS complex reaction {name}: {compEvents}, rate = {rate}', 3)
if self._parent._isSurfaceReac():
ilhs = simpLHS._GetStepsElems(Location.IN)
slhs = simpLHS._GetStepsElems(Location.SURF)
olhs = simpLHS._GetStepsElems(Location.OUT)
irhs = simpRHS._GetStepsElems(Location.IN)
srhs = simpRHS._GetStepsElems(Location.SURF)
orhs = simpRHS._GetStepsElems(Location.OUT)
if isinstance(rate, VDepRate):
raise NotImplementedError('Cannot use voltage dependent rates with STEPS complex reactions')
icompEvs = [ce._stepsObj for ce in compEvents if ce.loc == Location.IN]
scompEvs = [ce._stepsObj for ce in compEvents if ce.loc == Location.SURF]
ocompEvs = [ce._stepsObj for ce in compEvents if ce.loc == Location.OUT]
stepsReac = stepslib._py_ComplexSReac(
name,
surfSys.stepsSys,
ilhs=ilhs,
slhs=slhs,
olhs=olhs,
irhs=irhs,
srhs=srhs,
orhs=orhs,
icompEvs=icompEvs,
scompEvs=scompEvs,
ocompEvs=ocompEvs,
kcst=rate,
)
else:
lhs = simpLHS._GetStepsElems()
rhs = simpRHS._GetStepsElems()
if isinstance(rate, CompDepRate):
raise NotImplementedError(
f'{self._decl()}: Cannot use CompDepRate with STEPS complex '
f'reactions.'
)
compEvs = [ev._stepsObj for ev in compEvents]
stepsReac = stepslib._py_ComplexReac(name, volSys.stepsSys, lhs=lhs, rhs=rhs, compEvs=compEvs, kcst=rate)
return _SubReaction(self, stepsReac, rateMult=rateMult, name=name)
def _declareSimpleReac(self, name, lhs, rhs, mdl, volSys, surfSys, rate, rateMult, deps, antideps,
immobilization, max_dist):
"""Declare and return the steps reactions from the final sides lhs and rhs."""
# Compute rates and other parameters in case of complex dependencies
if isinstance(rate, CompDepRate):
rate = rate(lhs)
if isinstance(deps, CompDepFunc):
deps = deps(lhs)
if isinstance(immobilization, CompDepFunc):
immobilization = immobilization(lhs)
if isinstance(max_dist, CompDepFunc):
max_dist = max_dist(lhs)
nutils._print(f'\tAdding STEPS reaction {name}: {lhs} -> {rhs}, rate = {rate}, rateMult = {rateMult}', 3)
if self._parent._isSurfaceReac():
# Surface system reaction
ilhs = lhs._GetStepsElems(Location.IN)
slhs = lhs._GetStepsElems(Location.SURF)
olhs = lhs._GetStepsElems(Location.OUT)
irhs = rhs._GetStepsElems(Location.IN)
srhs = rhs._GetStepsElems(Location.SURF)
orhs = rhs._GetStepsElems(Location.OUT)
if isinstance(surfSys, VesicleSurfaceSystem):
# Vesicle surface reaction
if len(ilhs) > 0:
raise Exception(
f'{self._decl()}: Vesicle reactions cannot involve species inside the vesicle on'
f' the LHS.'
)
# On vesicle surface
_vlhs = lhs._GetStepsElems(Location.VESSURF)
_vrhs = rhs._GetStepsElems(Location.VESSURF)
vlhs = [se for se in _vlhs if isinstance(se, stepslib._py_Spec)]
llhs = [se for se in _vlhs if isinstance(se, stepslib._py_LinkSpec)]
vrhs = [se for se in _vrhs if isinstance(se, stepslib._py_Spec)]
lrhs = [se for se in _vrhs if isinstance(se, stepslib._py_LinkSpec)]
# Check that no forbidden locations are present
for elem in lhs:
if elem.loc not in [Location.IN, Location.OUT, Location.SURF, Location.VESSURF]:
raise Exception(
f'{self._decl()}: Element {elem} cannot appear in a vesicle surface reaction.'
)
# Dependencies
vdeps = []
if deps is not None:
deps = deps._toReactionSide()
if any(elem.loc != Location.VESSURF for elem in deps):
raise Exception(
f'{self._decl()}: Dependencies of reactions declared in a '
f'VesicleSurfaceSystem need to be explicitely declared on the surface of the '
f'vesicle (e.g. "spec.v"). The dependencies given do not satisfy this '
f'constraint: {deps}.'
)
vdeps = deps._GetStepsElems(Location.VESSURF)
# Immobilization
if immobilization not in [IMMOBILIZING, NO_EFFECT, MOBILIZING]:
raise ValueError(
f'{self._decl()}: Unsupported immobilization flag. IMMOBILIZING, NO_EFFECT'
f', or MOBILIZING.')
# Max distance
if max_dist is None:
max_dist = -1
elif max_dist < 0:
raise ValueError(
f'{self._decl()}: Maximum distance parameter should be a positive number, got '
f'{max_dist} instead.'
)
stepsReac = stepslib._py_VesSReac(
name,
surfSys._getStepsObjects()[0],
olhs=olhs,
slhs=slhs,
vlhs=vlhs,
llhs=llhs,
lrhs=lrhs,
vrhs=vrhs,
srhs=srhs,
orhs=orhs,
irhs=irhs,
vdeps=vdeps,
kcst=rate,
immobilization=immobilization,
max_dist=max_dist,
)
elif isinstance(surfSys, RaftSurfaceSystem):
# Raft surface reaction
# On vesicle surface
rlhs = lhs._GetStepsElems(Location.RAFTSURF)
rrhs = rhs._GetStepsElems(Location.RAFTSURF)
# Check that no forbidden locations are present
for elem in lhs:
if elem.loc not in [Location.IN, Location.OUT, Location.SURF, Location.RAFTSURF]:
raise Exception(
f'{self._decl()}: Element {elem} cannot appear in a raft surface reaction.'
)
# Dependencies
rdeps = []
if deps is not None:
deps = deps._toReactionSide()
if any(elem.loc != Location.RAFTSURF for elem in deps):
raise Exception(
f'{self._decl()}: Dependencies of reactions declared in a '
f'RaftSurfaceSystem need to be explicitely declared on the surface of the '
f'raft (e.g. "spec.r"). The dependencies given do not satisfy this '
f'constraint: {deps}.'
)
rdeps = deps._GetStepsElems(Location.RAFTSURF)
# Anti Dependencies
anti_rdeps = []
if antideps is not None:
antideps = antideps._toReactionSide()
if any(elem.loc != Location.RAFTSURF for elem in antideps):
raise Exception(
f'{self._decl()}: Anti-dependencies of reactions declared in a '
f'RaftSurfaceSystem need to be explicitely declared on the surface of the '
f'raft (e.g. "spec.r"). The anti-dependencies given do not satisfy this '
f'constraint: {antideps}.'
)
anti_rdeps = antideps._GetStepsElems(Location.RAFTSURF)
# Immobilization
if immobilization not in [IMMOBILIZING, NO_EFFECT, MOBILIZING]:
raise ValueError(
f'{self._decl()}: Unsupported immobilization flag. IMMOBILIZING, NO_EFFECT'
f', or MOBILIZING.')
stepsReac = stepslib._py_RaftSReac(
name,
surfSys._getStepsObjects()[0],
ilhs=ilhs,
olhs=olhs,
slhs=slhs,
rlhs=rlhs,
rrhs=rrhs,
srhs=srhs,
orhs=orhs,
irhs=irhs,
rdeps=rdeps,
anti_rdeps=anti_rdeps,
kcst=rate,
immobilization=immobilization,
)
else:
# Normal surface reaction
if isinstance(rate, VDepRate):
# Voltage dependent surface reaction
stepsReac = stepslib._py_VDepSReac(
name,
surfSys.stepsSys,
ilhs=ilhs,
slhs=slhs,
olhs=olhs,
irhs=irhs,
srhs=srhs,
orhs=orhs,
**rate._getVDepSReacParams(self._getRateUnits())
)
else:
stepsReac = stepslib._py_SReac(
name,
surfSys.stepsSys,
ilhs=ilhs,
slhs=slhs,
olhs=olhs,
irhs=irhs,
srhs=srhs,
orhs=orhs,
kcst=rate,
)
else:
# Volume system reaction
elhs = lhs._GetStepsElems()
erhs = rhs._GetStepsElems()
if isinstance(rate, VDepRate):
# Undeclare the parent reaction to avoid further exception raising when exiting context
# managers
# self._parent._declared = False
raise NotImplementedError(
f'{self._decl()}: Cannot declare a voltage-dependent volume '
f'reaction.'
)
else:
stepsReac = stepslib._py_Reac(name, volSys.stepsSys, lhs=elhs, rhs=erhs, kcst=rate)
# Check that there were no silently unused Parameters
if not isinstance(stepsReac, (stepslib._py_VesSReac, stepslib._py_RaftSReac)):
if deps is not None:
raise Exception(
f'{self._decl()}: Dependencies can only be used with vesicle and raft surface reactions.'
)
if immobilization != NO_EFFECT:
raise Exception(
f'{self._decl()}: Immobilization flags can only be used with vesicle and raft surface '
f'reactions.'
)
if not isinstance(stepsReac, stepslib._py_VesSReac):
if max_dist is not None:
raise Exception(
f'{self._decl()}: Maximum distances can only be used with vesicle surface reactions.'
)
return _SubReaction(self, stepsReac, lhs=lhs, rateMult=rateMult, name=name)
class _ComplexEvent:
"""Utility class for storing complex events"""
def __init__(self, comp, loc):
self.comp = comp
self.loc = loc
self._stepsObj = self._createStepsObj()
def getLocations(self):
return [self.loc]
class _ComplexCreateEvent(_ComplexEvent):
def __init__(self, comp, loc, init):
self.init = init
super().__init__(comp, loc)
def affectedSides(self):
return 'R'
def _createStepsObj(self):
return stepslib._py_ComplexCreateEvent(self.comp.name, self.init)
def __repr__(self):
return f'DeleteEv(comp={self.comp}, loc={self.loc}, init={self.init})'
class _ComplexLHSEvent(_ComplexEvent):
def __init__(self, comp, loc, filt):
self.filt = filt
super().__init__(comp, loc)
class _ComplexDeleteEvent(_ComplexLHSEvent):
def affectedSides(self):
return 'L'
def _createStepsObj(self):
return stepslib._py_ComplexDeleteEvent(self.comp.name, self.filt)
def __repr__(self):
return f'DeleteEv(comp={self.comp}, loc={self.loc}, filt={self.filt})'
class _ComplexUpdateEvent(_ComplexLHSEvent):
def __init__(self, comp, loc, filt, reac, upd, destLoc=None):
self.reac = reac
self.upd = upd
self.destLoc = destLoc
super().__init__(comp, loc, filt)
def affectedSides(self):
return 'LR'
def getLocations(self):
return [self.loc, self.destLoc]
def _createStepsObj(self):
return stepslib._py_ComplexUpdateEvent(self.comp.name, self.filt, self.reac, self.upd, self.destLoc)
def __repr__(self):
return f'UpdateEv(comp={self.comp}, loc={self.loc}, filt={self.filt}, reac={self.reac}, '\
f'upd={self.upd}, dest={self.destLoc})'
[docs]@nutils.FreezeAfterInit
class Reaction(
nutils.UsingObjects((Model, ComplexSelector), (VolumeSystem, SurfaceSystem)),
nutils.StepsWrapperObject,
nutils.ParameterizedObject,
):
"""A (possibly reversible) reaction.
A reaction needs to be declared within a :py:class:`Model` and within a
:py:class:`VolumeSystem` or a :py:class:`SurfaceSystem` (i.e. inside a ``with`` block,
see e.g. :py:class:`VolumeSystem`)
A :py:class:`Reaction` is defined by a left hand side (lhs) and a right hand side (rhs).
Both sides are instances of :py:class:`ReactionSide` and the overall reaction is declared
with a custom syntax designed to look like standard chemical reactions notation.
When :py:class:`Complex` are involved, a :py:class:`Reaction` object can represent several
elementary reactions, see user guide for details.
.. note::
Reactions should always be instantiated from a :py:class:`ReactionManager`.
Usage::
r = ReactionManager()
with vsys:
# Simple forward reaction
SA >r[1]> SB
r[1].K = 1 # Set forward rate
# Bidirectional reaction
SA + SB <r[1]> SC + SD
r[1].K = 2, 1 # Set forward and backward rates
# Chained reactions
SA <r[1]> SB >r[2]> SC
r[1].K = 2, 1
r[2].K = 3
.. warning::
The corresponding reactions are only truly added to steps after the `with vsys:` block
is exited.
"""
_FwdSpecifier = 'fwd'
_BkwSpecifier = 'bkw'
def __init__(self, _managerInfo=None, *args, **kwargs):
super().__init__(*args, **kwargs)
self.lhs = None
self.rhs = None
self._subReactions = collections.OrderedDict()
self._subReactions[Reaction._FwdSpecifier] = None
self._subReactions[Reaction._BkwSpecifier] = None
self._declared = False
self._added = False
self._bidir = None
self._declarationInfo = None
self._managerInfo = _managerInfo
self._usedObjects = list(self._getUsedObjects())
@classmethod
def _FromStepsObject(cls, obj, mdl):
"""Create the interface object from a STEPS object."""
reac = cls(name=obj.getID())
# Reaction sides
lhs, rhs = ReactionSide([]), ReactionSide([])
if isinstance(obj, stepslib._py_Reac):
for lhsSpec in obj.getLHS():
lhs += getattr(mdl, lhsSpec.getID())
for rhsSpec in obj.getRHS():
rhs += getattr(mdl, rhsSpec.getID())
elif isinstance(obj, (
stepslib._py_SReac, stepslib._py_VDepSReac,
stepslib._py_VesSReac, stepslib._py_RaftSReac)):
if not isinstance(obj, stepslib._py_VesSReac):
for lhsSpec in obj.getILHS():
lhs += getattr(mdl, lhsSpec.getID()).i
for lhsSpec in obj.getSLHS():
lhs += getattr(mdl, lhsSpec.getID()).s
for lhsSpec in obj.getOLHS():
lhs += getattr(mdl, lhsSpec.getID()).o
for rhsSpec in obj.getIRHS():
rhs += getattr(mdl, rhsSpec.getID()).i
for rhsSpec in obj.getSRHS():
rhs += getattr(mdl, rhsSpec.getID()).s
for rhsSpec in obj.getORHS():
rhs += getattr(mdl, rhsSpec.getID()).o
if isinstance(obj, stepslib._py_VesSReac):
for lhsSpec in obj.getVLHS() + obj.getLLHS():
lhs += getattr(mdl, lhsSpec.getID()).v
for rhsSpec in obj.getVRHS() + obj.getLRHS():
rhs += getattr(mdl, rhsSpec.getID()).v
if isinstance(obj, stepslib._py_RaftSReac):
for lhsSpec in obj.getRLHS():
lhs += getattr(mdl, lhsSpec.getID()).r
for rhsSpec in obj.getRRHS():
rhs += getattr(mdl, rhsSpec.getID()).r
else:
raise NotImplementedError(f'Cannot import from STEPS reaction {obj}')
lhs > reac > rhs
subreacLst = _SubReactionList(reac, obj.getID(), [(lhs, rhs, 1)], None)
subreacLst.append(_SubReaction(subreacLst, obj, lhs=lhs, rateMult=1.0, name=obj.getID()))
reac._subReactions[Reaction._FwdSpecifier] = subreacLst
# Rate
if not isinstance(obj, stepslib._py_VDepSReac):
# TODO not urgent: set K for VDepSReacs
reac.K = obj.getKcst()
# Immobilization
if isinstance(obj, (stepslib._py_VesSReac, stepslib._py_RaftSReac)):
reac.Immobilization = obj.getImmobilization()
# MaxDistance
if isinstance(obj, stepslib._py_VesSReac) and obj.getMaxDistance() != -1.0:
reac.MaxDistance = obj.getMaxDistance()
# Dependencies
if isinstance(obj, (stepslib._py_VesSReac, stepslib._py_RaftSReac)):
deps = ReactionSide([])
if isinstance(obj, stepslib._py_VesSReac):
for depSpec in obj.getVDeps():
deps += getattr(mdl, depSpec.getID()).v
else:
for depSpec in obj.getRDeps():
deps += getattr(mdl, depSpec.getID()).r
reac.Dependencies = deps
# Anti Dependencies
if isinstance(obj, (stepslib._py_RaftSReac, )):
antideps = ReactionSide([])
for adepSpec in obj.getAntiRDeps():
deps += getattr(mdl, adepSpec.getID()).r
reac.AntiDependencies = deps
reac._added = True
for loc in ALL_LOCATIONS:
subreacLst._reacElems['L'][loc] |= lhs._getAllElems(loc)
subreacLst._reacElems['R'][loc] |= rhs._getAllElems(loc)
return reac
@staticmethod
def _findMatchingComplexInd(csl, rhs):
"""Return the index of the complexSelector corresponding to csl in rhs."""
if csl.loc is None:
return rhs.index(csl) if csl in rhs else None
else:
# Handle the case in which complex selectors have different location
csrInds = [i for i, csr in enumerate(rhs) if (csl._elem, csl.stoich) == (csr._elem, csr.stoich)]
return csrInds[0] if len(csrInds) > 0 else None
def _expandSimpleReactions(self, lhs, rhs):
"""
Expand ComplexSelectors whose complex has '_statesAsSpecies = True' to a list of
lhs and rhs pairs.
"""
cre = _ComplexReactionElement
simpleCompLhs = [cs for cs in lhs._getElemsOfType(cre) if cs._elem._areStatesAsSpecies()]
simpleCompRhs = [cs for cs in rhs._getElemsOfType(cre) if cs._elem._areStatesAsSpecies()]
specLhs = lhs._getElemsOfType(Species, LinkSpecies, ComplexState)
specRhs = rhs._getElemsOfType(
Species, LinkSpecies, ComplexState
)
# TODO Later release: Using complexState in reaction when using C++ complex reacs: is it working?
# Initialize the list of reaction pairs with species
lhsRhsPairs = [(ReactionSide(specLhs), ReactionSide(specRhs), 1)]
# Group complexReactionElements by complex
comp2csls = {}
for csl in simpleCompLhs:
cmplx = csl._elem._getReferenceObject()
comp2csls.setdefault(cmplx, []).append(csl)
# Compute all reaction pairs for lhs ComplexSelectors
for comp, csls in comp2csls.items():
if len(csls) > 2:
raise Exception(
f'Cannot use more than 2 instances of a complex in a reaction, '
f'the current reaction uses {comp} {len(csls)} times.'
)
allStatePairs = []
for csl in csls:
csrInd = Reaction._findMatchingComplexInd(csl, simpleCompRhs)
if csrInd is not None:
# Compute state pairs if there is a matching ComplexSelector in the rhs
csr = simpleCompRhs[csrInd]
statePairs = []
for nl, nr, rm in csl._elem._getMatchingStates(csr._elem):
statePairs.append(
(
ReactingElement(nl, csl.stoich, csl.loc),
ReactingElement(nr, csr.stoich, csr.loc),
rm,
)
)
del simpleCompRhs[csrInd]
else:
if isinstance(csl, _ComplexReactants):
raise Exception(
f'ComplexReactants {csl} has no matching ComplexReactants on '
f'the RHS of the reaction.'
)
# Otherwise, just add all the states of csl
allStates = csl._elem._getAllStates()
statePairs = [
(ReactingElement(nl, csl.stoich, csl.loc), ReactionSide([]), 1) for nl in allStates
]
allStatePairs.append(statePairs)
# If 2 instances of the same complex are used, update the rate multipliers accordingly
newASP = []
if len(csls) == 2:
sameReacs = (
all(isinstance(csl._elem, _ComplexReactants) for csl in csls)
and csls[0]._elem._hasSameReactants(csls[1]._elem)
)
allL1 = {l1: rm1 for l1, r1, rm1 in allStatePairs[0]}
allL2 = {l2: rm2 for l2, r2, rm2 in allStatePairs[1]}
for l1, r1, rm1 in allStatePairs[0]:
for l2, r2, rm2 in allStatePairs[1]:
if sameReacs and l1 in allL2 and l2 in allL1:
newASP.append((l1 + l2, r1 + r2, rm1 * rm2 / 2))
else:
newASP.append((l1 + l2, r1 + r2, rm1 * rm2))
else:
newASP = allStatePairs[0]
# Update lhsRhsPairs with the matching states
newLRP = []
for l, r, rm in lhsRhsPairs:
for nl, nr, nrm in newASP:
newLRP.append((l + nl, r + nr, rm * nrm))
lhsRhsPairs = newLRP
# Add reaction pairs from remaining rhs elements if they are completely defined
# Otherwise, raise an exception.
if len(simpleCompRhs) > 0:
newLRP = []
for l, r, rm in lhsRhsPairs:
for csr in simpleCompRhs:
states = csr._elem._getAllStates()
if len(states) > 1:
raise Exception(
f'Complex selector {csr._elem} is used in the right hand side of a '
f'reaction but is not matching anything in the left hand side and is '
f'not fully defined. The reaction is ambiguous.'
)
newLRP.append((l, r + ReactingElement(states.pop(), csr.stoich, csr.loc), rm))
lhsRhsPairs = newLRP
return lhsRhsPairs
def _expandComplexReactions(self, lhs, rhs):
"""
Expand ComplexSelectors whose complex has '_statesAsSpecies = False' to a list of
complex events.
"""
cre = _ComplexReactionElement
clhs = [cs for cs in lhs._getElemsOfType(cre) if not cs._elem._areStatesAsSpecies()]
crhs = [cs for cs in rhs._getElemsOfType(cre) if not cs._elem._areStatesAsSpecies()]
# Initialize the list of complex reaction events
compEvents = [[]]
for cl in clhs:
comp = cl._elem._complex
crInd = Reaction._findMatchingComplexInd(cl, crhs)
if crInd is not None:
cr = crhs[crInd]._elem
events = [
_ComplexUpdateEvent(comp, cl.loc, filt, reac, upd, destLoc=crhs[crInd].loc)
for filt, reac, upd in cl._elem._getUpdateVects(cr)
]
del crhs[crInd]
else:
if isinstance(cl, _ComplexReactants):
raise Exception(
f'ComplexReactants {cl} has no matching ComplexReactants on '
f'the RHS of the reaction.'
)
events = [_ComplexDeleteEvent(comp, cl.loc, cl._elem._toUnorderedFilter())]
# Update compEvents with the new events
newCE = []
for evs in compEvents:
for e in events:
newCE.append(evs + [e])
compEvents = newCE
# Add creation events
for cr in crhs:
state = cr._elem._getFullyDeterminedState()
if state is None:
raise Exception(
f'Complex selector {cr._elem} is used on the RHS of the reaction and '
f'should lead to a new complex creation but is not fully defined. '
f'The reaction is ambiguous.'
)
for evs in compEvents:
evs.append(_ComplexCreateEvent(state._comp, cr.loc, state._toUnorderedVect()))
return compEvents
def _parseReac(self):
"""
Parse the reaction, expand complex reactions and setup data structures for declaring the
reactions in steps.
"""
mdl, globCompSels, volSys, surfSys = self._usedObjects
# Check that no SubUnitStates are wrongly used
if (self.lhs._isEmpty() and self.rhs._hasElemsOfType(SubUnitSelector)) or (
self.rhs._isEmpty() and self.lhs._hasElemsOfType(SubUnitSelector)
):
raise Exception('Cannot mix SubUnitState or SubUnitSelectors with empty reaction side.')
# Check that the correct spatial systems are being used
if self._isSurfaceReac() and surfSys is None:
raise Exception(
f'Reaction {self.name} is a surface reaction but is being declared '
f'outside of a surface system.'
)
if not self._isSurfaceReac() and volSys is None:
raise Exception(
f'Reaction {self.name} is a volume reaction but is being declared '
f'outside of a volume system.'
)
if not self._isSurfaceReac():
# Check that if reacElems have a location, it's the same for all reacElems and both sides
if self.lhs._hasLocation() or self.rhs._hasLocation():
if self.lhs._hasSingleLocation() and self.rhs._hasSingleLocation():
if self.lhs._hasVolumeLocation() and self.rhs._hasVolumeLocation():
# Add a constraint to the model that will be checked later,
# to make sure that the volume system matches the location of the surface system
loc = self.lhs._getLocation() if self.lhs._hasLocation() else self.rhs._getLocation()
mdl._addVolSysConstraint(self, volSys, surfSys, loc)
else:
raise Exception(
f'Reaction {self.name} is a volume reaction but has a side with '
f'surface located component.'
)
else:
raise Exception(
f'Reaction {self.name} is a volume reaction but has component in different locations.'
)
# Transform SubUnitSelectors to ComplexSelectors or ComplexReactants
lhs = self.lhs._getFinal(globCompSels)
rhs = self.rhs._getFinal(globCompSels)
# Transform ComplexSelectors whose complex has "statesAsSpecies = True" to species
fwdLRP = self._expandSimpleReactions(lhs, rhs)
fwdCompEvs = self._expandComplexReactions(lhs, rhs)
nameFwd = self.name + Reaction._FwdSpecifier if self._bidir else self.name
self._subReactions[Reaction._FwdSpecifier] = _SubReactionList(self, nameFwd, fwdLRP, fwdCompEvs, True)
if self._bidir:
bkwLRP = self._expandSimpleReactions(rhs, lhs)
bkwCompEvs = self._expandComplexReactions(rhs, lhs)
nameBkw = self.name + Reaction._BkwSpecifier
self._subReactions[Reaction._BkwSpecifier] = _SubReactionList(self, nameBkw, bkwLRP, bkwCompEvs, False)
# Retrieve filename and line number for easier debugging
for i, fi in enumerate(inspect.stack()):
if fi.function in ['__lt__', '__gt__']:
i += 1
break
frameInfo = inspect.stack()[i]
self._declarationInfo = (frameInfo.filename, frameInfo.lineno)
self._declared = True
def _comparison(self, other, gt):
if self._declared:
if self._managerInfo is None:
raise Exception('This reaction has already been declared.')
# Return the new reaction
# TODO not urgent: Manager info being a tuple hurts readability, make it a namedtuple
p = self._managerInfo[0]._popNewReacArgs()
if p is None:
self._managerInfo[0]._setNewReacArgs(other, gt)
# python exhibits special behavior for chained comparison operator
# see https://docs.python.org/3/reference/expressions.html#comparisons
# The first comparison needs to return True
return True
else:
r = self._managerInfo[0]._getNew(self._managerInfo[1])
r._comparison(*p)
return r._comparison(other, gt)
if isinstance(other, ReactionElement):
other = other._toReactionSide()
elif other is None:
other = ReactionSide([])
elif not isinstance(other, ReactionSide):
raise Exception(f'{other} cannot be one of the sides of reaction {self}.')
if self.lhs is None:
self.lhs = other._makeCopy()
self._bidir = gt
return self
else:
if not gt:
raise Exception('Used "<" instead of ">" after the reaction.')
self.rhs = other._makeCopy()
self._parseReac()
return other
def __gt__(self, other):
return self._comparison(other, True)
def __lt__(self, other):
return self._comparison(other, False)
def _getAllElems(self, loc, sides='LR'):
if self._subReactions[Reaction._FwdSpecifier] is not None:
res = self._subReactions[Reaction._FwdSpecifier]._getAllElems(loc, sides)
if self._subReactions[Reaction._BkwSpecifier] is not None:
res |= self._subReactions[Reaction._BkwSpecifier]._getAllElems(loc, sides)
return res
return set()
def _isSurfaceReac(self):
if self.lhs is None or self.rhs is None:
return None
lhsS, rhsS = self.lhs._isSurfaceSide(), self.rhs._isSurfaceSide()
return (lhsS and rhsS) or (None in [lhsS, rhsS] and True in [lhsS, rhsS])
def _isVesicleSurfReac(self):
return isinstance(self._usedObjects[3], VesicleSurfaceSystem)
def _isRaftSurfReac(self):
return isinstance(self._usedObjects[3], RaftSurfaceSystem)
[docs] def __getitem__(self, key):
"""Get sub-reactions from a reaction
Square bracket notation can be used on a reaction to retrieve either the forward
reaction(s), or the backward one(s). This should never be needed during model declaration
but it is very important for simulation control and data saving (see
:py:class:`steps.API_2.sim.SimPath`).
:param key: Use ``'{{model.Reaction._FwdSpecifier}}'`` to access the forward
reaction(s) or ``'{{model.Reaction._BkwSpecifier}}'`` to access the backward
reaction(s).
:type key: str
:returns: The querried reaction(s).
See :ref:`/API_2/STEPS_Tutorial_wm.ipynb` for usage example.
:meta public:
"""
if key not in self._subReactions:
raise KeyError(f'No reaction can be accessed with key {key}.')
return self._subReactions[key]
def __repr__(self):
return '{} {}-> {}'.format(self.lhs, '-' if not self._bidir else '<', self.rhs)
def _decl(self):
"""
Return a string that describes the declaration of the object. It returns its name by default
but it can contain e.g. the file and line at which it was declared.
"""
return f'{self.name} ({self._declarationInfo[0]}: {self._declarationInfo[1]})'
def _getStepsObjects(self):
"""Return a list of the steps objects that this named object holds."""
res = []
for specif, srl in self._subReactions.items():
if srl is not None:
res += srl._getStepsObjects()
return res
def _simPathWalkExpand(self):
"""Return an iterable of the elements that should be part of ResultSelector paths."""
return itertools.chain(
*[srl._simPathWalkExpand() for specif, srl in self._subReactions.items() if srl is not None]
)
def _getSubParameterizedObjects(self):
"""Return all subobjects that can hold parameters."""
return [srl for specif, srl in self._subReactions.items() if srl is not None]
def _includeinParamTables(self):
"""Whether the object should be included in parameter tables"""
# Only directed reactions should be included in parameter tables
return False
def _exitCallback(self, parent):
"""
Method to be called the first time we get out of a context manager in which the object as been
declared.
"""
# Only declare the reaction if we are exiting the volume / surface system
if isinstance(parent, SpaceSystem):
if self._declared and not self._added:
nutils._print(f'Adding pysteps reaction {self.name}: {self}', 2)
self._added = True
self._subReactions[Reaction._FwdSpecifier]._addToSteps()
if self._bidir:
self._subReactions[Reaction._BkwSpecifier]._addToSteps()
@staticmethod
def _getProperty(prop):
"""Build properties for the Reaction class based on the property prop."""
def getter(self):
if not self._declared:
raise Exception(f'Cannot get properties of a reaction that was not yet declared.')
if self._bidir:
return (
prop.fget(self._subReactions[Reaction._FwdSpecifier]),
prop.fget(self._subReactions[Reaction._BkwSpecifier]),
)
else:
return prop.fget(self._subReactions[Reaction._FwdSpecifier])
def setter(self, val):
if not self._declared:
raise Exception('Cannot set properties of a reaction that was not yet declared.')
if self._added:
raise Exception('Cannot set properties of a reaction that was already added to STEPS.')
if self._bidir:
if not (isinstance(val, tuple) and len(val) == 2):
raise Exception(
f'The reaction is reversible, two properties should be set, got {val} instead.'
)
prop.fset(self._subReactions[Reaction._FwdSpecifier], val[0])
prop.fset(self._subReactions[Reaction._BkwSpecifier], val[1])
else:
if not isinstance(val, tuple):
val = (val, )
if len(val) > 1:
raise Exception(f'The reaction is unidirectional but more than one property was set')
prop.fset(self._subReactions[Reaction._FwdSpecifier], val[0])
return property(getter, setter)
# Map _SubReactionList properties to Reaction
for propName in dir(_SubReactionList):
prop = getattr(_SubReactionList, propName)
if isinstance(prop, property):
newProp = Reaction._getProperty(prop)
newProp.__doc__ = prop.__doc__
setattr(Reaction, propName, newProp)
@nutils.FreezeAfterInit
class AnonymousReaction(Reaction):
"""
:meta private:
"""
def __init__(self, *args, **kwargs):
super().__init__(*args, **kwargs)
[docs]@nutils.FreezeAfterInit
class ReactionManager:
"""Class used to instantiate reactions
A :py:class:`ReactionManager` object allows to easily declare reactions by associating them
to an identifier using the square bracket syntax::
r = ReactionManager() # Only one ReactionManager is needed
S1 + S2 >r['R01']> S3 # String identifier for permanently naming the reaction
S1 + S2 >r[1]> S3 # Integer identifier for reactions that do not need explicit
# naming
String identifiers are used when the reaction will need to be accessed during simulation (see
:py:class:`steps.API_2.sim.SimPath` for usage), they cannot be reused for different reactions.
Integer identifiers are temporary identifiers and can be reused for declaring several un-named
reactions::
S1 + S2 >r[1]> S3 # First use of temporary identifier 1. r[1] now refers to this
r[1].K = ... # reaction until a new reaction is defined with r[1].
S1 + S3 <r[1]> S4 # r[1] now refers to this new reaction.
r[1].K = ...
In general, only one :py:class:`ReactionManager` should be instantiated and given a short
variable name (``r`` in all examples) for convenience. In some specific cases however, one might
need several independent reaction managers so this class is not implemented as a singleton.
"""
def __init__(self, *args, **kwargs):
super().__init__(*args, **kwargs)
self.reactions = {}
self._newReacArgs = None
[docs] def __getitem__(self, i):
"""Square brackets syntax for instantiating or accessing a reaction
:param i: Identifier of the reaction. Integer for a temporary id and string for a
permanent one.
:type i: Union[int, str]
:returns: The reaction corresponding to the given identifier
:rtype: :py:class:`Reaction`
:meta public:
"""
if i not in self.reactions:
if isinstance(i, str):
self.reactions[i] = Reaction(name=i, _managerInfo=(self, i))
else:
self.reactions[i] = AnonymousReaction(_managerInfo=(self, i))
return self.reactions[i]
def _getNew(self, i):
del self.reactions[i]
return self.__getitem__(i)
def _setNewReacArgs(self, *args):
self._newReacArgs = args
def _popNewReacArgs(self):
p = self._newReacArgs
self._newReacArgs = None
return p
###################################################################################################
# State-dependent rates
[docs]@nutils.FreezeAfterInit
class XDepFunc(nutils.NamedObject, nutils.ParameterizedObject):
r"""Base class for functions depending on some simulation or object state
This is a wrapper around a function that adds some functionalities.
:param func: lambda function or standard function. There is not support for variable named
arguments (\*\*kwargs).
:type func: Callable[[Any, ...], Number]
"""
_DEPENDS_ON_NAME = 'Depends on'
def __init__(self, func, *args, **kwargs):
super().__init__(*args, **kwargs)
if not hasattr(func, '__call__'):
raise TypeError(f'Expected a callable, got {func} instead.')
self._func = func
self._args = []
self._mult = 1
self._setParameter('Function', self._func)
self._subParamObjs = set()
def _getSubParameterizedObjects(self):
"""Return all subobjects that can hold parameters."""
return sorted(self._subParamObjs, key=lambda x: x.name)
[docs] def __mul__(self, m):
"""Multiplication with a constant or another function of the same class
:param m: Constant number or another :py:class:`XDepFunc` of the same type
:type m: Union[Number]
:returns: A new function, which corresponds to the multiplication of both operands
:rtype: :py:class:`XDepFunc`
:meta public:
"""
if isinstance(m, numbers.Number):
ret = copy.copy(self)
ret._mult *= m
return ret
else:
raise TypeError(f'Cannot multiply {self} with {m}.')
[docs] def __rmul__(self, m):
"""See :py:func:`__mul__`
:meta public:
"""
return self.__mul__(m)
# TODO Not urgent: add more operations
@property
def func(self):
"""Access the function itself
:type: Callable, read-only
"""
return self._func
[docs] def __call__(self, *lst, _recParams=True):
"""Call the wrapped function
:meta public:
"""
if _recParams:
nutils.Parameter._startRecording()
ret = self._mult * self._func(*lst)
if _recParams:
# Retrieve the parameters used when the function was called
allParams = [param for param in nutils.Parameter._endRecordingAndGetUsage() if param._isNamed() and param._isUserDefined()]
if len(allParams) > 0:
depon = self._getParameter(XDepFunc._DEPENDS_ON_NAME)
names = set([param.name for param in allParams])
if depon is not None:
names |= set(depon.value.lst)
allParams = list(set(depon._dependencies) | set(allParams))
nameList = nutils._ParamList(list(names))
self._setParameter(XDepFunc._DEPENDS_ON_NAME, nutils.Parameter(nameList, _dependencies=allParams))
if isinstance(ret, XDepFunc):
self._subParamObjs.add(ret)
return ret
[docs]@nutils.FreezeAfterInit
class VDepFunc(XDepFunc):
"""Hold a function and a voltage range for computing e.g. voltage-dependent reaction rates
For efficiency reasons, the values of these functions are pre-computed for specific voltage
values (by default between -150mV and 100mV with a 0.1mV step) and these values might be
changed by providing the optional ``vrange`` parameter.
:param func: lambda function or standard function that takes a membrane potential (in Volts)
as parameter and returns a single numerical value.
:type func: Callable[[float], Number]
:param vrange: Voltage range at which the function should be computed. Should be a 3-tuple
in the form ``(min_voltage, max_voltage, step_voltage)``, all expressed in Volts.
:type vrange: Tuple[float, float, float]
See :py:class:`VDepRate` for an example.
"""
_DEFAULT_VMIN = -150e-3
_DEFAULT_VMAX = 100e-3
_DEFAULT_DV = 1e-4
def __init__(self, func, vrange=None, *args, **kwargs):
super().__init__(func=func, *args, **kwargs)
if vrange is None:
vrange = [VDepFunc._DEFAULT_VMIN, VDepFunc._DEFAULT_VMAX, VDepFunc._DEFAULT_DV]
if vrange[0] >= vrange[1]:
raise Exception(f'The maximum voltage of the range is lower than the minimum.')
if vrange[2] > vrange[1] - vrange[0]:
raise Exception(f'The voltage step is greater than the voltage range.')
self.vrange = vrange
self._args = [self.vrange]
@property
@nutils.ParameterizedObject.RegisterGetter(units=nutils.Units('V'))
def vrange(self):
"""The associated voltage range
:type: Tuple[float, float, float]
"""
pass
@vrange.setter
@nutils.ParameterizedObject.RegisterSetter(units=nutils.Units('V'))
def vrange(self, val):
pass
[docs]@nutils.FreezeAfterInit
class CompDepFunc(XDepFunc):
"""Hold a function that depends on the fully specified state of a complex
:param func: lambda function or standard function that takes one or several fully specified
state(s) of complex(es) as parameter and returns a single numerical value.
:type func: Callable[[:py:class:`ComplexState`, ...], Number]
:param compsOrCompSels: A list of complexes or complex selectors in the same order as they
should be provided to the wrapped function.
:type compsOrCompSels: List[Union[:py:class:`Complex`, :py:class:`ComplexSelector`]]
See :py:class:`CompDepRate` for an example.
"""
def __init__(self, func, compsOrCompSels, *args, **kwargs):
super().__init__(func=func, *args, **kwargs)
if not hasattr(compsOrCompSels, '__iter__'):
raise TypeError(f'Expected an iterable containing complexes, got {compsOrCompSels} instead.')
self._complexes = []
for cs in compsOrCompSels:
if isinstance(cs, Complex):
self._complexes.append(cs)
elif isinstance(cs, ComplexSelector):
self._complexes.append(cs._complex)
elif isinstance(cs, ComplexState):
self._complexes.append(cs._comp)
else:
raise TypeError(f'Expected a complex or a complex selector, got {cs} instead.')
self._compsOrCompSels = compsOrCompSels
self._args = [self._compsOrCompSels]
self._setParameter('Complex', self._compsOrCompSels)
def _getStatesFromReacSide(self, lhs):
lhsCompStates = list(lhs._getElemsOfType(ComplexState))
lhsCompSels = list(lhs._getElemsOfType(ComplexSelector))
states = []
for cs in self._compsOrCompSels:
candidates = []
if isinstance(cs, ComplexSelector):
for re in lhsCompStates:
if (
isinstance(re._elem._parentCompSel, _ComplexReactants)
and re._elem._parentCompSel._filter == cs
) or re._elem._parentCompSel == cs:
candidates.append(re._elem)
for re in lhsCompSels:
if re._elem == cs:
for state in re._elem._getAllStates():
candidates.append(state)
elif isinstance(cs, Complex):
for re in lhsCompStates:
if re._elem._comp is cs:
candidates.append(re._elem)
for re in lhsCompSels:
if re._elem._complex is cs:
for state in re._elem._getAllStates():
candidates.append(state)
else:
raise TypeError(f'Expected a Complex or a ComplexSelector, got {cs} instead.')
if len(candidates) == 1:
states.append(candidates[0])
elif len(candidates) > 1:
raise Exception(
f'More than one complex state in {lhs} corresponded to {cs} for '
f'determining a rate: {candidates}'
)
else:
raise Exception(f'No complex states in {lhs} matched {cs}.')
return states
def __call__(self, *lst):
if all(isinstance(e, ComplexState) for e in lst):
return super().__call__(*lst)
elif len(lst) == 1 and isinstance(lst[0], ReactionSide):
return super().__call__(*self._getStatesFromReacSide(lst[0]))
else:
raise NotImplementedError()
[docs]@nutils.FreezeAfterInit
class VDepRate(VDepFunc):
"""Voltage dependent reaction rate
See base class :py:class:`VDepFunc` for how to create.
Usage::
with CaPchan[...]:
# VDepRate declaration
# alpha_cap and beta_cap are normal functions taking the membrane potential in Volts
# as parameter and returning the corresponding rate.
# with custom voltage range:
a_cap = VDepRate(alpha_cap, vrange = [-100.0e-3, 50e-3, 1e-4])
# with default voltage range
b_cap = VDepRate(beta_cap)
m0.s <r[1]> m1.s <r[2]> m2.s <r[3]> m3.s
r[1].K = 3 * a_cap, 1 * b_cap
r[2].K = 2 * a_cap, 2 * b_cap
r[3].K = 1 * a_cap, 3 * b_cap
# The transitions from m0 to m3 have rates that are multiples of a_cap and b_cap
# it is thus easier to first define these two VDepRate functions and reuse them
# when declaring reactions.
"""
def _getVDepSReacParams(self, units):
"""Return the parameters that are actually needed by the cython constructor"""
vmin, vmax, dv = self.vrange
tablesize = int((vmax - vmin) / dv) + 1
# Retrieve the parameters used when the function is called with vmin
# This can possibly miss some parameters if they are only involved in the computation
# for some specific values of the membrane potential
self(vmin)
ktab = []
for i in range(tablesize):
ret = self(vmin + i * dv, _recParams=False)
ktab.append(ret.valueIn(units) if isinstance(ret, nutils.Parameter) else ret)
return dict(ktab=ktab, vmin=vmin, vmax=vmax, dv=dv, tablesize=tablesize)
[docs]@nutils.FreezeAfterInit
class CompDepRate(CompDepFunc):
"""Complex state dependent reaction rate
See base class :py:class:`CompDepFunc` for how to create.
Usage::
with Comp[...]:
S1A >r[1]> S1B
rate = CompDepRate(lambda s: k1 + k2 * s.Count(S1B), [Comp])
r[1].K = rate
# The above reaction specifies that the subunit state S1A changes to S1B with
# a basal rate of k1 and this rate is increased by k2 for each subunit in the complex
# being in state S1B.
with Comp1[...] as C1, Comp2[...] as C2:
S1A[C1] + S3A[C2] >r[1]> S1B[C1] + S3B[C2]
rate = CompDepRate(lambda s1, s2: k3 + k4 * s1.Count(S2A) * s2.Count(S4A), [C1, C2])
r[1].K = rate
# The above reaction takes place between a subunit of Comp1 in state S1A and a
# subunit of Comp2 in state S3A. It takes place with a basal rate k3 and this rate
# is increased by a value that depends on the number of subunits in state S2A in
# Comp1 and the number of subunits in state S4A in Comp2.
"""
pass
[docs]@nutils.FreezeAfterInit
class CompDepCond(CompDepFunc):
"""Complex state dependent conductance
See base class :py:class:`CompDepFunc` for how to create.
"""
pass
[docs]@nutils.FreezeAfterInit
class CompDepP(CompDepFunc):
"""Complex state dependent permeability
See base class :py:class:`CompDepFunc` for how to create.
"""
pass
[docs]@nutils.FreezeAfterInit
class CompDepDcst(CompDepFunc):
"""Complex state dependent diffusion constant
See base class :py:class:`CompDepFunc` for how to create.
"""
pass
###################################################################################################
# Diffusion
class _SubDiffusion(nutils.NamedObject):
"""Wrapper for a single steps diffusion object"""
_DiffStr = 'Diff'
_SDiffStr = 'SDiff'
_DiffDirecStr = 'direction_tet'
_SDiffDirecStr = 'direction_tri'
def __init__(self, stepsDiff, parent, keyElem=None, *args, **kwargs):
super().__init__(*args, **kwargs)
self._stepsDiff = stepsDiff
self._elemStr = _SubDiffusion._SDiffStr if parent._isSurfaceDiff() else _SubDiffusion._DiffStr
self._keyElem = keyElem
self._direcStr = (
_SubDiffusion._SDiffDirecStr if parent._isSurfaceDiff() else _SubDiffusion._DiffDirecStr
)
self._direc = None
def _getStepsObjects(self):
"""Return a list of the steps objects that this named object holds."""
return [self._stepsDiff]
def _solverStr(self):
"""Return the string that is used as part of method names for this specific object."""
return self._elemStr
def _solverKeywordParams(self):
"""Return the additional keyword parameters that should be passed to the solver"""
if self._direc is not None:
return {self._direcStr: self._direc._solverId()[0]}
else:
return {}
def __call__(self, direc=None):
"""Return a specialized version of the subdiffusion"""
sd = copy.copy(self)
sd._direc = direc
return sd
def _solverSetValue(self, valName, v):
"""
Return the value that should actually be set in the solver when value 'v' is given by
the user.
"""
if valName == 'D' and isinstance(v, CompDepDcst):
return v(self._keyElem)
else:
return v
class _SubDiffusionList(nutils.SolverPathObject):
"""List of _SubDiffusion objects"""
def __init__(self, lst, *args, **kwargs):
super().__init__(*args, **kwargs)
self._diffLst = lst
def __iter__(self):
return iter(self._diffLst)
def _simPathWalkExpand(self):
"""Return an iterable of the elements that should be part of ResultSelector paths."""
return self.__iter__()
def __call__(self, direc=None):
return _SubDiffusionList([d.__call__(direc) for d in self._diffLst])
[docs]@nutils.FreezeAfterInit
class Diffusion(
nutils.UsingObjects(Model, (VolumeSystem, SurfaceSystem)),
nutils.StepsWrapperObject, nutils.ParameterizedObject
):
"""A diffusion rule for a chemical species in a volume or on a surface
:param elem: The element that subject to diffusion
:type elem: Union[:py:class:`Species`, :py:class:`ComplexState`, :py:class:`Complex`,
:py:class:`ComplexSelector`]
:param Dcst: The diffusion constant in S.I. units
:type Dcst: Union[float, :py:class:`CompDepDcst`]
Usage::
with vsys:
# Volume diffusion of Species S1
# The diffusion rule is named 'diffS1' for later access during simulation
diffS1 = Diffusion.Create(S1, 0.2e-9)
with ssys:
# Surface diffusion of Species S2
# The diffusion rule is not named and will thus not be accessible during simulation
Diffusion(S2, 0.08e-12)
"""
def __init__(self, elem, Dcst, *args, _createObj=True, **kwargs):
super().__init__(*args, **kwargs)
mdl, volSys, surfSys = self._getUsedObjects()
if volSys is not None and surfSys is not None:
raise Exception(
'A diffusion rule should only apply to a volume or a surface system, but not both.'
)
if not isinstance(elem, nutils.NamedObject):
raise TypeError(f'Expected a named ligand, got {elem} instead.')
self.sys = volSys if surfSys is None else surfSys
self._stepsDiffs = {}
self._elem = elem._getReferenceObject()
if _createObj:
if isinstance(self.sys, VesicleSurfaceSystem):
cls = stepslib._py_VesSDiff
elif self.sys.__class__ in [VolumeSystem, SurfaceSystem]:
cls = stepslib._py_Diff
else:
raise Exception(f'Cannot declare a diffusion rule in {self.sys}')
if (
isinstance(elem, (Complex, ComplexSelector))
) and elem._areStatesAsSpecies():
for state in elem:
name = f'{self.name}_{state.name}'
self._stepsDiffs[state] = _SubDiffusion(
cls(name, self.sys._getStepsObjects()[0], state._getStepsObjects()[0]),
self,
keyElem=state,
name=name,
)
elif isinstance(elem, (Species, ComplexState)):
self._stepsDiffs[elem] = _SubDiffusion(
cls(self.name, self.sys._getStepsObjects()[0], elem._getStepsObjects()[0]),
self,
keyElem=elem,
name=self.name,
)
else:
raise TypeError(f'Cannot declare a diffusion rule for {elem}.')
# Setting diffusion constant
self.Dcst = Dcst
def _getStepsObjects(self):
"""Return a list of the steps objects that this named object holds."""
return [sd._getStepsObjects()[0] for _, sd in self._stepsDiffs.items()]
@classmethod
def _FromStepsObject(cls, obj, mdl):
"""Create the interface object from a STEPS object."""
lig = getattr(mdl, obj.getLig().getID())
diff = cls(lig, obj.getDcst(), _createObj=False, name=obj.getID())
diff._stepsDiffs[lig] = _SubDiffusion(obj, diff, keyElem=lig, name=diff.name)
return diff
def _isSurfaceDiff(self):
return isinstance(self.sys, SurfaceSystem)
[docs] def __getitem__(self, key):
"""Get sub diffusion rules from a diffusion object
Square bracket notation can be used on a diffusion rule when it was defined for a complex
or a complex selector. This should never be needed during model declaration but it is very
important for simulation control and data saving (see :py:class:`steps.API_2.sim.SimPath`).
:param key: A sequence of one or more :py:class:`SubUnitState` or
:py:class:`SubUnitSelector` that corresponds to one or several :py:class:`ComplexState`.
:type key: Tuple[Union[:py:class:`SubUnitState`, :py:class:`SubUnitSelector`], ...]
:returns: The diffusion rule associated to ``key``
:meta public:
"""
if isinstance(self._elem, Complex):
e = self._elem.__getitem__(key)
if isinstance(e, ComplexState):
return self._stepsDiffs[e]
elif isinstance(e, ComplexSelector):
return _SubDiffusionList([self._stepsDiffs[s] for s in e])
raise NotImplementedError()
[docs] def __call__(self, direc=None):
"""Get a version of the diffusion rule with added information
Parentheses notation (function call) can be used on a diffusion rule to represent a
diffusion rule further specified with additional information. This should never be needed
during model declaration but can become handy for simulation control and data saving
(see :py:class:`steps.API_2.sim.SimPath`).
:param direc: A direction for the diffusion rule
:type key: :py:class:`Reference`
:returns: A version of the diffusion rule that is specific to direction ``direc``.
:meta public:
"""
return _SubDiffusionList([sd.__call__(direc=direc) for _, sd in self._stepsDiffs.items()])
def _getAllElems(self, loc):
if loc is None:
return [self._elem]
else:
if self._isSurfaceDiff():
return [self._elem] if loc == Location.SURF else []
else:
return []
def _simPathWalkExpand(self):
"""Return an iterable of the elements that should be part of ResultSelector paths."""
for _, sd in self._stepsDiffs.items():
yield sd
@property
@nutils.ParameterizedObject.RegisterGetter(units=nutils.Units('m^2 s^-1'))
def Dcst(self):
"""Diffusion constant (in m^2 s^-1)
:type: Union[float, :py:class:`steps.API_2.nutils.Parameter`]
"""
pass
@Dcst.setter
@nutils.ParameterizedObject.RegisterSetter(units=nutils.Units('m^2 s^-1'))
def Dcst(self, val):
if isinstance(val, numbers.Number):
for _, sd in self._stepsDiffs.items():
sd._stepsDiff.setDcst(val)
elif isinstance(val, CompDepFunc):
for cs, sd in self._stepsDiffs.items():
sd._stepsDiff.setDcst(val(cs))
else:
raise TypeError(f'{val} is not a valid diffusion constant.')
###################################################################################################
# Vesicle and Raft creation / deletion
[docs]@nutils.FreezeAfterInit
class Endocytosis(
nutils.UsingObjects(Model, SurfaceSystem),
nutils.StepsWrapperObject, nutils.ParameterizedObject
):
"""A specific type of interaction that models endocytosis of a vesicle
An endocytosis event models the process of vesicle endocytosis by creating
a new vesicle within a compartment. A vesicle will be created at a given
rate when an optional species signature is met within an endocytic zone of
a patch (see :py:class:`steps.API_2.geom.EndocyticZone`).
:param vesicle: Vesicle that should appear
:type vesicle: :py:class:`Vesicle`
:param kcst: Rate constant of the reaction in s^-1
:type kcst: float
:param deps: Species dependencies
:type deps: Union[None, :py:class:`ReactionElement`, :py:class:`ReactionSide`]
:param inside: If true, the vesicle should appear in the inside compartment, otherwise it should appear
in the outside compartment.
:type inside: bool
"""
def __init__(self, vesicle, kcst=0, deps=None, inside=True, _createObj=True, **kwargs):
super().__init__(**kwargs)
mdl, surfSys = self._getUsedObjects()
if isinstance(deps, ReactionElement):
deps = deps._toReactionSide()
elif deps is None:
deps = ReactionSide([])
elif not isinstance(deps, ReactionSide):
raise TypeError(f'Expected species dependency, got {deps} instead.')
if any(elem.loc not in [None, Location.SURF] for elem in deps):
raise Exception(f'Species dependencies cannot include species in volumes.')
if not isinstance(vesicle, Vesicle):
raise TypeError(f'Expected a vesicle, got {vesicle} instead.')
self._vesicle = vesicle
self._setParameter('Dependencies', deps)
self._setParameter('K', kcst, nutils.Units('s^-1'))
if _createObj:
deps = self.Dependencies._GetStepsElems()
if inside:
self.stepsEndo = stepslib._py_Endocytosis(
self.name, surfSys._getStepsObjects()[0], vesicle._getStepsObjects()[0], None,
deps, self.K
)
else:
self.stepsEndo = stepslib._py_Endocytosis(
self.name, surfSys._getStepsObjects()[0], None, vesicle._getStepsObjects()[0],
deps, self.K
)
else:
self.stepsEndo = None
def _getStepsObjects(self):
"""Return a list of the steps objects that this named object holds."""
return [self.stepsEndo]
def _getAllElems(self, loc=None):
"""
Return all elements that can directly be accessed by name from the parent physical
location.
"""
return self.Dependencies._getAllElems(loc)
def _simPathAutoMetaData(self):
"""Return a dictionary with string keys and string or numbers values."""
return {'obj_type': self._solverStr(), 'obj_id': self.name}
@classmethod
def _FromStepsObject(cls, obj, mdl):
"""Create the interface object from a STEPS object."""
ves = getattr(mdl, obj.getIRHS().getID())
deps = ReactionSide([])
for spec in obj.getSpecDeps():
deps += getattr(mdl, spec.getID())
endo = cls(ves, obj.getKcst(), deps=deps, inside=obj.getInner(), _createObj=False, name=obj.getID())
endo.stepsEndo = obj
return endo
@property
def Vesicle(self):
"""Vesicle associated with the endocytosis reaction
:type: :py:class:`Vesicle`, read-only
"""
return self._vesicle
@property
@nutils.ParameterizedObject.RegisterGetter(units=None)
def Dependencies(self):
"""Species dependencies for the endocytosis reaction
:type: :py:class:`ReactionSide`, read-only
"""
pass
@property
@nutils.ParameterizedObject.RegisterGetter(units=nutils.Units('s^-1'))
def K(self):
"""Endocytosis reaction constant (in s^-1)
:type: Union[float, :py:class:`steps.API_2.nutils.Parameter`]
"""
pass
@K.setter
@nutils.ParameterizedObject.RegisterSetter(units=nutils.Units('s^-1'))
def K(self, val):
if isinstance(val, numbers.Number):
self.stepsEndo.setKcst(val)
else:
raise TypeError(f'{val} is not a valid reaction constant.')
[docs]@nutils.FreezeAfterInit
class Exocytosis(
nutils.UsingObjects(Model, VesicleSurfaceSystem),
nutils.StepsWrapperObject, nutils.ParameterizedObject
):
"""A specific type of interaction that models exocytosis of a vesicle
An exocytosis event models the process of vesicle exocytosis. By default the vesicle
is destroyed, species in the vesicle surface are deposited in the patch at
the location at which exocytosis occurs, and species inside the vesicle lumen
are deposited in the opposite compartment.
:param kcst: Rate constant of the exocytosis reaction in s^-1
:type kcst: float
:param deps: Species dependencies on the vesicle surface
:type deps: Union[None, :py:class:`ReactionElement`, :py:class:`ReactionSide`]
:param raft: If this argument is given the exocytosis event will create a
raft on the patch at the location that exocytosis occurs. Any species on
the vesicle surface will be deposited in the raft instead of in the patch.
:type raft: :py:class:`Raft`
:param kissAndRun: This parameter can be set to True to model a kiss-and-run
exocytosis event. This does not result in collapse of the vesicle,
instead the vesicle is maintained in position after the release of
vesicle lumen contents into the opposite compartment. The vesicle
surface species are maintained on the vesicle.
:type kissAndRun: bool
:param knrSpecChanges: To aid kiss-and-run modeling, this argument can be a
list of tuples of two :py:class:`Species`; when the kiss-and-run event
happens, all vesicle surface species from the first element of the tuple
get changed to species of the second element. This feature can be used to
specify any species changes that take place on the vesicle surface upon a
kiss-and-run exocytosis event. This can be useful for modeling maturation
of a complex for which undocking may be dependent.
:type knrSpeChanges: List[Tuple[:py:class:`Species`, :py:class:`Species`]]
:param knrPartialRelease: To aid kiss-and-run modeling, this argument can
be used to specify partial release of any inner species that are released
upon a kiss-and-run exocytosis event. 1.0 means all species are released,
0.0 means none of the species are released. If not specified, it defaults
to 1.0.
:type knrPartialRelease: float
"""
def __init__(self, kcst=0, deps=None, raft=None, kissAndRun=False, knrSpecChanges=None,
knrPartialRelease=1, _createObj=True, **kwargs):
super().__init__(**kwargs)
mdl, surfSys = self._getUsedObjects()
if isinstance(deps, ReactionElement):
deps = deps._toReactionSide()
elif deps is None:
deps = ReactionSide([])
elif not isinstance(deps, ReactionSide):
raise TypeError(f'Expected species dependency, got {deps} instead.')
if any(elem.loc not in [None, Location.VESSURF] for elem in deps):
raise Exception(f'Species dependencies cannot include species in volumes.')
if raft is not None and not isinstance(raft, Raft):
raise TypeError(f'Expected a raft, got {raft} instad.')
self._raft = raft
self._setParameter('Dependencies', deps)
self._setParameter('KissAndRun', kissAndRun)
self._setParameter('KnRSpecChanges', knrSpecChanges)
self._setParameter('KnRPartialRelease', knrSpecChanges)
self._setParameter('K', kcst, nutils.Units('s^-1'))
if _createObj:
if self.KnRSpecChanges is not None:
knrSpecChanges = {
spec1._getStepsObjects()[0]: spec2._getStepsObjects()[0]
for spec1, spec2 in self.KnRSpecChanges
}
else:
knrSpecChanges = {}
deps = self.Dependencies._GetStepsElems()
self.stepsExo = stepslib._py_Exocytosis(
self.name, surfSys._getStepsObjects()[0], deps,
raft._getStepsObjects()[0] if raft is not None else None,
self.K, self.KissAndRun, knrSpecChanges, knrPartialRelease
)
else:
self.stepsExo = None
def _getStepsObjects(self):
"""Return a list of the steps objects that this named object holds."""
return [self.stepsExo]
def _simPathAutoMetaData(self):
"""Return a dictionary with string keys and string or numbers values."""
return {'obj_type': self._solverStr(), 'obj_id': self.name}
def _getAllElems(self, loc=None):
"""
Return all elements that can directly be accessed by name from the parent physical
location.
"""
return self.Dependencies._getAllElems(loc)
@classmethod
def _FromStepsObject(cls, obj, mdl):
"""Create the interface object from a STEPS object."""
if obj.getRaft() is not None:
raft = getattr(mdl, obj.getRaft().getID())
else:
raft = None
deps = ReactionSide([])
for spec in obj.getSpecDeps():
deps += getattr(mdl, spec.getID())
knrSpecChanges = list(obj.getKissAndRunSpecChanges().items())
exo = cls(
obj.getKcst(), deps=deps, raft=raft, kissAndRun=obj.getKissAndRun(),
knrSpecChanges=knrSpecChanges, knrPartialRelease=obj.getKissAndRunPartRelease(),
_createObj=False, name=obj.getID()
)
exo.stepsExo = obj
return exo
@property
def Raft(self):
"""Raft optionally associated with the exocytosis reaction
:type: Union[:py:class:`Raft`, None], read-only
"""
return self._raft
@property
@nutils.ParameterizedObject.RegisterGetter(units=None)
def Dependencies(self):
"""Species dependencies for the exocytosis reaction
:type: :py:class:`ReactionSide`, read-only
"""
pass
@property
@nutils.ParameterizedObject.RegisterGetter(units=None)
def KissAndRun(self):
"""Whether the exocytosis is a kiss-and-run exocytosis
:type: bool, read-only
"""
pass
@property
@nutils.ParameterizedObject.RegisterGetter(units=None)
def KnRSpecChanges(self):
"""Species changes upon kiss-and-run ecoxytosis reaction
:type: List[Tuple[:py:class:`Species`, :py:class:`Species`]], read-only
"""
pass
@property
@nutils.ParameterizedObject.RegisterGetter(units=None)
def KnRPartialRelease(self):
"""Fraction of inner species that are released upon a kiss-and-run exocytosis event
:type: float, read-only
"""
pass
@property
@nutils.ParameterizedObject.RegisterGetter(units=nutils.Units('s^-1'))
def K(self):
"""Exocytosys reaction constant (in s^-1)
:type: Union[float, :py:class:`steps.API_2.nutils.Parameter`]
"""
pass
@K.setter
@nutils.ParameterizedObject.RegisterSetter(units=nutils.Units('s^-1'))
def K(self, val):
if isinstance(val, numbers.Number):
self.stepsExo.setKcst(val)
else:
raise TypeError(f'{val} is not a valid reaction constant.')
[docs]@nutils.FreezeAfterInit
class RaftEndocytosis(
nutils.UsingObjects(Model, RaftSurfaceSystem),
nutils.StepsWrapperObject, nutils.ParameterizedObject
):
"""A specific type of interaction that models endocytosis of a raft to form a vesicle
A raft endocytosis event models the process of vesicle endocytosis by
creating a new vesicle within a compartment from a raft. A vesicle will be
created at a given rate when an optional species signature is met within
the raft.
:param vesicle: Vesicle that should appear
:type vesicle: :py:class:`Vesicle`
:param kcst: Rate constant of the reaction
:type kcst: float
:param deps: Species dependencies for endocytosis on the raft
:type deps: Union[None, :py:class:`ReactionElement`, :py:class:`ReactionSide`]
:param inside: If true, the vesicle should appear in the inside compartment, otherwise it should appear
in the outside compartment.
:type inside: bool
"""
def __init__(self, vesicle, kcst=0, deps=None, inside=True, _createObj=True, **kwargs):
super().__init__(**kwargs)
mdl, surfSys = self._getUsedObjects()
if isinstance(deps, ReactionElement):
deps = deps._toReactionSide()
elif deps is None:
deps = ReactionSide([])
elif not isinstance(deps, ReactionSide):
raise TypeError(f'Expected species dependency, got {deps} instead.')
if any(elem.loc not in [None, Location.RAFTSURF] for elem in deps):
raise Exception(f'Species dependencies cannot include species in volumes.')
if not isinstance(vesicle, Vesicle):
raise TypeError(f'Expected a vesicle, got {vesicle} instad.')
self._vesicle = vesicle
self._setParameter('Dependencies', deps)
self._setParameter('K', kcst, nutils.Units('s^-1'))
if _createObj:
deps = self.Dependencies._GetStepsElems()
if inside:
self.stepsRaftEndo = stepslib._py_RaftEndocytosis(
self.name, surfSys._getStepsObjects()[0], vesicle._getStepsObjects()[0], None,
deps, self.K
)
else:
self.stepsRaftEndo = stepslib._py_RaftEndocytosis(
self.name, surfSys._getStepsObjects()[0], None, vesicle._getStepsObjects()[0],
deps, self.K
)
else:
self.stepsRaftEndo = None
def _getAllElems(self, loc=None):
"""
Return all elements that can directly be accessed by name from the parent physical
location.
"""
return self.Dependencies._getAllElems(loc)
def _getStepsObjects(self):
"""Return a list of the steps objects that this named object holds."""
return [self.stepsRaftEndo]
def _simPathAutoMetaData(self):
"""Return a dictionary with string keys and string or numbers values."""
return {'obj_type': self._solverStr(), 'obj_id': self.name}
@classmethod
def _FromStepsObject(cls, obj, mdl):
"""Create the interface object from a STEPS object."""
ves = getattr(mdl, obj.getRHS().getID())
deps = ReactionSide([])
for spec in obj.getSpecDeps():
deps += getattr(mdl, spec.getID())
endo = cls(ves, obj.getKcst(), deps=deps, inside=obj.getInner(), _createObj=False, name=obj.getID())
endo.stepsRaftEndo = obj
return endo
@property
def Vesicle(self):
"""Vesicle associated with the raft endocytosis reaction
:type: :py:class:`Vesicle`, read-only
"""
return self._vesicle
@property
@nutils.ParameterizedObject.RegisterGetter(units=None)
def Dependencies(self):
"""Species dependencies for the raft endocytosis reaction
:type: :py:class:`ReactionSide`, read-only
"""
pass
@property
@nutils.ParameterizedObject.RegisterGetter(units=nutils.Units('s^-1'))
def K(self):
"""RaftEndocytosys reaction constant (in s^-1)
:type: Union[float, :py:class:`steps.API_2.nutils.Parameter`]
"""
pass
@K.setter
@nutils.ParameterizedObject.RegisterSetter(units=nutils.Units('s^-1'))
def K(self, val):
if isinstance(val, numbers.Number):
self.stepsRaftEndo.setKcst(val)
else:
raise TypeError(f'{val} is not a valid reaction constant.')
[docs]class RaftGen(
nutils.UsingObjects(SurfaceSystem),
nutils.StepsWrapperObject, nutils.ParameterizedObject
):
"""A specific type of interaction that models a raft genesis event
Generate a raft at a given rate when a defined species signature is met
within a patch triangle.
:param raft: Raft that should be created
:type raft: :py:class:`Raft`
:param kcst: Rate constant
:type kcst: float
:param deps: Species dependencies on the patch surface
:type deps: Union[None, :py:class:`ReactionElement`, :py:class:`ReactionSide`]
"""
def __init__(self, raft, kcst=0, deps=None, _createObj=True, **kwargs):
super().__init__(**kwargs)
surfSys, = self._getUsedObjects()
if isinstance(deps, ReactionElement):
deps = deps._toReactionSide()
elif deps is None:
deps = ReactionSide([])
elif not isinstance(deps, ReactionSide):
raise TypeError(f'Expected species dependency, got {deps} instead.')
if any(elem.loc not in [None, Location.SURF] for elem in deps):
raise Exception(f'Species dependencies cannot include species in volumes.')
if not isinstance(raft, Raft):
raise TypeError(f'Expected a raft, got {raft} instad.')
self._raft = raft
self._setParameter('Dependencies', deps)
self._setParameter('K', kcst, nutils.Units('s^-1'))
if _createObj:
deps = self.Dependencies._GetStepsElems()
self.stepsRaftGen = stepslib._py_RaftGen(
self.name, surfSys._getStepsObjects()[0], deps,
raft._getStepsObjects()[0], self.K
)
else:
self.stepsRaftGen = None
def _getStepsObjects(self):
"""Return a list of the steps objects that this named object holds."""
return [self.stepsRaftGen]
def _getAllElems(self, loc=None):
"""
Return all elements that can directly be accessed by name from the parent physical
location.
"""
return [self.Raft] + list(self.Dependencies._getAllElems(loc))
@classmethod
def _FromStepsObject(cls, obj, mdl):
"""Create the interface object from a STEPS object."""
if obj.getRaft() is not None:
raft = getattr(mdl, obj.getRaft().getID())
else:
raft = None
deps = ReactionSide([])
for spec in obj.getSpecDeps():
deps += getattr(mdl, spec.getID())
exo = cls(raft, obj.getKcst(), deps=deps, _createObj=False, name=obj.getID())
exo.stepsRaftGen = obj
return exo
@property
def Raft(self):
"""Raft associated with the genesis reaction
:type: :py:class:`Raft`, read-only
"""
return self._raft
@property
@nutils.ParameterizedObject.RegisterGetter(units=None)
def Dependencies(self):
"""Species dependencies for the raft genesis reaction
:type: :py:class:`ReactionSide`, read-only
"""
pass
@property
@nutils.ParameterizedObject.RegisterGetter(units=nutils.Units('s^-1'))
def K(self):
"""Raft genesis reaction constant (in s^-1)
:type: Union[float, :py:class:`steps.API_2.nutils.Parameter`]
"""
pass
@K.setter
@nutils.ParameterizedObject.RegisterSetter(units=nutils.Units('s^-1'))
def K(self, val):
if isinstance(val, numbers.Number):
self.stepsRaftGen.setKcst(val)
else:
raise TypeError(f'{val} is not a valid reaction constant.')
[docs]@nutils.FreezeAfterInit
class RaftDis(
nutils.UsingObjects(RaftSurfaceSystem),
nutils.StepsWrapperObject, nutils.ParameterizedObject
):
"""A specific type of interaction that models a raft dissolution event
A raft dissolution event results in removal of a raft at a given rate when
the population of species within a raft are **at or below** an
anti-dependencies signature. Any remaining species in the raft are
inserted into patch triangles.
:param kcst: Rate constant in s^-1
:type kcst: float
:param antideps: Species anti-dependencies on the raft surface. The raft
dissolution cannot happen if species count on the raft are **strictly
bigger** than at least one of the anti-dependencies.
:type deps: Union[None, :py:class:`ReactionElement`, :py:class:`ReactionSide`]
"""
def __init__(self, kcst=0, antideps=None, _createObj=True, **kwargs):
super().__init__(**kwargs)
surfSys, = self._getUsedObjects()
if isinstance(antideps, ReactionElement):
antideps = antideps._toReactionSide()
elif antideps is None:
antideps = ReactionSide([])
elif not isinstance(antideps, ReactionSide):
raise TypeError(f'Expected species dependency, got {antideps} instead.')
if any(elem.loc not in [None, Location.RAFTSURF] for elem in antideps):
raise Exception(f'Species dependencies cannot include species in volumes.')
self._setParameter('AntiDependencies', antideps)
self._setParameter('K', kcst, nutils.Units('s^-1'))
if _createObj:
antideps = self.AntiDependencies._GetStepsElems()
self.stepsRaftDis = stepslib._py_RaftDis(
self.name, surfSys._getStepsObjects()[0], antideps, self.K
)
else:
self.stepsRaftDis = None
def _getStepsObjects(self):
"""Return a list of the steps objects that this named object holds."""
return [self.stepsRaftDis]
def _getAllElems(self, loc=None):
"""
Return all elements that can directly be accessed by name from the parent physical
location.
"""
return self.AntiDependencies._getAllElems(loc)
@classmethod
def _FromStepsObject(cls, obj, mdl):
"""Create the interface object from a STEPS object."""
antideps = ReactionSide([])
for spec in obj.getSpecAntiDeps():
antideps += getattr(mdl, spec.getID())
exo = cls(obj.getKcst(), antideps=antideps, _createObj=False, name=obj.getID())
exo.stepsRaftDis = obj
return exo
@property
@nutils.ParameterizedObject.RegisterGetter(units=None)
def AntiDependencies(self):
"""Species antidependencies for the raft dissolution reaction
:type: :py:class:`ReactionSide`, read-only
"""
pass
@property
@nutils.ParameterizedObject.RegisterGetter(units=nutils.Units('s^-1'))
def K(self):
"""Raft dissolution reaction constant (in s^-1)
:type: Union[float, :py:class:`steps.API_2.nutils.Parameter`]
"""
pass
@K.setter
@nutils.ParameterizedObject.RegisterSetter(units=nutils.Units('s^-1'))
def K(self, val):
if isinstance(val, numbers.Number):
self.stepsRaftDis.setKcst(val)
else:
raise TypeError(f'{val} is not a valid reaction constant.')
[docs]@nutils.FreezeAfterInit
class VesicleBind(
nutils.UsingObjects(VolumeSystem),
nutils.StepsWrapperObject, nutils.ParameterizedObject
):
"""Vesicle binding reaction in a volume system
A vesicle binding event binds two vesicles by creating link species between them.
The types of vesicle are arbitrarily termed 1 and 2, and may be the same.
A reactant on vesicle surface of vesicle 1 binds with a reactant on vesicle
2, and link species are created, one on vesicle 1 and one on vesicle 2.
Internally these two link species will be associated with each other and
exist within a defined length bound.
:param vesicles: The two vesicles involved in the binding reaction (in a tuple)
:type vesicles: Tuple[:py:class:`Vesicle`, :py:class:`Vesicle`]
:param reactants: The species involved in the binding reaction for each vesicle (in a tuple,
the order should match the vesicles parameter)
:type reactants: Tuple[:py:class:`Species`, :py:class:`Species`]
:param linkProducts: The link species resulting from the binding for each vesicle (in a tuple,
the order should match the vesicles parameter)
:type linkProducts: Tuple[:py:class:`LinkSpecies`, :py:class:`LinkSpecies`]
:param lenMin: Minimum length of the link between the two vesicles, in m
:type lenMin: float
:param lenMax: Minimum length of the link between the two vesicles, in m
:type lenMax: float
:param kcst: Rate constant in M^-1 s^-1
:type kcst: float
:param immobilization: immobilization flag for the binding reaction (see
:py:attr:`VesicleBind.Immobilization`)
:type immobilization: Union[:py:data:`steps.API_2.model.IMMOBILIZING`,
:py:data:`steps.API_2.model.NO_EFFECT`, :py:data:`steps.API_2.model.MOBILIZING`]
:param deps: Dependencies for the binding reaction (in a tuple, the order should match the vesicle
parameter)
:type deps: Tuple[Union[None, :py:class:`ReactionElement`, :py:class:`ReactionSide`],
Union[None, :py:class:`ReactionElement`, :py:class:`ReactionSide`]]
"""
def __init__(
self, vesicles, reactants, linkProducts, lenMin, lenMax,
kcst=0, immobilization=NO_EFFECT, deps=(None, None), **kwargs
):
super().__init__(**kwargs)
self._volSys, = self._getUsedObjects()
self._stepsVesBind = None
self._declared = False
self._added = False
# Check that the arguments are of correct types
for val, cls in [(vesicles, Vesicle), (reactants, Species), (linkProducts, LinkSpecies)]:
if not (isinstance(val, tuple) and len(val) == 2 and all(isinstance(v, cls) for v in val)):
raise TypeError(f'Expected a 2-tuple of {cls.__name__}, got {val} instead.')
self._vesicles = vesicles
self._reactants = reactants
self._linkProducts = linkProducts
self.LengthMin = lenMin
self.LengthMax = lenMax
self.K = kcst
self.Immobilization = immobilization
self.Dependencies = deps
self._declared = True
def _getStepsObjects(self):
"""Return a list of the steps objects that this named object holds."""
return [self._stepsVesBind] if self._stepsVesBind is not None else []
def _getAllElems(self, loc=None):
"""
Return all elements that can directly be accessed by name from the parent physical
location.
"""
return list(self._vesicles)
def _getVesLinkSpecPairs(self):
"""Return pairs of (vesicle, link species) that are involved in the reaction"""
return zip(self._vesicles, self._linkProducts)
def _checkNotAdded(self):
"""Raise an exception if the vesbind reaction was already added to STEPS"""
if self._added:
raise Exception(
'Cannot set properties of a vesicle bind reaction that was already added to STEPS.'
)
def _exitCallback(self, parent):
"""
Method to be called the first time we get out of a context manager in which the object as been
declared.
"""
if not self._declared or self._added:
return
allDeps = []
for deps in self.Dependencies:
vdeps = []
ldeps = []
if deps is not None:
deps = deps._toReactionSide()
if any(elem.loc != Location.VESSURF for elem in deps):
raise Exception(
f'{self._decl()}: Dependencies should all be located on the vesicle surface '
f'(spec.v): {deps}'
)
for elem in deps._GetStepsElems():
if isinstance(elem, stepslib._py_Spec):
vdeps.append(elem)
elif isinstance(elem, stepslib._py_LinkSpec):
ldeps.append(elem)
else:
raise TypeError(f'{self._decl()}: Unexpected element in dependencies: {elem}')
allDeps.append((vdeps, ldeps))
self._stepsVesBind = stepslib._py_VesBind(
self.name, self._volSys._getStepsObjects()[0],
self._vesicles[0]._getStepsObjects()[0], self._reactants[0]._getStepsObjects()[0],
self._vesicles[1]._getStepsObjects()[0], self._reactants[1]._getStepsObjects()[0],
self._linkProducts[0]._getStepsObjects()[0], self._linkProducts[1]._getStepsObjects()[0],
self.LengthMax, self.LengthMin,
allDeps[0][0], allDeps[1][0],
allDeps[0][1], allDeps[1][1],
self.K, self.Immobilization
)
self._added = True
@classmethod
def _FromStepsObject(cls, obj, mdl):
"""Create the interface object from a STEPS object."""
ves1 = getattr(mdl, obj.getVesicle1().getID())
ves2 = getattr(mdl, obj.getVesicle2().getID())
r1 = getattr(mdl, obj.getReactant1().getID())
r2 = getattr(mdl, obj.getReactant2().getID())
l1 = getattr(mdl, obj.getProduct1().getID())
l2 = getattr(mdl, obj.getProduct2().getID())
dep1, dep2 = ReactionSide([]), ReactionSide([])
for spec in obj.getVDeps1():
dep1 += getattr(mdl, spec.getID()).v
for spec in obj.getVDeps2():
dep2 += getattr(mdl, spec.getID()).v
for spec in obj.getLDeps1():
dep1 += getattr(mdl, spec.getID()).v
for spec in obj.getLDeps2():
dep2 += getattr(mdl, spec.getID()).v
vb = cls(
(ves1, ves2), (r1, r2), (l1, l2),
obj.getLengthMin(), obj.getLengthMax(), obj.getKcst(), obj.getImmobilization(),
(dep1, dep2), name=obj.getID()
)
vb._stepsVesBind = obj
vb._added = True
return vb
@property
def Vesicles(self):
"""Vesicles associated with the binding reaction
:type: Tuple[:py:class:`Vesicle`, :py:class:`Vesicle`], read-only
"""
return self._vesicles
@property
@nutils.ParameterizedObject.RegisterGetter(units=nutils.Units('m'))
def LengthMin(self):
"""Minimum length of the link (in m)
:type: float
"""
pass
@LengthMin.setter
@nutils.ParameterizedObject.RegisterSetter(units=nutils.Units('m'))
def LengthMin(self, val):
self._checkNotAdded()
if not isinstance(val, numbers.Number):
raise TypeError(f'{val} is not a valid minimum link length.')
@property
@nutils.ParameterizedObject.RegisterGetter(units=nutils.Units('m'))
def LengthMax(self):
"""Maximum length of the link (in m)
:type: float
"""
pass
@LengthMax.setter
@nutils.ParameterizedObject.RegisterSetter(units=nutils.Units('m'))
def LengthMax(self, val):
self._checkNotAdded()
if not isinstance(val, numbers.Number):
raise TypeError(f'{val} is not a valid maximum link length.')
@property
@nutils.ParameterizedObject.RegisterGetter(units=nutils.Units('M^-1 s^-1'))
def K(self):
"""Binding reaction constant (in M^-1 s^-1)
:type: float
"""
return 0
@K.setter
@nutils.ParameterizedObject.RegisterSetter(units=nutils.Units('M^-1 s^-1'))
def K(self, val):
self._checkNotAdded()
if not isinstance(val, numbers.Number):
raise TypeError(f'{val} is not a valid reaction constant.')
@property
@nutils.ParameterizedObject.RegisterGetter(units=None)
def Immobilization(self):
"""Immobilization flag (:py:data:`IMMOBILIZING`, :py:data:`NO_EFFECT` or :py:data:`MOBILIZING`).
Defaults to :py:data:`NO_EFFECT`.
:type: Union[:py:data:`steps.API_2.model.IMMOBILIZING`, :py:data:`steps.API_2.model.NO_EFFECT`,
:py:data:`steps.API_2.model.MOBILIZING`]
"""
return NO_EFFECT
@Immobilization.setter
@nutils.ParameterizedObject.RegisterSetter(units=None)
def Immobilization(self, val):
if not isinstance(val, Immobilization):
raise TypeError(f'Expected an immobilization flag, got {val} instead.')
self._checkNotAdded()
@property
@nutils.ParameterizedObject.RegisterGetter(units=None)
def Dependencies(self):
"""Vesicle binding dependencies
A set of species or link species that need to be present in order for the binding reaction to happen.
The dependencies should be given in the shape of a reaction side, for example: `SA.v + SC.v`
The dependencies are empty by default.
:type: Tuple[Union[None, :py:class:`ReactionElement`, :py:class:`ReactionSide`],
Union[None, :py:class:`ReactionElement`, :py:class:`ReactionSide`]]
"""
# Return a default value if the property was not set
return None, None
@Dependencies.setter
@nutils.ParameterizedObject.RegisterSetter(units=None)
def Dependencies(self, deps):
self._checkNotAdded()
if (not isinstance(deps, tuple) or
not all(d is None or isinstance(d, (ReactionElement, ReactionSide)) for d in deps)):
raise TypeError(f'{deps} cannot be used as a dependency.')
[docs]@nutils.FreezeAfterInit
class VesicleUnbind(
nutils.UsingObjects(VolumeSystem),
nutils.StepsWrapperObject, nutils.ParameterizedObject
):
"""Vesicle unbinding reaction in a volume system
A vesicle unbinding event unbinds two vesicles bound by two link species.
The types of vesicle are arbitrarily termed 1 and 2, and may be the same.
Upon application of this reaction, the link species on vesicle 1 becomes a
species on the vesicle 1 surface and the link species on vesicle 2 becomes a
species on the vesicle 2 surface.
:param vesicles: The two vesicles involved in the binding reaction (in a tuple)
:type vesicles: Tuple[:py:class:`Vesicle`, :py:class:`Vesicle`]
:param linkReactants: The link species involved in the link for each vesicle (in a tuple,
the order should match the vesicles parameter)
:type linkReactants: Tuple[:py:class:`LinkSpecies`, :py:class:`LinkSpecies`]
:param products: The species resulting from the unbinding for each vesicle (in a tuple,
the order should match the vesicles parameter)
:type products: Tuple[:py:class:`Species`, :py:class:`Species`]
:param kcst: Rate constant in s^-1
:type kcst: float
:param immobilization: immobilization flag for the binding reaction (see
:py:attr:`VesicleUnbind.Immobilization`)
:type immobilization: Union[:py:data:`steps.API_2.model.IMMOBILIZING`,
:py:data:`steps.API_2.model.NO_EFFECT`, :py:data:`steps.API_2.model.MOBILIZING`]
"""
def __init__(self, vesicles, linkReactants, products, kcst=0, immobilization=NO_EFFECT, **kwargs):
super().__init__(**kwargs)
self._volSys, = self._getUsedObjects()
self._stepsVesUnbind = None
self._declared = False
self._added = False
# Check that the arguments are of correct types
for val, cls in [(vesicles, Vesicle), (products, Species), (linkReactants, LinkSpecies)]:
if not (isinstance(val, tuple) and len(val) == 2 and all(isinstance(v, cls) for v in val)):
raise TypeError(f'Expected a 2-tuple of {cls.__name__}, got {val} instead.')
self._vesicles = vesicles
self._linkReactants = linkReactants
self._products = products
self.K = kcst
self.Immobilization = immobilization
self._declared = True
def _getStepsObjects(self):
"""Return a list of the steps objects that this named object holds."""
return [self._stepsVesUnbind] if self._stepsVesUnbind is not None else []
def _getAllElems(self, loc=None):
"""
Return all elements that can directly be accessed by name from the parent physical
location.
"""
return list(self._vesicles)
def _getVesLinkSpecPairs(self):
"""Return pairs of (vesicle, link species) that are involved in the reaction"""
return zip(self._vesicles, self._linkReactants)
def _checkNotAdded(self):
"""Raise an exception if the vesunbind reaction was already added to STEPS"""
if self._added:
raise Exception(
'Cannot set properties of a vesicle unbind reaction that was already added to STEPS.'
)
def _exitCallback(self, parent):
"""
Method to be called the first time we get out of a context manager in which the object as been
declared.
"""
if not self._declared or self._added:
return
self._stepsVesUnbind = stepslib._py_VesUnbind(
self.name, self._volSys._getStepsObjects()[0],
self._linkReactants[0]._getStepsObjects()[0], self._linkReactants[1]._getStepsObjects()[0],
self._vesicles[0]._getStepsObjects()[0], self._products[0]._getStepsObjects()[0],
self._vesicles[1]._getStepsObjects()[0], self._products[1]._getStepsObjects()[0],
self.K, self.Immobilization
)
self._added = True
@classmethod
def _FromStepsObject(cls, obj, mdl):
"""Create the interface object from a STEPS object."""
ves1 = getattr(mdl, obj.getVesicle1().getID())
ves2 = getattr(mdl, obj.getVesicle2().getID())
l1 = getattr(mdl, obj.getLink1().getID())
l2 = getattr(mdl, obj.getLink2().getID())
p1 = getattr(mdl, obj.getProduct1().getID())
p2 = getattr(mdl, obj.getProduct2().getID())
vb = cls((ves1, ves2), (l1, l2), (p1, p2), obj.getKcst(), obj.getImmobilization(), name=obj.getID())
vb._stepsVesUnbind = obj
vb._added = True
return vb
@property
def Vesicles(self):
"""Vesicles associated with the binding reaction
:type: Tuple[:py:class:`Vesicle`, :py:class:`Vesicle`], read-only
"""
return self._vesicles
@property
@nutils.ParameterizedObject.RegisterGetter(units=nutils.Units('s^-1'))
def K(self):
"""Binding reaction constant (in s^-1)
:type: float
"""
return 0
@K.setter
@nutils.ParameterizedObject.RegisterSetter(units=nutils.Units('s^-1'))
def K(self, val):
self._checkNotAdded()
if not isinstance(val, numbers.Number):
raise TypeError(f'{val} is not a valid reaction constant.')
@property
@nutils.ParameterizedObject.RegisterGetter(units=None)
def Immobilization(self):
"""Immobilization flag (:py:data:`IMMOBILIZING`, :py:data:`NO_EFFECT` or :py:data:`MOBILIZING`)
Defaults to :py:data:`NO_EFFECT`.
:type: Union[:py:data:`steps.API_2.model.IMMOBILIZING`, :py:data:`steps.API_2.model.NO_EFFECT`,
:py:data:`steps.API_2.model.MOBILIZING`]
"""
return NO_EFFECT
@Immobilization.setter
@nutils.ParameterizedObject.RegisterSetter(units=None)
def Immobilization(self, val):
if not isinstance(val, Immobilization):
raise TypeError(f'Expected an immobilization flag, got {val} instead.')
self._checkNotAdded()
###################################################################################################
# Currents
class _SubCurrent(nutils.NamedObject):
"""Wrapper class for the steps OhmicCurr or GHKCurr object"""
def __init__(self, parent, obj, *args, **kwargs):
super().__init__(*args, name=obj.getID(), **kwargs)
self._parent = parent
self.stepsCurrent = obj
def _getStepsObjects(self):
"""Return a list of the steps objects that this named object holds."""
return [self.stepsCurrent]
def _solverStr(self):
"""Return the string that is used as part of method names for this specific object."""
return self._parent.__class__._currStr
def _simPathAutoMetaData(self):
"""Return a dictionary with string keys and string or numbers values."""
mtdt = {'obj_type': self._parent.__class__._currStr, 'obj_id': self.name}
# Add information about parent comp or patch
for key, val in self._parent._simPathAutoMetaData().items():
mtdt['parent_' + key] = val
return mtdt
class _SubCurrentList:
"""List of _SubCurrent objects"""
def __init__(self, lst, *args, **kwargs):
super().__init__(*args, **kwargs)
self._currentLst = lst
def __iter__(self):
return iter(self._currentLst)
def _simPathWalkExpand(self):
"""Return an iterable of the elements that should be part of ResultSelector paths."""
return self.__iter__()
def _simPathCombinerClass(self):
"""Return the class that needs to be used to combine expanded elements."""
return nutils.SumSimPathComb
[docs]@nutils.FreezeAfterInit
class Current(nutils.UsingObjects(SurfaceSystem), nutils.StepsWrapperObject, nutils.ParameterizedObject):
"""Base class representing a current associated with one or more channel states
All :py:class:`Current` objects need to be declared inside a :py:class:`SurfaceSystem`.
.. note::
A :py:class:`Current` object should not be directly instantiated by the user,
only its subclasses should be instantiated. It is only documented for clarity.
"""
def __init__(self, *args, **kwargs):
super().__init__(*args, **kwargs)
self._currents = {}
self._complex = None
self._added = False
self._declared = False
def _getAllElems(self, loc):
return (
list(set(state._comp for state, curr in self._currents.items()))
+ [state for state, curr in self._currents.items()]
+ [curr for state, curr in self._currents.items()]
)
def _getComplexStates(self, states):
self._setParameter('Opened state', states)
if isinstance(states, Complex):
return states, states[...]._getAllStates()
if isinstance(states, ComplexSelector):
return states._complex, states._getAllStates()
elif isinstance(states, ComplexState):
return states._comp, [states]
else:
raise TypeError(f'Expected a ComplexSelector or a ComplexState, got {states} instead.')
def _getCurrentFromState(self, state):
if state not in self._currents:
raise Exception(f'ComplexState {state} is not associated to a _SubCurrent in {self}.')
return self._currents[state]
[docs] def __getitem__(self, key):
"""Get sub current from a current object
Square bracket notation can be used on a current to retrieve the specific current object
corresponding to a specific channel state. This should never be needed during model
declaration but it is very important for simulation control and data saving (see
:py:class:`steps.API_2.sim.SimPath`). If a :py:class:`ComplexSelector` covering several complex
states is used, an iterable list of subcurrents is returned.
:param key: A sequence of one or more :py:class:`SubUnitState` or
:py:class:`SubUnitSelector` that corresponds to one or several :py:class:`ComplexState`.
:type key: Tuple[Union[:py:class:`SubUnitState`, :py:class:`SubUnitSelector`], ...]
:returns: The sub current(s) associated to ``key``
:meta public:
"""
state = self._complex.__getitem__(key)
if isinstance(state, ComplexState):
return self._getCurrentFromState(state)
elif isinstance(state, ComplexSelector):
return _SubCurrentList([self._getCurrentFromState(s) for s in state._getAllStates()])
else:
raise NotImplementedError()
def __iter__(self):
return iter([curr for state, curr in self._currents.items()])
def _simPathWalkExpand(self):
"""Return an iterable of the elements that should be part of ResultSelector paths."""
return self.__iter__()
def _simPathCombinerClass(self):
"""Return the class that needs to be used to combine expanded elements."""
if self._complex._statesAsSpecies:
return nutils.SumSimPathComb
else:
return None
def _getStepsObjects(self):
"""Return a list of the steps objects that this named object holds."""
return [obj.stepsCurrent for state, obj in self._currents.items()]
def _createStepsObj(self):
"""Create the actual STEPS objects and populate self._currents"""
raise NotImplementedError()
def _exitCallback(self, parent):
"""
Method to be called the first time we get out of a context manager in which the object as been
declared.
"""
self._added = True
self._createStepsObj()
def _checkNotAdded(self):
"""Raise an exception if the currents were already added"""
if self._added:
raise Exception(f'Cannot modify the properties of a current that was already added to STEPS.')
def _simPathAutoMetaData(self):
"""Return a dictionary with string keys and string or numbers values."""
return {'obj_type': self.__class__._currStr, 'obj_id': self.name}
[docs]@nutils.FreezeAfterInit
class OhmicCurr(Current):
"""An ohmic current through a channel in one or more conducting states
:param states: The conducting state(s) of the channel
:type states: Union[:py:class:`ComplexState`, :py:class:`ComplexSelector`]
:param conduct: A constant conductance in Siemens or a :py:class:`CompDepCond` function that
returns a conductance in Siemens as a function of a complex state. Defaults to 0 Siemens.
:type conduct: Union[float, :py:class:`CompDepCond`]
:param rev_pot: Reversal potential in Volts. Defaults to 0 Volts.
:type rev_pot: float
An ohmic current object is a simple current based on a (possibly state specific) fixed value
for single-channel conductance and a constant reversal potential. An ohmic current
does not result in movement of ions between compartments, but simply
contributes a continuous current to the EField solver for every conducting channel
on a membrane surface that contains the ohmic current.
Usage::
with ssys:
# Constant conductance:
OC_SK = OhmicCurr.Create(SKchan[sko1|sko2], SK_G, SK_rev)
# The above line declares an ohmic current named OC_SK for channel SKchan when its
# single subunit is in states sko1 or sko2.
# Complex dependent conductance:
cond = CompDepCond(lambda s: k1 + k2 * s.Count(S1A), [Chan1])
OC_Chan1 = OhmicCurr.Create(Chan1[:, :, S2A], cond, e_rev)
# The above 2 lines declare an ohmic current from channel Chan1 when its 3rd subunit
# is in state S2A. Thanks to the CompDepCond function, the conductance changes
# depending on the specific states of the other subunits. In this case, there is a
# basal conductance k1 and this conductance is increased by k2 for each subunit in
# state S1A.
.. warning::
The corresponding currents are only truly added to steps after the `with ssys:` block
is exited.
"""
_currStr = 'Ohmic'
def __init__(self, states, conduct=0, rev_pot=0, *args, **kwargs):
super().__init__(*args, **kwargs)
self._complex, self._states = self._getComplexStates(states)
# Properties
self.G = conduct
self.ERev = rev_pot
self._declared = True
def _createStepsObj(self):
"""Create the actual STEPS objects and populate self._currents"""
if not self._declared:
return
surfSys = self._getParentOfType(SurfaceSystem)
conduct = self.G
rev_pot = self.ERev
if isinstance(conduct, CompDepFunc):
if len(conduct._complexes) != 1 or conduct._complexes[0] is not self._complex:
raise Exception('The CompDepFunc used as conductance is not compatible with the channel.')
condStates = [(s, conduct(s)) for s in self._states]
else:
condStates = [(s, conduct) for s in self._states]
for state, g in condStates:
chanState = self._complex._compStates[state]
self._currents[state] = _SubCurrent(
self,
stepslib._py_OhmicCurr(
f'{self.name}_{state.name}', surfSys.stepsSys, chanstate=chanState, g=g, erev=rev_pot
),
)
@property
@nutils.ParameterizedObject.RegisterGetter(units=nutils.Units('V'))
def ERev(self):
"""Reversal potential in Volts
Can only be set in the `with ssys:` block in which the current has been declared.
:type: Union[float, :py:class:`steps.API_2.nutils.Parameter`]
"""
pass
@ERev.setter
@nutils.ParameterizedObject.RegisterSetter(units=nutils.Units('V'))
def ERev(self, val):
self._checkNotAdded()
if not isinstance(val, numbers.Number):
raise TypeError(f'Expected a float, got {val} instead')
@property
@nutils.ParameterizedObject.RegisterGetter(units=nutils.Units('S'))
def G(self):
"""Conductance in Siemens
Can only be set in the `with ssys:` block in which the current has been declared.
:type: Union[float, :py:class:`CompDepCond`]
"""
pass
@G.setter
@nutils.ParameterizedObject.RegisterSetter(units=nutils.Units('S'))
def G(self, val):
self._checkNotAdded()
if not isinstance(val, (numbers.Number, CompDepCond)):
raise TypeError(f'Expected a float or a CompDepCond object, got {val} instead')
[docs]@nutils.FreezeAfterInit
class GHKCurr(Current):
"""A Goldman-Hodgkin-Katz current through a channel in one or more conducting states
:param states: The conducting state(s) of the channel.
:type states: Union[:py:class:`ComplexState`, :py:class:`ComplexSelector`]
:param spec: The species involved in the current, must have a non zero
:py:func:`Species.valence`.
:type spec: :py:class:`Species`
:param P: A single channel permeability in m^3/s, or a :py:class:`PInfo` object, or a
:py:class:`CompDepP` function that returns a permeability as a function of a complex state.
:type P: Union[float, :py:class:`PInfo`, :py:class:`CompDepP`]
:param computeflux: If True, then the current will result in movement of ions between
compartments, if False the current will be calculated but will not correspond to a real
ion flux.
:type computeflux: bool
:param virtual_oconc: A 'virtual outer concentration' that can be specified so that the outer
compartment does not have to be explicitly simulated (if it retains default `None`
value then the outer concentration of ion must be simulated).
:type virtual_oconc: float
:param vshift:
:type vshift: float
Each GHK current in the simulation is solved within the SSA with a rate determined
from the simulation state, i.e. membrane potential, 'outer' and 'inner'
concentration of the ion and temperature (which is fixed), and constant parameters
ion valence and permeability. The user can supply a call to :py:func:`PInfo` instead of a
permeability, STEPS will compute the corresponding permeability.
A GHK current involves, optionally, a real transfer of ions between comparments separated
by a membrane in the STEPS simulation. It is important that these ions implement diffusion
if compartments are not well-mixed.
Usage::
with mdl:
Ca = Species.Create()
Ca.valence = 2
with ssys:
# Constant permeability:
GC_CaP = GHKCurr.Create(CaPchan[m3], Ca, CaP_P, computeflux = True)
# The above line declares GHK current named OC_CaP for Ca2+ ions and for
# channel CaPchan when its single subunit is in state m3.
# Complex dependent permeability:
P = CompDepP(lambda s: p1 + p2 * s.Count(S1A), [Chan2])
GC_Chan2 = GHKCurr.Create(Chan2[:, :, S2A], Ca, P, computeflux = True)
# The above 2 lines declare GHK current for Ca2+ ions and for channel Chan2
# when its 3rd subunit is in state S2A. Thanks to the CompDepP function,
# the permeability changes depending on the specific states of the other
# subunits. In this case, there is a basal permeability p1 and this
# permeability is increased by p2 for each subunit in state S1A.
"""
_currStr = 'GHK'
class _GHKCurrPInfo(nutils.ParameterizedObject):
"""Simple wrapper for providing Pinfo instead of permeability"""
def __init__(self, g, V, T, oconc, iconc, *args, **kwargs):
super().__init__(*args, **kwargs)
self.g = g
self.V = V
self.T = T
self.oconc = oconc
self.iconc = iconc
self._valence = None
def _setValence(self, valence):
self._valence = valence
@property
def value(self):
if self._valence is None:
raise Exception(
f'The PInfo object was not associated with any GHK current, cannot compute permeability.'
)
return stepslib._py_permeability(self.g, self.V, self._valence, self.T, self.iconc, self.oconc)
_GHKCurrPInfo.RegisterParameter('g', nutils.Units('S'))
_GHKCurrPInfo.RegisterParameter('V', nutils.Units('V'))
_GHKCurrPInfo.RegisterParameter('T', nutils.Units('K'))
_GHKCurrPInfo.RegisterParameter('oconc', nutils.Units('M'))
_GHKCurrPInfo.RegisterParameter('iconc', nutils.Units('M'))
def __init__(self, states, spec, P, computeflux=True, virtual_oconc=None, vshift=0, *args, **kwargs):
super().__init__(*args, **kwargs)
if not isinstance(spec, Species):
raise TypeError(f'Expected a Species, got {spec} instead.')
self._spec = spec
self._complex, self._states = self._getComplexStates(states)
self._computeflux = computeflux
self._vshift = vshift
# Properties
self.P = copy.copy(P)
if isinstance(P, GHKCurr._GHKCurrPInfo):
self.P._setValence(spec.valence)
self.VOConc = virtual_oconc
self._declared = True
def _createStepsObj(self):
"""Create the actual STEPS objects and populate self._currents"""
if not self._declared:
return
surfSys = self._getParentOfType(SurfaceSystem)
P = self.P
virtual_oconc = self.VOConc if self.VOConc is not None else -1
if isinstance(P, CompDepFunc):
if len(P._complexes) != 1 or P._complexes[0] is not self._complex:
raise Exception('The CompDepFunc used as permeability is not compatible with the channel.')
condStates = [(s, P(s)) for s in self._states]
else:
condStates = [(s, P) for s in self._states]
for state, perm in condStates:
chanState = self._complex._compStates[state]
self._currents[state] = _SubCurrent(
self,
stepslib._py_GHKcurr(
f'{self.name}_{state.name}',
surfSys.stepsSys,
chanState,
self._spec.stepsSpecies,
computeflux=self._computeflux,
virtual_oconc=virtual_oconc,
vshift=self._vshift,
),
)
if isinstance(perm, numbers.Number):
self._currents[state].stepsCurrent.setP(perm)
elif isinstance(perm, GHKCurr._GHKCurrPInfo):
self._currents[state].stepsCurrent.setPInfo(
perm.g, perm.V, perm.T, perm.oconc, perm.iconc
)
else:
raise TypeError(f'Expected a permeability, got {perm} instead.')
[docs] @classmethod
def PInfo(cls, g, V, T, oconc, iconc):
"""Provide information from channel measurement instead of permeability
:param g: A measured single-channel conductance in Siemens
:type g: Union[float, :py:class:`steps.API_2.nutils.Parameter`]
:param V: The potential in Volts
:type V: Union[float, :py:class:`steps.API_2.nutils.Parameter`]
:param T: The temperature in Kelvin (may be different from the STEPS simulation
temperature)
:type T: Union[float, :py:class:`steps.API_2.nutils.Parameter`]
:param oconc: The 'outer' concentration of the ion in molar units
:type oconc: Union[float, :py:class:`steps.API_2.nutils.Parameter`]
:param iconc: The 'inner' concentration of the ion in molar units
:type iconc: Union[float, :py:class:`steps.API_2.nutils.Parameter`]
:returns: An aggregate object. Does not return the permeability.
Usage::
Pinfos = GHKCurr.PInfo(g=20e-12, V=-22e-3, T=293.15, oconc=4e-3, iconc=155e-3)
GC_CaP = GHKCurr.Create(CaPchan[m3], Ca, Pinfos, computeflux = True)
The actual permeability can be accessed with `GC_CaP.P.value`.
"""
return GHKCurr._GHKCurrPInfo(g, V, T, oconc, iconc)
def _getAllElems(self, loc):
if loc != Location.SURF:
return super()._getAllElems(loc) + [self._spec]
else:
return super()._getAllElems(loc)
@property
@nutils.ParameterizedObject.RegisterGetter(units=nutils.Units('M'))
def VOConc(self):
"""Virtual outer concentration, in Molars
Can only be set in the `with ssys:` block in which the current has been declared.
:type: Union[float, :py:class:`steps.API_2.nutils.Parameter`]
"""
pass
@VOConc.setter
@nutils.ParameterizedObject.RegisterSetter(units=nutils.Units('M'))
def VOConc(self, val):
self._checkNotAdded()
@property
@nutils.ParameterizedObject.RegisterGetter(units=nutils.Units('m^3 s^-1'))
def P(self):
"""Single channel permeability in m^3 s^-1
Can only be set in the `with ssys:` block in which the current has been declared.
:type: Union[float, :py:class:`steps.API_2.nutils.Parameter`, :py:class:`PInfo`, :py:class:`CompDepP`]
"""
pass
@P.setter
@nutils.ParameterizedObject.RegisterSetter(units=nutils.Units('m^3 s^-1'))
def P(self, val):
self._checkNotAdded()