37#include "NodesOnlyMesh.hpp"
38#include "ChasteCuboid.hpp"
40template<
unsigned SPACE_DIM>
43 mMaximumInteractionDistance(1.0),
45 mMinimumNodeDomainBoundarySeparation(1.0),
46 mMaxAddedNodeIndex(0u),
47 mpBoxCollection(nullptr),
48 mCalculateNodeNeighbours(true)
52template<
unsigned SPACE_DIM>
59template<
unsigned SPACE_DIM>
62 assert(maxInteractionDistance > 0.0 && maxInteractionDistance < DBL_MAX);
63 mMaximumInteractionDistance = maxInteractionDistance;
65 mMinimumNodeDomainBoundarySeparation = mMaximumInteractionDistance;
69 SetUpBoxCollection(rNodes);
71 mLocalInitialNodes.resize(rNodes.size(),
false);
73 for (
unsigned i=0; i<rNodes.size(); i++)
75 if (mpBoxCollection->IsOwned(rNodes[i]))
77 assert(!rNodes[i]->IsDeleted());
79 mLocalInitialNodes[i] =
true;
82 c_vector<double, SPACE_DIM> location = rNodes[i]->rGetLocation();
88 if (rNodes[i]->HasNodeAttributes())
93 this->mNodes.push_back(p_node_copy);
96 mNodesMapping[p_node_copy->
GetIndex()] = this->mNodes.size()-1;
101template<
unsigned SPACE_DIM>
105 std::vector<Node<SPACE_DIM>*> temp_nodes(rNodes.size());
106 for (
unsigned idx = 0; idx < rNodes.size(); idx++)
108 temp_nodes[idx] = rNodes[idx].get();
111 ConstructNodesWithoutMesh(temp_nodes, maxInteractionDistance);
114template<
unsigned SPACE_DIM>
117 ConstructNodesWithoutMesh(rGeneratingMesh.
mNodes, maxInteractionDistance);
120template<
unsigned SPACE_DIM>
123 return mLocalInitialNodes;
126template<
unsigned SPACE_DIM>
129 std::map<unsigned, unsigned>::const_iterator node_position = mNodesMapping.find(index);
131 if (node_position == mNodesMapping.end())
136 return node_position->second;
139template<
unsigned SPACE_DIM>
146 mNodesMapping.clear();
151template<
unsigned SPACE_DIM>
154 return mpBoxCollection;
157template <
unsigned SPACE_DIM>
162 std::map<unsigned, unsigned>::const_iterator node_position = mHaloNodesMapping.find(index);
164 if (node_position != mHaloNodesMapping.end())
166 p_node = mHaloNodes[node_position->second].get();
170 p_node = this->GetNode(index);
173 assert(p_node !=
nullptr);
178template<
unsigned SPACE_DIM>
181 return mpBoxCollection->IsOwned(location);
184template<
unsigned SPACE_DIM>
187 return this->mNodes.size() - this->mDeletedNodeIndices.size();
190template<
unsigned SPACE_DIM>
196template<
unsigned SPACE_DIM>
199 mMaximumInteractionDistance = maxDistance;
202template<
unsigned SPACE_DIM>
205 return mMaximumInteractionDistance;
208template<
unsigned SPACE_DIM>
219template<
unsigned SPACE_DIM>
222 mCalculateNodeNeighbours = calculateNodeNeighbours;
225template<
unsigned SPACE_DIM>
228 assert(mpBoxCollection);
230 mpBoxCollection->CalculateInteriorNodePairs(this->mNodes, rNodePairs);
233template<
unsigned SPACE_DIM>
236 assert(mpBoxCollection);
238 mpBoxCollection->CalculateBoundaryNodePairs(this->mNodes, rNodePairs);
241template<
unsigned SPACE_DIM>
246 RemoveDeletedNodes(map);
248 this->mDeletedNodeIndices.clear();
249 this->mAddedNodes =
false;
253 this->SetMeshHasChangedSinceLoading();
256template<
unsigned SPACE_DIM>
259 typename std::vector<Node<SPACE_DIM>* >::iterator node_iter = this->mNodes.begin();
260 while (node_iter != this->mNodes.end())
262 if ((*node_iter)->IsDeleted())
266 mNodesMapping.erase((*node_iter)->GetIndex());
270 node_iter = this->mNodes.erase(node_iter);
279template<
unsigned SPACE_DIM>
282 mNodesMapping.clear();
283 for (
unsigned location_in_vector=0; location_in_vector < this->mNodes.size(); location_in_vector++)
285 unsigned global_index = this->mNodes[location_in_vector]->GetIndex();
286 mNodesMapping[global_index] = location_in_vector;
290template<
unsigned SPACE_DIM>
293 mNodesToSendRight.clear();
294 mNodesToSendLeft.clear();
297 node_iter != this->GetNodeIteratorEnd();
300 unsigned owning_process = mpBoxCollection->GetProcessOwningNode(&(*node_iter));
307 mNodesToSendRight.push_back(node_iter->GetIndex());
311 mNodesToSendLeft.push_back(node_iter->GetIndex());
323 mNodesToSendLeft.push_back(node_iter->GetIndex());
325 else if ( owning_process == 0 )
328 mNodesToSendRight.push_back(node_iter->GetIndex());
334template<
unsigned SPACE_DIM>
337 return mNodesToSendLeft;
340template<
unsigned SPACE_DIM>
343 return mNodesToSendRight;
346template<
unsigned SPACE_DIM>
349 return mpBoxCollection->rGetHaloNodesRight();
352template<
unsigned SPACE_DIM>
355 return mpBoxCollection->rGetHaloNodesLeft();
358template<
unsigned SPACE_DIM>
361 unsigned location_in_nodes_vector = 0;
363 if (this->mDeletedNodeIndices.empty())
365 this->mNodes.push_back(pNewNode);
366 location_in_nodes_vector = this->mNodes.size() - 1;
370 location_in_nodes_vector = this->mDeletedNodeIndices.back();
371 this->mDeletedNodeIndices.pop_back();
372 delete this->mNodes[location_in_nodes_vector];
373 this->mNodes[location_in_nodes_vector] = pNewNode;
376 this->mAddedNodes =
true;
378 mMaxAddedNodeIndex = (pNewNode->
GetIndex() > mMaxAddedNodeIndex) ? pNewNode->
GetIndex() : mMaxAddedNodeIndex;
381 mNodesMapping[pNewNode->
GetIndex()] = location_in_nodes_vector;
387template<
unsigned SPACE_DIM>
390 mHaloNodes.push_back(pNewNode);
391 mHaloNodesMapping[pNewNode->
GetIndex()] = mHaloNodes.size() - 1;
394template<
unsigned SPACE_DIM>
399 mHaloNodesMapping.clear();
402template<
unsigned SPACE_DIM>
405 unsigned fresh_global_index = GetNextAvailableIndex();
406 pNewNode->
SetIndex(fresh_global_index);
408 AddNodeWithFixedIndex(pNewNode);
410 return fresh_global_index;
413template<
unsigned SPACE_DIM>
417 assert(!concreteMove);
420 this->GetNode(nodeIndex)->SetPoint(point);
423template<
unsigned SPACE_DIM>
427 unsigned index = pMovedNode->GetIndex();
428 c_vector<double, SPACE_DIM> location = pMovedNode->rGetLocation();
432 if (pMovedNode->HasNodeAttributes())
434 double radius = pMovedNode->GetRadius();
437 unsigned region = pMovedNode->GetRegion();
440 bool is_particle = pMovedNode->IsParticle();
443 for (
unsigned i=0; i<pMovedNode->GetNumNodeAttributes(); i++)
445 double attribute = pMovedNode->rGetNodeAttributes()[i];
450 AddNodeWithFixedIndex(p_node);
453template<
unsigned SPACE_DIM>
456 if (this->GetNode(index)->IsDeleted())
458 EXCEPTION(
"Trying to delete a deleted node");
461 unsigned local_index = SolveNodeMapping(index);
463 this->mNodes[local_index]->MarkAsDeleted();
464 this->mDeletedNodeIndices.push_back(local_index);
465 mDeletedGlobalNodeIndices.push_back(index);
468template<
unsigned SPACE_DIM>
474 mDeletedGlobalNodeIndices.pop_back();
477template<
unsigned SPACE_DIM>
480 assert(!(separation < 0.0));
482 mMinimumNodeDomainBoundarySeparation = separation;
485template<
unsigned SPACE_DIM>
490 if (!this->mDeletedGlobalNodeIndices.empty())
492 index = this->mDeletedGlobalNodeIndices.back();
493 this->mDeletedGlobalNodeIndices.pop_back();
497 unsigned counter = mIndexCounter;
505template <
unsigned SPACE_DIM>
508 assert(mpBoxCollection);
510 int num_local_rows = mpBoxCollection->GetNumLocalRows();
513 c_vector<double, 2 * SPACE_DIM> current_domain_size = mpBoxCollection->rGetDomainSize();
514 c_vector<double, 2 * SPACE_DIM> new_domain_size;
515 new_domain_size = current_domain_size;
517 double fudge = 1e-14;
518 c_vector<bool, SPACE_DIM> is_periodic = mpBoxCollection->GetIsPeriodicAllDims();
519 for (
unsigned d = 0; d < SPACE_DIM; d++)
524 new_domain_size[2 * d] = current_domain_size[2 * d] - (mMaximumInteractionDistance - fudge);
525 new_domain_size[2 * d + 1] = current_domain_size[2 * d + 1] + (mMaximumInteractionDistance - fudge);
528 SetUpBoxCollection(mMaximumInteractionDistance, new_domain_size, new_local_rows);
531template<
unsigned SPACE_DIM>
534 assert(mpBoxCollection);
536 int is_local_node_close = 0;
537 c_vector<double, 2*SPACE_DIM> domain_boundary = mpBoxCollection->rGetDomainSize();
540 node_iter != this->GetNodeIteratorEnd();
543 c_vector<double, SPACE_DIM> location = node_iter->rGetLocation();
546 c_vector<bool, SPACE_DIM> is_periodic = mpBoxCollection->GetIsPeriodicAllDims();
547 for (
unsigned d=0; d<SPACE_DIM; d++)
549 if ( !is_periodic(d) &&
550 ( location[d] < (domain_boundary[2*d] + mMinimumNodeDomainBoundarySeparation) || location[d] > (domain_boundary[2*d+1] - mMinimumNodeDomainBoundarySeparation) ) )
552 is_local_node_close = 1;
556 if (is_local_node_close)
563 int is_any_node_close = 0;
564 MPI_Allreduce(&is_local_node_close, &is_any_node_close, 1, MPI_INT, MPI_SUM,
PetscTools::GetWorld());
566 return (is_any_node_close > 0);
569template<
unsigned SPACE_DIM>
574 delete mpBoxCollection;
576 mpBoxCollection =
nullptr;
579template<
unsigned SPACE_DIM>
582 this->SetUpBoxCollection(maxInteractionDistance, domainSize);
585template<
unsigned SPACE_DIM>
588 ClearBoxCollection();
597 swell_factor = (1 + swell_factor) * 1e-14;
599 c_vector<double, 2*SPACE_DIM> domain_size;
600 for (
unsigned i=0; i < SPACE_DIM; i++)
605 SetUpBoxCollection(mMaximumInteractionDistance, domain_size);
608template<
unsigned SPACE_DIM>
611 ClearBoxCollection();
613 bool isPeriodicInX =
false;
614 bool isPeriodicInY =
false;
615 bool isPeriodicInZ =
false;
617 for (
unsigned i=0; i<SPACE_DIM; i++ )
621 isPeriodicInX = isDimPeriodic(i);
625 isPeriodicInY = isDimPeriodic(i);
629 isPeriodicInZ = isDimPeriodic(i);
633 mpBoxCollection->SetupLocalBoxesHalfOnly();
634 mpBoxCollection->SetCalculateNodeNeighbours(mCalculateNodeNeighbours);
637template<
unsigned SPACE_DIM>
642 node_iter != this->GetNodeIteratorEnd();
645 unsigned box_index = mpBoxCollection->CalculateContainingBox(&(*node_iter));
646 mpBoxCollection->rGetBox(box_index).AddNode(&(*node_iter));
650template<
unsigned SPACE_DIM>
654 for (
typename std::vector<boost::shared_ptr<
Node<SPACE_DIM> > >::iterator halo_node_iter = mHaloNodes.begin();
655 halo_node_iter != mHaloNodes.end();
658 unsigned box_index = mpBoxCollection->CalculateContainingBox((*halo_node_iter).get());
659 mpBoxCollection->rGetHaloBox(box_index).AddNode((*halo_node_iter).get());
663template<
unsigned SPACE_DIM>
666 assert(mpBoxCollection);
669 mpBoxCollection->EmptyBoxes();
673 mpBoxCollection->UpdateHaloBoxes();
676template<
unsigned SPACE_DIM>
679 if (!mpBoxCollection)
681 SetUpBoxCollection(this->mNodes);
684 while (IsANodeCloseToDomainBoundary())
686 EnlargeBoxCollection();
690template<
unsigned SPACE_DIM>
693 return mpBoxCollection->GetIsPeriodicAcrossProcs();
696template<
unsigned SPACE_DIM>
699 std::vector<int> local_node_distribution = mpBoxCollection->CalculateNumberOfNodesInEachStrip();
701 unsigned new_rows = mpBoxCollection->LoadBalance(local_node_distribution);
703 c_vector<double, 2*SPACE_DIM> current_domain_size = mpBoxCollection->rGetDomainSize();
706 double fudge = 1e-14;
707 for (
unsigned d=0; d < SPACE_DIM; d++)
709 current_domain_size[2*d] = current_domain_size[2*d] + fudge;
710 current_domain_size[2*d+1] = current_domain_size[2*d+1] - fudge;
712 SetUpBoxCollection(mMaximumInteractionDistance, current_domain_size, new_rows);
715template<
unsigned SPACE_DIM>
721 for (
unsigned i=0; i<this->mNodes.size(); i++)
723 this->mNodes[i]->SetIndex(GetNextAvailableIndex());
727template<
unsigned SPACE_DIM>
730 std::vector<unsigned> indices(GetNumNodes());
731 unsigned live_index=0;
732 for (
unsigned i=0; i<this->mNodes.size(); i++)
735 if (!this->mNodes[i]->IsDeleted())
737 indices[live_index] = this->mNodes[i]->GetIndex();
#define EXCEPTION(message)
#define EXPORT_TEMPLATE_CLASS_SAME_DIMS(CLASS)
virtual double GetWidth(const unsigned &rDimension) const
std::vector< Node< SPACE_DIM > * > mNodes
const ChastePoint< SPACE_DIM > & rGetUpperCorner() const
const ChastePoint< SPACE_DIM > & rGetLowerCorner() const
unsigned GetLongestAxis() const
double GetWidth(unsigned rDimension) const
void SetDeleted(unsigned index)
void SetRegion(unsigned region)
void AddNodeAttribute(double attribute)
std::vector< double > & rGetNodeAttributes()
void SetIndex(unsigned index)
void SetRadius(double radius)
unsigned GetIndex() const
void SetIsParticle(bool isParticle)
void ClearBoxCollection()
void AddMovedNode(boost::shared_ptr< Node< SPACE_DIM > > pMovedNode)
void ConstructNodesWithoutMesh(const std::vector< Node< SPACE_DIM > * > &rNodes, double maxInteractionDistance)
double GetMaximumInteractionDistance()
void EnlargeBoxCollection()
unsigned GetNumNodes() const
virtual unsigned GetMaximumNodeIndex()
void AddNodeWithFixedIndex(Node< SPACE_DIM > *pNewNode)
void SetMaximumInteractionDistance(double maxDistance)
unsigned AddNode(Node< SPACE_DIM > *pNewNode)
void DeleteNode(unsigned index)
unsigned GetNextAvailableIndex()
void UpdateBoxCollection()
std::vector< unsigned > GetAllNodeIndices() const
void CalculateNodesOutsideLocalDomain()
void SetInitialBoxCollection(const c_vector< double, 2 *SPACE_DIM > domainSize, double maxInteractionDistance)
DistributedBoxCollection< SPACE_DIM > * GetBoxCollection()
std::vector< bool > & rGetInitiallyOwnedNodes()
void SetUpBoxCollection(const std::vector< Node< SPACE_DIM > * > &rNodes)
void SetCalculateNodeNeighbours(bool calculateNodeNeighbours)
void ConstructFromMeshReader(AbstractMeshReader< SPACE_DIM, SPACE_DIM > &rMeshReader)
void SetNode(unsigned nodeIndex, ChastePoint< SPACE_DIM > point, bool concreteMove=false)
Node< SPACE_DIM > * GetNodeOrHaloNode(unsigned index) const
bool IsANodeCloseToDomainBoundary()
std::vector< unsigned > & rGetHaloNodesToSendRight()
double GetWidth(const unsigned &rDimension) const
void RemoveDeletedNodes(NodeMap &map)
std::vector< unsigned > & rGetNodesToSendLeft()
unsigned SolveNodeMapping(unsigned index) const
void CalculateBoundaryNodePairs(std::vector< std::pair< Node< SPACE_DIM > *, Node< SPACE_DIM > * > > &rNodePairs)
void ResizeBoxCollection()
void SetMinimumNodeDomainBoundarySeparation(double separation)
void AddHaloNode(boost::shared_ptr< Node< SPACE_DIM > > pNewNode)
std::vector< unsigned > & rGetHaloNodesToSendLeft()
void AddHaloNodesToBoxes()
void DeleteMovedNode(unsigned index)
void CalculateInteriorNodePairs(std::vector< std::pair< Node< SPACE_DIM > *, Node< SPACE_DIM > * > > &rNodePairs)
bool GetIsPeriodicAcrossProcsFromBoxCollection() const
bool IsOwned(c_vector< double, SPACE_DIM > &location)
std::vector< unsigned > & rGetNodesToSendRight()
void ConstructFromMeshReader(AbstractMeshReader< ELEMENT_DIM, SPACE_DIM > &rMeshReader)