Documentation for Release 2024.2

Cell Based Demo

This tutorial is automatically generated from TestCellBasedDemoTutorial.hpp at revision 409e06cb314b. Note that the code is given in full at the bottom of the page.

Examples showing how to create, run and cell-based simulations in Chaste

Introduction

This tutorial is designed to give you a quick introduction to running cell-based simulations in Chaste. Full details are postponed until later tutorials.

We begin with a simple monolayer simulation and see how to:

change the cell-level model;
how to impose boundaries;
how to impose periodic conditions;
how to specify how to remove cells; and
how to change cell-cycle models.

The test

We begin by including the necessary header files. These will be described in detail in subsequent cell-based tutorials.

#include <cxxtest/TestSuite.h>
#include "CellBasedSimulationArchiver.hpp"
#include "SmartPointers.hpp"
#include "AbstractCellBasedTestSuite.hpp"
#include "AdhesionPottsUpdateRule.hpp"
#include "CellsGenerator.hpp"
#include "CylindricalHoneycombMeshGenerator.hpp"
#include "GeneralisedLinearSpringForce.hpp"
#include "HoneycombMeshGenerator.hpp"
#include "HoneycombVertexMeshGenerator.hpp"
#include "MeshBasedCellPopulationWithGhostNodes.hpp"
#include "NagaiHondaForce.hpp"
#include "SimpleTargetAreaModifier.hpp"
#include "NodeBasedCellPopulation.hpp"
#include "OffLatticeSimulation.hpp"
#include "OnLatticeSimulation.hpp"
#include "PlaneBoundaryCondition.hpp"
#include "PottsBasedCellPopulation.hpp"
#include "PottsMeshGenerator.hpp"
#include "RandomCellKiller.hpp"
#include "RepulsionForce.hpp"
#include "UniformG1GenerationalCellCycleModel.hpp"
#include "SurfaceAreaConstraintPottsUpdateRule.hpp"
#include "TysonNovakCellCycleModel.hpp"
#include "VertexBasedCellPopulation.hpp"
#include "VolumeConstraintPottsUpdateRule.hpp"
#include "VoronoiDataWriter.hpp"

#include "FakePetscSetup.hpp"

Next, we define the test class which inherits from AbstractCellBasedTestSuite. We inherit from AbstractCellBasedTestSuite rather than CxxTest::TestSuite directly because this class sets up and destroys some singleton objects for us. Singletons are objects that we want to exist only once in each simulation and will be covered in detail in later tutorials. Since we are using AbstractCellBasedTestSuite the singleton SimulationTime is initialised to zero at the beginning of the test and destroyed at the end of the test; RandomNumberGenerator is re-seeded with zero at the beginning and destroyed at the end of the test; and CellPropertyRegistry (which stores CellProperties, you learn about these in a later tutorial CreatingAndUsingANewCellProperty) is cleared at the beginning of the test. This makes for cleaner code.

class TestCellBasedDemoTutorial : public AbstractCellBasedTestSuite
{
public:

Test 1 - a basic vertex-based simulation

In the first test, we run a simple vertex-based simulation of an epithelial monolayer. Each cell in the simulation is assigned a simple stochastic cell-cycle model, the cells will divide randomly and never stop proliferating.

