This tutorial was generated from the file projects/Harvey2015/test/TestUnitValidationLiteratePaper.hpp at revision r23521. Note that the code is given in full at the bottom of the page.
Validate a simple three cell simulation in parallel and serial (Section 3.1)
On this wiki page we describe in detail some code that is used to compare the results of a simple three cell simulation in parallel and serial. The results should be the same (to machine output precision) as described in Section 3.1.
Use
Both tests in this file are designed to be run twice:
## in serial
scons build=GccOptNative projects/Harvey2015/test/TestUnitValidationLiteratePaper.hpp
## In parallel
scons build=GccOptNative_2 projects/Harvey2015/test/TestUnitValidationLiteratePaper.hpp
After this the positional output may be checked to machine output precision:
cd /tmp/$USER/testoutput
## first test comparison
diff ValidateTwoCells_1_Procs/results_from_time_0/results.viznodes ValidateTwoCells_2_Procs/results_from_time_0/results.viznodes
## second test
diff ValidateTwoCellsOneProc_1_Procs/results_from_time_0/results.viznodes ValidateTwoCellsOneProc_2_Procs/results_from_time_0/results.viznodes
The results.viznode
file show one line per timestep with the x,y,z coordinates of each cell listed in order.
Note that the output is not given to machine precision (only C++ stdio
precision) unless you amend the trunk with setprecision(..)
in
the NodeLocationWriter
(see wiki:PaperTutorials/Harvey2015/ValidateSimulation).
VTK files will contain full machine precision position information, together with process ownership.
paraview --data=/tmp/$USER/testoutput/ValidateTwoCells_2_Procs/results_from_time_0/results.pvd
## View the cells by adding the Glyph filter and rotating the z-axis
Code overview
The first thing to do is to include the necessary header files.
// For any extra output
#include <iostream>
// The testing framework we use
#include <cxxtest/TestSuite.h>
#include "AbstractCellBasedTestSuite.hpp"
// The main implementation used in this publication
#include "NodeBasedCellPopulation.hpp"
// Other Chaste cell_based code
#include "OffLatticeSimulation.hpp"
#include "GeneralisedLinearSpringForce.hpp"
#include "CellsGenerator.hpp"
#include "FixedDurationGenerationBasedCellCycleModel.hpp"
// Required to run the simulation in parallel
#include "PetscSetupAndFinalize.hpp"
The test suite
class TestUnitValidation : public AbstractCellBasedTestSuite
{
public:
First unit test
This test places the three cells at
- z=0.0
- z=0.4
- z=0.6
The natural strip size of the simulation is 0.5 so, when run on 2 or more process, the start configuration is
- z=0.0 on Process 0
- z=0.4 on Process 0
- z=0.6 on Process 1 All other processes are assigned no cells.
The forces between cells are such that the cells are in repulsion
void TestTwoCellOnTwoProcesses() throw (Exception)
{
Here we make define the node positions. When run in parallel this is done on every process.
std::vector<Node<3>* > nodes;
nodes.push_back(new Node<3>(0, false, 0.0, 0.0, 0.0));
nodes.push_back(new Node<3>(0, false, 0.0, 0.0, 0.4));
nodes.push_back(new Node<3>(0, false, 0.0, 0.0, 0.6));
The nodes are used to construct the NodesOnlyMesh
object with a “maximum interaction distance” of 0.5.
This parameter informs strip size for the parallelisation. The bounding geometric region
NodesOnlyMesh<3> mesh;
mesh.ConstructNodesWithoutMesh(nodes, 0.5);
The following is to confirm that when run parallel in the distributed mesh
- process 0 owns the first two nodes
- process 1 owns the 3rd node
- other processes do not take part
std::cout << "Node/cell ownership on process "<<PetscTools::GetMyRank();
for (unsigned i=0; i<3; i++)
{
std::cout<<"\t"<<mesh.rGetInitiallyOwnedNodes()[i];
}
std::cout<<"\n";
Cells are constructed for each of the 3 nodes
std::vector<CellPtr> cells;
CellsGenerator<FixedDurationGenerationBasedCellCycleModel, 3> cells_generator;
cells_generator.GenerateBasicRandom(cells, mesh.GetNumNodes());
A population is constructed. The absolute movement threshold is there to alter us if a cell moves further than expected (since such a cell might move over multiple strips).
NodeBasedCellPopulation<3> population(mesh, cells);
population.SetAbsoluteMovementThreshold(0.5);
Name the output folder based on how many processes we are running.
This folder is to be found relative to CHASTE_TEST_OUTPUT
which by default is
/tmp/$USER/testoutput
std::ostringstream num_procs;
num_procs << PetscTools::GetNumProcs();
std::string output_directory = "ValidateTwoCells_" + num_procs.str() + "_Procs";
A simulation is set up to run for one hour (of simulated time) in 120 discrete steps.
OffLatticeSimulation<3> simulator(population);
simulator.SetOutputDirectory(output_directory);
simulator.SetEndTime(1.0);
simulator.SetDt(1.0/120.0);
Create a force law and pass it to the OffLatticeSimulation
. The cells have been positioned
such that the force will repel them from each other.
MAKE_PTR(GeneralisedLinearSpringForce<3>, p_force);
p_force->SetCutOffLength(0.5);
simulator.AddForce(p_force);
Run the simulation
simulator.Solve();
}
Second unit test
This test places the three cells at
- z=0.0
- z=0.4
- z=0.45
The natural strip size of the simulation is 0.5 so, when run on any number of processes, the start configuration has
- All cells on Process 0
- No cells on Process 1, 2, …
The forces between cells are such that the cells are in repulsion. Cells will migrate during the simulation.
void TestTwoCellOnOneProcesses() throw (Exception)
{
std::vector<Node<3>* > nodes;
nodes.push_back(new Node<3>(0, false, 0.0, 0.0, 0.0));
nodes.push_back(new Node<3>(0, false, 0.0, 0.0, 0.4));
nodes.push_back(new Node<3>(0, false, 0.0, 0.0, 0.45));
NodesOnlyMesh<3> mesh;
mesh.ConstructNodesWithoutMesh(nodes, 0.5);
std::cout << "Node/cell ownership on process "<<PetscTools::GetMyRank();
for (unsigned i=0; i<3; i++)
{
std::cout<<"\t"<<mesh.rGetInitiallyOwnedNodes()[i];
}
std::cout<<"\n";
// Cells
std::vector<CellPtr> cells;
CellsGenerator<FixedDurationGenerationBasedCellCycleModel, 3> cells_generator;
cells_generator.GenerateBasicRandom(cells, mesh.GetNumNodes());
// Population
NodeBasedCellPopulation<3> population(mesh, cells);
population.SetAbsoluteMovementThreshold(0.5);
// Name the output file based on how many procs we are running.
std::ostringstream num_procs;
num_procs << PetscTools::GetNumProcs();
std::string output_directory = "ValidateTwoCellsOneProc_" + num_procs.str() + "_Procs";
// Simulation
OffLatticeSimulation<3> simulator(population);
simulator.SetOutputDirectory(output_directory);
simulator.SetEndTime(1.0);
simulator.SetDt(1.0/120.0);
// Create a force law and pass it to the OffLatticeSimulation
MAKE_PTR(GeneralisedLinearSpringForce<3>, p_force);
p_force->SetCutOffLength(0.5);
simulator.AddForce(p_force);
// Run
simulator.Solve();
}
};
Code
The full code is given below
File name TestUnitValidationLiteratePaper.hpp
// For any extra output
#include <iostream>
// The testing framework we use
#include <cxxtest/TestSuite.h>
#include "AbstractCellBasedTestSuite.hpp"
// The main implementation used in this publication
#include "NodeBasedCellPopulation.hpp"
// Other Chaste cell_based code
#include "OffLatticeSimulation.hpp"
#include "GeneralisedLinearSpringForce.hpp"
#include "CellsGenerator.hpp"
#include "FixedDurationGenerationBasedCellCycleModel.hpp"
// Required to run the simulation in parallel
#include "PetscSetupAndFinalize.hpp"
class TestUnitValidation : public AbstractCellBasedTestSuite
{
public:
void TestTwoCellOnTwoProcesses() throw (Exception)
{
std::vector<Node<3>* > nodes;
nodes.push_back(new Node<3>(0, false, 0.0, 0.0, 0.0));
nodes.push_back(new Node<3>(0, false, 0.0, 0.0, 0.4));
nodes.push_back(new Node<3>(0, false, 0.0, 0.0, 0.6));
NodesOnlyMesh<3> mesh;
mesh.ConstructNodesWithoutMesh(nodes, 0.5);
std::cout << "Node/cell ownership on process "<<PetscTools::GetMyRank();
for (unsigned i=0; i<3; i++)
{
std::cout<<"\t"<<mesh.rGetInitiallyOwnedNodes()[i];
}
std::cout<<"\n";
std::vector<CellPtr> cells;
CellsGenerator<FixedDurationGenerationBasedCellCycleModel, 3> cells_generator;
cells_generator.GenerateBasicRandom(cells, mesh.GetNumNodes());
NodeBasedCellPopulation<3> population(mesh, cells);
population.SetAbsoluteMovementThreshold(0.5);
std::ostringstream num_procs;
num_procs << PetscTools::GetNumProcs();
std::string output_directory = "ValidateTwoCells_" + num_procs.str() + "_Procs";
OffLatticeSimulation<3> simulator(population);
simulator.SetOutputDirectory(output_directory);
simulator.SetEndTime(1.0);
simulator.SetDt(1.0/120.0);
MAKE_PTR(GeneralisedLinearSpringForce<3>, p_force);
p_force->SetCutOffLength(0.5);
simulator.AddForce(p_force);
simulator.Solve();
}
void TestTwoCellOnOneProcesses() throw (Exception)
{
std::vector<Node<3>* > nodes;
nodes.push_back(new Node<3>(0, false, 0.0, 0.0, 0.0));
nodes.push_back(new Node<3>(0, false, 0.0, 0.0, 0.4));
nodes.push_back(new Node<3>(0, false, 0.0, 0.0, 0.45));
NodesOnlyMesh<3> mesh;
mesh.ConstructNodesWithoutMesh(nodes, 0.5);
std::cout << "Node/cell ownership on process "<<PetscTools::GetMyRank();
for (unsigned i=0; i<3; i++)
{
std::cout<<"\t"<<mesh.rGetInitiallyOwnedNodes()[i];
}
std::cout<<"\n";
// Cells
std::vector<CellPtr> cells;
CellsGenerator<FixedDurationGenerationBasedCellCycleModel, 3> cells_generator;
cells_generator.GenerateBasicRandom(cells, mesh.GetNumNodes());
// Population
NodeBasedCellPopulation<3> population(mesh, cells);
population.SetAbsoluteMovementThreshold(0.5);
// Name the output file based on how many procs we are running.
std::ostringstream num_procs;
num_procs << PetscTools::GetNumProcs();
std::string output_directory = "ValidateTwoCellsOneProc_" + num_procs.str() + "_Procs";
// Simulation
OffLatticeSimulation<3> simulator(population);
simulator.SetOutputDirectory(output_directory);
simulator.SetEndTime(1.0);
simulator.SetDt(1.0/120.0);
// Create a force law and pass it to the OffLatticeSimulation
MAKE_PTR(GeneralisedLinearSpringForce<3>, p_force);
p_force->SetCutOffLength(0.5);
simulator.AddForce(p_force);
// Run
simulator.Solve();
}
};