Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
NodalPrismEvenlySpaced.cpp
Go to the documentation of this file.
1 ///////////////////////////////////////////////////////////////////////////////
2 //
3 // File NodalPrismEvenlySpaced.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: 3D Nodal Prism Evenly Spaced Point Definitions
33 //
34 ///////////////////////////////////////////////////////////////////////////////
35 
36 
41 #include <vector>
42 
43 namespace Nektar
44 {
45  namespace LibUtilities
46  {
47  namespace
48  {
49  bool isVertex(int x, int y, int z, int npts)
50  {
51  return (x == 0 && y == 0 && z == 0 ) ||
52  (x == (npts-1) && y == 0 && z == 0 ) ||
53  (x == (npts-1) && y == (npts-1) && z == 0 ) ||
54  (x == 0 && y == (npts-1) && z == 0 ) ||
55  (x == 0 && y == 0 && z == (npts-1)) ||
56  (x == 0 && y == (npts-1) && z == (npts-1));
57  }
58 
59  bool isEdge_01(int x, int y, int z, int npts)
60  {
61  return y == 0 && z == 0;
62  }
63 
64  bool isEdge_12(int x, int y, int z, int npts)
65  {
66  return x == (npts-1) && z == 0;
67  }
68 
69  bool isEdge_23(int x, int y, int z, int npts)
70  {
71  return y == (npts-1) && z == 0;
72  }
73 
74  bool isEdge_30(int x, int y, int z, int npts)
75  {
76  return x == 0 && z == 0;
77  }
78 
79  bool isEdge_04(int x, int y, int z, int npts)
80  {
81  return x == 0 && y == 0;
82  }
83 
84  bool isEdge_14(int x, int y, int z, int npts)
85  {
86  return x + z == (npts-1) && y == 0;
87  }
88 
89  bool isEdge_25(int x, int y, int z, int npts)
90  {
91  return x + z == (npts-1) && y == (npts-1);
92  }
93 
94  bool isEdge_35(int x, int y, int z, int npts)
95  {
96  return x == 0 && y == (npts-1);
97  }
98 
99  bool isEdge_45(int x, int y, int z, int npts)
100  {
101  return x == 0 && z == (npts-1);
102  }
103 
104  bool isEdge(int x, int y, int z, int npts){
105  return isEdge_01(x,y,z,npts) || isEdge_12(x,y,z,npts) ||
106  isEdge_23(x,y,z,npts) || isEdge_30(x,y,z,npts) ||
107  isEdge_04(x,y,z,npts) || isEdge_14(x,y,z,npts) ||
108  isEdge_25(x,y,z,npts) || isEdge_35(x,y,z,npts) ||
109  isEdge_45(x,y,z,npts);
110  }
111 
112  bool isFace_0123(int x, int y, int z, int npts)
113  {
114  return z == 0;
115  }
116 
117  bool isFace_014(int x, int y, int z, int npts)
118  {
119  return y == 0;
120  }
121 
122  bool isFace_1254(int x, int y, int z, int npts)
123  {
124  return x + z == npts-1;
125  }
126 
127  bool isFace_325(int x, int y, int z, int npts)
128  {
129  return y == (npts-1);
130  }
131 
132  bool isFace_0354(int x, int y, int z, int npts)
133  {
134  return x == 0;
135  }
136 
137  bool isFace(int x, int y, int z, int npts){
138  return isFace_0123(x,y,z,npts) || isFace_014(x,y,z,npts) ||
139  isFace_1254(x,y,z,npts) || isFace_325(x,y,z,npts) ||
140  isFace_0354(x,y,z,npts);
141  }
142  }
143 
144  // Calculate evenly spaced number of points
146  {
147  // Allocate the storage for points
149 
150  // Populate m_points
151  unsigned int npts = GetNumPoints();
152  NekDouble delta = 2.0/(npts - 1.0);
153  for(unsigned int z=0, index=0; z<npts; ++z){
154  for(int y=0; y<npts; ++y){
155  for(int x=0; x<npts-z; ++x, ++index){
156  NekDouble xi = -1.0 + x*delta;
157  NekDouble yi = -1.0 + y*delta;
158  NekDouble zi = -1.0 + z*delta;
159 
160  m_points[0][index] = xi;
161  m_points[1][index] = yi;
162  m_points[2][index] = zi;
163  }
164  }
165  }
166 
168  }
169 
171  {
172  unsigned int npts = GetNumPoints();
173  using std::vector;
174  vector<int> vertex;
175  vector<int> iEdge_01; // interior edge 0
176  vector<int> iEdge_12; // interior edge 1
177  vector<int> iEdge_23; // interior edge 2
178  vector<int> iEdge_30; // interior edge 3
179  vector<int> iEdge_04; // interior edge 4
180  vector<int> iEdge_14; // interior edge 5
181  vector<int> iEdge_25; // interior edge 6
182  vector<int> iEdge_35; // interior edge 7
183  vector<int> iEdge_45; // interior edge 8
184  vector<int> iFace_0123; // interior face 0
185  vector<int> iFace_014; // interior face 1
186  vector<int> iFace_1254; // interior face 2
187  vector<int> iFace_325; // interior face 3
188  vector<int> iFace_0354; // interior face 4
189  vector<int> interiorVolumePoints; // interior volume points
190  vector<int> map;
191 
192  // Build the lattice prism left to right - bottom to top
193  for(int z=0, index=0; z<npts; ++z){
194  for(int y=0; y<npts; ++y){
195  for(int x=0; x<npts-z; ++x, ++index){
196  if (isVertex(x,y,z,npts))
197  {
198  vertex.push_back(index);
199  }
200  else if (isEdge(x,y,z,npts))
201  {
202  if (isEdge_01(x,y,z,npts))
203  {
204  iEdge_01.push_back(index);
205  }
206  else if (isEdge_12(x,y,z,npts))
207  {
208  iEdge_12.push_back(index);
209  }
210  else if (isEdge_23(x,y,z,npts))
211  {
212  iEdge_23.push_back(index);
213  }
214  else if (isEdge_30(x,y,z,npts))
215  {
216  iEdge_30.push_back(index);
217  }
218  else if (isEdge_04(x,y,z,npts))
219  {
220  iEdge_04.push_back(index);
221  }
222  else if (isEdge_14(x,y,z,npts))
223  {
224  iEdge_14.push_back(index);
225  }
226  else if (isEdge_25(x,y,z,npts))
227  {
228  iEdge_25.push_back(index);
229  }
230  else if (isEdge_35(x,y,z,npts))
231  {
232  iEdge_35.push_back(index);
233  }
234  else if (isEdge_45(x,y,z,npts))
235  {
236  iEdge_45.push_back(index);
237  }
238  }
239  else if (isFace(x,y,z,npts))
240  {
241  if (isFace_0123(x,y,z,npts))
242  {
243  iFace_0123.push_back(index);
244  }
245  else if (isFace_014(x,y,z,npts))
246  {
247  iFace_014.push_back(index);
248  }
249  else if (isFace_1254(x,y,z,npts))
250  {
251  iFace_1254.push_back(index);
252  }
253  else if (isFace_325(x,y,z,npts))
254  {
255  iFace_325.push_back(index);
256  }
257  else if (isFace_0354(x,y,z,npts))
258  {
259  iFace_0354.push_back(index);
260  }
261  }
262  else
263  {
264  interiorVolumePoints.push_back(index);
265  }
266  }
267  }
268  }
269 
270  for (unsigned int n=0; n<vertex.size(); ++n)
271  {
272  map.push_back(vertex[n]);
273  }
274 
275  for (unsigned int n=0; n<iEdge_01.size(); ++n)
276  {
277  map.push_back(iEdge_01[n]);
278  }
279 
280  for (unsigned int n=0; n<iEdge_12.size(); ++n)
281  {
282  map.push_back(iEdge_12[n]);
283  }
284 
285  for (unsigned int n=0; n<iEdge_23.size(); ++n)
286  {
287  map.push_back(iEdge_23[n]);
288  }
289 
290  for (unsigned int n=0; n<iEdge_30.size(); ++n)
291  {
292  map.push_back(iEdge_30[n]);
293  }
294 
295  for (unsigned int n=0; n<iEdge_04.size(); ++n)
296  {
297  map.push_back(iEdge_04[n]);
298  }
299 
300  for (unsigned int n=0; n<iEdge_14.size(); ++n)
301  {
302  map.push_back(iEdge_14[n]);
303  }
304 
305  for (unsigned int n=0; n<iEdge_25.size(); ++n)
306  {
307  map.push_back(iEdge_25[n]);
308  }
309 
310  for (unsigned int n=0; n<iEdge_35.size(); ++n)
311  {
312  map.push_back(iEdge_35[n]);
313  }
314 
315  for (unsigned int n=0; n<iEdge_45.size(); ++n)
316  {
317  map.push_back(iEdge_45[n]);
318  }
319 
320  for (unsigned int n=0; n<iFace_0123.size(); ++n)
321  {
322  map.push_back(iFace_0123[n]);
323  }
324 
325  for (unsigned int n=0; n<iFace_014.size(); ++n)
326  {
327  map.push_back(iFace_014[n]);
328  }
329 
330  for(unsigned int n=0; n<iFace_1254.size(); ++n)
331  {
332  map.push_back(iFace_1254[n]);
333  }
334 
335  for(unsigned int n=0; n<iFace_325.size(); ++n)
336  {
337  map.push_back(iFace_325[n]);
338  }
339 
340  for(unsigned int n=0; n<iFace_0354.size(); ++n)
341  {
342  map.push_back(iFace_0354[n]);
343  }
344 
345  for(unsigned int n=0; n<interiorVolumePoints.size(); ++n)
346  {
347  map.push_back(interiorVolumePoints[n]);
348  }
349 
350  Array<OneD, NekDouble> points[3];
354 
355  for(unsigned int index=0; index<map.size(); ++index)
356  {
357  points[0][index] = m_points[0][index];
358  points[1][index] = m_points[1][index];
359  points[2][index] = m_points[2][index];
360  }
361 
362  for(unsigned int index=0; index<map.size(); ++index)
363  {
364  m_points[0][index] = points[0][map[index]];
365  m_points[1][index] = points[1][map[index]];
366  m_points[2][index] = points[2][map[index]];
367  }
368  }
369 
371  {
372 
373  }
374 
375 
376  // ////////////////////////////////////////
377  // CalculateInterpMatrix()
379  const Array<OneD, const NekDouble>& yia,
380  const Array<OneD, const NekDouble>& zia,
381  Array<OneD, NekDouble>& interp)
382  {
383  ASSERTL0(false, "Not yet implemented");
384  }
385 
386  // ////////////////////////////////////////
387  // CalculateDerivMatrix()
389  {
390 
391  }
392 
393  boost::shared_ptr<PointsBaseType> NodalPrismEvenlySpaced::Create(const PointsKey &key)
394  {
395  boost::shared_ptr<PointsBaseType> returnval(MemoryManager<NodalPrismEvenlySpaced>::AllocateSharedPtr(key));
396 
397  returnval->Initialize();
398 
399  return returnval;
400  }
401  } // end of namespace
402 } // end of namespace
403 
void CalculateInterpMatrix(const Array< OneD, const NekDouble > &xi, const Array< OneD, const NekDouble > &yi, const Array< OneD, const NekDouble > &zi, Array< OneD, NekDouble > &interp)
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
General purpose memory allocation routines with the ability to allocate from thread specific memory p...
Array< OneD, DataType > m_points[3]
Definition: Points.h:349
unsigned int GetNumPoints() const
Definition: Points.h:246
static std::string npts
Definition: InputFld.cpp:43
Defines a specification for a set of points.
Definition: Points.h:58
double NekDouble
static boost::shared_ptr< PointsBaseType > Create(const PointsKey &key)
unsigned int GetTotNumPoints() const
Definition: Points.h:251