Numerical Solvers

The Virtual Cell converts the biological description defined in a BioModel into a corresponding mathematical system of equations (e.g. ordinary and/or partial differential equations), VCell will then solve the equations by applying numerical solvers to perform and analyze simulations.

In Virtual Cell, numerical solvers can be divided into the following 4 types.

Deterministic Compartmental Solvers The compartmental simulation executes a single point simulation based on the defined physiological model, the geometric assumptions, and structure sizes or surface to volume ratios and volume fractions. The ordinary differential equations (ODE), for single point approximations, representing the reactions kinetics are generated and passed to an interpreted ODE solver. The set of nonlinear ODE equations are typically solved in seconds. It allows an interactive modification of parameters and a quick determination of the effect over time. In addition, the Virtual Cell software computes the local sensitivity (Sensitivity Analysis) of any species concentration to any parameter as a function of time evaluated at the nominal solution.

• IDA (Variable Order, Variable Time Step, ODE/DAE): IDA addresses systems of differential-algebraic equations (DAEs), and uses Backward Differentiation Formula methods. ODEs are a subset of DAEs, therefore IDA may be used for solving ODEs.
• Forward Euler (First Order, Fixed Time Step): Forward Euler method is a fixed time step mehtod to solve ordinary differential equations.
• Combined Stiff Solver (IDA/CVODE) This chooses between IDA and CVODE depending on the problem to be solved.
• CVODE (Variable Order, Variable Time Step): CVODE is used for solving initial value problems for ordinary differential equations. It solves both stiff and nonstiff systems, using variable-coefficient Adams and BDF methods. In the stiff case, options for treating the Jacobian of the system include dense and band matrix solvers, and a preconditioned Krylov (iterative) solver. In the highly modular organization of CVODE, the core integrator module is independent of the linear system solvers, and all operations on N-vectors are isolated in a module of vector kernels.
• Runge-Kutta (Fourth Order, Fixed Time Step) and Runge-Kutta (Second Order, Fixed Time Step): The Runge Kutta methods provide further systematic improvement in the spirit of the modified Euler method. Second order fixed time step method is also called the midpoint method.
• Adams-Moulton (Fifth Order, Fixed Time Step): The methods such as Foward Euler, Runge-Kutta etc. are called single-step methods because they use only the information from one previous point to compute the successive point. Adams-Moulton methods are explicit linear multistep methods that depend on multiple previous solution points to generate a new approximate solution point. It is a fixed time step method.
• Runge-Kutta-Fehlberg (Fifth Order, Variable Time Step): The Runge-Kutta-Fehlberg integrator is primarily designed to solve non-stiff and mildly stiff differential equations when derivative evaluations are inexpensive. It should generally not be used when the user is demanding high accuracy. It is a variable time step method.

Stochastic Compartmental Solvers

Different numerical methods are required to simulate different biological systems. Deterministic simulation is valid as long as the concentrations of simulated populations are high. However, the concentrations of these factors could be very low in living cells. In these latter cases, the biology is more accurately simulated by discrete approaches. As a result, stochastic simulation has been implemented in Virtual Cell to allow users to describe the discrete nature of changes in cell systems.

• Hybrid (Gibson + Euler-Maruyama Method): This is a hybrid stochastic method. It partitions the system into subsets of fast and slow reactions and approximates the fast reactions as a continuous Markov process, using a chemical Langevin equation, and accurately describes the slow dynamics using the Gibson algorithm. Fixed time step Euler-Maruyama is used for approximate numerical solution of CLE.
• Hybrid (Gibson + Milstein Method): This is a hybrid stochastic method. It partitions the system into subsets of fast and slow reactions and approximates the fast reactions as a continuous Markov process, using a chemical Langevin equation, and accurately describes the slow dynamics using the Gibson algorithm. Fixed time step Milstein is used for approximate numerical solution of CLE.
• Hybrid (Adaptive Gibson + Milstein Method): This is a hybrid stochastic method. It partitions the system into subsets of fast and slow reactions and approximates the fast reactions as a continuous Markov process, using a chemical Langevin equation, and accurately describes the slow dynamics using the Gibson algorithm. Adaptive time step Milstein is used for approximate numerical solution of CLE.
• Gibson (Next Reaction Stochastic Method): Gibson-Bruck is an improved exact stochastic method based on Gillespie's SSA. It uses only a single random number per simulation event and takes time proportional to the logarithm of the number of reactions. Better performance is also acheived by utilizing a dependency graph and an indexed priority queue.