    void TestVertexBasedMonolayer()
    {

The first thing we define is a 2D (specified by the <2,2>) mesh which holds the spatial information of the simulation. To do this we use one of a number of MeshGenerators.

        HoneycombVertexMeshGenerator generator(2, 2);
        boost::shared_ptr<MutableVertexMesh<2,2> > p_mesh = generator.GetMesh();

We now generate a collection of cells. We do this by using a CellsGenerator and we specify the proliferative behaviour of the cell by choosing a CellCycleModel, here we choose a UniformG1GenerationalCellCycleModel where each cell is given a division time, drawn from a uniform distribution when it is created. (Note that here we need to use a phase based cell cycle model so that we can use the target area modifiers which are needed by the vertex based simulations). For a vertex simulation we need as may cells as elements in the mesh.

        std::vector<CellPtr> cells;
        MAKE_PTR(TransitCellProliferativeType, p_transit_type);
        CellsGenerator<UniformG1GenerationalCellCycleModel, 2> cells_generator;
        cells_generator.GenerateBasicRandom(cells, p_mesh->GetNumElements(), p_transit_type);

We now create a CellPopulation object (passing in the mesh and cells) to connect the mesh and the cells together. Here that is a VertexBasedCellPopulation and the dimension is <2>.

        VertexBasedCellPopulation<2> cell_population(*p_mesh, cells);

We now create an OffLatticeSimulation object and pass in the CellPopulation. We also set some options on the simulation like output directory, output multiple (so we don’t visualize every timestep), and end time.

        OffLatticeSimulation<2> simulator(cell_population);
        simulator.SetOutputDirectory("CellBasedDemo1");
        simulator.SetSamplingTimestepMultiple(200);
        simulator.SetEndTime(20.0);

To specify how cells move around, we create a “shared pointer” to a Force object and pass it to the OffLatticeSimulation. This is done using the MAKE_PTR macro as follows.

        MAKE_PTR(NagaiHondaForce<2>, p_force);
        simulator.AddForce(p_force);

We now pass a “target area modifier” to the simulation. This modifies the target area of individual cells and thus alters the relative forces between neighbouring cells.

        MAKE_PTR(SimpleTargetAreaModifier<2>, p_growth_modifier);
        simulator.AddSimulationModifier(p_growth_modifier);

Finally we call the Solve method on the simulation to run the simulation.

        simulator.Solve();

The next two lines are for test purposes only and are not part of this tutorial. We are checking that we reached the end time of the simulation with the correct number of cells. If different simulation input parameters are being explored the lines should be removed.

        TS_ASSERT_EQUALS(cell_population.GetNumRealCells(), 16u);
        TS_ASSERT_DELTA(SimulationTime::Instance()->GetTime(), 20.0, 1e-10);
    }

To visualize the results, open a new terminal, cd to the Chaste directory, then cd to anim. Then do: java Visualize2dVertexCells $CHASTE_TEST_OUTPUT/CellBasedDemo1/results_from_time_0. We may have to do: javac Visualize2dVertexCells.java beforehand to create the java executable.

The make_a_movie script can be used to generate a video based on the results of your simulation. To do this, first visualize the results using Visualize2dVertexCells as described above. Click on the box marked “Output” and play through the whole simulation to generate a sequence of .png images, one for each time step. Next, still in the anim folder, do: ./make_a_movie. This reads in the .png files and creates a video file called simulation.mpeg.

Results can also be visualized using Paraview. See the Visualizing With Paraview tutorial for more information.

Test 2 - basic node-based simulation

We next show how to modify the previous test to implement a ’node-based’ simulation, in which cells are represented by overlapping spheres (actually circles, since we’re in 2D).

    void TestNodeBasedMonolayer()
    {

We now need to create a NodesOnlyMesh we do this by first creating a MutableMesh and passing this to a helper method ConstructNodesWithoutMesh along with a interaction cut off length that defines the connectivity in the mesh.

        HoneycombMeshGenerator generator(2, 2); //**Changed**//
        boost::shared_ptr<MutableMesh<2,2> > p_generating_mesh = generator.GetMesh(); //**Changed**//
        NodesOnlyMesh<2> mesh; //**Changed**//
        mesh.ConstructNodesWithoutMesh(*p_generating_mesh, 1.5); //**Changed**//

We create the cells as before, only this time we need one cell per node.

        std::vector<CellPtr> cells;
        MAKE_PTR(TransitCellProliferativeType, p_transit_type);
        CellsGenerator<UniformG1GenerationalCellCycleModel, 2> cells_generator;
        cells_generator.GenerateBasicRandom(cells, mesh.GetNumNodes(), p_transit_type); //**Changed**//

This time we create a NodeBasedCellPopulation as we are using a NodesOnlyMesh.

        NodeBasedCellPopulation<2> cell_population(mesh, cells);//**Changed**//

We create an OffLatticeSimulation object as before, all we change is the output directory and output results more often as a larger default timestep is used for these simulations.

        OffLatticeSimulation<2> simulator(cell_population);
        simulator.SetOutputDirectory("CellBasedDemo2"); //**Changed**//
        simulator.SetSamplingTimestepMultiple(12); //**Changed**//
        simulator.SetEndTime(20.0);

We use a different Force which is suitable for node based simulations.

        MAKE_PTR(RepulsionForce<2>, p_force); //**Changed**//
        simulator.AddForce(p_force);

In all types of simulation you may specify how cells are removed from the simulation by specifying a CellKiller. You create these in the same was as the Force and pass them to the CellBasedSimulation. Note that here the constructor for RandomCellKiller requires some arguments to be passed to it, therefore we use the MAKE_PTR_ARGS macro.

        MAKE_PTR_ARGS(RandomCellKiller<2>, p_cell_killer, (&cell_population, 0.01)); //**Changed**//
        simulator.AddCellKiller(p_cell_killer);

Again we call the Solve method on the simulation to run the simulation.

        simulator.Solve();

The next two lines are for test purposes only and are not part of this tutorial. Again, we are checking that we reached the end time of the simulation with the correct number of cells.

        TS_ASSERT_EQUALS(cell_population.GetNumRealCells(), 7u);
        TS_ASSERT_DELTA(SimulationTime::Instance()->GetTime(), 20.0, 1e-10);
    }

To visualize the results, open a new terminal, cd to the Chaste directory, then cd to anim. Then do: java Visualize2dCentreCells $CHASTE_TEST_OUTPUT/CellBasedDemo2/results_from_time_0. We may have to do: javac Visualize2dCentreCells.java beforehand to create the java executable.

As described above, the make_a_movie script can be used to generate a video based on the results of your simulation.

Test 3 - basic mesh-based simulation

We next show how to modify the previous test to implement a ‘mesh-based’ simulation, in which cells are represented by their centres and a Voronoi tessellation is used to find nearest neighbours.

