jLM.RDME#

Spatial model

Classes

File(fname[, replicate, latticeType, ...])

Load a previously defined simulation

Sim(name, filename, dimensions, ...[, ...])

Define and run an RDME simulation

SpatialModel(name, filename, dimensions, ...)

Base RDME model

partial

partial(func, *args, **keywords) - new function with partial application of the given arguments and keywords.

tqdm(*_, **__)

Decorate an iterable object, returning an iterator which acts exactly like the original iterable, but prints a dynamically updating progressbar every time a value is requested.

class jLM.RDME.File(fname, replicate=1, latticeType=None, sample_frame=False, max_frames=100)[source]#

Bases: TrajAnalysisMixin, FileJupyterMixin, LmReadMixin, SpatialModel

Load a previously defined simulation

Load a RDME simulation file

Parameters:
  • fname (str) – LM data filename

  • replicate (int) – Replicate to load initially

  • latticeType (str) – “Byte” / “Int”

  • sample_frame (bool) – Sample frames or not

  • max_frames (int) – Number of frames after sampling

Note

jLM includes some metadata in the Lattice Microbes HDF5 file which is necessary to reload the model.

assignReaction(reaction, region)#

Assign a reaction to a region

Parameters:
close_hdf_file()[source]#
construct()#

Track newly created model objects and display in Notebook.

Context manager which tracks new species, reactions, etc., and displays a HTML summary when used in Jupyter

diffusionConst(rate, value, **kwargs)#

Lookup/create diffusion constant instance

Parameters:
  • rate (str) – Name of diffusion constant

  • value (float) – Value of diffusion constant

  • texRepr (str) – TeX math mode representation of the rate constant

  • annotation (str) – Description text

Returns:

DiffusionConst instance

Return type:

DiffusionConst

displayGeometry(**kwargs)#
distributeConcentration(sp, reg, conc)#

Distribute a concentration of particles uniformly through a region

Parameters:
  • sp (Species]) – Species type

  • reg (Region]) – Region type

  • conc (float) – Concentration of particles

distributeNumber(sp, reg, count)#

Distribute a number of particles uniformly through a region

Parameters:
  • sp (Species]) – Species type

  • reg (Region]) – Region type

  • count (int) – Number of particles

getLatticeHistogram(regex=None, species=None, integrated='xyz', replicate=None, frameStart=None, frameEnd=None, timeStart=None, timeEnd=None)#

Calculate the spatial distribution of particles over an interval

Species can be selected individually with the species key. If a list of species are given instead the sum of the particle numbers for those species types will be returned. A regular expression can also be used to select the species to return. Species types or species names in the form of a string are both acceptible to use in the species key. The regex and species options are mutually exclusive. The interval of time in which the histogram is computed can be selected either through the frame numbers (frameStart, frameEnd) or through the simulation time (timeStart, timeEnd). If the full 3-D lattice is not needed, any of the directions can be integrated out.

Parameters:
  • regex (str) – Regular expression matching species names. Multiple matches will be summed.

  • species ([Species]) – Species

  • integrated (str) – Combination of ‘x’, ‘y’, ‘z’ specifying the directions to integrate out

  • replicate (int) – Replicate to return, default is current replicate.

  • frameStart (int) – Starting frame

  • frameEnd (int) – Ending frame

  • timeStart (int) – Starting time

  • timeEnd (int) – Ending Time

Returns:

The average particle counts

Return type:

ndarray(dtype=float)

Note

Since this calculation requires reading the entire particle lattice trajectory, it can be slow.

getLatticeTrajectory(regex=None, species=None, integrated='', replicate=None, frameStart=None, frameEnd=None, timeStart=None, timeEnd=None)#

Calculate the spatial distribution of particles versus time

Species can be selected individually with the species key. If a list of species are given instead the sum of the particle numbers for those species types will be returned. A regular expression can also be used to select the species to return. Species types or species names in the form of a string are both acceptible to use in the species key. The regex and species options are mutually exclusive. Instead of returning the entire trajectory, a segment can be selected either through the frame numbers (frameStart, frameEnd) or through the simulation time (timeStart, timeEnd). If the full 3-D lattice is not needed, any of the directions can be integrated out.

