43#include "Toroidal2dMesh.hpp"
45#include "VtkMeshWriter.hpp"
68 for (
unsigned index=0; index<nodes.size(); index++)
70 Node<2>* p_temp_node = nodes[index];
73 assert( 0 <= x && x < width);
76 assert( 0 <= y && y < depth);
77 mNodes.push_back(p_temp_node);
111 c_vector<double, 2> location;
112 unsigned this_node_index = node_iter->GetIndex();
115 assert(0.0 <= node_iter->rGetLocation()[1]);
116 assert(node_iter->rGetLocation()[1] <=
mHeight);
126 c_vector<double, 2> location;
128 location[1] = location[1] +
mHeight;
138 c_vector<double, 2> location;
140 location[1] = location[1] -
mHeight;
157 c_vector<double, 2> location;
158 unsigned this_node_index = node_iter->GetIndex();
161 assert(0.0 <= node_iter->rGetLocation()[0]);
162 assert(node_iter->rGetLocation()[0] <=
mWidth);
171 c_vector<double, 2> location;
173 location[0] = location[0] +
mWidth;
183 c_vector<double, 2> location;
185 location[0] = location[0] -
mWidth;
203 rMap.
Resize(old_num_all_nodes);
207 for (
unsigned i=0; i<old_num_all_nodes; i++)
209 if (
mNodes[i]->IsDeleted())
242 double left_boundary = -0.5*
mWidth;
243 double right_boundary = 1.5*
mWidth;
244 double bottom_boundary = -0.5*
mHeight;
245 double top_boundary = 1.5*
mHeight;
251 unsigned num_nodes_outside = 0;
252 for (
unsigned j=0; j<3; j++)
254 Node<2>* p_node = this->
GetNode(elem_iter->GetNodeGlobalIndex(j));
256 c_vector<double, 2> location;
258 double this_node_x_location = location[0];
259 double this_node_y_location = location[1];
261 if ((this_node_x_location < left_boundary) ||
262 (this_node_x_location > right_boundary) ||
263 (this_node_y_location < bottom_boundary) ||
264 (this_node_y_location > top_boundary))
269 if (num_nodes_outside==3)
271 elem_iter->MarkAsDeleted();
280 unsigned num_nodes_outside = 0;
281 for (
unsigned j=0; j<2; j++)
283 Node<2>* p_node = this->
GetNode((*elem_iter)->GetNodeGlobalIndex(j));
285 c_vector<double, 2> location;
287 double this_node_x_location = location[0];
288 double this_node_y_location = location[1];
290 if ((this_node_x_location < left_boundary) ||
291 (this_node_x_location > right_boundary) ||
292 (this_node_y_location < bottom_boundary) ||
293 (this_node_y_location > top_boundary))
298 if (num_nodes_outside==2)
300 (*elem_iter)->MarkAsDeleted();
392 bool boundary_element_made =
false;
393 unsigned elem_index = 0;
395 while (elem_index<num_elements && !boundary_element_made)
400 boundary_element_made =
true;
401 std::vector<Node<2>*> nodes;
402 nodes.push_back(p_element->
GetNode(0));
403 nodes.push_back(p_element->
GetNode(1));
422 for (
unsigned i=0; i<rMap.
GetSize(); i++)
462 unsigned number_of_left_image_nodes = 0;
463 unsigned number_of_right_image_nodes = 0;
465 for (
unsigned i=0; i<3; i++)
467 unsigned this_node_index = elem_iter->GetNodeGlobalIndex(i);
471 number_of_left_image_nodes++;
475 number_of_right_image_nodes++;
480 if (number_of_right_image_nodes >= 1)
482 elem_iter->MarkAsDeleted();
487 if (number_of_left_image_nodes == 3)
489 elem_iter->MarkAsDeleted();
498 if (number_of_left_image_nodes == 1 || number_of_left_image_nodes == 2)
500 for (
unsigned i=0; i<3; i++)
502 unsigned this_node_index = elem_iter->GetNodeGlobalIndex(i);
506 elem_iter->ReplaceNode(
mNodes[this_node_index],
mNodes[it->second]);
599 unsigned number_of_bottom_image_nodes = 0;
600 unsigned number_of_top_image_nodes = 0;
602 for (
unsigned i=0; i<3; i++)
604 unsigned this_node_index = elem_iter->GetNodeGlobalIndex(i);
608 number_of_bottom_image_nodes++;
612 number_of_top_image_nodes++;
617 if (number_of_top_image_nodes >= 1)
619 elem_iter->MarkAsDeleted();
624 if (number_of_bottom_image_nodes == 3)
626 elem_iter->MarkAsDeleted();
635 if (number_of_bottom_image_nodes == 1 || number_of_bottom_image_nodes == 2)
637 for (
unsigned i=0; i<3; i++)
639 unsigned this_node_index = elem_iter->GetNodeGlobalIndex(i);
643 elem_iter->ReplaceNode(
mNodes[this_node_index],
mNodes[it->second]);
727 c_vector<double, 2> vector = rLocation2 - rLocation1;
728 vector[0] = fmod(vector[0],
mWidth);
729 vector[1] = fmod(vector[1],
mHeight);
735 if (vector[0] > 0.5*
mWidth)
739 else if (vector[0] < -0.5*
mWidth)
747 else if (vector[1] < -0.5*
mHeight)
788 assert(rDimension==0 || rDimension==1);
806 SetNode(node_index, new_node_point,
false);
835 unsigned elem_index = *left_iter;
843 c_vector<unsigned,3> original_element_node_indices;
844 c_vector<unsigned,3> corresponding_element_node_indices;
845 for (
unsigned i=0; i<3; i++)
856 unsigned corresponding_elem_index = *right_iter;
860 bool is_corresponding_node =
true;
862 for (
unsigned i=0; i<3; i++)
864 if ((corresponding_element_node_indices[i] != p_corresponding_element->
GetNodeGlobalIndex(0)) &&
865 (corresponding_element_node_indices[i] != p_corresponding_element->
GetNodeGlobalIndex(1)) &&
866 (corresponding_element_node_indices[i] != p_corresponding_element->
GetNodeGlobalIndex(2)) )
868 is_corresponding_node =
false;
873 if (is_corresponding_node)
897 if (temp_left_hand_side_elements.empty() || temp_right_hand_side_elements.empty())
905 if (temp_right_hand_side_elements.size() == 2 && temp_left_hand_side_elements.size() == 2)
941 unsigned elem_index = *bottom_iter;
949 c_vector<unsigned,3> original_element_node_indices;
950 c_vector<unsigned,3> corresponding_element_node_indices;
951 for (
unsigned i=0; i<3; i++)
962 unsigned corresponding_elem_index = *top_iter;
966 bool is_corresponding_node =
true;
968 for (
unsigned i=0; i<3; i++)
970 if ((corresponding_element_node_indices[i] != p_corresponding_element->
GetNodeGlobalIndex(0)) &&
971 (corresponding_element_node_indices[i] != p_corresponding_element->
GetNodeGlobalIndex(1)) &&
972 (corresponding_element_node_indices[i] != p_corresponding_element->
GetNodeGlobalIndex(2)) )
974 is_corresponding_node =
false;
979 if (is_corresponding_node)
987 temp_top_hand_side_elements.erase(corresponding_elem_index);
1004 if (temp_bottom_hand_side_elements.empty() || temp_top_hand_side_elements.empty())
1012 if (temp_top_hand_side_elements.size() == 2 && temp_bottom_hand_side_elements.size() == 2)
1029 unsigned incidences_of_zero_left_image_nodes = 0;
1030 unsigned incidences_of_zero_right_image_nodes = 0;
1037 unsigned number_of_left_image_nodes = 0;
1038 unsigned number_of_right_image_nodes = 0;
1040 for (
unsigned i=0; i<3; i++)
1042 unsigned this_node_index = elem_iter->GetNodeGlobalIndex(i);
1046 number_of_left_image_nodes++;
1050 number_of_right_image_nodes++;
1054 if ((number_of_left_image_nodes == 0) && (number_of_right_image_nodes == 1 || number_of_right_image_nodes == 2) )
1056 incidences_of_zero_left_image_nodes++;
1058 if ((number_of_right_image_nodes == 0) && (number_of_left_image_nodes == 1 || number_of_left_image_nodes == 2) )
1060 incidences_of_zero_right_image_nodes++;
1074 if (number_of_right_image_nodes == 1 || number_of_right_image_nodes == 2)
1080 if (number_of_left_image_nodes == 1 || number_of_left_image_nodes == 2)
1122 unsigned incidences_of_zero_bottom_image_nodes = 0;
1123 unsigned incidences_of_zero_top_image_nodes = 0;
1130 unsigned number_of_bottom_image_nodes = 0;
1131 unsigned number_of_top_image_nodes = 0;
1133 for (
unsigned i=0; i<3; i++)
1135 unsigned this_node_index = elem_iter->GetNodeGlobalIndex(i);
1139 number_of_bottom_image_nodes++;
1143 number_of_top_image_nodes++;
1147 if ((number_of_bottom_image_nodes == 0) && (number_of_top_image_nodes == 1 || number_of_top_image_nodes == 2) )
1149 incidences_of_zero_bottom_image_nodes++;
1151 if ((number_of_top_image_nodes == 0) && (number_of_bottom_image_nodes == 1 || number_of_bottom_image_nodes == 2) )
1153 incidences_of_zero_top_image_nodes++;
1157 if (number_of_top_image_nodes == 1 || number_of_top_image_nodes == 2)
1163 if (number_of_bottom_image_nodes == 1 || number_of_bottom_image_nodes == 2)
1182 unsigned corresponding_node_index = UINT_MAX;
1224 assert(corresponding_node_index != UINT_MAX);
1225 return corresponding_node_index;
1231 unsigned corresponding_node_index = UINT_MAX;
1242 std::vector<unsigned>::iterator top_im_iter = std::find(
mTopImages.begin(),
mTopImages.end(), nodeIndex);
1273 assert(corresponding_node_index != UINT_MAX);
1274 return corresponding_node_index;
1281 unsigned num_nodes =
mNodes.size();
1282 for (
unsigned i=0; i<num_nodes; i++)
1284 double& x_location = (
mNodes[i]->rGetModifiableLocation())[0];
1285 if (x_location < 0.0)
1289 else if (x_location >=
mWidth)
1291 x_location = fmod(x_location,
mWidth);
1293 double& y_location = (
mNodes[i]->rGetModifiableLocation())[1];
1294 if (y_location < 0.0)
1298 else if (y_location >=
mHeight)
1300 y_location = fmod(y_location,
mHeight);
#define CHASTE_CLASS_EXPORT(T)
Node< SPACE_DIM > * GetNode(unsigned localIndex) const
unsigned GetNodeGlobalIndex(unsigned localIndex) const
virtual unsigned GetNumAllNodes() const
NodeIterator GetNodeIteratorEnd()
Node< SPACE_DIM > * GetNode(unsigned index) const
std::vector< Node< SPACE_DIM > * > mNodes
NodeIterator GetNodeIteratorBegin(bool skipDeletedNodes=true)
ElementIterator GetElementIteratorBegin(bool skipDeletedElements=true)
std::vector< BoundaryElement< ELEMENT_DIM-1, SPACE_DIM > * >::const_iterator BoundaryElementIterator
ElementIterator GetElementIteratorEnd()
Element< ELEMENT_DIM, SPACE_DIM > * GetElement(unsigned index) const
BoundaryElementIterator GetBoundaryElementIteratorBegin() const
unsigned GetNumAllElements() const
BoundaryElementIterator GetBoundaryElementIteratorEnd() const
std::vector< BoundaryElement< ELEMENT_DIM-1, SPACE_DIM > * > mBoundaryElements
c_vector< double, DIM > & rGetLocation()
void SetCoordinate(unsigned i, double value)
void ReIndex(NodeMap &map)
virtual void SetNode(unsigned index, ChastePoint< SPACE_DIM > point, bool concreteMove=true)
std::vector< unsigned > mDeletedBoundaryElementIndices
std::vector< unsigned > mDeletedElementIndices
std::vector< unsigned > mDeletedNodeIndices
virtual unsigned AddNode(Node< SPACE_DIM > *pNewNode)
void SetDeleted(unsigned index)
unsigned GetNewIndex(unsigned oldIndex) const
void SetNewIndex(unsigned oldIndex, unsigned newIndex)
void Resize(unsigned size)
bool IsDeleted(unsigned index)
const c_vector< double, SPACE_DIM > & rGetLocation() const
ChastePoint< SPACE_DIM > GetPoint() const
std::vector< double > mBoundaryElementJacobianDeterminants
std::vector< c_vector< double, SPACE_DIM > > mBoundaryElementWeightedDirections
std::vector< unsigned > mBottomOriginals
std::map< unsigned, unsigned > mImageToRightOriginalNodeMap
void GenerateVectorsOfElementsStraddlingToroidalPeriodicBoundaries()
std::vector< unsigned > mLeftImages
std::vector< unsigned > mRightOriginals
Toroidal2dMesh(double width, double depth)
std::set< unsigned > mBottomPeriodicBoundaryElementIndices
void ReconstructToroidalMesh()
std::map< unsigned, unsigned > mImageToTopOriginalNodeMap
std::set< unsigned > mTopPeriodicBoundaryElementIndices
std::vector< unsigned > mLeftOriginals
unsigned GetCorrespondingToroidalNodeIndex(unsigned nodeIndex)
void SetNode(unsigned index, ChastePoint< 2 > point, bool concreteMove)
void CorrectToroidalNonPeriodicMesh()
std::set< unsigned > mLeftPeriodicBoundaryElementIndices
std::map< unsigned, unsigned > mImageToBottomOriginalNodeMap
c_vector< double, 2 > GetVectorFromAtoB(const c_vector< double, 2 > &rLocation1, const c_vector< double, 2 > &rLocation2)
std::vector< unsigned > mRightImages
std::set< unsigned > mRightPeriodicBoundaryElementIndices
std::vector< unsigned > mTopOriginals
unsigned GetCorrespondingCylindricalNodeIndex(unsigned nodeIndex)
void ReconstructCylindricalMesh()
double GetWidth(const unsigned &rDimension) const
std::vector< unsigned > mBottomImages
std::map< unsigned, unsigned > mImageToLeftOriginalNodeMap
unsigned AddNode(Node< 2 > *pNewNode)
void CorrectCylindricalNonPeriodicMesh()
std::vector< unsigned > mTopImages
void GenerateVectorsOfElementsStraddlingCylindricalPeriodicBoundaries()