#include <VoxelPools.h>
Public Member Functions | |
VoxelPools () | |
virtual | ~VoxelPools () |
void | reinit (double dt) |
void | setStoich (Stoich *stoich, const OdeSystem *ode) |
void | advance (const ProcInfo *p) |
Do the numerical integration. Advance the simulation. | |
void | setInitDt (double dt) |
Set initial timestep to use by the solver. | |
void | setVolumeAndDependencies (double vol) |
void | updateAllRateTerms (const vector< RateTerm * > &rates, unsigned int numCoreRates) |
Updates all the rate constants from the reference rates vector. | |
void | updateRateTerms (const vector< RateTerm * > &rates, unsigned int numCoreRates, unsigned int index) |
void | updateRates (const double *s, double *yprime) const |
void | updateReacVelocities (const double *s, vector< double > &v) const |
void | print () const |
Used for debugging. | |
Static Public Member Functions | |
static int | gslFunc (double t, const double *y, double *dydt, void *params) |
This is the function which evaluates the rates. |
This is the class for handling reac-diff voxels used for deterministic computations.
VoxelPools::VoxelPools | ( | ) |
VoxelPools::~VoxelPools | ( | ) | [virtual] |
References VoxelPoolsBase::rates_.
void VoxelPools::advance | ( | const ProcInfo * | p | ) |
Do the numerical integration. Advance the simulation.
References ProcInfo::currTime, ProcInfo::dt, and VoxelPoolsBase::varS().
int VoxelPools::gslFunc | ( | double | t, | |
const double * | y, | |||
double * | dydt, | |||
void * | params | |||
) | [static] |
This is the function which evaluates the rates.
References q, VoxelPoolsBase::stoichPtr_, Stoich::updateFuncs(), and updateRates().
Referenced by Ksolve::setStoich().
void VoxelPools::print | ( | ) | const |
Used for debugging.
For debugging: Print contents of voxel pool.
Changes cross rate terms to zero if there is no junction void filterCrossRateTerms( const vector< pair< Id, Id > >& vec );
Reimplemented from VoxelPoolsBase.
References Stoich::getNumCoreRates(), VoxelPoolsBase::rates_, and VoxelPoolsBase::stoichPtr_.
void VoxelPools::reinit | ( | double | dt | ) |
void VoxelPools::setInitDt | ( | double | dt | ) |
Set initial timestep to use by the solver.
setStoich: Assigns the ODE system and Stoich. Stoich may be modified to create a new RateTerm vector in case the volume of this VoxelPools is new.
References OdeSystem::epsAbs, OdeSystem::epsRel, OdeSystem::initStepSize, VoxelPoolsBase::reinit(), and VoxelPoolsBase::stoichPtr_.
Referenced by SteadyState::setStoich().
void VoxelPools::setVolumeAndDependencies | ( | double | vol | ) | [virtual] |
Handles volume change and subsequent cascading updates.
Handle volume updates. Inherited Virtual func.
Reimplemented from VoxelPoolsBase.
References Stoich::getNumCoreRates(), Stoich::getRateTerms(), Stoich::setupCrossSolverReacVols(), VoxelPoolsBase::stoichPtr_, and updateAllRateTerms().
void VoxelPools::updateAllRateTerms | ( | const vector< RateTerm * > & | rates, | |
unsigned int | numCoreRates | |||
) | [virtual] |
Updates all the rate constants from the reference rates vector.
Implements VoxelPoolsBase.
Referenced by SteadyState::setStoich(), and setVolumeAndDependencies().
void VoxelPools::updateRates | ( | const double * | s, | |
double * | yprime | |||
) | const |
Core computation function. Updates the reaction velocities vector yprime given the current mol 'n' vector s.
References KinSparseMatrix::computeRowRate(), Stoich::getNumAllPools(), Stoich::getNumBufPools(), Stoich::getNumProxyPools(), Stoich::getNumVarPools(), Stoich::getStoichiometryMatrix(), N, SparseMatrix< T >::nColumns(), SparseMatrix< T >::nRows(), VoxelPoolsBase::rates_, and VoxelPoolsBase::stoichPtr_.
Referenced by SteadyState::classifyState(), and gslFunc().
void VoxelPools::updateRateTerms | ( | const vector< RateTerm * > & | rates, | |
unsigned int | numCoreRates, | |||
unsigned int | index | |||
) | [virtual] |
updateRateTerms updates the rate consts of a belonging to the specified index on the rates vector. It does recaling and assigning using values from the internal rates vector.
Implements VoxelPoolsBase.
void VoxelPools::updateReacVelocities | ( | const double * | s, | |
vector< double > & | v | |||
) | const |
updateReacVelocities computes the velocity *v* of each reaction from the vector *s* of pool s. This is a utility function for programs like SteadyState that need to analyze velocity.