Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
BLMesh.cpp
Go to the documentation of this file.
1 ////////////////////////////////////////////////////////////////////////////////
2 //
3 // File: TetMesh.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: tet meshing methods
33 //
34 ////////////////////////////////////////////////////////////////////////////////
35 
38 
39 using namespace std;
40 namespace Nektar
41 {
42 namespace NekMeshUtils
43 {
44 
45 void BLMesh::Mesh()
46 {
47  // At this stage the surface mesh is complete and the elements know their
48  // neighbours through element links in the edges, this includes quads.
49 
50  // here elements are made for the boundary layer they will need to know
51  // links (maybe facelinks), so that the tetmeshing module can extract the
52  // surface upon which it needs to mesh (top of the bl and the rest of the
53  // surface).
54 
55  vector<ElementSharedPtr> quad;
56  vector<ElementSharedPtr> ptri; // triangles to grow prisms onto
57 
58  for (int i = 0; i < m_mesh->m_element[2].size(); i++)
59  {
60  bool onblsurf = false;
61  for (int j = 0; j < m_blsurfs.size(); j++)
62  {
63  if (m_mesh->m_element[2][i]->CADSurfId == m_blsurfs[j])
64  {
65  onblsurf = true;
66  break;
67  }
68  }
69  if (m_mesh->m_element[2][i]->GetConf().m_e ==
71  {
72  quad.push_back(m_mesh->m_element[2][i]);
73  }
74  else if (onblsurf)
75  {
76  ptri.push_back(m_mesh->m_element[2][i]);
77  }
78  }
79 
80  map<int, NodeSharedPtr> blpair;
81  for (int i = 0; i < quad.size(); i++)
82  {
83  vector<EdgeSharedPtr> e = quad[i]->GetEdgeList();
84  for (int j = 0; j < e.size(); j++)
85  {
86  // if both or none are on curve skip
87  if ((e[j]->m_n1->GetNumCadCurve() > 0 &&
88  e[j]->m_n2->GetNumCadCurve() > 0) ||
89  (!(e[j]->m_n1->GetNumCadCurve() > 0) &&
90  !(e[j]->m_n2->GetNumCadCurve() > 0)))
91  {
92  continue;
93  }
94 
95  if (e[j]->m_n1->GetNumCadCurve() > 0)
96  {
97  blpair[e[j]->m_n1->m_id] = e[j]->m_n2;
98  }
99  else if (e[j]->m_n2->GetNumCadCurve() > 0)
100  {
101  blpair[e[j]->m_n2->m_id] = e[j]->m_n1;
102  }
103  else
104  {
105  ASSERTL0(false, "that failed");
106  }
107  }
108  }
109 
110  map<int, int> nm;
111  nm[0] = 0;
112  nm[1] = 3;
113  nm[2] = 4;
114  nm[3] = 5;
115  nm[4] = 1;
116  nm[5] = 2;
117 
118  for (int i = 0; i < ptri.size(); i++)
119  {
120  vector<NodeSharedPtr> pn(6); // all prism nodes
121  vector<NodeSharedPtr> n = ptri[i]->GetVertexList();
122 
123  vector<pair<int, CADSurfSharedPtr> > tmpss = n[0]->GetCADSurfs();
124  CADSurfSharedPtr tmps;
125 
126  for (int j = 0; j < tmpss.size(); j++)
127  {
128  if (tmpss[j].first == ptri[i]->CADSurfId)
129  {
130  tmps = tmpss[j].second;
131  break;
132  }
133  }
134 
135  if (tmps->IsReversedNormal())
136  {
137  nm[0] = 0;
138  nm[1] = 3;
139  nm[2] = 1;
140  nm[3] = 2;
141  nm[4] = 4;
142  nm[5] = 5;
143  }
144 
145  for (int j = 0; j < n.size(); j++)
146  {
147  pn[nm[j * 2]] = n[j];
148 
150  it = blpair.find(n[j]->m_id);
151  if (it != blpair.end())
152  {
153  pn[nm[j * 2 + 1]] = blpair[n[j]->m_id];
154  }
155  else
156  {
158 
159  vector<pair<int, CADSurfSharedPtr> > surfs =
160  n[j]->GetCADSurfs();
161 
162  for (int s = 0; s < surfs.size(); s++)
163  {
164  Array<OneD, NekDouble> N = surfs[s].second->N(
165  n[j]->GetCADSurfInfo(surfs[s].first));
166  for (int k = 0; k < 3; k++)
167  {
168  AN[k] += N[k];
169  }
170  }
171 
172  NekDouble mag =
173  sqrt(AN[0] * AN[0] + AN[1] * AN[1] + AN[2] * AN[2]);
174 
175  for (int k = 0; k < 3; k++)
176  {
177  AN[k] /= mag;
178  }
179 
180  Array<OneD, NekDouble> loc = n[j]->GetLoc();
182  for (int k = 0; k < 3; k++)
183  np[k] = loc[k] + AN[k] * m_bl;
184  NodeSharedPtr nn = boost::shared_ptr<Node>(
185  new Node(m_mesh->m_numNodes++, np[0], np[1], np[2]));
186  pn[nm[j * 2 + 1]] = nn;
187  blpair[n[j]->m_id] = nn;
188  }
189  }
190 
191  ElmtConfig conf(LibUtilities::ePrism, 1, false, false);
192  vector<int> tags;
193  tags.push_back(1); // all prisms are comp 1
195  LibUtilities::ePrism, conf, pn, tags);
196 
197  m_mesh->m_element[3].push_back(E);
198 
199  // need to give the surface element some information about
200  // which prism is above it
201  // so that tetmesh can infer the pseudo surface
202  vector<NodeSharedPtr> faceNodes;
203  vector<EdgeSharedPtr> edgeList = ptri[i]->GetEdgeList();
204  FaceSharedPtr F = boost::shared_ptr<Face>(new Face(
205  n, faceNodes, edgeList, ptri[i]->GetConf().m_faceCurveType));
206  vector<FaceSharedPtr> f = E->GetFaceList();
207  for (int j = 0; j < f.size(); j++)
208  {
209  if (f[j]->m_vertexList.size() != 3) // quad
210  continue;
211 
212  // only two triangle faces so if its not this one, this is the
213  // pseudo surfaces
214  if (!(F == f[j]))
215  {
216  m_surftopriface[ptri[i]->GetId()] = f[j];
217  }
218  }
219  }
220 }
221 }
222 }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
Basic information about an element.
Definition: Element.h:58
tBaseSharedPtr CreateInstance(tKey idKey BOOST_PP_COMMA_IF(MAX_PARAM) BOOST_PP_ENUM_BINARY_PARAMS(MAX_PARAM, tParam, x))
Create an instance of the class referred to by idKey.
Definition: NekFactory.hpp:162
Represents a face comprised of three or more edges.
Definition: Face.h:60
STL namespace.
ElementFactory & GetElementFactory()
Definition: Element.cpp:47
Represents a point in the domain.
Definition: Node.h:60
boost::shared_ptr< Node > NodeSharedPtr
Definition: Node.h:50
double NekDouble
boost::shared_ptr< CADSurf > CADSurfSharedPtr
Definition: CADSurf.h:217
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
boost::shared_ptr< Element > ElementSharedPtr
Definition: Edge.h:52
boost::shared_ptr< Face > FaceSharedPtr
Shared pointer to a face.
Definition: Face.h:378