    void TestMeshBasedMonolayer()
    {

This time we just create a MutableMesh and use that to specify the spatial locations of cells.

        HoneycombMeshGenerator generator(2, 2);
        boost::shared_ptr<MutableMesh<2,2> > p_mesh = generator.GetMesh();  //**Changed**//

We create the same number of cells as the previous test.

        std::vector<CellPtr> cells;
        MAKE_PTR(TransitCellProliferativeType, p_transit_type);
        CellsGenerator<UniformG1GenerationalCellCycleModel, 2> cells_generator;
        cells_generator.GenerateBasicRandom(cells, p_mesh->GetNumNodes(), p_transit_type);

This time we create a MeshBasedCellPopulation as we are using a MutableMesh.

        MeshBasedCellPopulation<2> cell_population(*p_mesh, cells); //**Changed**//

To view the results of this and the subsequent mesh based tutorials in Paraview it is necessary to explicitly generate the required .vtu files. This is detailed in the VisualizingWithParaview tutorial. Note that the results in Paraview may appear different to those in the java based visualizer. This is related to the different methods used to generate voronoi tesselations in each and is resolved through the use of ‘ghost nodes’, as shown in the next test.

        cell_population.AddPopulationWriter<VoronoiDataWriter>();

We create an OffLatticeSimulation object as before, all we change is the output directory.

        OffLatticeSimulation<2> simulator(cell_population);
        simulator.SetOutputDirectory("CellBasedDemo3"); //**Changed**//
        simulator.SetSamplingTimestepMultiple(12);
        simulator.SetEndTime(20.0);

We use a different Force which is suitable for mesh based simulations.

        MAKE_PTR(GeneralisedLinearSpringForce<2>, p_force); //**Changed**//
        simulator.AddForce(p_force);

Again we call the Solve method on the simulation to run the simulation.

        simulator.Solve();

The next two lines are for test purposes only and are not part of this tutorial.

        TS_ASSERT_EQUALS(cell_population.GetNumRealCells(), 16u);
        TS_ASSERT_DELTA(SimulationTime::Instance()->GetTime(), 20.0, 1e-10);
    }

The results may be visualized using Visualize2dCentreCells as described in the previous test, with the results directory changed from CellBasedDemo2 to CellBasedDemo3.

Test 4 - basic mesh-based simulation with ghost nodes

We next show how to modify the previous test to include ‘ghost nodes’, which do not correspond to cells but are sometimes needed when using a Voronoi tessellation. We will discuss ghost nodes in more detail in subsequent cell-based tutorials.