Parameters:
  • regex (str) – Regular expression matching species names. Multiple matches will be summed.

  • species ([Species]) – Species

  • integrated (str) – Combination of ‘x’, ‘y’, ‘z’ specifying the directions to integrate out

  • replicate (int) – Replicate to return, default is current replicate.

  • frameStart (int) – Starting frame

  • frameEnd (int) – Ending frame

  • timeStart (int) – Starting time

  • timeEnd (int) – Ending Time

Returns:

The evaluation times and particle counts

Return type:

(ndarray(shape=(nt,), dtype=float), ndarray(dtype=float))

Note

Since this calculation requires reading the entire particle lattice trajectory, it can be slow.

getNumberTrajectory(regex=None, species=None, replicate=None, frameStart=None, frameEnd=None, timeStart=None, timeEnd=None)#

Calculate particle number trajectories

The time course of the particle numbers in the simulation can be queried with this function. Species can be selected individually with the species key. If a list of species are given instead the sum of the particle numbers for those species types will be returned. A regular expression can also be used to select the species to return. Species types or species names in the form of a string are both acceptible to use in the species key. The regex and species options are mutually exclusive. Instead of returning the entire trajectory, a segment can be selected either through the frame numbers (frameStart, frameEnd) or through the simulation time (timeStart, timeEnd).

Parameters:
  • regex (str) – Regular expression matching species names. Multiple matches will be summed.

  • species ([Species]) – Species

  • replicate (int) – Replicate to return, default is current replicate.

  • frameStart (int) – Starting frame

  • frameEnd (int) – Ending frame

  • timeStart (int) – Starting time

  • timeEnd (int) – Ending Time

Returns:

The evaluation times and particle counts

Return type:

(ndarray(shape=(nt,), dtype=float), ndarray(shape=(nt,), dtype=float))

getNumberTrajectoryFromRegion(spRegex=None, species=None, regRegex=None, region=None, replicate=None, frameStart=None, frameEnd=None, frameDownScale=None, timeStart=None, timeEnd=None, timeDownScale=None)#

Calculate particle number trajectories for specific regions

The time course of the particle numbers in the simulation within specific regions can be queried with this function. Species can be selected individually with the species key. If a list of species are given instead the sum of the particle numbers for those species types will be returned. A regular expression can also be used to select the species to return. Species types or species names in the form of a string are both acceptible to use in the species key. The spRegex and species options are mutually exclusive. Instead of returning the entire trajectory, a segment can be selected either through the frame numbers (frameStart, frameEnd) or through the simulation time (timeStart, timeEnd).If both given, the frame numbers will be used. The regions to compute particle numbers over are selected similar to the species through the options regRegex and region.

Parameters:
  • spRegex (str) – Regular expression matching species names. Multiple matches will be summed.

  • species ([Species]) – Species

  • regRegex (str) – Regular expression matching region names. Multiple matches will be summed.

  • region ([Region]) – Region

  • replicate (int) – Replicate to return, default is current replicate.

  • frameStart (int) – Starting frame

  • frameEnd (int) – Ending frame

  • timeStart (int) – Starting time

  • timeEnd (int) – Ending Time

Returns:

The evaluation times and particle counts

Return type:

(ndarray(shape=(nt,), dtype=float), ndarray(shape=(nt,), dtype=float))

Note

Since this calculation requires reading the entire particle lattice trajectory, it can be slow.

init_worker(file_name)[source]#
maxDiffusionRate(latticeSpacing=None, dt=None)#

Compute max allowed diffusion constant for the simulation

Parameters:
  • latticeSpacing (float) – Lattice spacing in meters

  • dt (float) – Timestep in seconds

Returns:

Maximum diffusion constant in m^/s

Return type:

float

particleStatistics(particleLattice=None, siteLattice=None)#

Compute properties of the particle lattice

Dictionary keys

countBySpeciesRegion

Particle counts for each species type for each region.

concBySpeciesRegion

Concentration for each species type for each region.

countByRegion

Total particle counts for each region.

concByRegion

Total concentration for each region.

countBySpecies

Total particle counts for each species.

concBySpecies

Total concentration of each species, averaged over simulation volume

