Chaste Commit::baa90ac2819b962188b7562f2326be23c47859a7
AbstractCardiacTissue.hpp
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#ifndef ABSTRACTCARDIACTISSUE_HPP_
36#define ABSTRACTCARDIACTISSUE_HPP_
37
38#include <set>
39#include <vector>
40#include <boost/shared_ptr.hpp>
41
43
45#include "ClassIsAbstract.hpp"
47#include <boost/serialization/base_object.hpp>
48#include <boost/serialization/shared_ptr.hpp>
49#include <boost/serialization/vector.hpp>
50#include <boost/serialization/string.hpp>
51#include <boost/serialization/split_member.hpp>
52
53#include "AbstractCardiacCellInterface.hpp"
54#include "FakeBathCell.hpp"
55#include "AbstractCardiacCellFactory.hpp"
56#include "AbstractConductivityTensors.hpp"
57#include "AbstractPurkinjeCellFactory.hpp"
58#include "ReplicatableVector.hpp"
59#include "HeartConfig.hpp"
60#include "ArchiveLocationInfo.hpp"
61#include "AbstractDynamicallyLoadableEntity.hpp"
62#include "DynamicModelLoaderRegistry.hpp"
63#include "AbstractConductivityModifier.hpp"
64
80#ifdef DOXYGEN_CHASTE_ISSUE_199 // See https://github.com/Chaste/Chaste/issues/199
81template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
82#else
83template<unsigned ELEMENT_DIM, unsigned SPACE_DIM = ELEMENT_DIM>
84#endif
85class AbstractCardiacTissue : private boost::noncopyable
86{
87private:
88
91 friend class TestMonodomainTissue;
92
99 template<class Archive>
100 void save(Archive & archive, const unsigned int version) const
101 {
102 if (version >= 3)
103 {
104 archive & mHasPurkinje;
105 }
106 if (version >= 2)
107 {
108 archive & mExchangeHalos;
109 }
110 // Don't use the std::vector serialization for cardiac cells, so that we can load them
111 // more cleverly when migrating checkpoints.
113
114 // archive & mpMesh; Archived in save/load_constructs at the bottom of Mono/BidomainTissue.hpp
115 // archive & mpIntracellularConductivityTensors; Loaded from HeartConfig every time constructor is called
116 if (HeartConfig::Instance()->IsMeshProvided() && HeartConfig::Instance()->GetLoadMesh())
117 {
118 switch (HeartConfig::Instance()->GetConductivityMedia())
119 {
120 case cp::media_type::Orthotropic:
121 {
123 assert(source_file.Exists());
125
126 TRY_IF_MASTER(source_file.CopyTo(dest_file));
127 break;
128 }
129
130 case cp::media_type::Axisymmetric:
131 {
133 assert(source_file.Exists());
136
137 TRY_IF_MASTER(source_file.CopyTo(dest_file));
138 break;
139 }
140
141 case cp::media_type::NoFibreOrientation:
142 break;
143
144 default:
146 }
147 }
148
149 // archive & mIionicCacheReplicated; // will be regenerated
150 // archive & mIntracellularStimulusCacheReplicated; // will be regenerated
151 archive & mDoCacheReplication;
152 // archive & mMeshUnarchived; Not archived since set to true when archiving constructor is called.
153
155
156 // Paranoia: check we agree with the mesh on who owns what
157 assert(mpDistributedVectorFactory == mpMesh->GetDistributedVectorFactory());
158 assert(mpDistributedVectorFactory->GetLow()==mpMesh->GetDistributedVectorFactory()->GetLow());
159 assert(mpDistributedVectorFactory->GetLocalOwnership()==mpMesh->GetDistributedVectorFactory()->GetLocalOwnership());
160 }
161
168 template<class Archive>
169 void load(Archive & archive, const unsigned int version)
170 {
171 // archive & mpMesh; Archived in save/load_constructs at the bottom of Mono/BidomainTissue.hpp
172 // archive & mpIntracellularConductivityTensors; Loaded from HeartConfig every time constructor is called
173
174 if (version >= 3)
175 {
176 archive & mHasPurkinje;
177 if (mHasPurkinje)
178 {
180 }
181 }
182 if (version >= 2)
183 {
184 archive & mExchangeHalos;
185 if (mExchangeHalos)
186 {
189 unsigned num_halo_nodes = mHaloNodes.size();
190 mHaloCellsDistributed.resize( num_halo_nodes );
191 for (unsigned local_index = 0; local_index < num_halo_nodes; local_index++)
192 {
193 unsigned global_index = mHaloNodes[local_index];
194 mHaloGlobalToLocalIndexMap[global_index] = local_index;
195 }
196 }
197 }
198
199 // mCellsDistributed & mHaloCellsDistributed:
201
202 // archive & mIionicCacheReplicated; // will be regenerated
203 // archive & mIntracellularStimulusCacheReplicated; // will be regenerated
204 archive & mDoCacheReplication;
205
206 // we no longer have a bool mDoOneCacheReplication, but to maintain backwards compatibility
207 // we archive something if version==0
208 if (version==0)
209 {
210 bool do_one_cache_replication = true;
211 archive & do_one_cache_replication;
212 }
213
215
216 // Paranoia: check we agree with the mesh on who owns what
217 assert(mpDistributedVectorFactory == mpMesh->GetDistributedVectorFactory());
218 assert(mpDistributedVectorFactory->GetLow()==mpMesh->GetDistributedVectorFactory()->GetLow());
219 assert(mpDistributedVectorFactory->GetLocalOwnership()==mpMesh->GetDistributedVectorFactory()->GetLocalOwnership());
220 // archive & mMeshUnarchived; Not archived since set to true when archiving constructor is called.
221
222 // not archiving mpConductivityModifier for the time being (mechanics simulations are only use-case at the moment, and they
223 // do not get archived...). mpConductivityModifier has to be reset to NULL upon load.
225 }
226 BOOST_SERIALIZATION_SPLIT_MEMBER()
227
228
232
233protected:
234
236 AbstractTetrahedralMesh<ELEMENT_DIM,SPACE_DIM>* mpMesh;
237
241
244
248
254
260
266
272
275
285
290
297
300
311
316
322
324 std::vector<unsigned> mHaloNodes;
325
328
330 std::map<unsigned, unsigned> mHaloGlobalToLocalIndexMap;
331
337 std::vector<std::vector<unsigned> > mNodesToSendPerProcess;
338
343 std::vector<std::vector<unsigned> > mNodesToReceivePerProcess;
344
351
358 void SetUpHaloCells(AbstractCardiacCellFactory<ELEMENT_DIM,SPACE_DIM>* pCellFactory);
359
360public:
371 AbstractCardiacTissue(AbstractCardiacCellFactory<ELEMENT_DIM,SPACE_DIM>* pCellFactory, bool exchangeHalos=false);
372
378 AbstractCardiacTissue(AbstractTetrahedralMesh<ELEMENT_DIM,SPACE_DIM>* pMesh);
379
381 virtual ~AbstractCardiacTissue();
382
384 bool HasPurkinje();
385
392 void SetCacheReplication(bool doCacheReplication);
393
400
404 const c_matrix<double, SPACE_DIM, SPACE_DIM>& rGetIntracellularConductivityTensor(unsigned elementIndex);
405
413 virtual const c_matrix<double, SPACE_DIM, SPACE_DIM>& rGetExtracellularConductivityTensor(unsigned elementIndex);
414
423 AbstractCardiacCellInterface* GetCardiacCell( unsigned globalIndex );
424
433 AbstractCardiacCellInterface* GetPurkinjeCell( unsigned globalIndex );
434
444
454 virtual void SolveCellSystems(Vec existingSolution, double time, double nextTime, bool updateVoltage=false);
455
458
461
464
467
468
476 void UpdateCaches(unsigned globalIndex, unsigned localIndex, double nextTime);
477
485 void UpdatePurkinjeCaches(unsigned globalIndex, unsigned localIndex, double nextTime);
486
490 void ReplicateCaches();
491
495 const std::vector<AbstractCardiacCellInterface*>& rGetCellsDistributed() const;
496
501
507 const AbstractTetrahedralMesh<ELEMENT_DIM,SPACE_DIM>* pGetMesh() const;
508
514 void SetConductivityModifier(AbstractConductivityModifier<ELEMENT_DIM,SPACE_DIM>* pModifier);
515
527 template<class Archive>
528 void SaveCardiacCells(Archive & archive, const unsigned int version) const
529 {
530 const std::vector<AbstractCardiacCellInterface*> & r_cells_distributed = rGetCellsDistributed();
531 assert(mpDistributedVectorFactory == this->mpMesh->GetDistributedVectorFactory());
532 archive & mpDistributedVectorFactory; // Needed when loading
533 const unsigned num_cells = r_cells_distributed.size();
534 archive & num_cells;
535 for (unsigned i=0; i<num_cells; i++)
536 {
537 AbstractDynamicallyLoadableEntity* p_entity = dynamic_cast<AbstractDynamicallyLoadableEntity*>(r_cells_distributed[i]);
538 bool is_dynamic = (p_entity != NULL);
539 archive & is_dynamic;
540 if (is_dynamic)
541 {
542#ifdef CHASTE_CAN_CHECKPOINT_DLLS
543 archive & p_entity->GetLoader()->GetLoadableModulePath();
544#else
545 // We should have thrown an exception before this point
547#endif // CHASTE_CAN_CHECKPOINT_DLLS
548 }
549 archive & r_cells_distributed[i];
550 if (mHasPurkinje)
551 {
552 archive & rGetPurkinjeCellsDistributed()[i];
553 }
554 }
555 }
556
569 template<class Archive>
570 void LoadCardiacCells(Archive & archive, const unsigned int version)
571 {
572 // Note that p_factory loaded from this archive might not be the same as our mesh's factory,
573 // since we're loading a process-specific archive that could have been written by any process.
574 // We therefore need to use p_mesh_factory to determine the partitioning in use for the resumed
575 // simulation, and p_factory to determine what the original partitioning was when the simulation
576 // was saved.
577 DistributedVectorFactory* p_factory;
578 DistributedVectorFactory* p_mesh_factory = this->mpMesh->GetDistributedVectorFactory();
579 archive & p_factory;
580 unsigned num_cells;
581 archive & num_cells;
582 if (mCellsDistributed.empty())
583 {
584 mCellsDistributed.resize(p_mesh_factory->GetLocalOwnership());
585#ifndef NDEBUG
586 // Paranoia
587 for (unsigned i=0; i<mCellsDistributed.size(); i++)
588 {
589 assert(mCellsDistributed[i] == NULL);
590 }
591#endif
592 }
593 else
594 {
595 assert(mCellsDistributed.size() == p_mesh_factory->GetLocalOwnership());
596 }
597 if (mHasPurkinje)
598 {
599 if (mPurkinjeCellsDistributed.empty())
600 {
601 mPurkinjeCellsDistributed.resize(p_mesh_factory->GetLocalOwnership());
602 }
603 else
604 {
605 assert(mPurkinjeCellsDistributed.size() == p_mesh_factory->GetLocalOwnership());
606 }
607 }
608
609 // We don't store a cell index in the archive, so need to work out what global index this tissue starts at.
610 // If we have an original factory we use the original low index; otherwise we use the current low index.
611 unsigned index_low = p_factory->GetOriginalFactory() ? p_factory->GetOriginalFactory()->GetLow() : p_mesh_factory->GetLow();
612
613 // Track fake cells (which might have multiple pointers to the same object) to make sure we only delete non-local ones
614 std::set<FakeBathCell*> fake_cells_non_local, fake_cells_local;
615
616 /*
617 * Historical note:
618 *
619 * We always do a dumb partition when we unarchive.
620 *
621 * When unarchive was first implemented in parallel (migration #1199) it was thought that we might want to repartition the mesh. This would be feasible and would give
622 * better partitions when we move to different numbers of processes. However it would require us to apply a new permutation to entire data structure.
623 *
624 * In the case where the original mesh was permuted and *copied* into the archive, we need to apply the stored permutation to the mesh but not to the archive (cells). That
625 * is, any permutation stored with the mesh can be safely ignored. (If we also had to repartition/permute the archive, then we would be applying a double permutation to the
626 * mesh and a single permutation to archive.)
627 *
628 */
629 for (unsigned local_index=0; local_index<num_cells; local_index++)
630 {
631 // Figure out where this cell goes
632 unsigned global_index = index_low + local_index;
633 bool local = p_mesh_factory->IsGlobalIndexLocal(global_index);
634
635 // Check if this will be a halo cell
636 std::map<unsigned, unsigned>::const_iterator halo_position;
637 bool halo = ((halo_position=mHaloGlobalToLocalIndexMap.find(global_index)) != mHaloGlobalToLocalIndexMap.end());
638 // halo_position->second is local halo index
639
640 bool is_dynamic;
641 archive & is_dynamic;
642 if (is_dynamic)
643 {
644#ifdef CHASTE_CAN_CHECKPOINT_DLLS
645 // Ensure the shared object file for this cell model is loaded.
646 // We need to do this here, rather than in the class' serialization code,
647 // because that code won't be available until this is done...
648 std::string shared_object_path;
649 archive & shared_object_path;
650 DynamicModelLoaderRegistry::Instance()->GetLoader(shared_object_path);
651#else
652 // Could only happen on Mac OS X, and will probably be trapped earlier.
654#endif // CHASTE_CAN_CHECKPOINT_DLLS
655 }
657 archive & p_cell;
658 AbstractCardiacCellInterface* p_purkinje_cell = NULL;
659 if (mHasPurkinje)
660 {
661 archive & p_purkinje_cell;
662 }
663 // Check if it's a fake cell
664 FakeBathCell* p_fake = dynamic_cast<FakeBathCell*>(p_cell);
665 if (p_fake)
666 {
667 if (halo || local)
668 {
669 fake_cells_local.insert(p_fake);
670 }
671 else
672 {
673 fake_cells_non_local.insert(p_fake);
674 }
675 }
676 FakeBathCell* p_fake_purkinje = dynamic_cast<FakeBathCell*>(p_purkinje_cell);
677 if (p_fake_purkinje)
678 {
679 if (halo || local)
680 {
681 fake_cells_local.insert(p_fake_purkinje);
682 }
683 else
684 {
685 fake_cells_non_local.insert(p_fake_purkinje);
686 }
687 }
688 // Add real cells to the local or halo vectors
689 if (local)
690 {
691 // Note that the original local_index was relative to the archived partition (distributed vector)
692 // The new_local_index is local relative to the new partition in memory
693 unsigned new_local_index = global_index - p_mesh_factory->GetLow();
694 assert(mCellsDistributed[new_local_index] == NULL);
695 mCellsDistributed[new_local_index] = p_cell;
696 if (mHasPurkinje)
697 {
698 assert(mPurkinjeCellsDistributed[new_local_index] == NULL);
699 mPurkinjeCellsDistributed[new_local_index] = p_purkinje_cell;
700 }
701 }
702 else if (halo)
703 {
704 assert(mHaloCellsDistributed[halo_position->second] == NULL);
705 mHaloCellsDistributed[halo_position->second] = p_cell;
706 }
707 else
708 {
709 if (!p_fake)
710 {
711 // Non-local real cell, so free the memory.
712 delete p_cell;
713 }
714 if (!p_fake_purkinje)
715 {
716 // This will be NULL if there's no Purkinje, so a delete is OK.
717 delete p_purkinje_cell;
718 }
719 }
720 }
721
722 // Delete any unused fake cells
723 for (std::set<FakeBathCell*>::iterator it = fake_cells_non_local.begin();
724 it != fake_cells_non_local.end();
725 ++it)
726 {
727 if (fake_cells_local.find(*it) == fake_cells_local.end())
728 {
729 delete (*it);
730 }
731 }
732 }
733};
734
736
737namespace boost {
738namespace serialization {
745template <unsigned ELEMENT_DIM, unsigned SPACE_DIM>
746struct version<AbstractCardiacTissue<ELEMENT_DIM, SPACE_DIM> >
747{
750};
751} // namespace serialization
752} // namespace boost
753
754#endif /*ABSTRACTCARDIACTISSUE_HPP_*/
755
gcov doesn't like this file...
#define TEMPLATED_CLASS_IS_ABSTRACT_2_UNSIGNED(T)
#define NEVER_REACHED
#define TRY_IF_MASTER(method)
const std::vector< AbstractCardiacCellInterface * > & rGetPurkinjeCellsDistributed() const
ReplicatableVector mPurkinjeIntracellularStimulusCacheReplicated
ReplicatableVector & rGetIionicCacheReplicated()
ReplicatableVector mIntracellularStimulusCacheReplicated
const AbstractTetrahedralMesh< ELEMENT_DIM, SPACE_DIM > * pGetMesh() const
AbstractCardiacCellInterface * GetCardiacCellOrHaloCell(unsigned globalIndex)
void UpdateCaches(unsigned globalIndex, unsigned localIndex, double nextTime)
void UpdatePurkinjeCaches(unsigned globalIndex, unsigned localIndex, double nextTime)
const std::vector< AbstractCardiacCellInterface * > & rGetCellsDistributed() const
std::vector< std::vector< unsigned > > mNodesToReceivePerProcess
ReplicatableVector & rGetPurkinjeIntracellularStimulusCacheReplicated()
virtual void SolveCellSystems(Vec existingSolution, double time, double nextTime, bool updateVoltage=false)
ReplicatableVector & rGetIntracellularStimulusCacheReplicated()
DistributedVectorFactory * mpDistributedVectorFactory
void save(Archive &archive, const unsigned int version) const
std::vector< AbstractCardiacCellInterface * > mHaloCellsDistributed
void SaveCardiacCells(Archive &archive, const unsigned int version) const
AbstractConductivityTensors< ELEMENT_DIM, SPACE_DIM > * mpIntracellularConductivityTensors
void SetCacheReplication(bool doCacheReplication)
ReplicatableVector & rGetPurkinjeIionicCacheReplicated()
std::vector< AbstractCardiacCellInterface * > mPurkinjeCellsDistributed
AbstractCardiacCellInterface * GetPurkinjeCell(unsigned globalIndex)
std::vector< std::vector< unsigned > > mNodesToSendPerProcess
void SetConductivityModifier(AbstractConductivityModifier< ELEMENT_DIM, SPACE_DIM > *pModifier)
std::vector< AbstractCardiacCellInterface * > mCellsDistributed
friend class boost::serialization::access
void LoadCardiacCells(Archive &archive, const unsigned int version)
ReplicatableVector mPurkinjeIionicCacheReplicated
void SetUpHaloCells(AbstractCardiacCellFactory< ELEMENT_DIM, SPACE_DIM > *pCellFactory)
void load(Archive &archive, const unsigned int version)
AbstractConductivityModifier< ELEMENT_DIM, SPACE_DIM > * mpConductivityModifier
AbstractCardiacCellInterface * GetCardiacCell(unsigned globalIndex)
virtual const c_matrix< double, SPACE_DIM, SPACE_DIM > & rGetExtracellularConductivityTensor(unsigned elementIndex)
const c_matrix< double, SPACE_DIM, SPACE_DIM > & rGetIntracellularConductivityTensor(unsigned elementIndex)
std::map< unsigned, unsigned > mHaloGlobalToLocalIndexMap
AbstractTetrahedralMesh< ELEMENT_DIM, SPACE_DIM > * mpMesh
std::vector< unsigned > mHaloNodes
ReplicatableVector mIionicCacheReplicated
Forward declaration which is going to be used for friendship.
const DynamicCellModelLoaderPtr GetLoader() const
static std::string GetMeshFilename()
static std::string GetArchiveRelativePath()
DistributedVectorFactory * GetOriginalFactory()
bool IsGlobalIndexLocal(unsigned globalIndex)
DynamicCellModelLoaderPtr GetLoader(const std::string &rPath)
static DynamicModelLoaderRegistry * Instance()
bool Exists() const
FileFinder CopyTo(const FileFinder &rDest) const
static HeartConfig * Instance()
void Resize(unsigned size)
CHASTE_VERSION_CONTENT(3)
Macro to set the version number of templated archive in known versions of Boost.