    void TestMeshBasedMonolayerWithGhostNodes()
    {

This time we just create a MutableMesh and use that to specify the spatial locations of cells. Here we pass an extra argument to the HoneycombMeshGenerator which adds another 2 rows of nodes round the mesh, known as ghost nodes.

        HoneycombMeshGenerator generator(2, 2, 2); //**Changed**//
        boost::shared_ptr<MutableMesh<2,2> > p_mesh = generator.GetMesh();

We only want to create cells for non ghost nodes. To find these we get them from the HoneycombMeshGenerator using the method GetCellLocationIndices. We also use a different CellCycleModel. Here we use a TysonNovakCellCycleModel which solves a coupled set of ODEs for each cell to calculate when each cell divides.

        std::vector<unsigned> location_indices = generator.GetCellLocationIndices();//**Changed**//
        std::vector<CellPtr> cells;
        MAKE_PTR(TransitCellProliferativeType, p_transit_type);
        CellsGenerator<TysonNovakCellCycleModel, 2> cells_generator; //**Changed**//
        cells_generator.GenerateBasicRandom(cells, location_indices.size(), p_transit_type); //**Changed**//

This time we create a MeshBasedCellPopulation as we are using a MutableMesh and have ghost nodes. We also need to pass the indices of non ghost nodes as an extra argument.

        MeshBasedCellPopulationWithGhostNodes<2> cell_population(*p_mesh, cells, location_indices); //**Changed**//

Again Paraview output is explicitly requested.

        cell_population.AddPopulationWriter<VoronoiDataWriter>();

We create an OffLatticeSimulation object as before, all we change is the output directory and the end time. The Tyson Novak model is for yeast cells and therefore cells proliferate much more often and so we run the simulation for less time to keep cell numbers relatively small for this demo.

        OffLatticeSimulation<2> simulator(cell_population);
        simulator.SetOutputDirectory("CellBasedDemo4"); //**Changed**//
        simulator.SetSamplingTimestepMultiple(12);
        simulator.SetEndTime(2.0); //**Changed**//

We use the same Force as before and run the simulation in the same way.

        MAKE_PTR(GeneralisedLinearSpringForce<2>, p_force);
        simulator.AddForce(p_force);
        simulator.Solve();

The next two lines are for test purposes only and are not part of this tutorial.

        TS_ASSERT_EQUALS(cell_population.GetNumRealCells(), 16u);
        TS_ASSERT_DELTA(SimulationTime::Instance()->GetTime(), 2.0, 1e-10);
    }

The results may be visualized using Visualize2dCentreCells as described in the previous test, with the results directory changed from CellBasedDemo3 to CellBasedDemo4.

Test 5 - basic periodic mesh-based simulation

We next show how to modify the previous test to implement a periodic boundary to the left and right of the domain.

    void TestMeshBasedMonolayerPeriodic()
    {

We now want to impose periodic boundaries on the domain. To do this we create a Cylindrical2dMesh using a CylindricalHoneycombMeshGenerator.

        CylindricalHoneycombMeshGenerator generator(5, 2, 2); //**Changed**//
        boost::shared_ptr<Cylindrical2dMesh> p_mesh = generator.GetCylindricalMesh(); //**Changed**//

Again we create one cell for each non ghost node. Note that we have changed back to using a UniformG1GenerationalCellCycleModel.

        std::vector<unsigned> location_indices = generator.GetCellLocationIndices();
        std::vector<CellPtr> cells;
        MAKE_PTR(TransitCellProliferativeType, p_transit_type);
        CellsGenerator<UniformG1GenerationalCellCycleModel, 2> cells_generator; //**Changed**//
        cells_generator.GenerateBasicRandom(cells, location_indices.size(), p_transit_type);

We use the same CellPopulation, CellBasedSimulation (only changing the output directory and end time) and Force as before and run the simulation.

        MeshBasedCellPopulationWithGhostNodes<2> cell_population(*p_mesh, cells, location_indices);

Again Paraview output is explicitly requested.

        cell_population.AddPopulationWriter<VoronoiDataWriter>();

        OffLatticeSimulation<2> simulator(cell_population);
        simulator.SetOutputDirectory("CellBasedDemo5"); //**Changed**//
        simulator.SetSamplingTimestepMultiple(12);
        simulator.SetEndTime(20.0); //**Changed**//

        MAKE_PTR(GeneralisedLinearSpringForce<2>, p_force);
        simulator.AddForce(p_force);

        simulator.Solve();

The next two lines are for test purposes only and are not part of this tutorial.

        TS_ASSERT_EQUALS(cell_population.GetNumRealCells(), 29u);
        TS_ASSERT_DELTA(SimulationTime::Instance()->GetTime(), 20.0, 1e-10);
    }

The results may be visualized using Visualize2dCentreCells as described in the previous test, with the results directory changed from CellBasedDemo4 to CellBasedDemo5.

Test 6 - basic periodic mesh-based simulation with obstructions

We next show how to modify the previous test to include one or more ‘obstructions’ within the domain.

    void TestMeshBasedMonolayerPeriodicSolidBottomBoundary()
    {

We make the same Mesh, Cells, CellPopulation, CellBasedSimulation and forces as before, all we change is the output directory.

        CylindricalHoneycombMeshGenerator generator(5, 2, 2);
        boost::shared_ptr<Cylindrical2dMesh> p_mesh = generator.GetCylindricalMesh();

        std::vector<unsigned> location_indices = generator.GetCellLocationIndices();
        std::vector<CellPtr> cells;
        MAKE_PTR(StemCellProliferativeType, p_stem_type);
        CellsGenerator<UniformG1GenerationalCellCycleModel, 2> cells_generator;
        cells_generator.GenerateBasicRandom(cells, location_indices.size(), p_stem_type); //**Changed**//

        MeshBasedCellPopulationWithGhostNodes<2> cell_population(*p_mesh, cells, location_indices);
        cell_population.AddPopulationWriter<VoronoiDataWriter>();

        OffLatticeSimulation<2> simulator(cell_population);
        simulator.SetOutputDirectory("CellBasedDemo6"); //**Changed**//
        simulator.SetSamplingTimestepMultiple(50);
        simulator.SetEndTime(20.0);

        MAKE_PTR(GeneralisedLinearSpringForce<2>, p_force);
        simulator.AddForce(p_force);

We now want to impose the condition y>0 on the cells. To do this we create a “shared pointer” to a PlaneBoundaryCondition. Much like the RandomCellKiller earlier we pass arguments to the constructor (a point (0,0) on the plane (line in 2D) and an outward pointing normal to the plane (0,-1) ) using the MAKE_PTR_ARGS macro.

        c_vector<double,2> point = zero_vector<double>(2);
        c_vector<double,2> normal = zero_vector<double>(2);
        normal(1) = -1.0;
        MAKE_PTR_ARGS(PlaneBoundaryCondition<2>, p_bc, (&cell_population, point, normal));
        simulator.AddCellPopulationBoundaryCondition(p_bc);

Finally we call the Solve method as in all other simulations.

        simulator.Solve();

The next two lines are for test purposes only and are not part of this tutorial.

        TS_ASSERT_EQUALS(cell_population.GetNumRealCells(), 23u);
        TS_ASSERT_DELTA(SimulationTime::Instance()->GetTime(), 20.0, 1e-10);
    }

The results may be visualized using Visualize2dCentreCells as described in the previous test, with the results directory changed from CellBasedDemo5 to CellBasedDemo6.

Test 7 - basic Potts-based simulation

In the final test we show how to modify the earlier tests (using off lattice models) to implement a ‘Potts-based’ simulation, in which cells are represented by collections of sites on a fixed lattice.

    void TestPottsBasedMonolayer()
    {

In common with the off lattice simulations we begin by creating a mesh. Here we use the PottsMeshGenerator class to generate a PottsMesh each element in the mesh is a collection of lattice sites (represented by nodes at their centres). All the connectivity between lattice sites is defined by the PottsMeshGenerator, and there are arguments to make the domains periodic.

        PottsMeshGenerator<2> generator(20, 2, 4, 20, 2, 4); //**Changed**//
        boost::shared_ptr<PottsMesh<2> > p_mesh = generator.GetMesh(); //**Changed**//

We generate one cell for each element as in vertex based simulations.

        std::vector<CellPtr> cells;
        MAKE_PTR(TransitCellProliferativeType, p_transit_type);
        CellsGenerator<UniformG1GenerationalCellCycleModel, 2> cells_generator;
        cells_generator.GenerateBasicRandom(cells, p_mesh->GetNumElements(), p_transit_type);

As we have a PottsMesh we use a PottsBasedCellPopulation. Note here we also change the “temperature” of the Potts simulation to make cells more motile.

        PottsBasedCellPopulation<2> cell_population(*p_mesh, cells);//**Changed**//
        cell_population.SetTemperature(1.0);

As a Potts simulation is restricted to a lattice we create a OnSimulation object and pass in the CellPopulation in much the same way as an OffLatticeSimulation in the above examples. We also set some options on the simulation like output directory and end time.

        OnLatticeSimulation<2> simulator(cell_population);//**Changed**//
        simulator.SetOutputDirectory("CellBasedDemo7"); //**Changed**//
        simulator.SetEndTime(20.0);

In order to specify how cells move around we create “shared pointers” to UpdateRule objects and pass them to the OnLatticeSimulation. This is analogous to Forces in earlier examples.

        MAKE_PTR(VolumeConstraintPottsUpdateRule<2>, p_volume_constraint_update_rule); //**Changed**//
        simulator.AddUpdateRule(p_volume_constraint_update_rule); //**Changed**//
        MAKE_PTR(SurfaceAreaConstraintPottsUpdateRule<2>, p_surface_area_update_rule); //**Changed**//
        simulator.AddUpdateRule(p_surface_area_update_rule); //**Changed**//
        MAKE_PTR(AdhesionPottsUpdateRule<2>, p_adhesion_update_rule); //**Changed**//
        simulator.AddUpdateRule(p_adhesion_update_rule); //**Changed**//

We can add CellKillers as before.

        MAKE_PTR_ARGS(RandomCellKiller<2>, p_cell_killer, (&cell_population, 0.01));
        simulator.AddCellKiller(p_cell_killer);

Again we run the simulation by calling the Solve method.

        simulator.Solve();

The next two lines are for test purposes only and are not part of this tutorial.

        TS_ASSERT_EQUALS(cell_population.GetNumRealCells(), 16u);
        TS_ASSERT_DELTA(SimulationTime::Instance()->GetTime(), 20.0, 1e-10);
    }
};

To visualize the results, open a new terminal, cd to the Chaste directory, then cd to anim. Then do: java Visualize2dVertexCells $CHASTE_TEST_OUTPUT/CellBasedDemo7/results_from_time_0. We may have to do: javac Visualize2dVertexCells.java beforehand to create the java executable.

Full code

#include <cxxtest/TestSuite.h>
#include "CellBasedSimulationArchiver.hpp"
#include "SmartPointers.hpp"
#include "AbstractCellBasedTestSuite.hpp"
#include "AdhesionPottsUpdateRule.hpp"
#include "CellsGenerator.hpp"
#include "CylindricalHoneycombMeshGenerator.hpp"
#include "GeneralisedLinearSpringForce.hpp"
#include "HoneycombMeshGenerator.hpp"
#include "HoneycombVertexMeshGenerator.hpp"
#include "MeshBasedCellPopulationWithGhostNodes.hpp"
#include "NagaiHondaForce.hpp"
#include "SimpleTargetAreaModifier.hpp"
#include "NodeBasedCellPopulation.hpp"
#include "OffLatticeSimulation.hpp"
#include "OnLatticeSimulation.hpp"
#include "PlaneBoundaryCondition.hpp"
#include "PottsBasedCellPopulation.hpp"
#include "PottsMeshGenerator.hpp"
#include "RandomCellKiller.hpp"
#include "RepulsionForce.hpp"
#include "UniformG1GenerationalCellCycleModel.hpp"
#include "SurfaceAreaConstraintPottsUpdateRule.hpp"
#include "TysonNovakCellCycleModel.hpp"
#include "VertexBasedCellPopulation.hpp"
#include "VolumeConstraintPottsUpdateRule.hpp"
#include "VoronoiDataWriter.hpp"

#include "FakePetscSetup.hpp"

class TestCellBasedDemoTutorial : public AbstractCellBasedTestSuite
{
public:
    void TestVertexBasedMonolayer()
    {
        HoneycombVertexMeshGenerator generator(2, 2);
        boost::shared_ptr<MutableVertexMesh<2,2> > p_mesh = generator.GetMesh();

        std::vector<CellPtr> cells;
        MAKE_PTR(TransitCellProliferativeType, p_transit_type);
        CellsGenerator<UniformG1GenerationalCellCycleModel, 2> cells_generator;
        cells_generator.GenerateBasicRandom(cells, p_mesh->GetNumElements(), p_transit_type);

        VertexBasedCellPopulation<2> cell_population(*p_mesh, cells);

        OffLatticeSimulation<2> simulator(cell_population);
        simulator.SetOutputDirectory("CellBasedDemo1");
        simulator.SetSamplingTimestepMultiple(200);
        simulator.SetEndTime(20.0);

        MAKE_PTR(NagaiHondaForce<2>, p_force);
        simulator.AddForce(p_force);

        MAKE_PTR(SimpleTargetAreaModifier<2>, p_growth_modifier);
        simulator.AddSimulationModifier(p_growth_modifier);

        simulator.Solve();

        TS_ASSERT_EQUALS(cell_population.GetNumRealCells(), 16u);
        TS_ASSERT_DELTA(SimulationTime::Instance()->GetTime(), 20.0, 1e-10);
    }

    void TestNodeBasedMonolayer()
    {
        HoneycombMeshGenerator generator(2, 2); //**Changed**//
        boost::shared_ptr<MutableMesh<2,2> > p_generating_mesh = generator.GetMesh(); //**Changed**//
        NodesOnlyMesh<2> mesh; //**Changed**//
        mesh.ConstructNodesWithoutMesh(*p_generating_mesh, 1.5); //**Changed**//

        std::vector<CellPtr> cells;
        MAKE_PTR(TransitCellProliferativeType, p_transit_type);
        CellsGenerator<UniformG1GenerationalCellCycleModel, 2> cells_generator;
        cells_generator.GenerateBasicRandom(cells, mesh.GetNumNodes(), p_transit_type); //**Changed**//

        NodeBasedCellPopulation<2> cell_population(mesh, cells);//**Changed**//

        OffLatticeSimulation<2> simulator(cell_population);
        simulator.SetOutputDirectory("CellBasedDemo2"); //**Changed**//
        simulator.SetSamplingTimestepMultiple(12); //**Changed**//
        simulator.SetEndTime(20.0);

        MAKE_PTR(RepulsionForce<2>, p_force); //**Changed**//
        simulator.AddForce(p_force);

        MAKE_PTR_ARGS(RandomCellKiller<2>, p_cell_killer, (&cell_population, 0.01)); //**Changed**//
        simulator.AddCellKiller(p_cell_killer);

        simulator.Solve();

        TS_ASSERT_EQUALS(cell_population.GetNumRealCells(), 7u);
        TS_ASSERT_DELTA(SimulationTime::Instance()->GetTime(), 20.0, 1e-10);
    }

    void TestMeshBasedMonolayer()
    {
        HoneycombMeshGenerator generator(2, 2);
        boost::shared_ptr<MutableMesh<2,2> > p_mesh = generator.GetMesh();  //**Changed**//

        std::vector<CellPtr> cells;
        MAKE_PTR(TransitCellProliferativeType, p_transit_type);
        CellsGenerator<UniformG1GenerationalCellCycleModel, 2> cells_generator;
        cells_generator.GenerateBasicRandom(cells, p_mesh->GetNumNodes(), p_transit_type);

        MeshBasedCellPopulation<2> cell_population(*p_mesh, cells); //**Changed**//

        cell_population.AddPopulationWriter<VoronoiDataWriter>();

        OffLatticeSimulation<2> simulator(cell_population);
        simulator.SetOutputDirectory("CellBasedDemo3"); //**Changed**//
        simulator.SetSamplingTimestepMultiple(12);
        simulator.SetEndTime(20.0);

        MAKE_PTR(GeneralisedLinearSpringForce<2>, p_force); //**Changed**//
        simulator.AddForce(p_force);

        simulator.Solve();

        TS_ASSERT_EQUALS(cell_population.GetNumRealCells(), 16u);
        TS_ASSERT_DELTA(SimulationTime::Instance()->GetTime(), 20.0, 1e-10);
    }

    void TestMeshBasedMonolayerWithGhostNodes()
    {
        HoneycombMeshGenerator generator(2, 2, 2); //**Changed**//
        boost::shared_ptr<MutableMesh<2,2> > p_mesh = generator.GetMesh();

        std::vector<unsigned> location_indices = generator.GetCellLocationIndices();//**Changed**//
        std::vector<CellPtr> cells;
        MAKE_PTR(TransitCellProliferativeType, p_transit_type);
        CellsGenerator<TysonNovakCellCycleModel, 2> cells_generator; //**Changed**//
        cells_generator.GenerateBasicRandom(cells, location_indices.size(), p_transit_type); //**Changed**//

        MeshBasedCellPopulationWithGhostNodes<2> cell_population(*p_mesh, cells, location_indices); //**Changed**//

        cell_population.AddPopulationWriter<VoronoiDataWriter>();

        OffLatticeSimulation<2> simulator(cell_population);
        simulator.SetOutputDirectory("CellBasedDemo4"); //**Changed**//
        simulator.SetSamplingTimestepMultiple(12);
        simulator.SetEndTime(2.0); //**Changed**//

        MAKE_PTR(GeneralisedLinearSpringForce<2>, p_force);
        simulator.AddForce(p_force);
        simulator.Solve();

        TS_ASSERT_EQUALS(cell_population.GetNumRealCells(), 16u);
        TS_ASSERT_DELTA(SimulationTime::Instance()->GetTime(), 2.0, 1e-10);
    }

    void TestMeshBasedMonolayerPeriodic()
    {
        CylindricalHoneycombMeshGenerator generator(5, 2, 2); //**Changed**//
        boost::shared_ptr<Cylindrical2dMesh> p_mesh = generator.GetCylindricalMesh(); //**Changed**//

        std::vector<unsigned> location_indices = generator.GetCellLocationIndices();
        std::vector<CellPtr> cells;
        MAKE_PTR(TransitCellProliferativeType, p_transit_type);
        CellsGenerator<UniformG1GenerationalCellCycleModel, 2> cells_generator; //**Changed**//
        cells_generator.GenerateBasicRandom(cells, location_indices.size(), p_transit_type);

        MeshBasedCellPopulationWithGhostNodes<2> cell_population(*p_mesh, cells, location_indices);

        cell_population.AddPopulationWriter<VoronoiDataWriter>();

        OffLatticeSimulation<2> simulator(cell_population);
        simulator.SetOutputDirectory("CellBasedDemo5"); //**Changed**//
        simulator.SetSamplingTimestepMultiple(12);
        simulator.SetEndTime(20.0); //**Changed**//

        MAKE_PTR(GeneralisedLinearSpringForce<2>, p_force);
        simulator.AddForce(p_force);

        simulator.Solve();

        TS_ASSERT_EQUALS(cell_population.GetNumRealCells(), 29u);
        TS_ASSERT_DELTA(SimulationTime::Instance()->GetTime(), 20.0, 1e-10);
    }

