Chaste Commit::333233612e3dcc941d260418acd7b5ccdc30390e
AbstractCellBasedSimulation.cpp
1/*
2
3Copyright (c) 2005-2024, University of Oxford.
4All rights reserved.
5
6University of Oxford means the Chancellor, Masters and Scholars of the
7University of Oxford, having an administrative office at Wellington
8Square, Oxford OX1 2JD, UK.
9
10This file is part of Chaste.
11
12Redistribution and use in source and binary forms, with or without
13modification, are permitted provided that the following conditions are met:
14 * Redistributions of source code must retain the above copyright notice,
15 this list of conditions and the following disclaimer.
16 * Redistributions in binary form must reproduce the above copyright notice,
17 this list of conditions and the following disclaimer in the documentation
18 and/or other materials provided with the distribution.
19 * Neither the name of the University of Oxford nor the names of its
20 contributors may be used to endorse or promote products derived from this
21 software without specific prior written permission.
22
23THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
24AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
25IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
26ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE
27LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
28CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE
29GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
30HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
31LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT
32OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
33
34*/
35
36#include <cmath>
37#include <iostream>
38#include <fstream>
39#include <set>
40#include "AbstractCellBasedSimulation.hpp"
41#include "CellBasedEventHandler.hpp"
42#include "LogFile.hpp"
43#include "ExecutableSupport.hpp"
44#include "AbstractPdeModifier.hpp"
45#include "CellDivisionLocationsWriter.hpp"
46#include "CellRemovalLocationsWriter.hpp"
47
48template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
50 bool deleteCellPopulationInDestructor,
51 bool initialiseCells)
52 : mDt(DOUBLE_UNSET),
53 mEndTime(DOUBLE_UNSET), // hours - this is set later on
54 mrCellPopulation(rCellPopulation),
55 mDeleteCellPopulationInDestructor(deleteCellPopulationInDestructor),
56 mInitialiseCells(initialiseCells),
57 mNoBirth(false),
58 mUpdateCellPopulation(true),
59 mOutputDirectory(""),
60 mSimulationOutputDirectory(mOutputDirectory),
61 mNumBirths(0),
62 mNumDeaths(0),
63 mOutputDivisionLocations(false),
64 mOutputCellVelocities(false),
65 mSamplingTimestepMultiple(1),
66 mUpdatingTimestepMultiple(1)
67{
68 // Set a random seed of 0 if it wasn't specified earlier
70
72 {
73 mrCellPopulation.InitialiseCells();
74 }
75
76 /*
77 * Specify a default time step to use, which may depend on the cell population type.
78 * Note that the time step may be reset using SetDt().
79 */
80 mDt = rCellPopulation.GetDefaultTimeStep();
81}
82
83template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
85{
86 if (mDeleteCellPopulationInDestructor)
87 {
88 delete &mrCellPopulation;
89 }
90}
91
92template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
94{
95 if (mNoBirth)
96 {
97 return 0;
98 }
99
100 unsigned num_births_this_step = 0;
101
102 // Iterate over all cells, seeing if each one can be divided
103 for (typename AbstractCellPopulation<ELEMENT_DIM,SPACE_DIM>::Iterator cell_iter = mrCellPopulation.Begin();
104 cell_iter != mrCellPopulation.End();
105 ++cell_iter)
106 {
107 // Check if this cell is ready to divide
108 double cell_age = cell_iter->GetAge();
109 if (cell_age > 0.0)
110 {
111 if (cell_iter->ReadyToDivide())
112 {
113 // Check if there is room into which the cell may divide
114 if (mrCellPopulation.IsRoomToDivide(*cell_iter))
115 {
116 // Store parent ID for output if required
117 unsigned parent_cell_id = cell_iter->GetCellId();
118
119 // Create a new cell
120 CellPtr p_new_cell = cell_iter->Divide();
121
128 if (mrCellPopulation.template HasWriter<CellDivisionLocationsWriter>())
129 {
130 c_vector<double, SPACE_DIM> cell_location = mrCellPopulation.GetLocationOfCellCentre(*cell_iter);
131
132 std::stringstream division_info;
133 division_info << SimulationTime::Instance()->GetTime() << "\t";
134 for (unsigned i = 0; i < SPACE_DIM; i++)
135 {
136 division_info << cell_location[i] << "\t";
137 }
138 division_info << "\t" << cell_age << "\t" << parent_cell_id << "\t" << cell_iter->GetCellId() << "\t" << p_new_cell->GetCellId() << "\t";
139
140 mrCellPopulation.AddDivisionInformation(division_info.str());
141 }
142
143 // Add the new cell to the cell population
144 mrCellPopulation.AddCell(p_new_cell, *cell_iter);
145
146 // Update counter
147 num_births_this_step++;
148 }
149 }
150 }
151 }
152 return num_births_this_step;
153}
154
155template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
157{
158 unsigned num_deaths_this_step = 0;
159
160 /*
161 * This labels cells as dead or apoptosing. It does not actually remove the cells,
162 * mrCellPopulation.RemoveDeadCells() needs to be called for this.
163 */
164 for (typename std::vector<boost::shared_ptr<AbstractCellKiller<SPACE_DIM> > >::iterator killer_iter = mCellKillers.begin();
165 killer_iter != mCellKillers.end();
166 ++killer_iter)
167 {
168 (*killer_iter)->CheckAndLabelCellsForApoptosisOrDeath();
169 }
170
171 num_deaths_this_step += mrCellPopulation.RemoveDeadCells();
172
173 return num_deaths_this_step;
174}
175
176template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
178{
179 assert(dt > 0);
180 mDt = dt;
181}
182
183template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
188
189template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
194
195template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
200
201template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
203{
204 assert(endTime > 0);
205 mEndTime = endTime;
206}
207
208template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
210{
211 mOutputDirectory = outputDirectory;
212 mSimulationOutputDirectory = mOutputDirectory;
213}
214
215template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
220
221template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
223{
224 assert(samplingTimestepMultiple > 0);
225 mSamplingTimestepMultiple = samplingTimestepMultiple;
226}
227
228template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
230{
231 assert(updatingTimestepMultiple > 0);
232 mUpdatingTimestepMultiple = updatingTimestepMultiple;
233}
235template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
240
241template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
246
247template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
249{
250 mUpdateCellPopulation = updateCellPopulation;
251}
252
253template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
258
259template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
261{
262 mNoBirth = noBirth;
263}
264
265template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
267{
268 mCellKillers.push_back(pCellKiller);
269}
270
271template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
276
277template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
279{
280 mSimulationModifiers.push_back(pSimulationModifier);
281}
283template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
284std::vector<boost::shared_ptr<AbstractCellBasedSimulationModifier<ELEMENT_DIM, SPACE_DIM> > >* AbstractCellBasedSimulation<ELEMENT_DIM,SPACE_DIM>::GetSimulationModifiers()
285{
286 return &mSimulationModifiers;
287}
288
289template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
292 mTopologyUpdateSimulationModifiers.push_back(pSimulationModifier);
293}
294
295template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
296std::vector<boost::shared_ptr<AbstractCellBasedSimulationModifier<ELEMENT_DIM, SPACE_DIM> > >* AbstractCellBasedSimulation<ELEMENT_DIM,SPACE_DIM>::GetTopologyUpdateSimulationModifiers()
297{
298 return &mTopologyUpdateSimulationModifiers;
300
301template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
302std::vector<double> AbstractCellBasedSimulation<ELEMENT_DIM,SPACE_DIM>::GetNodeLocation(const unsigned& rNodeIndex)
303{
304 std::vector<double> location;
305 for (unsigned i=0; i<SPACE_DIM; i++)
306 {
307 location.push_back(mrCellPopulation.GetNode(rNodeIndex)->rGetLocation()[i]);
308 }
309 return location;
310}
311
312template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
315 CellBasedEventHandler::BeginEvent(CellBasedEventHandler::EVERYTHING);
316 CellBasedEventHandler::BeginEvent(CellBasedEventHandler::SETUP);
317
318 // Set up the simulation time
320 double current_time = p_simulation_time->GetTime();
321
322 assert(mDt != DOUBLE_UNSET); //Subclass constructors take care of this
323
324 if (mEndTime == DOUBLE_UNSET)
325 {
326 EXCEPTION("SetEndTime has not yet been called.");
327 }
328
329 /*
330 * Note that mDt is used here for "ideal time step". If this step doesn't divide the time remaining
331 * then a *different* time step will be taken by the time-stepper. The real time-step (used in the
332 * SimulationTime singleton) is currently not available to this class.
334 * \todo Should we over-write the value of mDt, or change this behaviour? (see #2159)
335 */
336 unsigned num_time_steps = (unsigned) ((mEndTime-current_time)/mDt+0.5);
337 if (current_time > 0) // use the reset function if necessary
338 {
339 p_simulation_time->ResetEndTimeAndNumberOfTimeSteps(mEndTime, num_time_steps);
341 else
342 {
343 if (p_simulation_time->IsEndTimeAndNumberOfTimeStepsSetUp())
344 {
345 EXCEPTION("End time and number of timesteps already setup. You should not use SimulationTime::SetEndTimeAndNumberOfTimeSteps in cell-based tests.");
346 }
347 else
349 p_simulation_time->SetEndTimeAndNumberOfTimeSteps(mEndTime, num_time_steps);
350 }
351 }
352
353 if (mOutputDirectory == "")
354 {
355 EXCEPTION("OutputDirectory not set");
357
358 double time_now = p_simulation_time->GetTime();
359 std::ostringstream time_string;
360 time_string << time_now;
361
362 std::string results_directory = mOutputDirectory +"/results_from_time_" + time_string.str();
363 mSimulationOutputDirectory = results_directory;
364
365 // Set up simulation
366
367 // Create output files for the visualizer
368 OutputFileHandler output_file_handler(results_directory + "/", true);
369
370 if (mOutputDivisionLocations)
371 {
372 mrCellPopulation.template AddCellPopulationEventWriter<CellDivisionLocationsWriter>();
373 mrCellPopulation.template AddCellPopulationEventWriter<CellRemovalLocationsWriter>();
374 }
375 if (mOutputCellVelocities)
376 {
377 OutputFileHandler output_file_handler2(this->mSimulationOutputDirectory + "/", false);
378 mpCellVelocitiesFile = output_file_handler2.OpenOutputFile("cellvelocities.dat");
379 }
380
381 mrCellPopulation.OpenWritersFiles(output_file_handler);
382
385 mpVizSetupFile = output_file_handler.OpenOutputFile("results.vizsetup");
386
387 for (typename std::vector<boost::shared_ptr<AbstractCellBasedSimulationModifier<ELEMENT_DIM, SPACE_DIM> > >::iterator iter = mSimulationModifiers.begin();
388 iter != mSimulationModifiers.end();
389 ++iter)
390 {
391 if (boost::dynamic_pointer_cast<AbstractPdeModifier<SPACE_DIM> >(*iter))
392 {
393 *this->mpVizSetupFile << "PDE \n";
394 }
395 }
397
398 this->mrCellPopulation.SimulationSetupHook(this);
399
400 SetupSolve();
402 // Call SetupSolve() on each modifier
403 for (typename std::vector<boost::shared_ptr<AbstractCellBasedSimulationModifier<ELEMENT_DIM, SPACE_DIM> > >::iterator iter = mSimulationModifiers.begin();
404 iter != mSimulationModifiers.end();
405 ++iter)
406 {
407 (*iter)->SetupSolve(this->mrCellPopulation,this->mSimulationOutputDirectory);
408 }
410 // Call SetupSolve() on each topology update modifier
411 for (typename std::vector<boost::shared_ptr<AbstractCellBasedSimulationModifier<ELEMENT_DIM, SPACE_DIM> > >::iterator iter = mTopologyUpdateSimulationModifiers.begin();
412 iter != mTopologyUpdateSimulationModifiers.end();
413 ++iter)
415 (*iter)->SetupSolve(this->mrCellPopulation,this->mSimulationOutputDirectory);
416 }
417
418 /*
419 * Age the cells to the correct time. Note that cells are created with
420 * negative birth times so that some are initially almost ready to divide.
421 */
422 LOG(1, "Setting up cells...");
423 for (typename AbstractCellPopulation<ELEMENT_DIM,SPACE_DIM>::Iterator cell_iter = mrCellPopulation.Begin();
424 cell_iter != mrCellPopulation.End();
425 ++cell_iter)
426 {
427 /*
428 * We don't use the result; this call is just to force the cells to age
429 * to the current time running their cell-cycle models to get there.
430 */
431 cell_iter->ReadyToDivide();
432 }
433 LOG(1, "\tdone\n");
434
435 // Write initial conditions to file for the visualizer
436 WriteVisualizerSetupFile();
437
439 {
440 *mpVizSetupFile << std::flush;
441 }
442
443 mrCellPopulation.WriteResultsToFiles(results_directory+"/");
444
445 OutputSimulationSetup();
446 CellBasedEventHandler::EndEvent(CellBasedEventHandler::SETUP);
447
448 // Enter main time loop
449 while (!( p_simulation_time->IsFinished() || StoppingEventHasOccurred() ) )
450 {
451 LOG(1, "--TIME = " << p_simulation_time->GetTime() << "\n");
452
453
454 if (p_simulation_time->GetTimeStepsElapsed()%mUpdatingTimestepMultiple == 0)
455 {
456 // This function calls DoCellRemoval(), DoCellBirth() and CellPopulation::Update()
457 UpdateCellPopulation();
458 }
459
460 /*
461 * Store whether we are sampling results at the current timestep note we have +1 here as we caculate this
462 * before the timestep is increaced.
463 */
465 bool at_sampling_timestep = ((p_time->GetTimeStepsElapsed()+1)%this->mSamplingTimestepMultiple == 0);
466
467 /*
468 * If required, store the current locations of cell centres. Note that we need to
469 * use a std::map between cells and locations, rather than (say) a std::vector with
470 * location indices corresponding to cells, since once we call UpdateCellLocations()
471 * the location index of each cell may change. This is especially true in the case
472 * of a CaBasedCellPopulation.
473 */
474 std::map<CellPtr, c_vector<double, SPACE_DIM> > old_cell_locations;
475 if (mOutputCellVelocities && at_sampling_timestep)
476 {
477 for (typename AbstractCellPopulation<ELEMENT_DIM, SPACE_DIM>::Iterator cell_iter = mrCellPopulation.Begin();
478 cell_iter != mrCellPopulation.End();
479 ++cell_iter)
481 old_cell_locations[*cell_iter] = mrCellPopulation.GetLocationOfCellCentre(*cell_iter);
482 }
483 }
484
485 /* Call UpdateAtEndOfTimeStep() on each topology update modifier */
486 CellBasedEventHandler::BeginEvent(CellBasedEventHandler::UPDATESIMULATION);
487 for (typename std::vector<boost::shared_ptr<AbstractCellBasedSimulationModifier<ELEMENT_DIM, SPACE_DIM> > >::iterator iter = mTopologyUpdateSimulationModifiers.begin();
488 iter != mTopologyUpdateSimulationModifiers.end();
489 ++iter)
490 {
491 (*iter)->UpdateAtEndOfTimeStep(this->mrCellPopulation);
493 CellBasedEventHandler::EndEvent(CellBasedEventHandler::UPDATESIMULATION);
494
495 // Update cell locations and topology
496 UpdateCellLocationsAndTopology();
498 // Now write cell velocities to file if required
499 if (mOutputCellVelocities && at_sampling_timestep)
500 {
501 // Offset as doing this before we increase time by mDt
502 *mpCellVelocitiesFile << p_time->GetTime() + mDt<< "\t";
503
504 for (typename AbstractCellPopulation<ELEMENT_DIM, SPACE_DIM>::Iterator cell_iter = mrCellPopulation.Begin();
505 cell_iter != mrCellPopulation.End();
506 ++cell_iter)
507 {
508 unsigned index = mrCellPopulation.GetLocationIndexUsingCell(*cell_iter);
509 const c_vector<double,SPACE_DIM>& position = mrCellPopulation.GetLocationOfCellCentre(*cell_iter);
510
511 c_vector<double, SPACE_DIM> velocity; // Two lines for profile build
512 velocity = (position - old_cell_locations[*cell_iter])/mDt;
513
514 *mpCellVelocitiesFile << index << " ";
515 for (unsigned i=0; i<SPACE_DIM; i++)
516 {
517 *mpCellVelocitiesFile << position[i] << " ";
518 }
519
520 for (unsigned i=0; i<SPACE_DIM; i++)
521 {
522 *mpCellVelocitiesFile << velocity[i] << " ";
523 }
524 }
525 *mpCellVelocitiesFile << "\n";
526 }
527
528 // Update the assignment of cells to processes.
529 mrCellPopulation.UpdateCellProcessLocation();
530
531 // Increment simulation time here, so results files look sensible
532 p_simulation_time->IncrementTimeOneStep();
533
534 /* Call UpdateAtEndOfTimeStep() on each modifier */
535 CellBasedEventHandler::BeginEvent(CellBasedEventHandler::UPDATESIMULATION);
536 for (typename std::vector<boost::shared_ptr<AbstractCellBasedSimulationModifier<ELEMENT_DIM, SPACE_DIM> > >::iterator iter = mSimulationModifiers.begin();
537 iter != mSimulationModifiers.end();
538 ++iter)
539 {
540 (*iter)->UpdateAtEndOfTimeStep(this->mrCellPopulation);
541 }
542 CellBasedEventHandler::EndEvent(CellBasedEventHandler::UPDATESIMULATION);
543
544 // Output current results to file
545 CellBasedEventHandler::BeginEvent(CellBasedEventHandler::OUTPUT);
546 if (p_simulation_time->GetTimeStepsElapsed()%mSamplingTimestepMultiple == 0)// should be at_sampling_timestep !
547 {
548 mrCellPopulation.WriteResultsToFiles(results_directory+"/");
549
550 // Call UpdateAtEndOfOutputTimeStep() on each modifier
551 for (typename std::vector<boost::shared_ptr<AbstractCellBasedSimulationModifier<ELEMENT_DIM, SPACE_DIM> > >::iterator iter = mSimulationModifiers.begin();
552 iter != mSimulationModifiers.end();
553 ++iter)
554 {
555 (*iter)->UpdateAtEndOfOutputTimeStep(this->mrCellPopulation);
556 }
557 }
558 CellBasedEventHandler::EndEvent(CellBasedEventHandler::OUTPUT);
559 }
560
561 LOG(1, "--END TIME = " << p_simulation_time->GetTime() << "\n");
562
563 /*
564 * Carry out a final update so that cell population is coherent with new cell positions.
565 * Note that cell birth and death still need to be checked because they may be spatially
566 * dependent.
567 */
568 UpdateCellPopulation();
569
570 CellBasedEventHandler::BeginEvent(CellBasedEventHandler::UPDATESIMULATION);
571 // Call UpdateAtEndOfSolve(), on each modifier
572 for (typename std::vector<boost::shared_ptr<AbstractCellBasedSimulationModifier<ELEMENT_DIM, SPACE_DIM> > >::iterator iter = mSimulationModifiers.begin();
573 iter != mSimulationModifiers.end();
574 ++iter)
575 {
576 (*iter)->UpdateAtEndOfSolve(this->mrCellPopulation);
577 }
578 CellBasedEventHandler::EndEvent(CellBasedEventHandler::UPDATESIMULATION);
579
580 CellBasedEventHandler::BeginEvent(CellBasedEventHandler::OUTPUT);
581
582 mrCellPopulation.CloseWritersFiles();
583
584 if (mOutputCellVelocities)
585 {
586 mpCellVelocitiesFile->close();
587 }
588
590 {
591 *mpVizSetupFile << "Complete\n";
592 mpVizSetupFile->close();
593 }
594
595 CellBasedEventHandler::EndEvent(CellBasedEventHandler::OUTPUT);
596 CellBasedEventHandler::EndEvent(CellBasedEventHandler::EVERYTHING);
597}
598
599template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
604
605template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
607{
608 // Remove dead cells
609 CellBasedEventHandler::BeginEvent(CellBasedEventHandler::DEATH);
610 unsigned deaths_this_step = DoCellRemoval();
611 mNumDeaths += deaths_this_step;
612 LOG(1, "\tNum deaths = " << mNumDeaths << "\n");
613 CellBasedEventHandler::EndEvent(CellBasedEventHandler::DEATH);
614
615 // Divide cells
616 CellBasedEventHandler::BeginEvent(CellBasedEventHandler::BIRTH);
617 unsigned births_this_step = DoCellBirth();
618 mNumBirths += births_this_step;
619 LOG(1, "\tNum births = " << mNumBirths << "\n");
620 CellBasedEventHandler::EndEvent(CellBasedEventHandler::BIRTH);
621
622 // This allows NodeBasedCellPopulation::Update() to do the minimum amount of work
623 bool births_or_death_occurred = ((births_this_step>0) || (deaths_this_step>0));
624
625 // Update topology of cell population
626 CellBasedEventHandler::BeginEvent(CellBasedEventHandler::UPDATECELLPOPULATION);
628 if (mUpdateCellPopulation && (p_time->GetTimeStepsElapsed() % mUpdatingTimestepMultiple == 0) )
629 {
630 LOG(1, "\tUpdating cell population...");
631 mrCellPopulation.Update(births_or_death_occurred);
632 LOG(1, "\tdone.\n");
633 }
634 else if (births_or_death_occurred)
635 {
636 if (!mUpdateCellPopulation)
637 {
638 EXCEPTION("CellPopulation has had births or deaths but mUpdateCellPopulation is set to false, please set it to true.");
639 }
640 else if ((p_time->GetTimeStepsElapsed() % mUpdatingTimestepMultiple != 0))
641 {
642 EXCEPTION("CellPopulation has had births or deaths but you were on a non update step, make sure your cell cycle model and killer only operate on update steps.");
643 }
644 else
645 {
647 }
648 }
649 CellBasedEventHandler::EndEvent(CellBasedEventHandler::UPDATECELLPOPULATION);
650}
651
652template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
654{
655 return mOutputDivisionLocations;
656}
657
658template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
660{
661 mOutputDivisionLocations = outputDivisionLocations;
662}
663
664template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
669
670template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
672{
673 mOutputCellVelocities = outputCellVelocities;
674}
675
676template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
678{
679 OutputFileHandler output_file_handler(this->mSimulationOutputDirectory + "/", false);
680
681 // Output machine information
684
686 {
687 // Output Chaste provenance information
688 out_stream build_info_file = output_file_handler.OpenOutputFile("build.info");
689 std::string build_info;
691 *build_info_file << build_info;
692 build_info_file->close();
693
694 // Output simulation parameter and setup details
695 out_stream parameter_file = output_file_handler.OpenOutputFile("results.parameters");
696
697 // Output simulation details
698 std::string simulation_type = GetIdentifier();
699
700 *parameter_file << "<Chaste>\n";
701 *parameter_file << "\n\t<" << simulation_type << ">\n";
702 OutputSimulationParameters(parameter_file);
703 *parameter_file << "\t</" << simulation_type << ">\n";
704 *parameter_file << "\n";
705
706 // Output cell population details (includes cell-cycle model details)
707 mrCellPopulation.OutputCellPopulationInfo(parameter_file);
708
709 // Loop over cell killers
710 *parameter_file << "\n\t<CellKillers>\n";
711 for (typename std::vector<boost::shared_ptr<AbstractCellKiller<SPACE_DIM> > >::iterator iter = mCellKillers.begin();
712 iter != mCellKillers.end();
713 ++iter)
714 {
715 // Output cell killer details
716 (*iter)->OutputCellKillerInfo(parameter_file);
717 }
718 *parameter_file << "\t</CellKillers>\n";
719
720 // Iterate over simulationmodifiers
721 *parameter_file << "\n\t<SimulationModifiers>\n";
722 for (typename std::vector<boost::shared_ptr<AbstractCellBasedSimulationModifier<ELEMENT_DIM, SPACE_DIM> > >::iterator iter = mSimulationModifiers.begin();
723 iter != mSimulationModifiers.end();
724 ++iter)
725 {
726 // Output simulation modifier details
727 (*iter)->OutputSimulationModifierInfo(parameter_file);
728 }
729 *parameter_file << "\t</SimulationModifiers>\n";
730
731 // This is used to output information about subclasses
732 OutputAdditionalSimulationSetup(parameter_file);
733
734 *parameter_file << "\n</Chaste>\n";
735 parameter_file->close();
736 }
737}
738
739template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
741{
742 *rParamsFile << "\t\t<Dt>" << mDt << "</Dt>\n";
743 *rParamsFile << "\t\t<EndTime>" << mEndTime << "</EndTime>\n";
744 *rParamsFile << "\t\t<UpdateCellPopulation>" << mUpdateCellPopulation << "</UpdateCellPopulation>\n";
745 *rParamsFile << "\t\t<SamplingTimestepMultiple>" << mSamplingTimestepMultiple << "</SamplingTimestepMultiple>\n";
746 *rParamsFile << "\t\t<UpdatingTimestepMultiple>" << mUpdatingTimestepMultiple << "</UpdatingTimestepMultiple>\n";
747 *rParamsFile << "\t\t<OutputDivisionLocations>" << mOutputDivisionLocations << "</OutputDivisionLocations>\n";
748 *rParamsFile << "\t\t<OutputCellVelocities>" << mOutputCellVelocities << "</OutputCellVelocities>\n";
749}
750
751// Explicit instantiation
const double DOUBLE_UNSET
Definition Exception.hpp:57
#define EXCEPTION(message)
#define NEVER_REACHED
AbstractCellPopulation< ELEMENT_DIM, SPACE_DIM > & mrCellPopulation
void SetUpdateCellPopulationRule(bool updateCellPopulation)
void SetOutputCellVelocities(bool outputCellVelocities)
void SetUpdatingTimestepMultiple(unsigned updatingTimestepMultiple)
AbstractCellPopulation< ELEMENT_DIM, SPACE_DIM > & rGetCellPopulation()
std::vector< boost::shared_ptr< AbstractCellBasedSimulationModifier< ELEMENT_DIM, SPACE_DIM > > > * GetSimulationModifiers()
void AddSimulationModifier(boost::shared_ptr< AbstractCellBasedSimulationModifier< ELEMENT_DIM, SPACE_DIM > > pSimulationModifier)
void AddTopologyUpdateSimulationModifier(boost::shared_ptr< AbstractCellBasedSimulationModifier< ELEMENT_DIM, SPACE_DIM > > pSimulationModifier)
virtual void OutputSimulationParameters(out_stream &rParamsFile)=0
void SetOutputDirectory(std::string outputDirectory)
void AddCellKiller(boost::shared_ptr< AbstractCellKiller< SPACE_DIM > > pCellKiller)
std::vector< boost::shared_ptr< AbstractCellBasedSimulationModifier< ELEMENT_DIM, SPACE_DIM > > > * GetTopologyUpdateSimulationModifiers()
AbstractCellBasedSimulation(AbstractCellPopulation< ELEMENT_DIM, SPACE_DIM > &rCellPopulation, bool deleteCellPopulationInDestructor=false, bool initialiseCells=true)
std::vector< double > GetNodeLocation(const unsigned &rNodeIndex)
void SetSamplingTimestepMultiple(unsigned samplingTimestepMultiple)
void SetOutputDivisionLocations(bool outputDivisionLocations)
virtual double GetDefaultTimeStep()=0
static void SetOutputDirectory(const std::string &rOutputDirectory)
static void GetBuildInfo(std::string &rInfo)
static void WriteMachineInfoFile(std::string fileBaseName)
std::string GetOutputDirectoryFullPath() const
out_stream OpenOutputFile(const std::string &rFileName, std::ios_base::openmode mode=std::ios::out|std::ios::trunc) const
static bool AmMaster()
static RandomNumberGenerator * Instance()
bool IsFinished() const
double GetTime() const
void SetEndTimeAndNumberOfTimeSteps(double endTime, unsigned totalTimeStepsInSimulation)
bool IsEndTimeAndNumberOfTimeStepsSetUp() const
void IncrementTimeOneStep()
static SimulationTime * Instance()
void ResetEndTimeAndNumberOfTimeSteps(const double &rEndTime, const unsigned &rNumberOfTimeStepsInThisRun)
unsigned GetTimeStepsElapsed() const