176{
178 using std::vector;
179 vector<int> vertex;
180 vector<int> iEdge_01;
181 vector<int> iEdge_12;
182 vector<int> iEdge_23;
183 vector<int> iEdge_30;
184 vector<int> iEdge_04;
185 vector<int> iEdge_14;
186 vector<int> iEdge_25;
187 vector<int> iEdge_35;
188 vector<int> iEdge_45;
189 vector<int> iFace_0123;
190 vector<int> iFace_014;
191 vector<int> iFace_1254;
192 vector<int> iFace_325;
193 vector<int> iFace_0354;
194 vector<int> interiorVolumePoints;
195 vector<int> map;
196
197
198 for (size_t y = 0, index = 0; y < npts; y++)
199 {
200 for (size_t t = 0; t < npts * (npts + 1) / 2; t++, index++)
201 {
202 if (isVertex(t, y, npts))
203 {
204 vertex.push_back(index);
205 }
206 else if (isEdge(t, y, npts))
207 {
208 if (isEdge_01(t, y, npts))
209 {
210 iEdge_01.push_back(index);
211 }
212 else if (isEdge_12(t, y, npts))
213 {
214 iEdge_12.push_back(index);
215 }
216 else if (isEdge_23(t, y, npts))
217 {
218 iEdge_23.push_back(index);
219 }
220 else if (isEdge_30(t, y, npts))
221 {
222 iEdge_30.push_back(index);
223 }
224 else if (isEdge_04(t, y, npts))
225 {
226 iEdge_04.push_back(index);
227 }
228 else if (isEdge_14(t, y, npts))
229 {
230 iEdge_14.push_back(index);
231 }
232 else if (isEdge_25(t, y, npts))
233 {
234 iEdge_25.push_back(index);
235 }
236 else if (isEdge_35(t, y, npts))
237 {
238 iEdge_35.push_back(index);
239 }
240 else if (isEdge_45(t, y, npts))
241 {
242 iEdge_45.push_back(index);
243 }
244 }
245 else if (isFace(t, y, npts))
246 {
247 if (isFace_0123(t, y, npts))
248 {
249 iFace_0123.push_back(index);
250 }
251 else if (isFace_014(t, y, npts))
252 {
253 iFace_014.push_back(index);
254 }
255 else if (isFace_1254(t, y, npts))
256 {
257 iFace_1254.push_back(index);
258 }
259 else if (isFace_325(t, y, npts))
260 {
261 iFace_325.push_back(index);
262 }
263 else if (isFace_0354(t, y, npts))
264 {
265 iFace_0354.push_back(index);
266 }
267 }
268 }
269 }
270
271
272 std::swap(vertex[2], vertex[4]);
273
274 std::reverse(iEdge_23.begin(), iEdge_23.end());
275 std::reverse(iEdge_30.begin(), iEdge_30.end());
276 std::reverse(iEdge_04.begin(), iEdge_04.end());
277 std::reverse(iEdge_35.begin(), iEdge_35.end());
278
279
280 for (size_t i = 0; i < npts - 2; i++)
281 {
282 for (size_t j = i + 1; j < npts - 2; j++)
283 {
284 std::swap(iFace_1254[i * (npts - 2) + j],
285 iFace_1254[j * (npts - 2) + i]);
286 }
287 }
288 for (size_t i = 0; i < npts - 2; i++)
289 {
290 std::reverse(iFace_0354.begin() + (i * (npts - 2)),
291 iFace_0354.begin() + (i * (npts - 2) + npts - 2));
292 }
293 for (size_t i = 0; i < npts - 2; i++)
294 {
295 for (size_t j = i + 1; j < npts - 2; j++)
296 {
297 std::swap(iFace_0354[i * (npts - 2) + j],
298 iFace_0354[j * (npts - 2) + i]);
299 }
300 }
301
302
303
304
305 size_t index = npts * (npts + 1) / 2 + 3 * (npts - 1);
306
307 for (size_t k = 0; k < npts - 2; k++)
308 {
309 size_t index1 = index;
310 for (size_t j = 0; j < npts - 2; j++)
311 {
312 for (size_t i = 0; i < npts - 3 - k; i++)
313 {
314 interiorVolumePoints.push_back(index1++);
315 }
316 index1 += npts * (npts + 1) / 2 - (npts - 3 - k);
317 }
318 index += npts - 3 - k;
319 }
320
321 for (size_t n = 0; n < vertex.size(); ++n)
322 {
323 map.push_back(vertex[n]);
324 }
325
326 for (size_t n = 0; n < iEdge_01.size(); ++n)
327 {
328 map.push_back(iEdge_01[n]);
329 }
330
331 for (size_t n = 0; n < iEdge_12.size(); ++n)
332 {
333 map.push_back(iEdge_12[n]);
334 }
335
336 for (size_t n = 0; n < iEdge_23.size(); ++n)
337 {
338 map.push_back(iEdge_23[n]);
339 }
340
341 for (size_t n = 0; n < iEdge_30.size(); ++n)
342 {
343 map.push_back(iEdge_30[n]);
344 }
345
346 for (size_t n = 0; n < iEdge_04.size(); ++n)
347 {
348 map.push_back(iEdge_04[n]);
349 }
350
351 for (size_t n = 0; n < iEdge_14.size(); ++n)
352 {
353 map.push_back(iEdge_14[n]);
354 }
355
356 for (size_t n = 0; n < iEdge_25.size(); ++n)
357 {
358 map.push_back(iEdge_25[n]);
359 }
360
361 for (size_t n = 0; n < iEdge_35.size(); ++n)
362 {
363 map.push_back(iEdge_35[n]);
364 }
365
366 for (size_t n = 0; n < iEdge_45.size(); ++n)
367 {
368 map.push_back(iEdge_45[n]);
369 }
370
371 for (size_t n = 0; n < iFace_0123.size(); ++n)
372 {
373 map.push_back(iFace_0123[n]);
374 }
375
376 for (size_t n = 0; n < iFace_014.size(); ++n)
377 {
378 map.push_back(iFace_014[n]);
379 }
380
381 for (size_t n = 0; n < iFace_1254.size(); ++n)
382 {
383 map.push_back(iFace_1254[n]);
384 }
385
386 for (size_t n = 0; n < iFace_325.size(); ++n)
387 {
388 map.push_back(iFace_325[n]);
389 }
390
391 for (size_t n = 0; n < iFace_0354.size(); ++n)
392 {
393 map.push_back(iFace_0354[n]);
394 }
395
396 for (size_t n = 0; n < interiorVolumePoints.size(); ++n)
397 {
398 map.push_back(interiorVolumePoints[n]);
399 }
400
401 Array<OneD, NekDouble> points[3];
405
406 for (size_t index = 0; index < map.size(); ++index)
407 {
408 points[0][index] =
m_points[0][index];
409 points[1][index] =
m_points[1][index];
410 points[2][index] =
m_points[2][index];
411 }
412
413 for (size_t index = 0; index < map.size(); ++index)
414 {
415 m_points[0][index] = points[0][map[index]];
416 m_points[1][index] = points[1][map[index]];
417 m_points[2][index] = points[2][map[index]];
418 }
419}
Array< OneD, DataType > m_points[3]
Storage for the point locations, allowing for up to a 3D points storage.
size_t GetNumPoints() const
size_t GetTotNumPoints() const