count

Total particle count

conc

Total concentration, averaged over simulation volume

vol

Total simulation volume

siteCount

Total number of subvolumes

regionVol

Volume of each region

regionCounts

Number of subvolumes in each region

Returns:

Particle lattice statistics

Return type:

dict

placeNumber(sp, x, y, z, n)#

Place a number of particles at a specific subvolume

Parameters:
  • sp (Species) – Species type

  • x (int) – x-coordinate

  • y (int) – y-coordinate

  • z (int) – z-coordinate

  • n (int) – Number of particles

process_frame(frame_index, siteTypesLattice, numberSpecies, numberSiteTypes)[source]#
rateConst(rate, value, order, **kwargs)#

Lookup/create reaction rate constant instance

Parameters:
  • rate (str) – Name of rate constant

  • value (float) – Value of rate constant

  • order (int) – Order of reaction

  • texRepr (str) – TeX math mode representation of the rate constant

  • annotation (str) – Description text

Returns:

RateConst instance

Return type:

RateConst

reaction(reactants, products, rate, value=None, regions=None, **kwargs)#

Lookup/create reaction instance

The product/reactants can be specified as a string, Species instance, list of strings, list of Species instances, or the empty list, denoting no species. The reaction rate can be created in place if a string is given for rate and the value of the rate constant is specified by value. Rates must be provided in standard chemical kinetic units, e.g. for a reaction

\[\mathrm{A} + \mathrm{B} \overset{k}{\to} \mathrm{C},\]

the rate coefficent, \(k\) has units of \(\mathrm{M}^{-1}\mathrm{s}^{-1}\)

Parameters:
  • reactants ([Species]) – Species, or list of species acting as reactants

  • products ([Species]) – Species, or list of species acting as products

  • rate ([RateConst]) – Rate constant

  • value (float) – Optional value of rate constant

  • regions (py:class:~jLM.Types.Region]) – List of applicable regions

Returns:

Reaction instance

Return type:

Reaction

reduce_trajectory_file(file_name, max_frames)[source]#
region(name, **kwargs)#

Lookup/create region instance

Parameters:
  • name (str) – Name of the region

  • texRepr (str) – TeX math mode representation of the region

  • annotation (str) – Description text

Returns:

Region instance

Return type:

Region

resizeLattice(dimensions, latticeSpacing, latticeType=None)#
run(solver=None, replicate=1, seed=None, cudaDevices=None, checkpointInterval=0, sample_frame=False, max_frames=100)#

Run the RDME simulation

Parameters:
  • solver (lm.RDMESolver) – Rdme solver

  • replicate (int) – Index of replicate

  • seed (int) – RNG seed

  • cudaDevices ([int]) – List of CUDA device indexes

  • checkpointInterval (int) – Number of seconds between checkpoints

  • sample_frame (bool) – Sample frames or not

  • max_frames (int) – Number of frames after sampling

Returns:

Simulation result

Return type:

jLM.File

setMaximumTimestep()#

Set the simulation timestep using the fastest diffusion rate

setReplicate(replicate=1)#

Select replicate from HDF5 file

Parameters:

replicate (int) – Replicate index (1-based)

showAllParameters(**kwargs)#
showAllSpecies(**kwargs)#
showDiffusionConstants(**kwargs)#
showRateConstants(**kwargs)#
showReactions(**kwargs)#
showRegion(**kwargs)#
showRegionStack(**kwargs)#
showSpecies(**kwargs)#
species(name, **kwargs)#

Lookup/create species instance

Parameters:
  • name (str) – Name of the species

  • texRepr (str) – TeX math mode representation of the species

  • annotation (str) – Description text

Returns:

Species instance

Return type:

Species

transitionRate(sp, rFrom, rTo, rate, value=None)#

Define the diffusive transition rate between regions

This allows for fine grained detail over the diffusion matrix. The transition is defined for the species sp from region rFrom to rTo. If None is given for any of these parameters, then the entire axis of the matrix is affected, i.e. None is a wildcard. The rate is specified with the rate parameter, either as a DiffusionConst instance, or as a string to lookup the previously created DiffusionConst. If a string is provided and a DiffusionConst instance has not been created, providing the parameter value will create the new constant.

Parameters:
  • sp (Species]) – Species type or None

  • rFrom (Region]) – Region type or None

  • rTo (Region]) – Region type or None

  • rate (DiffusionConst]) – Diffusion rate

  • value (float) – Value of diffusion constant if new

vmd_precomp(file_name)[source]#
write_to_h5(file_name, maxParticleCounts, maxSiteCounts)[source]#
NA#

Avogadro’s number

Type:

float

property diffusionFast#

Lookup/create a maximum-valued diffusion constant instance

Returns:

DiffusionConst instance

Return type:

DiffusionConst

property diffusionZero#

Lookup/create a zero-valued diffusion constant instance

Returns:

DiffusionConst instance

Return type:

DiffusionConst

filename#

LM data filename

Type:

str

property hookInterval#

Simulation time between calls to hookSimulation

latticeType#

Lattice type (Int or Byte)

Type:

str

property latticeWriteInterval#

Simulation time between writing the lattice to the HDF5 file

name#

Simulation name

Type:

str

property perfPrintInterval#

Real elapsed time between writing simulation progress to stdout

property simulationTime#

Total simulated time

property speciesWriteInterval#

Simulation time between writing the species counts to the HDF5 file

class jLM.RDME.Sim(name, filename, dimensions, latticeSpacing, regionName, latticeType=None, dt=None)[source]#

Bases: LmWriteMixin, SpatialModel

Define and run an RDME simulation

Create new RDME object

Parameters:
  • name (str) – Name of simulation

  • filename (str) – LM data filename

  • dimensions ((nx,ny,nz)) – Dimensions of lattice

  • latticeSpacing (float) – Lattice spacing in meters

  • regionName (str) – Name of the default region

  • latticeType (str) – “Byte” / “Int”

  • dt (float) – Timestep

assignReaction(reaction, region)#

Assign a reaction to a region

Parameters:
construct()#

Track newly created model objects and display in Notebook.

Context manager which tracks new species, reactions, etc., and displays a HTML summary when used in Jupyter

copyParticleLattice(fname, diffusion_model=False, replicate=1, frame=-1)#

Copy particle lattice from existing simulation

Parameters:

fname (str) – LM simulation H5 file

Keyword Arguments:
  • diffusion_model (bool) – Use the initial condition particle lattice

  • replicate (int) – LM simulation replicate

  • frame (int) – LM simulation particle lattice frame. Negative values index from end.

diffusionConst(rate, value, **kwargs)#

Lookup/create diffusion constant instance

Parameters:
  • rate (str) – Name of diffusion constant

  • value (float) – Value of diffusion constant

  • texRepr (str) – TeX math mode representation of the rate constant

  • annotation (str) – Description text

Returns:

DiffusionConst instance

Return type:

DiffusionConst

displayGeometry(**kwargs)#
distributeConcentration(sp, reg, conc)#

Distribute a concentration of particles uniformly through a region

Parameters:
  • sp (Species]) – Species type

  • reg (Region]) – Region type

  • conc (float) – Concentration of particles

distributeNumber(sp, reg, count)#

Distribute a number of particles uniformly through a region

Parameters:
  • sp (Species]) – Species type

  • reg (Region]) – Region type

  • count (int) – Number of particles

finalize()#
h5()#
maxDiffusionRate(latticeSpacing=None, dt=None)#

Compute max allowed diffusion constant for the simulation

Parameters:
  • latticeSpacing (float) – Lattice spacing in meters

  • dt (float) – Timestep in seconds

Returns:

Maximum diffusion constant in m^/s

Return type:

float

particleStatistics(particleLattice=None, siteLattice=None)[source]#

Compute properties of the particle lattice

Dictionary keys

countBySpeciesRegion

Particle counts for each species type for each region.

concBySpeciesRegion

Concentration for each species type for each region.

countByRegion

Total particle counts for each region.

concByRegion

Total concentration for each region.

countBySpecies

Total particle counts for each species.

concBySpecies

Total concentration of each species, averaged over simulation volume

count

Total particle count

conc

Total concentration, averaged over simulation volume

vol

Total simulation volume

siteCount

Total number of subvolumes

regionVol

Volume of each region

regionCounts

Number of subvolumes in each region

Returns:

Particle lattice statistics

Return type:

dict

placeNumber(sp, x, y, z, n)#

Place a number of particles at a specific subvolume

Parameters:
  • sp (Species) – Species type

  • x (int) – x-coordinate

  • y (int) – y-coordinate

  • z (int) – z-coordinate

  • n (int) – Number of particles

rateConst(rate, value, order, **kwargs)#

Lookup/create reaction rate constant instance

Parameters:
  • rate (str) – Name of rate constant

  • value (float) – Value of rate constant

  • order (int) – Order of reaction

  • texRepr (str) – TeX math mode representation of the rate constant

  • annotation (str) – Description text

Returns:

RateConst instance

Return type:

RateConst

reaction(reactants, products, rate, value=None, regions=None, **kwargs)#

Lookup/create reaction instance

The product/reactants can be specified as a string, Species instance, list of strings, list of Species instances, or the empty list, denoting no species. The reaction rate can be created in place if a string is given for rate and the value of the rate constant is specified by value. Rates must be provided in standard chemical kinetic units, e.g. for a reaction

\[\mathrm{A} + \mathrm{B} \overset{k}{\to} \mathrm{C},\]

the rate coefficent, \(k\) has units of \(\mathrm{M}^{-1}\mathrm{s}^{-1}\)

Parameters:
  • reactants ([Species]) – Species, or list of species acting as reactants

  • products ([Species]) – Species, or list of species acting as products

  • rate ([RateConst]) – Rate constant

  • value (float) – Optional value of rate constant

  • regions (py:class:~jLM.Types.Region]) – List of applicable regions

Returns:

Reaction instance

Return type:

Reaction

region(name, **kwargs)#

Lookup/create region instance

Parameters:
  • name (str) – Name of the region

  • texRepr (str) – TeX math mode representation of the region

  • annotation (str) – Description text

Returns:

Region instance

Return type:

Region

resizeLattice(dimensions, latticeSpacing, latticeType=None)#
run(solver=None, replicate=1, seed=None, cudaDevices=None, checkpointInterval=0, sample_frame=False, max_frames=100)#

Run the RDME simulation

Parameters:
  • solver (lm.RDMESolver) – Rdme solver

  • replicate (int) – Index of replicate

  • seed (int) – RNG seed

  • cudaDevices ([int]) – List of CUDA device indexes

  • checkpointInterval (int) – Number of seconds between checkpoints

  • sample_frame (bool) – Sample frames or not

  • max_frames (int) – Number of frames after sampling

Returns:

Simulation result

Return type:

jLM.File

setMaximumTimestep()#

Set the simulation timestep using the fastest diffusion rate

showAllParameters(**kwargs)#
showAllSpecies(**kwargs)#
showDiffusionConstants(**kwargs)#
showRateConstants(**kwargs)#
showReactions(**kwargs)#
showRegion(**kwargs)#
showRegionStack(**kwargs)#
showSpecies(**kwargs)#
species(name, **kwargs)#

Lookup/create species instance

Parameters:
  • name (str) – Name of the species

  • texRepr (str) – TeX math mode representation of the species

  • annotation (str) – Description text

Returns:

Species instance

Return type:

Species

transitionRate(sp, rFrom, rTo, rate, value=None)#

Define the diffusive transition rate between regions

This allows for fine grained detail over the diffusion matrix. The transition is defined for the species sp from region rFrom to rTo. If None is given for any of these parameters, then the entire axis of the matrix is affected, i.e. None is a wildcard. The rate is specified with the rate parameter, either as a DiffusionConst instance, or as a string to lookup the previously created DiffusionConst. If a string is provided and a DiffusionConst instance has not been created, providing the parameter value will create the new constant.

Parameters:
  • sp (Species]) – Species type or None

  • rFrom (Region]) – Region type or None

  • rTo (Region]) – Region type or None

  • rate (DiffusionConst]) – Diffusion rate

  • value (float) – Value of diffusion constant if new

NA#

Avogadro’s number

Type:

float

dc#

Convienient diffusion constant access

Type:

Namespace

diffRateList#

List of diffusion constants

Type:

SimObjs

property diffusionFast#

Lookup/create a maximum-valued diffusion constant instance

Returns:

DiffusionConst instance

Return type:

DiffusionConst

property diffusionZero#

Lookup/create a zero-valued diffusion constant instance

Returns:

DiffusionConst instance

Return type:

DiffusionConst

filename#

LM data filename

Type:

str

property hookInterval#

Simulation time between calls to hookSimulation

latticeType#

Lattice type (Int or Byte)

Type:

str

property latticeWriteInterval#

Simulation time between writing the lattice to the HDF5 file

name#

Simulation name

Type:

str

property perfPrintInterval#

Real elapsed time between writing simulation progress to stdout

rc#

Convienient reaction rate access

Type:

Namespace

reactionList#

List of reactions

Type:

SimObjs

reg#

Convienient region access

Type:

Namespace

regionList#

List of regions

Type:

SimObjs

rxnRateList#

List of reaction rates

Type:

SimObjs

property simulationTime#

Total simulated time

sp#

Convienient species access

Type:

Namespace

speciesList#

List of species types

Type:

SimObjs

property speciesWriteInterval#

Simulation time between writing the species counts to the HDF5 file

class jLM.RDME.SpatialModel(name, filename, dimensions, latticeSpacing, latticeType=None)[source]#

Bases: LatticeAnalysisMixin, JupyterDisplayMixin

Base RDME model

assignReaction(reaction, region)[source]#

Assign a reaction to a region

Parameters:
construct()#

Track newly created model objects and display in Notebook.

Context manager which tracks new species, reactions, etc., and displays a HTML summary when used in Jupyter

diffusionConst(rate, value, **kwargs)[source]#

Lookup/create diffusion constant instance

Parameters:
  • rate (str) – Name of diffusion constant

  • value (float) – Value of diffusion constant

  • texRepr (str) – TeX math mode representation of the rate constant

  • annotation (str) – Description text

Returns:

DiffusionConst instance

Return type:

DiffusionConst

displayGeometry(**kwargs)#
distributeConcentration(sp, reg, conc)[source]#

Distribute a concentration of particles uniformly through a region

Parameters:
  • sp (Species]) – Species type

  • reg (Region]) – Region type

  • conc (float) – Concentration of particles

distributeNumber(sp, reg, count)[source]#

Distribute a number of particles uniformly through a region

Parameters:
  • sp (Species]) – Species type

  • reg (Region]) – Region type

  • count (int) – Number of particles

maxDiffusionRate(latticeSpacing=None, dt=None)[source]#

Compute max allowed diffusion constant for the simulation

Parameters:
  • latticeSpacing (float) – Lattice spacing in meters

  • dt (float) – Timestep in seconds

Returns:

Maximum diffusion constant in m^/s

Return type:

float

particleStatistics(particleLattice=None, siteLattice=None)#

Compute properties of the particle lattice

Dictionary keys

countBySpeciesRegion

Particle counts for each species type for each region.

concBySpeciesRegion

Concentration for each species type for each region.

countByRegion

Total particle counts for each region.

concByRegion

Total concentration for each region.

countBySpecies

Total particle counts for each species.

concBySpecies

Total concentration of each species, averaged over simulation volume

count

Total particle count

conc

Total concentration, averaged over simulation volume

vol

Total simulation volume

siteCount

Total number of subvolumes

regionVol

Volume of each region

regionCounts

Number of subvolumes in each region

Returns:

Particle lattice statistics

Return type:

dict

placeNumber(sp, x, y, z, n)[source]#

Place a number of particles at a specific subvolume

Parameters:
  • sp (Species) – Species type

  • x (int) – x-coordinate

  • y (int) – y-coordinate

  • z (int) – z-coordinate

  • n (int) – Number of particles

rateConst(rate, value, order, **kwargs)[source]#

Lookup/create reaction rate constant instance

Parameters:
  • rate (str) – Name of rate constant

  • value (float) – Value of rate constant

  • order (int) – Order of reaction

  • texRepr (str) – TeX math mode representation of the rate constant

  • annotation (str) – Description text

Returns:

RateConst instance

Return type:

RateConst

reaction(reactants, products, rate, value=None, regions=None, **kwargs)[source]#

Lookup/create reaction instance

The product/reactants can be specified as a string, Species instance, list of strings, list of Species instances, or the empty list, denoting no species. The reaction rate can be created in place if a string is given for rate and the value of the rate constant is specified by value. Rates must be provided in standard chemical kinetic units, e.g. for a reaction

\[\mathrm{A} + \mathrm{B} \overset{k}{\to} \mathrm{C},\]

the rate coefficent, \(k\) has units of \(\mathrm{M}^{-1}\mathrm{s}^{-1}\)

Parameters:
  • reactants ([Species]) – Species, or list of species acting as reactants

  • products ([Species]) – Species, or list of species acting as products

  • rate ([RateConst]) – Rate constant

  • value (float) – Optional value of rate constant

  • regions (py:class:~jLM.Types.Region]) – List of applicable regions

Returns:

Reaction instance

Return type:

Reaction

region(name, **kwargs)[source]#

Lookup/create region instance

Parameters:
  • name (str) – Name of the region

  • texRepr (str) – TeX math mode representation of the region

  • annotation (str) – Description text

Returns:

Region instance

Return type:

Region

resizeLattice(dimensions, latticeSpacing, latticeType=None)[source]#
run(solver=None, replicate=1, seed=None, cudaDevices=None, checkpointInterval=0, sample_frame=False, max_frames=100)[source]#

Run the RDME simulation

Parameters:
  • solver (lm.RDMESolver) – Rdme solver

  • replicate (int) – Index of replicate

  • seed (int) – RNG seed

  • cudaDevices ([int]) – List of CUDA device indexes

  • checkpointInterval (int) – Number of seconds between checkpoints

  • sample_frame (bool) – Sample frames or not

  • max_frames (int) – Number of frames after sampling

Returns:

Simulation result

Return type:

jLM.File

setMaximumTimestep()[source]#

Set the simulation timestep using the fastest diffusion rate

showAllParameters(**kwargs)#
showAllSpecies(**kwargs)#
showDiffusionConstants(**kwargs)#
showRateConstants(**kwargs)#
showReactions(**kwargs)#
showRegion(**kwargs)#
showRegionStack(**kwargs)#
showSpecies(**kwargs)#
species(name, **kwargs)[source]#

Lookup/create species instance

Parameters:
  • name (str) – Name of the species

  • texRepr (str) – TeX math mode representation of the species

  • annotation (str) – Description text

Returns:

Species instance

Return type:

Species

transitionRate(sp, rFrom, rTo, rate, value=None)[source]#

Define the diffusive transition rate between regions

This allows for fine grained detail over the diffusion matrix. The transition is defined for the species sp from region rFrom to rTo. If None is given for any of these parameters, then the entire axis of the matrix is affected, i.e. None is a wildcard. The rate is specified with the rate parameter, either as a DiffusionConst instance, or as a string to lookup the previously created DiffusionConst. If a string is provided and a DiffusionConst instance has not been created, providing the parameter value will create the new constant.

Parameters:
  • sp (Species]) – Species type or None

  • rFrom (Region]) – Region type or None

  • rTo (Region]) – Region type or None

  • rate (DiffusionConst]) – Diffusion rate

  • value (float) – Value of diffusion constant if new

NA#

Avogadro’s number

Type:

float

property diffusionFast#

Lookup/create a maximum-valued diffusion constant instance

Returns:

DiffusionConst instance

Return type:

DiffusionConst

property diffusionZero#

Lookup/create a zero-valued diffusion constant instance

Returns:

DiffusionConst instance

Return type:

DiffusionConst

filename#

LM data filename

Type:

str

property hookInterval#

Simulation time between calls to hookSimulation

latticeType#

Lattice type (Int or Byte)

Type:

str

property latticeWriteInterval#

Simulation time between writing the lattice to the HDF5 file

name#

Simulation name

Type:

str

property perfPrintInterval#

Real elapsed time between writing simulation progress to stdout

property simulationTime#

Total simulated time

property speciesWriteInterval#

Simulation time between writing the species counts to the HDF5 file