Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
Quadrilateral.cpp
Go to the documentation of this file.
1 ////////////////////////////////////////////////////////////////////////////////
2 //
3 // File: Quadrilateral.cpp
4 //
5 // For more information, please see: http://www.nektar.info/
6 //
7 // The MIT License
8 //
9 // Copyright (c) 2006 Division of Applied Mathematics, Brown University (USA),
10 // Department of Aeronautics, Imperial College London (UK), and Scientific
11 // Computing and Imaging Institute, University of Utah (USA).
12 //
13 // License for the specific language governing rights and limitations under
14 // Permission is hereby granted, free of charge, to any person obtaining a
15 // copy of this software and associated documentation files (the "Software"),
16 // to deal in the Software without restriction, including without limitation
17 // the rights to use, copy, modify, merge, publish, distribute, sublicense,
18 // and/or sell copies of the Software, and to permit persons to whom the
19 // Software is furnished to do so, subject to the following conditions:
20 //
21 // The above copyright notice and this permission notice shall be included
22 // in all copies or substantial portions of the Software.
23 //
24 // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
25 // OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
26 // FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
27 // THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
28 // LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
29 // FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
30 // DEALINGS IN THE SOFTWARE.
31 //
32 // Description: Mesh quad object.
33 //
34 ////////////////////////////////////////////////////////////////////////////////
35 
36 #include <LocalRegions/QuadExp.h>
38 
40 
41 using namespace std;
42 
43 namespace Nektar
44 {
45 namespace NekMeshUtils
46 {
47 
48 LibUtilities::ShapeType Quadrilateral::m_type =
50  LibUtilities::eQuadrilateral, Quadrilateral::create, "Quadrilateral");
51 
52 /**
53  * @brief Create a quadrilateral element.
54  */
55 Quadrilateral::Quadrilateral(ElmtConfig pConf,
56  vector<NodeSharedPtr> pNodeList,
57  vector<int> pTagList)
58  : Element(pConf, GetNumNodes(pConf), pNodeList.size())
59 {
60  m_tag = "Q";
61  m_dim = 2;
62  m_taglist = pTagList;
63  int n = m_conf.m_order - 1;
64 
65  // Create a map to relate edge nodes to a pair of vertices
66  // defining an edge. This is based on the ordering produced by
67  // gmsh.
68  map<pair<int, int>, int> edgeNodeMap;
69  map<pair<int, int>, int>::iterator it;
70  edgeNodeMap[pair<int, int>(1, 2)] = 5;
71  edgeNodeMap[pair<int, int>(2, 3)] = 5 + n;
72  edgeNodeMap[pair<int, int>(3, 4)] = 5 + 2 * n;
73  edgeNodeMap[pair<int, int>(4, 1)] = 5 + 3 * n;
74 
75  // Add vertices. This logic will determine (in 2D) whether the
76  // element is clockwise (sum > 0) or counter-clockwise (sum < 0).
77  NekDouble sum = 0.0;
78  for (int i = 0; i < 4; ++i)
79  {
80  int o = (i + 1) % 4;
81  m_vertex.push_back(pNodeList[i]);
82  sum += (pNodeList[o]->m_x - pNodeList[i]->m_x) *
83  (pNodeList[o]->m_y + pNodeList[i]->m_y);
84  }
85 
86  // Create edges (with corresponding set of edge points)
87  for (it = edgeNodeMap.begin(); it != edgeNodeMap.end(); ++it)
88  {
89  vector<NodeSharedPtr> edgeNodes;
90  if (m_conf.m_order > 1)
91  {
92  for (int j = it->second; j < it->second + n; ++j)
93  {
94  edgeNodes.push_back(pNodeList[j - 1]);
95  }
96  }
97  m_edge.push_back(EdgeSharedPtr(new Edge(pNodeList[it->first.first - 1],
98  pNodeList[it->first.second - 1],
99  edgeNodes,
101  }
102 
103  if (pConf.m_reorient)
104  {
105  if (sum > 0.0)
106  {
107  reverse(m_edge.begin(), m_edge.end());
108  }
109  }
110 
111  if (m_conf.m_faceNodes)
112  {
113  m_volumeNodes.insert(m_volumeNodes.begin(),
114  pNodeList.begin() + 4 * m_conf.m_order,
115  pNodeList.end());
116  }
117 }
118 
120  int edgeId, EdgeSharedPtr edge)
121 {
122  int locVert = edgeId;
123  if (edge->m_n1 == m_vertex[locVert])
124  {
125  return StdRegions::eForwards;
126  }
127  else if (edge->m_n2 == m_vertex[locVert])
128  {
129  return StdRegions::eBackwards;
130  }
131  else
132  {
133  ASSERTL1(false, "Edge is not connected to this quadrilateral.");
134  }
135 
137 }
138 
142  int coordDim,
143  int &id,
144  bool justConfig)
145 {
146  m_conf.m_order = order;
147  m_curveType = pType;
148  m_conf.m_volumeNodes = false;
149  m_volumeNodes.clear();
150 
151  // Quadrilaterals of order == 1 have no interior volume points.
152  if (order == 1)
153  {
154  m_conf.m_faceNodes = false;
155  return;
156  }
157 
158  m_conf.m_faceNodes = true;
159 
160  if (justConfig)
161  {
162  return;
163  }
164 
165  int nPoints = order + 1;
166  StdRegions::StdExpansionSharedPtr xmap = geom->GetXmap();
167 
169  LibUtilities::PointsKey pKey(nPoints, pType);
170  ASSERTL1(pKey.GetPointsDim() == 1, "Points distribution must be 1D");
171  LibUtilities::PointsManager()[pKey]->GetPoints(px);
172 
173  Array<OneD, Array<OneD, NekDouble> > phys(coordDim);
174 
175  for (int i = 0; i < coordDim; ++i)
176  {
177  phys[i] = Array<OneD, NekDouble>(xmap->GetTotPoints());
178  xmap->BwdTrans(geom->GetCoeffs(i), phys[i]);
179  }
180 
181  int nQuadIntPts = (nPoints - 2) * (nPoints - 2);
182  m_volumeNodes.resize(nQuadIntPts);
183 
184  for (int i = 1, cnt = 0; i < nPoints-1; ++i)
185  {
186  for (int j = 1; j < nPoints-1; ++j, ++cnt)
187  {
189  xp[0] = px[j];
190  xp[1] = px[i];
191 
192  Array<OneD, NekDouble> x(3, 0.0);
193  for (int k = 0; k < coordDim; ++k)
194  {
195  x[k] = xmap->PhysEvaluate(xp, phys[k]);
196  }
197 
198  m_volumeNodes[cnt] = boost::shared_ptr<Node>(
199  new Node(id++, x[0], x[1], x[2]));
200  }
201  }
202 }
203 
205 {
209 
210  for (int i = 0; i < 4; ++i)
211  {
212  edges[i] = m_edge[i]->GetGeom(coordDim);
213  verts[i] = m_vertex[i]->GetGeom(coordDim);
214  }
215 
216  StdRegions::Orientation edgeorient[4] = {
217  SpatialDomains::SegGeom::GetEdgeOrientation(*edges[0], *edges[1]),
218  SpatialDomains::SegGeom::GetEdgeOrientation(*edges[1], *edges[2]),
219  SpatialDomains::SegGeom::GetEdgeOrientation(*edges[2], *edges[3]),
220  SpatialDomains::SegGeom::GetEdgeOrientation(*edges[3], *edges[0])};
221 
223  m_id, verts, edges, edgeorient);
224 
225  return ret;
226 }
227 
228 void Quadrilateral::GetCurvedNodes(std::vector<NodeSharedPtr> &nodeList) const
229 {
230  int n = m_edge[0]->GetNodeCount();
231  nodeList.resize(n * n);
232 
233  // Write vertices
234  nodeList[0] = m_vertex[0];
235  nodeList[n - 1] = m_vertex[1];
236  nodeList[n * n - 1] = m_vertex[2];
237  nodeList[n * (n - 1)] = m_vertex[3];
238 
239  // Write edge-interior
240  int skips[4][2] = {
241  {0, 1}, {n - 1, n}, {n * n - 1, -1}, {n * (n - 1), -n}};
242  for (int i = 0; i < 4; ++i)
243  {
244  bool reverseEdge = m_edge[i]->m_n1 == m_vertex[i];
245 
246  if (!reverseEdge)
247  {
248  for (int j = 1; j < n - 1; ++j)
249  {
250  nodeList[skips[i][0] + j * skips[i][1]] =
251  m_edge[i]->m_edgeNodes[n - 2 - j];
252  }
253  }
254  else
255  {
256  for (int j = 1; j < n - 1; ++j)
257  {
258  nodeList[skips[i][0] + j * skips[i][1]] =
259  m_edge[i]->m_edgeNodes[j - 1];
260  }
261  }
262  }
263 
264  // Write interior
265  for (int i = 1; i < n - 1; ++i)
266  {
267  for (int j = 1; j < n - 1; ++j)
268  {
269  nodeList[i * n + j] =
270  m_volumeNodes[(i - 1) * (n - 2) + (j - 1)];
271  }
272  }
273 }
274 
275 
276 /**
277  * @brief Return the number of nodes defining a quadrilateral.
278  */
280 {
281  int n = pConf.m_order;
282  if (!pConf.m_faceNodes)
283  return 4 * n;
284  else
285  return (n + 1) * (n + 1);
286 }
287 }
288 }
bool m_faceNodes
Denotes whether the element contains face nodes. For 2D elements, if this is true then the element co...
Definition: ElementConfig.h:80
Basic information about an element.
Definition: ElementConfig.h:50
Represents an edge which joins two points.
Definition: Edge.h:58
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
static NEKMESHUTILS_EXPORT unsigned int GetNumNodes(ElmtConfig pConf)
Return the number of nodes defining a quadrilateral.
virtual NEKMESHUTILS_EXPORT SpatialDomains::GeometrySharedPtr GetGeom(int coordDim)
Generate a Nektar++ geometry object for this element.
boost::shared_ptr< QuadGeom > QuadGeomSharedPtr
Definition: QuadGeom.h:54
STL namespace.
ElementFactory & GetElementFactory()
Definition: Element.cpp:47
ElmtConfig m_conf
Contains configuration of the element.
Definition: Element.h:380
static StdRegions::Orientation GetEdgeOrientation(const SegGeom &edge1, const SegGeom &edge2)
Get the orientation of edge1.
Definition: SegGeom.cpp:293
std::vector< int > m_taglist
List of integers specifying properties of the element.
Definition: Element.h:384
LibUtilities::PointsType m_edgeCurveType
Distribution of points in edges.
Definition: ElementConfig.h:92
unsigned int m_order
Order of the element.
Definition: ElementConfig.h:87
std::vector< NodeSharedPtr > m_vertex
List of element vertex nodes.
Definition: Element.h:386
unsigned int m_dim
Dimension of the element.
Definition: Element.h:378
bool m_volumeNodes
Denotes whether the element contains volume (i.e. interior) nodes. These are not supported by either ...
Definition: ElementConfig.h:85
boost::shared_ptr< SegGeom > SegGeomSharedPtr
Definition: Geometry2D.h:60
std::vector< EdgeSharedPtr > m_edge
List of element edges.
Definition: Element.h:388
unsigned int GetPointsDim() const
Definition: Points.h:149
PointsManagerT & PointsManager(void)
Defines a specification for a set of points.
Definition: Points.h:58
double NekDouble
std::vector< NodeSharedPtr > m_volumeNodes
List of element volume nodes.
Definition: Element.h:392
virtual NEKMESHUTILS_EXPORT StdRegions::Orientation GetEdgeOrient(int edgeId, EdgeSharedPtr edge)
Get the edge orientation of edge with respect to the local element, which lies at edge index edgeId...
std::string m_tag
Tag character describing the element.
Definition: Element.h:382
boost::shared_ptr< Edge > EdgeSharedPtr
Shared pointer to an edge.
Definition: Edge.h:135
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
unsigned int m_id
ID of the element.
Definition: Element.h:376
LibUtilities::PointsType m_curveType
Volume curve type.
Definition: Element.h:394
bool m_reorient
Denotes whether the element needs to be re-orientated for a spectral element framework.
Definition: ElementConfig.h:90
virtual NEKMESHUTILS_EXPORT void MakeOrder(int order, SpatialDomains::GeometrySharedPtr geom, LibUtilities::PointsType pType, int coordDim, int &id, bool justConfig=false)
Insert interior (i.e. volume) points into this element to make the geometry an order order representa...
boost::shared_ptr< StdExpansion > StdExpansionSharedPtr
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:228
boost::shared_ptr< Geometry > GeometrySharedPtr
Definition: Geometry.h:53
boost::shared_ptr< PointGeom > PointGeomSharedPtr
Definition: Geometry.h:60
Base class for element definitions.
Definition: Element.h:59
virtual NEKMESHUTILS_EXPORT void GetCurvedNodes(std::vector< NodeSharedPtr > &nodeList) const
get list of volume interior nodes
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, tDescription pDesc="")
Register a class with the factory.
Definition: NekFactory.hpp:215