CryptSimulation2d.cpp
00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029 #include "CryptSimulation2d.hpp"
00030 #include "WntConcentration.hpp"
00031 #include "VanLeeuwen2009WntSwatCellCycleModelHypothesisOne.hpp"
00032 #include "VanLeeuwen2009WntSwatCellCycleModelHypothesisTwo.hpp"
00033 #include "SmartPointers.hpp"
00034
00035 CryptSimulation2d::CryptSimulation2d(AbstractCellPopulation<2>& rCellPopulation,
00036 bool deleteCellPopulationInDestructor,
00037 bool initialiseCells)
00038 : OffLatticeSimulation<2>(rCellPopulation,
00039 deleteCellPopulationInDestructor,
00040 initialiseCells),
00041 mWriteBetaCatenin(false),
00042 mUsingMeshBasedCellPopulation(false)
00043 {
00044
00045
00046
00047
00048 if ( (dynamic_cast<VertexBasedCellPopulation<2>*>(&rCellPopulation) == NULL)
00049 &&(dynamic_cast<MeshBasedCellPopulation<2>*>(&rCellPopulation) == NULL) )
00050 {
00051 EXCEPTION("CryptSimulation2d is to be used with MeshBasedCellPopulation or VertexBasedCellPopulation (or subclasses) only");
00052 }
00053
00054 if (dynamic_cast<MeshBasedCellPopulation<2>*>(&mrCellPopulation))
00055 {
00056 mUsingMeshBasedCellPopulation = true;
00057 }
00058
00059
00060
00061
00062
00063
00064
00065 if (dynamic_cast<AbstractVanLeeuwen2009WntSwatCellCycleModel*>(mrCellPopulation.Begin()->GetCellCycleModel()))
00066 {
00067 mWriteBetaCatenin = true;
00068 }
00069
00070 if (!mDeleteCellPopulationInDestructor)
00071 {
00072
00073 MAKE_PTR_ARGS(CryptSimulationBoundaryCondition<2>, p_bc, (&rCellPopulation));
00074 AddCellPopulationBoundaryCondition(p_bc);
00075 }
00076 }
00077
00078 CryptSimulation2d::~CryptSimulation2d()
00079 {
00080 }
00081
00082 c_vector<double, 2> CryptSimulation2d::CalculateCellDivisionVector(CellPtr pParentCell)
00083 {
00084 if (mUsingMeshBasedCellPopulation)
00085 {
00086
00087 c_vector<double, 2> parent_coords = mrCellPopulation.GetLocationOfCellCentre(pParentCell);
00088 c_vector<double, 2> daughter_coords;
00089
00090
00091 double separation =
00092 static_cast<MeshBasedCellPopulation<2>*>(&mrCellPopulation)->GetMeinekeDivisionSeparation();
00093
00094
00095 c_vector<double, 2> random_vector;
00096
00097
00098
00099
00100
00101
00102 double random_angle = RandomNumberGenerator::Instance()->ranf();
00103 random_angle *= 2.0*M_PI;
00104
00105 random_vector(0) = 0.5*separation*cos(random_angle);
00106 random_vector(1) = 0.5*separation*sin(random_angle);
00107
00108 c_vector<double, 2> proposed_new_parent_coords = parent_coords - random_vector;
00109 c_vector<double, 2> proposed_new_daughter_coords = parent_coords + random_vector;
00110
00111 if ((proposed_new_parent_coords(1) >= 0.0) && (proposed_new_daughter_coords(1) >= 0.0))
00112 {
00113
00114 parent_coords = proposed_new_parent_coords;
00115 daughter_coords = proposed_new_daughter_coords;
00116 }
00117 else
00118 {
00119 proposed_new_daughter_coords = parent_coords + 2.0*random_vector;
00120 while (proposed_new_daughter_coords(1) < 0.0)
00121 {
00122 random_angle = RandomNumberGenerator::Instance()->ranf();
00123 random_angle *= 2.0*M_PI;
00124
00125 random_vector(0) = separation*cos(random_angle);
00126 random_vector(1) = separation*sin(random_angle);
00127 proposed_new_daughter_coords = parent_coords + random_vector;
00128 }
00129 daughter_coords = proposed_new_daughter_coords;
00130 }
00131
00132 assert(daughter_coords(1) >= 0.0);
00133 assert(parent_coords(1) >= 0.0);
00134
00135
00136 ChastePoint<2> parent_coords_point(parent_coords);
00137
00138 unsigned node_index = mrCellPopulation.GetLocationIndexUsingCell(pParentCell);
00139 mrCellPopulation.SetNode(node_index, parent_coords_point);
00140
00141 return daughter_coords;
00142 }
00143 else
00144 {
00145 c_vector<double, 2> axis_of_division = zero_vector<double>(2);
00146
00147
00148 bool is_wnt_included = WntConcentration<2>::Instance()->IsWntSetUp();
00149 if (!is_wnt_included)
00150 {
00151 WntConcentration<2>::Destroy();
00152 if (pParentCell->GetCellCycleModel()->GetCellProliferativeType() == STEM)
00153 {
00154 axis_of_division(0) = 1.0;
00155 axis_of_division(1) = 0.0;
00156 }
00157 }
00158 return axis_of_division;
00159 }
00160 }
00161
00162 void CryptSimulation2d::SetupWriteBetaCatenin()
00163 {
00164 OutputFileHandler output_file_handler(this->mSimulationOutputDirectory + "/", false);
00165 mVizBetaCateninResultsFile = output_file_handler.OpenOutputFile("results.vizbetacatenin");
00166 *mpVizSetupFile << "BetaCatenin\n";
00167 }
00168
00169 void CryptSimulation2d::WriteBetaCatenin(double time)
00170 {
00171 *mVizBetaCateninResultsFile << time << "\t";
00172
00173 for (AbstractCellPopulation<2>::Iterator cell_iter = mrCellPopulation.Begin();
00174 cell_iter != mrCellPopulation.End();
00175 ++cell_iter)
00176 {
00177 unsigned global_index = mrCellPopulation.GetLocationIndexUsingCell(*cell_iter);
00178 double x = mrCellPopulation.GetLocationOfCellCentre(*cell_iter)[0];
00179 double y = mrCellPopulation.GetLocationOfCellCentre(*cell_iter)[1];
00180
00181
00182 assert(mWriteBetaCatenin);
00183
00184 AbstractVanLeeuwen2009WntSwatCellCycleModel* p_model = dynamic_cast<AbstractVanLeeuwen2009WntSwatCellCycleModel*>(cell_iter->GetCellCycleModel());
00185 double b_cat_membrane = p_model->GetMembraneBoundBetaCateninLevel();
00186 double b_cat_cytoplasm = p_model->GetCytoplasmicBetaCateninLevel();
00187 double b_cat_nuclear = p_model->GetNuclearBetaCateninLevel();
00188
00189 *mVizBetaCateninResultsFile << global_index << " " << x << " " << y << " " << b_cat_membrane << " " << b_cat_cytoplasm << " " << b_cat_nuclear << " ";
00190 }
00191
00192 *mVizBetaCateninResultsFile << "\n";
00193 }
00194
00195 void CryptSimulation2d::SetupSolve()
00196 {
00197
00198 OffLatticeSimulation<2>::SetupSolve();
00199
00200
00201
00202
00203
00204
00205 bool any_cells_present = (mrCellPopulation.Begin() != mrCellPopulation.End());
00206 if (any_cells_present && mWriteBetaCatenin)
00207 {
00208 SetupWriteBetaCatenin();
00209 double current_time = SimulationTime::Instance()->GetTime();
00210 WriteBetaCatenin(current_time);
00211 }
00212 }
00213
00214 void CryptSimulation2d::PostSolve()
00215 {
00216
00217 OffLatticeSimulation<2>::PostSolve();
00218
00219 SimulationTime* p_time = SimulationTime::Instance();
00220
00221 if ((p_time->GetTimeStepsElapsed()+1)%mSamplingTimestepMultiple == 0)
00222 {
00223
00224
00225
00226
00227
00228 bool any_cells_present = (mrCellPopulation.Begin() != mrCellPopulation.End());
00229 if (any_cells_present && mWriteBetaCatenin)
00230 {
00231 double time_next_step = p_time->GetTime() + p_time->GetTimeStep();
00232 WriteBetaCatenin(time_next_step);
00233 }
00234 }
00235 }
00236
00237 void CryptSimulation2d::AfterSolve()
00238 {
00239
00240 OffLatticeSimulation<2>::AfterSolve();
00241
00242
00243
00244
00245
00246 bool any_cells_present = (mrCellPopulation.Begin() != mrCellPopulation.End());
00247 if (any_cells_present && mWriteBetaCatenin)
00248 {
00249 mVizBetaCateninResultsFile->close();
00250 }
00251 }
00252
00253 void CryptSimulation2d::UseJiggledBottomCells()
00254 {
00255
00256 boost::static_pointer_cast<CryptSimulationBoundaryCondition<2> >(mBoundaryConditions[0])->SetUseJiggledBottomCells(true);
00257 }
00258
00259 void CryptSimulation2d::SetBottomCellAncestors()
00260 {
00261
00262
00263
00264
00265
00266 double threshold_height = 1.0;
00267 if (mUsingMeshBasedCellPopulation)
00268 {
00269 threshold_height = 0.5;
00270 }
00271
00272 unsigned index = 0;
00273 for (AbstractCellPopulation<2>::Iterator cell_iter = mrCellPopulation.Begin();
00274 cell_iter != mrCellPopulation.End();
00275 ++cell_iter)
00276 {
00277 if (mrCellPopulation.GetLocationOfCellCentre(*cell_iter)[1] < threshold_height)
00278 {
00279 cell_iter->SetAncestor(index++);
00280 }
00281 }
00282 }
00283
00284 void CryptSimulation2d::OutputSimulationParameters(out_stream& rParamsFile)
00285 {
00286 double width = mrCellPopulation.GetWidth(0);
00287 bool use_jiggled_bottom_cells = boost::static_pointer_cast<CryptSimulationBoundaryCondition<2> >(mBoundaryConditions[0])->GetUseJiggledBottomCells();
00288
00289 *rParamsFile << "\t\t<CryptCircumference>" << width << "</CryptCircumference>\n";
00290 *rParamsFile << "\t\t<UseJiggledBottomCells>" << use_jiggled_bottom_cells << "</UseJiggledBottomCells>\n";
00291
00292
00293 OffLatticeSimulation<2>::OutputSimulationParameters(rParamsFile);
00294 }
00295
00296
00297 #include "SerializationExportWrapperForCpp.hpp"
00298 CHASTE_CLASS_EXPORT(CryptSimulation2d)