    void TestMeshBasedMonolayerPeriodicSolidBottomBoundary()
    {
        CylindricalHoneycombMeshGenerator generator(5, 2, 2);
        boost::shared_ptr<Cylindrical2dMesh> p_mesh = generator.GetCylindricalMesh();

        std::vector<unsigned> location_indices = generator.GetCellLocationIndices();
        std::vector<CellPtr> cells;
        MAKE_PTR(StemCellProliferativeType, p_stem_type);
        CellsGenerator<UniformG1GenerationalCellCycleModel, 2> cells_generator;
        cells_generator.GenerateBasicRandom(cells, location_indices.size(), p_stem_type); //**Changed**//

        MeshBasedCellPopulationWithGhostNodes<2> cell_population(*p_mesh, cells, location_indices);
        cell_population.AddPopulationWriter<VoronoiDataWriter>();

        OffLatticeSimulation<2> simulator(cell_population);
        simulator.SetOutputDirectory("CellBasedDemo6"); //**Changed**//
        simulator.SetSamplingTimestepMultiple(50);
        simulator.SetEndTime(20.0);

        MAKE_PTR(GeneralisedLinearSpringForce<2>, p_force);
        simulator.AddForce(p_force);

        c_vector<double,2> point = zero_vector<double>(2);
        c_vector<double,2> normal = zero_vector<double>(2);
        normal(1) = -1.0;
        MAKE_PTR_ARGS(PlaneBoundaryCondition<2>, p_bc, (&cell_population, point, normal));
        simulator.AddCellPopulationBoundaryCondition(p_bc);

        simulator.Solve();

        TS_ASSERT_EQUALS(cell_population.GetNumRealCells(), 23u);
        TS_ASSERT_DELTA(SimulationTime::Instance()->GetTime(), 20.0, 1e-10);
    }
    void TestPottsBasedMonolayer()
    {
        PottsMeshGenerator<2> generator(20, 2, 4, 20, 2, 4); //**Changed**//
        boost::shared_ptr<PottsMesh<2> > p_mesh = generator.GetMesh(); //**Changed**//

        std::vector<CellPtr> cells;
        MAKE_PTR(TransitCellProliferativeType, p_transit_type);
        CellsGenerator<UniformG1GenerationalCellCycleModel, 2> cells_generator;
        cells_generator.GenerateBasicRandom(cells, p_mesh->GetNumElements(), p_transit_type);

        PottsBasedCellPopulation<2> cell_population(*p_mesh, cells);//**Changed**//
        cell_population.SetTemperature(1.0);

        OnLatticeSimulation<2> simulator(cell_population);//**Changed**//
        simulator.SetOutputDirectory("CellBasedDemo7"); //**Changed**//
        simulator.SetEndTime(20.0);

        MAKE_PTR(VolumeConstraintPottsUpdateRule<2>, p_volume_constraint_update_rule); //**Changed**//
        simulator.AddUpdateRule(p_volume_constraint_update_rule); //**Changed**//
        MAKE_PTR(SurfaceAreaConstraintPottsUpdateRule<2>, p_surface_area_update_rule); //**Changed**//
        simulator.AddUpdateRule(p_surface_area_update_rule); //**Changed**//
        MAKE_PTR(AdhesionPottsUpdateRule<2>, p_adhesion_update_rule); //**Changed**//
        simulator.AddUpdateRule(p_adhesion_update_rule); //**Changed**//

        MAKE_PTR_ARGS(RandomCellKiller<2>, p_cell_killer, (&cell_population, 0.01));
        simulator.AddCellKiller(p_cell_killer);

        simulator.Solve();

        TS_ASSERT_EQUALS(cell_population.GetNumRealCells(), 16u);
        TS_ASSERT_DELTA(SimulationTime::Instance()->GetTime(), 20.0, 1e-10);
    }
};

Cell Based Demo

Cell-based Chaste Code Structure

Cell Based Demo

Examples showing how to create, run and cell-based simulations in Chaste#

Introduction#

The test#

Test 1 - a basic vertex-based simulation#

Test 2 - basic node-based simulation#

Test 3 - basic mesh-based simulation#

Test 4 - basic mesh-based simulation with ghost nodes#

Test 5 - basic periodic mesh-based simulation#

Test 6 - basic periodic mesh-based simulation with obstructions#

Test 7 - basic Potts-based simulation#

Full code#

Examples showing how to create, run and cell-based simulations in Chaste

Introduction

The test

Test 1 - a basic vertex-based simulation

Test 2 - basic node-based simulation

Test 3 - basic mesh-based simulation

Test 4 - basic mesh-based simulation with ghost nodes

Test 5 - basic periodic mesh-based simulation

Test 6 - basic periodic mesh-based simulation with obstructions

Test 7 - basic Potts-based simulation

Full code