Commit 8873d0a4 authored by Michael Turner's avatar Michael Turner
Browse files

split mesh elements code

parent 389d011c
SET(MESHUTILS_SOURCES
MeshElements/MeshElements.cpp
MeshElements/Mesh.cpp
MeshElements/Element.cpp
MeshElements/BooleanOperators.cpp
MeshElements/Point.cpp
MeshElements/Line.cpp
MeshElements/Triangle.cpp
MeshElements/Quadrilateral.cpp
MeshElements/Tetrahedron.cpp
MeshElements/Pyramid.cpp
MeshElements/Prism.cpp
MeshElements/Hexahedron.cpp
)
IF(NEKTAR_USE_MESH)
SET(MESHUTILS_SOURCES ${MESHUTILS_SOURCES}
......
////////////////////////////////////////////////////////////////////////////////
//
// File: MeshElements.cpp
//
// For more information, please see: http://www.nektar.info/
//
// The MIT License
//
// Copyright (c) 2006 Division of Applied Mathematics, Brown University (USA),
// Department of Aeronautics, Imperial College London (UK), and Scientific
// Computing and Imaging Institute, University of Utah (USA).
//
// License for the specific language governing rights and limitations under
// Permission is hereby granted, free of charge, to any person obtaining a
// copy of this software and associated documentation files (the "Software"),
// to deal in the Software without restriction, including without limitation
// the rights to use, copy, modify, merge, publish, distribute, sublicense,
// and/or sell copies of the Software, and to permit persons to whom the
// Software is furnished to do so, subject to the following conditions:
//
// The above copyright notice and this permission notice shall be included
// in all copies or substantial portions of the Software.
//
// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
// OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
// THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
// DEALINGS IN THE SOFTWARE.
//
// Description: Mesh manipulation objects.
//
////////////////////////////////////////////////////////////////////////////////
#include "MeshElements.h"
using namespace std;
namespace Nektar
{
namespace MeshUtils
{
bool operator==(ElmtConfig const &c1, ElmtConfig const &c2)
{
return (c1.m_e == c2.m_e && c1.m_order == c2.m_order);
}
/**
* @brief Compares two %HOSurf objects referred to as shared pointers.
*
* Two %HOSurf objects are defined to be equal if they contain identical
* vertex ids contained in HOSurf::vertId.
*/
bool operator==(HOSurfSharedPtr const &p1, HOSurfSharedPtr const &p2)
{
if (p1->vertId.size() != p2->vertId.size())
{
return false;
}
vector<int> ids1 = p1->vertId;
vector<int> ids2 = p2->vertId;
sort(ids1.begin(), ids1.end());
sort(ids2.begin(), ids2.end());
for (int i = 0; i < ids1.size(); ++i)
{
if (ids1[i] != ids2[i])
return false;
}
return true;
}
/**
* @brief Test equality of two conditions - i.e. compare types, fields
* and values but _not_ composite ids.
*/
bool operator==(ConditionSharedPtr const &c1, ConditionSharedPtr const &c2)
{
int i, n = c1->type.size();
if (n != c2->type.size())
{
return false;
}
for (i = 0; i < n; ++i)
{
if (c1->type[i] != c2->type[i])
{
return false;
}
if (c1->field[i] != c2->field[i] || c1->value[i] != c2->value[i])
{
return false;
}
}
return true;
}
/**
* @brief Defines equality between two #NodeSharedPtr objects.
*/
bool operator==(NodeSharedPtr const &p1, NodeSharedPtr const &p2)
{
return *p1 == *p2;
}
bool operator!=(NodeSharedPtr const &p1, NodeSharedPtr const &p2)
{
if(p1->m_id != p2->m_id)
{
return true;
}
else
{
return false;
}
}
/**
* @brief Defines ordering between two #NodeSharedPtr objects.
*/
bool operator< (NodeSharedPtr const &p1, NodeSharedPtr const &p2)
{
return *p1 < *p2;
}
std::ostream &operator<<(std::ostream &os, const NodeSharedPtr &n)
{
os << n->m_x << " " << n->m_y << " " << n->m_z;
return os;
}
/**
* @brief Defines equality of two edges (equal if IDs of end nodes
* match in either ordering).
*/
bool operator==(EdgeSharedPtr const &p1, EdgeSharedPtr const &p2)
{
return ( ((*(p1->m_n1) == *(p2->m_n1)) && (*(p1->m_n2) == *(p2->m_n2)))
|| ((*(p1->m_n2) == *(p2->m_n1)) && (*(p1->m_n1) == *(p2->m_n2))));
}
/**
* @brief Defines ordering between two edges (based on ID of edges).
*/
bool operator< (EdgeSharedPtr const &p1, EdgeSharedPtr const &p2)
{
return p1->m_id < p2->m_id;
}
/**
* @brief Defines equality of two faces (equal if IDs of vertices are
* the same.)
*/
bool operator==(FaceSharedPtr const &p1, FaceSharedPtr const &p2)
{
std::vector<NodeSharedPtr>::iterator it1, it2;
for (it1 = p1->m_vertexList.begin(); it1 != p1->m_vertexList.end(); ++it1)
{
if (find(p2->m_vertexList.begin(), p2->m_vertexList.end(), *it1)
== p2->m_vertexList.end())
{
return false;
}
}
return true;
}
/**
* @brief Defines ordering between two faces (depending on ID of
* faces).
*/
bool operator< (FaceSharedPtr const &p1, FaceSharedPtr const &p2)
{
return p1->m_id < p2->m_id;
}
}
}
......@@ -49,8 +49,62 @@ namespace MeshUtils
public:
Composite() : m_reorder(true) {}
/// Generate the list of IDs of elements within this composite.
std::string GetXmlString(bool doSort=true);
/**
* @brief Generate a Nektar++ string describing the composite.
*
* The list of composites may include individual element IDs or ranges
* of element IDs.
*/
std::string GetXmlString(bool doSort=true)
{
#if 0 // turn this option off since causes problem with InputNekpp.cpp
if (doSort)
{
element_id_less_than sortOperator;
sort(m_items.begin(), m_items.end(), sortOperator);
}
#endif
stringstream st;
vector<ElementSharedPtr>::iterator it;
bool range = false;
int vId = m_items[0]->GetId();
int prevId = vId;
st << " " << m_tag << "[" << vId;
for (it = m_items.begin()+1; it != m_items.end(); ++it){
// store previous element ID and get current one
prevId = vId;
vId = (*it)->GetId();
// continue an already started range
if (prevId > -1 && vId == prevId + 1)
{
range = true;
// if this is the last element, it's the end of a range, so write
if (*it == m_items.back())
{
st << "-" << vId;
}
continue;
}
// terminate a range, if present
if (range)
{
st << "-" << prevId;
range = false;
}
// write what will be either a single entry or start of new range
st << "," << vId;
}
// terminate
st << "] ";
return st.str();
}
/// ID of composite.
unsigned int m_id;
......
////////////////////////////////////////////////////////////////////////////////
//
// File: MeshElements.cpp
//
// For more information, please see: http://www.nektar.info/
//
// The MIT License
//
// Copyright (c) 2006 Division of Applied Mathematics, Brown University (USA),
// Department of Aeronautics, Imperial College London (UK), and Scientific
// Computing and Imaging Institute, University of Utah (USA).
//
// License for the specific language governing rights and limitations under
// Permission is hereby granted, free of charge, to any person obtaining a
// copy of this software and associated documentation files (the "Software"),
// to deal in the Software without restriction, including without limitation
// the rights to use, copy, modify, merge, publish, distribute, sublicense,
// and/or sell copies of the Software, and to permit persons to whom the
// Software is furnished to do so, subject to the following conditions:
//
// The above copyright notice and this permission notice shall be included
// in all copies or substantial portions of the Software.
//
// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
// OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
// THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
// DEALINGS IN THE SOFTWARE.
//
// Description: Mesh manipulation objects.
//
////////////////////////////////////////////////////////////////////////////////
#include <loki/Singleton.h>
#include "MeshElements.h"
using namespace std;
namespace Nektar
{
namespace MeshUtils
{
ElementFactory& GetElementFactory()
{
typedef Loki::SingletonHolder<ElementFactory,
Loki::CreateUsingNew,
Loki::NoDestroy > Type;
return Type::Instance();
}
Element::Element(ElmtConfig pConf,
unsigned int pNumNodes,
unsigned int pGotNodes)
: m_conf(pConf),
m_curveType(LibUtilities::ePolyEvenlySpaced),
m_geom()
{
if (pNumNodes != pGotNodes)
{
cerr << "Number of modes mismatch for type "
<< pConf.m_e << "! Should be " << pNumNodes
<< " but got " << pGotNodes << " nodes." << endl;
abort();
}
}
/**
* @brief Replace a vertex in the element.
*
* When a vertex is replaced, the element edges and faces are also
* searched and the corresponding edge/face nodes are updated to
* maintain consistency.
*
* @param p Index of the vertex to replace.
* @param pNew New vertex.
*/
void Element::SetVertex(unsigned int p, NodeSharedPtr pNew)
{
NodeSharedPtr vOld = m_vertex[p];
m_vertex[p] = pNew;
for (unsigned int i = 0; i < m_edge.size(); ++i)
{
if (m_edge[i]->m_n1 == vOld)
{
m_edge[i]->m_n1 = pNew;
}
else if (m_edge[i]->m_n2 == vOld)
{
m_edge[i]->m_n2 = pNew;
}
}
for (unsigned int i = 0; i < m_face.size(); ++i)
{
// Replace vertices in faces
for (unsigned int j = 0; j < m_face[i]->m_vertexList.size(); ++j)
{
if (m_face[i]->m_vertexList[j] == vOld)
{
m_face[i]->m_vertexList[j] = pNew;
}
}
for (unsigned int j = 0; j < m_face[i]->m_edgeList.size(); ++j)
{
if (m_face[i]->m_edgeList[j]->m_n1 == vOld)
{
m_face[i]->m_edgeList[j]->m_n1 = pNew;
}
else if (m_face[i]->m_edgeList[j]->m_n2 == vOld)
{
m_face[i]->m_edgeList[j]->m_n2 = pNew;
}
}
}
}
/**
* @brief Replace an edge in the element.
*
* When an edge is replaced, the element faces are also searched and
* the corresponding face edges are updated to maintain consistency.
*
* @param p Index of the edge to replace.
* @param pNew New edge.
*/
void Element::SetEdge(unsigned int p, EdgeSharedPtr pNew)
{
EdgeSharedPtr vOld = m_edge[p];
m_edge[p] = pNew;
for (unsigned int i = 0; i < m_face.size(); ++i)
{
for (unsigned int j = 0; j < m_face[i]->m_edgeList.size(); ++j)
{
if (m_face[i]->m_edgeList[j] == vOld)
{
m_face[i]->m_edgeList[j] = pNew;
}
}
}
}
/**
* @brief Replace a face in the element.
*
* When a face is replaced, no other consistency checks are required.
*
* @param p Index of the face to replace.
* @param pNew New face.
*/
void Element::SetFace(unsigned int p, FaceSharedPtr pNew)
{
m_face[p] = pNew;
}
/**
* @brief Obtain the order of an element by looking at edges.
*/
int Element::GetMaxOrder()
{
int i, ret = 1;
for (i = 0; i < m_edge.size(); ++i)
{
int edgeOrder = m_edge[i]->GetNodeCount()-1;
if (edgeOrder > ret)
{
ret = edgeOrder;
}
}
return ret;
}
}
}
////////////////////////////////////////////////////////////////////////////////
//
// File: MeshElements.cpp
//
// For more information, please see: http://www.nektar.info/
//
// The MIT License
//
// Copyright (c) 2006 Division of Applied Mathematics, Brown University (USA),
// Department of Aeronautics, Imperial College London (UK), and Scientific
// Computing and Imaging Institute, University of Utah (USA).
//
// License for the specific language governing rights and limitations under
// Permission is hereby granted, free of charge, to any person obtaining a
// copy of this software and associated documentation files (the "Software"),
// to deal in the Software without restriction, including without limitation
// the rights to use, copy, modify, merge, publish, distribute, sublicense,
// and/or sell copies of the Software, and to permit persons to whom the
// Software is furnished to do so, subject to the following conditions:
//
// The above copyright notice and this permission notice shall be included
// in all copies or substantial portions of the Software.
//
// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
// OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
// THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
// DEALINGS IN THE SOFTWARE.
//
// Description: Mesh manipulation objects.
//
////////////////////////////////////////////////////////////////////////////////
#include <SpatialDomains/HexGeom.h>
#include "MeshElements.h"
using namespace std;
namespace Nektar
{
namespace MeshUtils
{
LibUtilities::ShapeType Hexahedron::m_type = GetElementFactory().
RegisterCreatorFunction(LibUtilities::eHexahedron, Hexahedron::create, "Hexahedron");
/**
* @brief Create a hexahedral element.
*/
Hexahedron::Hexahedron(ElmtConfig pConf,
vector<NodeSharedPtr> pNodeList,
vector<int> pTagList)
: Element(pConf, GetNumNodes(pConf), pNodeList.size())
{
m_tag = "H";
m_dim = 3;
m_taglist = pTagList;
int n = m_conf.m_order-1;
// Create a map to relate edge nodes to a pair of vertices defining an edge
// This is based on the ordering produced by gmsh.
map<pair<int,int>, int> edgeNodeMap;
map<pair<int,int>, int>::iterator it;
edgeNodeMap[pair<int,int>(1,2)] = 9;
edgeNodeMap[pair<int,int>(2,3)] = 9 + n;
edgeNodeMap[pair<int,int>(3,4)] = 9 + 2*n;
edgeNodeMap[pair<int,int>(4,1)] = 9 + 3*n;
edgeNodeMap[pair<int,int>(1,5)] = 9 + 4*n;
edgeNodeMap[pair<int,int>(2,6)] = 9 + 5*n;
edgeNodeMap[pair<int,int>(3,7)] = 9 + 6*n;
edgeNodeMap[pair<int,int>(4,8)] = 9 + 7*n;
edgeNodeMap[pair<int,int>(5,6)] = 9 + 8*n;
edgeNodeMap[pair<int,int>(6,7)] = 9 + 9*n;
edgeNodeMap[pair<int,int>(7,8)] = 9 + 10*n;
edgeNodeMap[pair<int,int>(8,5)] = 9 + 11*n;
// Add vertices
for (int i = 0; i < 8; ++i) {
m_vertex.push_back(pNodeList[i]);
}
// Create edges (with corresponding set of edge points)
for (it = edgeNodeMap.begin(); it != edgeNodeMap.end(); ++it)
{
vector<NodeSharedPtr> edgeNodes;
if (m_conf.m_order > 1) {
for (int j = it->second; j < it->second + n; ++j) {
edgeNodes.push_back(pNodeList[j-1]);
}
}
m_edge.push_back(EdgeSharedPtr(new Edge(pNodeList[it->first.first-1],
pNodeList[it->first.second-1],
edgeNodes,
m_conf.m_edgeCurveType)));
}
// Create faces
int face_edges[6][4];
int face_ids[6][4] = {
{0,1,2,3},{0,1,5,4},{1,2,6,5},{3,2,6,7},{0,3,7,4},{4,5,6,7}};
for (int j = 0; j < 6; ++j)
{
vector<NodeSharedPtr> faceVertices;
vector<EdgeSharedPtr> faceEdges;
vector<NodeSharedPtr> faceNodes;
for (int k = 0; k < 4; ++k)
{
faceVertices.push_back(m_vertex[face_ids[j][k]]);
NodeSharedPtr a = m_vertex[face_ids[j][k]];
NodeSharedPtr b = m_vertex[face_ids[j][(k+1)%4]];
for (unsigned int i = 0; i < m_edge.size(); ++i)
{
if ( ((*(m_edge[i]->m_n1)==*a) && (*(m_edge[i]->m_n2)==*b))
|| ((*(m_edge[i]->m_n1)==*b) && (*(m_edge[i]->m_n2) == *a)) )
{
face_edges[j][k] = i;
faceEdges.push_back(m_edge[i]);
break;
}
}
}
if (m_conf.m_faceNodes)
{
int N = 8 + 12*n + j*n*n;
for (int i = 0; i < n*n; ++i)
{
faceNodes.push_back(pNodeList[N+i]);
}
}
m_face.push_back(FaceSharedPtr(
new Face(faceVertices, faceNodes, faceEdges, m_conf.m_faceCurveType)));
}
// Reorder edges to be consistent with Nektar++ ordering.
vector<EdgeSharedPtr> tmp(12);
tmp[ 0] = m_edge[face_edges[0][0]];
tmp[ 1] = m_edge[face_edges[0][1]];
tmp[ 2] = m_edge[face_edges[0][2]];
tmp[ 3] = m_edge[face_edges[0][3]];
tmp[ 4] = m_edge[face_edges[1][3]];
tmp[ 5] = m_edge[face_edges[1][1]];
tmp[ 6] = m_edge[face_edges[2][1]];
tmp[ 7] = m_edge[face_edges[4][1]];
tmp[ 8] = m_edge[face_edges[5][0]];
tmp[ 9] = m_edge[face_edges[5][1]];
tmp[10] = m_edge[face_edges[5][2]];
tmp[11] = m_edge[face_edges[5][3]];
m_edge = tmp;