35 #include "AbstractCardiacProblem.hpp" 37 #include "DistributedVector.hpp" 39 #include "GenericMeshReader.hpp" 40 #include "Hdf5ToCmguiConverter.hpp" 41 #include "Hdf5ToMeshalyzerConverter.hpp" 42 #include "Hdf5ToVtkConverter.hpp" 43 #include "HeartConfig.hpp" 44 #include "HeartEventHandler.hpp" 45 #include "LinearSystem.hpp" 47 #include "PostProcessingWriter.hpp" 48 #include "ProgressReporter.hpp" 49 #include "TimeStepper.hpp" 51 template <
unsigned ELEMENT_DIM,
unsigned SPACE_DIM,
unsigned PROBLEM_DIM>
55 mAllocatedMemoryForMesh(false),
58 mpCardiacTissue(NULL),
60 mpCellFactory(pCellFactory),
64 mpTimeAdaptivityController(NULL),
66 mUseHdf5DataWriterCache(false),
67 mHdf5DataWriterChunkSizeAndAlignment(0)
72 EXCEPTION(
"AbstractCardiacProblem: Please supply a cell factory pointer to your cardiac problem constructor.");
77 template <
unsigned ELEMENT_DIM,
unsigned SPACE_DIM,
unsigned PROBLEM_DIM>
101 template <
unsigned ELEMENT_DIM,
unsigned SPACE_DIM,
unsigned PROBLEM_DIM>
116 template <
unsigned ELEMENT_DIM,
unsigned SPACE_DIM,
unsigned PROBLEM_DIM>
124 WARNING(
"Using a non-distributed mesh in a parallel simulation is not a good idea.");
135 std::shared_ptr<AbstractMeshReader<ELEMENT_DIM, SPACE_DIM> > p_mesh_reader
137 mpMesh->ConstructFromMeshReader(*p_mesh_reader);
149 c_vector<double, 1> fibre_length;
151 mpMesh->ConstructRegularSlabMesh(inter_node_space, fibre_length[0]);
156 c_vector<double, 2> sheet_dimensions;
158 mpMesh->ConstructRegularSlabMesh(inter_node_space, sheet_dimensions[0], sheet_dimensions[1]);
163 c_vector<double, 3> slab_dimensions;
165 mpMesh->ConstructRegularSlabMesh(inter_node_space, slab_dimensions[0], slab_dimensions[1], slab_dimensions[2]);
181 EXCEPTION(std::string(
"No mesh given: define it in XML parameters file or call SetMesh()\n") + e.
GetShortMessage());
216 template <
unsigned ELEMENT_DIM,
unsigned SPACE_DIM,
unsigned PROBLEM_DIM>
222 template <
unsigned ELEMENT_DIM,
unsigned SPACE_DIM,
unsigned PROBLEM_DIM>
228 template <
unsigned ELEMENT_DIM,
unsigned SPACE_DIM,
unsigned PROBLEM_DIM>
233 EXCEPTION(
"Cardiac tissue is null, Initialise() probably hasn't been called");
237 EXCEPTION(
"End time should be in the future");
243 EXCEPTION(
"Either explicitly specify not to print output (call PrintOutput(false)) or specify the output directory and filename prefix");
258 if (fabs(end_time - pde_time * round(end_time / pde_time)) > 1e-10)
260 EXCEPTION(
"PDE timestep does not seem to divide end time - check parameters");
264 template <
unsigned ELEMENT_DIM,
unsigned SPACE_DIM,
unsigned PROBLEM_DIM>
268 Vec initial_condition = p_factory->
CreateVec(PROBLEM_DIM);
270 std::vector<DistributedVector::Stripe> stripe;
271 stripe.reserve(PROBLEM_DIM);
273 for (
unsigned i = 0; i < PROBLEM_DIM; i++)
282 stripe[0][index] =
mpCardiacTissue->GetCardiacCell(index.Global)->GetVoltage();
283 if (PROBLEM_DIM == 2)
285 stripe[1][index] = 0;
291 return initial_condition;
294 template <
unsigned ELEMENT_DIM,
unsigned SPACE_DIM,
unsigned PROBLEM_DIM>
302 assert(pMesh != NULL);
307 template <
unsigned ELEMENT_DIM,
unsigned SPACE_DIM,
unsigned PROBLEM_DIM>
313 template <
unsigned ELEMENT_DIM,
unsigned SPACE_DIM,
unsigned PROBLEM_DIM>
319 template <
unsigned ELEMENT_DIM,
unsigned SPACE_DIM,
unsigned PROBLEM_DIM>
325 template <
unsigned ELEMENT_DIM,
unsigned SPACE_DIM,
unsigned PROBLEM_DIM>
328 return mpMesh->GetDistributedVectorFactory()->CreateDistributedVector(
mSolution);
331 template <
unsigned ELEMENT_DIM,
unsigned SPACE_DIM,
unsigned PROBLEM_DIM>
337 template <
unsigned ELEMENT_DIM,
unsigned SPACE_DIM,
unsigned PROBLEM_DIM>
344 template <
unsigned ELEMENT_DIM,
unsigned SPACE_DIM,
unsigned PROBLEM_DIM>
349 EXCEPTION(
"Tissue not yet set up, you may need to call Initialise() before GetTissue().");
354 template <
unsigned ELEMENT_DIM,
unsigned SPACE_DIM,
unsigned PROBLEM_DIM>
370 template <
unsigned ELEMENT_DIM,
unsigned SPACE_DIM,
unsigned PROBLEM_DIM>
375 std::vector<double> additional_stopping_times;
382 additional_stopping_times);
391 for (
unsigned problem_index = 0; problem_index < PROBLEM_DIM; problem_index++)
402 Vec initial_condition;
412 std::string progress_reporter_dir;
417 bool extending_file =
false;
458 progress_reporter_dir =
"";
460 BOOST_FOREACH (boost::shared_ptr<AbstractOutputModifier> p_output_modifier,
mOutputModifiers)
462 p_output_modifier->InitialiseAtStart(this->
mpMesh->GetDistributedVectorFactory(), this->
mpMesh->rGetNodePermutation());
463 p_output_modifier->ProcessSolutionAtTimeStep(stepper.
GetTime(), initial_condition, PROBLEM_DIM);
486 mpSolver->SetInitialCondition(initial_condition);
552 BOOST_FOREACH (boost::shared_ptr<AbstractOutputModifier> p_output_modifier,
mOutputModifiers)
554 p_output_modifier->ProcessSolutionAtTimeStep(stepper.
GetTime(),
mSolution, PROBLEM_DIM);
578 BOOST_FOREACH (boost::shared_ptr<AbstractOutputModifier> p_output_modifier,
mOutputModifiers)
580 p_output_modifier->FinaliseAtEnd();
586 template <
unsigned ELEMENT_DIM,
unsigned SPACE_DIM,
unsigned PROBLEM_DIM>
683 template <
unsigned ELEMENT_DIM,
unsigned SPACE_DIM,
unsigned PROBLEM_DIM>
696 if (
mpMesh->rGetNodePermutation().size() > 0)
700 EXCEPTION(
"HeartConfig setting `GetOutputUsingOriginalNodeOrdering` is meaningless when outputting particular nodes in parallel. (Nodes are written with their original indices by default).");
702 std::vector<unsigned> nodes_to_output_permuted(
mNodesToOutput.size());
729 template <
unsigned ELEMENT_DIM,
unsigned SPACE_DIM,
unsigned PROBLEM_DIM>
737 std::vector<std::string> output_variables;
739 const unsigned num_vars = output_variables.size();
743 for (
unsigned var_index = 0; var_index < num_vars; var_index++)
746 std::string var_name = output_variables[var_index];
767 template <
unsigned ELEMENT_DIM,
unsigned SPACE_DIM,
unsigned PROBLEM_DIM>
771 std::vector<std::string> output_variables;
777 assert(output_variables.size() == num_vars);
780 for (
unsigned var_index = 0; var_index < num_vars; var_index++)
783 Vec variable_data = this->
mpMesh->GetDistributedVectorFactory()->CreateVec();
784 DistributedVector distributed_var_data = this->
mpMesh->GetDistributedVectorFactory()->CreateDistributedVector(variable_data);
788 index != distributed_var_data.
End();
796 distributed_var_data[index] = 0.0;
801 distributed_var_data[index] = this->
mpCardiacTissue->GetCardiacCell(index.Global)->GetAnyVariable(output_variables[var_index],
mCurrentTime);
804 distributed_var_data.
Restore();
813 template <
unsigned ELEMENT_DIM,
unsigned SPACE_DIM,
unsigned PROBLEM_DIM>
843 EXCEPTION(
"Attempting to extend " << h5_file.
GetAbsolutePath() <<
" with results from time = " <<
mCurrentTime <<
", but it already contains results up to time = " << times.back() <<
"." 844 " Calling HeartConfig::Instance()->SetOutputDirectory() before Solve() will direct results elsewhere.");
877 if (success ==
false)
892 template <
unsigned ELEMENT_DIM,
unsigned SPACE_DIM,
unsigned PROBLEM_DIM>
898 template <
unsigned ELEMENT_DIM,
unsigned SPACE_DIM,
unsigned PROBLEM_DIM>
904 template <
unsigned ELEMENT_DIM,
unsigned SPACE_DIM,
unsigned PROBLEM_DIM>
910 template <
unsigned ELEMENT_DIM,
unsigned SPACE_DIM,
unsigned PROBLEM_DIM>
915 EXCEPTION(
"Data reader invalid as data writer cannot be initialised");
920 template <
unsigned ELEMENT_DIM,
unsigned SPACE_DIM,
unsigned PROBLEM_DIM>
926 template <
unsigned ELEMENT_DIM,
unsigned SPACE_DIM,
unsigned PROBLEM_DIM>
virtual void WriteOneStep(double time, Vec voltageVec)=0
BccType mpDefaultBoundaryConditionsContainer
std::string GetAbsolutePath() const
bool ApplyPermutation(const std::vector< unsigned > &rPermutation, bool unsafeExtendingMode=false)
void WritePostProcessingFiles()
std::string GetMeshName() const
static bool IsRegionBath(HeartRegionType regionId)
Hdf5DataWriter * mpWriter
std::vector< boost::shared_ptr< AbstractOutputModifier > > mOutputModifiers
std::string mMeshFilename
double GetInterNodeSpace() const
unsigned mVoltageColumnId
void PrintOutput(bool rPrintOutput)
int GetVariableByName(const std::string &rVariableName)
DistributedTetrahedralMeshPartitionType::type GetMeshPartitioning() const
void Update(double currentTime)
void GetOutputVariables(std::vector< std::string > &rOutputVariables) const
DistributedVector CreateDistributedVector(Vec vec, bool readOnly=false)
void Write(bool useArchiveLocationInfo=false, std::string subfolderName="output")
hsize_t mHdf5DataWriterChunkSizeAndAlignment
#define EXCEPTION(message)
static void BeginEvent(unsigned event)
virtual Vec CreateInitialCondition()
virtual void CreateMeshFromHeartConfig()
bool GetOutputUsingOriginalNodeOrdering()
AbstractTetrahedralMesh< ELEMENT_DIM, SPACE_DIM > & rGetMesh()
double GetNextTime() const
void SetOutputUsingOriginalNodeOrdering(bool useOriginal)
std::vector< unsigned > mExtraVariablesId
void GetFibreLength(c_vector< double, 1 > &fibreLength) const
AbstractTimeAdaptivityController * mpTimeAdaptivityController
virtual void SetUpAdditionalStoppingTimes(std::vector< double > &rAdditionalStoppingTimes)
virtual void PreSolveChecks()
void SetTargetChunkSize(hsize_t targetSize)
void AdvanceOneTimeStep()
std::string GetOutputFilenamePrefix() const
void CloseFilesAndPostProcess()
AbstractTetrahedralMesh< ELEMENT_DIM, SPACE_DIM > * mpMesh
void SetHdf5DataWriterTargetChunkSizeAndAlignment(hsize_t size)
bool GetVisualizeWithMeshalyzer() const
void GetSlabDimensions(c_vector< double, 3 > &slabDimensions) const
void GetSheetDimensions(c_vector< double, 2 > &sheetDimensions) const
bool mAllocatedMemoryForMesh
Hdf5DataReader GetDataReader()
void SetUseTimeAdaptivityController(bool useAdaptivity, AbstractTimeAdaptivityController *pController=NULL)
virtual void WriteInfo(double time)=0
BccType mpBoundaryConditionsContainer
AbstractCardiacTissue< ELEMENT_DIM, SPACE_DIM > * GetTissue()
void AdvanceAlongUnlimitedDimension()
std::string GetOutputDirectory() const
std::vector< double > GetUnlimitedDimensionValues()
unsigned EstimateTimeSteps() const
int DefineVariable(const std::string &rVariableName, const std::string &rVariableUnits)
std::vector< unsigned > mNodesToOutput
virtual void OnEndOfTimestep(double time)
std::string GetSubdirectory()
void PutVector(int variableID, Vec petscVector)
void SetMesh(AbstractTetrahedralMesh< ELEMENT_DIM, SPACE_DIM > *pMesh)
void SetUseHdf5DataWriterCache(bool useCache=true)
AbstractCardiacCellFactory< ELEMENT_DIM, SPACE_DIM > * mpCellFactory
bool mUseHdf5DataWriterCache
virtual bool GetHasBath()
void DefineFixedDimension(long dimensionSize)
void SetOutputNodes(std::vector< unsigned > &rNodesToOutput)
double GetSimulationDuration() const
AbstractDynamicLinearPdeSolver< ELEMENT_DIM, SPACE_DIM, PROBLEM_DIM > * mpSolver
virtual void DefineWriterColumns(bool extending)
virtual ~AbstractCardiacProblem()
virtual AbstractCardiacTissue< ELEMENT_DIM, SPACE_DIM > * CreateCardiacTissue()=0
void SetBoundaryConditionsContainer(BccType pBcc)
virtual void SetElectrodes()
void SetAlignment(hsize_t alignment)
AbstractCardiacTissue< ELEMENT_DIM, SPACE_DIM > * mpCardiacTissue
static void EndEvent(unsigned event)
static std::string GetChasteTestOutputDirectory()
void SetWriteInfo(bool writeInfo=true)
virtual AbstractDynamicLinearPdeSolver< ELEMENT_DIM, SPACE_DIM, PROBLEM_DIM > * CreateSolver()=0
void WriteExtraVariablesOneStep()
DistributedVector GetSolutionDistributedVector()
double GetPdeTimeStep() const
std::string GetShortMessage() const
static HeartConfig * Instance()
virtual void AtBeginningOfTimestep(double time)
void DefineExtraVariablesWriterColumns(bool extending)
virtual void EndDefineMode()
void DefineUnlimitedDimension(const std::string &rVariableName, const std::string &rVariableUnits, unsigned estimatedLength=1)