Note: to exercise hybrid solvers you need to go to species specifications and check the checkbox for "Force Continuous" for at least one species. See Species Specifications.

Deterministic Spatial Solvers

Virtual Cell uses the finite volume method to deterministically solve spatial problems. The finite volume method is a method of spatial discretization of partial differential equations which exactly preserves conservation laws. Similar to the finite difference method, values are calculated at discrete places on a meshed geometry. The standalone version of the finite volume method is a little slower but gives better error messages.

• Fully-Implicit Finite Volume, Regular Grid (Variable Time Step): This is our fully implicit, adaptive time step finite volume method. The finite volume method represents partial differential equations as algebraic discretization equations which exactly preserves conservation laws. Similar to the finite difference method, values are calculated at discrete places on a meshed geometry. This method employs Sundials stiff solver CVODE for time stepping (method of lines). Please note that relative and absolute tolerances affect the accuracy of time descritization only, therefore spatial discritization is the only significant source of solution error.
• EBChombo, Semi-Implicit (Fixed Time Step, Experimental): Chombo provides a set of tools for implementing finite difference methods for the solution of partial differential equations on block-structured adaptively refined rectangular grids. Both elliptic and time-dependent modules are included. Chombo supports calculations in complex geometries with both embedded boundaries and mapped grids, and Chombo also supports particle methods. Most parallel platforms are supported, and cross-platform self-describing file formats are included.
• Semi-Implicit Finite Volume-Particle Hybrid, Regular Grid (Fixed Time Step): This was our interpreted standalone version of the finite volume method. The finite volume method is a method for representing and evaluating partial differential equations as algebraic discretization equations which exactly preserves conservation laws. Similar to the finite difference method, values were calculated at discrete places on a meshed geometry.
• Semi-Implicit Finite Volume Compiled, Regular Grid (Fixed Time Step)(DISABLED, No Longer Available): use Semi-Implicit Finite Volume-Particle Hybrid, Regular Grid (Fixed Time Step) instead). It is possible to view results run in prior versions of Virtual Cell under this solver, but new simulations can't be run, nor can models that reference this solver be edited and saved. The finite volume method is a method for representing and evaluating partial differential equations as algebraic discretization equations which exactly preserves conservation laws. Similar to the finite difference method, values are calculated at discrete places on a meshed geometry.

Stochastic Spatial Solvers

Virtual Cell incorporates the Smoldyn(Smoluchowski Dynamics) to stochastically solve spatial problems. Smoldyn is a computer program for cell-scale biochemical simulations. It simulates each molecule of interest individually to capture natural stochasticity and for nanometer-scale spatial resolution. It treats other molecules implicity, so it can simulate tens of thousands of molecules over several minutes of real time. Simulated molecules diffuse, react, are confined by surfaces, and bind to membranes much as they would in a real biological system.

• Smoldyn (Spatial Stochastic Simulator): Smoldyn is a Brownian dynamics simulator. It represents space as a 1-, 2-, or 3-dimensional continuum, as opposed to a lattice, and it steps through time using finite length time steps. Smoldyn represents molecules as individual point-like particles and membranes as infinitely thin surfaces. Smoldyn simulates molecular diffusion, chemical reactions between individual molecules, and a wide variety of molecule-surface interactions. So far, Smoldyn has been used primarily for either detailed biophysics research problems, such as on diffusion-influenced reaction dynamics, or for investigating the effects of spatial organization on simple biological systems, such as the Escherichia coli chemotaxis system.