Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
TriangleInterface.cpp
Go to the documentation of this file.
1 ////////////////////////////////////////////////////////////////////////////////
2 //
3 // File: TriangleInterface.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: Interface to triangle mesher
33 //
34 ////////////////////////////////////////////////////////////////////////////////
35 
37 
38 #include <sstream>
39 
40 using namespace std;
41 namespace Nektar
42 {
43 namespace NekMeshUtils
44 {
45 
46 void TriangleInterface::Mesh(bool Quiet, bool Quality)
47 {
48  SetUp();
49 
50  int numPoints = 0;
51  int numSeg = 0;
52  for (int i = 0; i < m_boundingloops.size(); i++)
53  {
54  numSeg += m_boundingloops[i].size();
55  }
56  numPoints = numSeg + m_stienerpoints.size();
57 
58  stringstream ss;
59  ss << "3 points required for triangulation, " << numPoints << " provided";
60 
61  ASSERTL0(numPoints > 2, ss.str());
62 
63  in.numberofpoints = numPoints;
64  in.numberofpointattributes = 0;
65  in.pointlist = new REAL[in.numberofpoints * 2];
66  in.pointmarkerlist = new int [in.numberofpoints];
67 
68  int pointc = 0;
69 
70  for (int i = 0; i < m_boundingloops.size(); i++)
71  {
72  for (int j = 0; j < m_boundingloops[i].size(); j++, pointc++)
73  {
74  nodemap[pointc] = m_boundingloops[i][j];
75 
77  m_boundingloops[i][j]->GetCADSurfInfo(sid);
78  in.pointlist[pointc * 2 + 0] = uv[0] * m_str;
79  in.pointlist[pointc * 2 + 1] = uv[1];
80  }
81  }
82 
83  for (int i = 0; i < m_stienerpoints.size(); i++, pointc++)
84  {
85  nodemap[pointc] = m_stienerpoints[i];
86 
87  Array<OneD, NekDouble> uv = m_stienerpoints[i]->GetCADSurfInfo(sid);
88  in.pointlist[pointc * 2 + 0] = uv[0] * m_str;
89  in.pointlist[pointc * 2 + 1] = uv[1];
90  }
91 
92  in.numberofsegments = numSeg;
93  in.segmentlist = new int[in.numberofsegments * 2];
94  pointc = 0;
95  for (int i = 0; i < m_boundingloops.size(); i++, pointc++)
96  {
97  int first = pointc;
98  for (int j = 0; j < m_boundingloops[i].size() - 1; j++, pointc++)
99  {
100  in.segmentlist[pointc * 2 + 0] = pointc;
101  in.segmentlist[pointc * 2 + 1] = pointc + 1;
102  }
103  in.segmentlist[pointc * 2 + 0] = pointc;
104  in.segmentlist[pointc * 2 + 1] = first;
105  }
106 
107  in.numberofregions = 0;
108  in.numberofholes = m_centers.size() - 1;
109  in.holelist = new REAL[in.numberofholes * 2];
110 
111  for (int i = 1; i < m_centers.size(); i++)
112  {
113  in.holelist[(i - 1) * 2 + 0] = m_centers[i][0] * m_str;
114  in.holelist[(i - 1) * 2 + 1] = m_centers[i][1];
115  }
116 
117  if (Quiet && Quality)
118  {
119  triangulate("pzenqQYY", &in, &out, NULL);
120  }
121  else if (Quiet && !Quality)
122  {
123  triangulate("pzenYYQ", &in, &out, NULL);
124  }
125  else if (!Quiet && Quality)
126  {
127  triangulate("pzenqYY", &in, &out, NULL);
128  }
129  else if (!Quiet && !Quality)
130  {
131  triangulate("pzenYY", &in, &out, NULL);
132  }
133 
134  // verify the mesh a bit
135  if (out.numberofpoints - out.numberofedges + out.numberoftriangles !=
136  2 - m_centers.size())
137  {
138  cout << endl << "epc wrong" << endl;
139  cout << out.numberofpoints - out.numberofedges + out.numberoftriangles
140  << " " << m_centers.size() << " " << sid << endl;
141  }
142 }
143 
144 void TriangleInterface::SetUp()
145 {
146  in.pointlist = (REAL *)NULL;
147  in.pointattributelist = (REAL *)NULL;
148  in.pointmarkerlist = (int *)NULL;
149  in.numberofpoints = 0;
150  in.numberofpointattributes = 0;
151  //
152  in.trianglelist = (int *)NULL;
153  in.triangleattributelist = (REAL *)NULL;
154  in.trianglearealist = (REAL *)NULL;
155  in.neighborlist = (int *)NULL;
156  in.numberoftriangles = 0;
157  in.numberofcorners = 0;
158  in.numberoftriangleattributes = 0;
159  //
160  in.segmentlist = (int *)NULL;
161  in.segmentmarkerlist = (int *)NULL;
162  in.numberofsegments = 0;
163  //
164  in.holelist = (REAL *)NULL;
165  in.numberofholes = 0;
166  //
167  in.regionlist = (REAL *)NULL;
168  in.numberofregions = 0;
169  //
170  in.edgelist = (int *)NULL;
171  in.edgemarkerlist = (int *)NULL;
172  in.normlist = (REAL *)NULL;
173  in.numberofedges = 0;
174  //
175  out.pointlist = (REAL *)NULL;
176  out.pointattributelist = (REAL *)NULL;
177  out.pointmarkerlist = (int *)NULL;
178  out.numberofpoints = 0;
179  out.numberofpointattributes = 0;
180  //
181  out.trianglelist = (int *)NULL;
182  out.triangleattributelist = (REAL *)NULL;
183  out.trianglearealist = (REAL *)NULL;
184  out.neighborlist = (int *)NULL;
185  out.numberoftriangles = 0;
186  out.numberofcorners = 0;
187  out.numberoftriangleattributes = 0;
188  //
189  out.segmentlist = (int *)NULL;
190  out.segmentmarkerlist = (int *)NULL;
191  out.numberofsegments = 0;
192  //
193  out.holelist = (REAL *)NULL;
194  out.numberofholes = 0;
195  //
196  out.regionlist = (REAL *)NULL;
197  out.numberofregions = 0;
198  //
199  out.edgelist = (int *)NULL;
200  out.edgemarkerlist = (int *)NULL;
201  out.normlist = (REAL *)NULL;
202  out.numberofedges = 0;
203 }
204 
205 void TriangleInterface::Extract(
206  std::vector<std::vector<NodeSharedPtr> > &Connec)
207 {
208  Connec.clear();
209  for (int i = 0; i < out.numberoftriangles; i++)
210  {
212  n1 = nodemap.find(out.trianglelist[i * 3 + 0]);
213  n2 = nodemap.find(out.trianglelist[i * 3 + 1]);
214  n3 = nodemap.find(out.trianglelist[i * 3 + 2]);
215 
216  ASSERTL0(n1 != nodemap.end() && n2 != nodemap.end() &&
217  n3 != nodemap.end(),
218  "node index error");
219 
220  vector<NodeSharedPtr> tri(3);
221  tri[0] = n1->second;
222  tri[1] = n2->second;
223  tri[2] = n3->second;
224  Connec.push_back(tri);
225  }
226 }
227 }
228 }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
STL namespace.
#define REAL
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator