Qspectral_elements.cc
Go to the documentation of this file.
1// LIC// ====================================================================
2// LIC// This file forms part of oomph-lib, the object-oriented,
3// LIC// multi-physics finite-element library, available
4// LIC// at http://www.oomph-lib.org.
5// LIC//
6// LIC// Copyright (C) 2006-2025 Matthias Heil and Andrew Hazel
7// LIC//
8// LIC// This library is free software; you can redistribute it and/or
9// LIC// modify it under the terms of the GNU Lesser General Public
10// LIC// License as published by the Free Software Foundation; either
11// LIC// version 2.1 of the License, or (at your option) any later version.
12// LIC//
13// LIC// This library is distributed in the hope that it will be useful,
14// LIC// but WITHOUT ANY WARRANTY; without even the implied warranty of
15// LIC// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
16// LIC// Lesser General Public License for more details.
17// LIC//
18// LIC// You should have received a copy of the GNU Lesser General Public
19// LIC// License along with this library; if not, write to the Free Software
20// LIC// Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
21// LIC// 02110-1301 USA.
22// LIC//
23// LIC// The authors may be contacted at oomph-lib@maths.man.ac.uk.
24// LIC//
25// LIC//====================================================================
26#include "Qspectral_elements.h"
27#include "Qelement_face_coordinate_translation_schemes.h"
28
29
30namespace oomph
31{
32 std::map<unsigned, Vector<double>> OneDLegendreShapeParam::z;
33
34
35 //////////////////////////////////////////////////////////////////////////
36 /// 1D QLegendreElements
37 //////////////////////////////////////////////////////////////////////////
38
39
40 //=======================================================================
41 /// Assign the static integral
42 //=======================================================================
43 template<unsigned NNODE_1D>
44 GaussLobattoLegendre<1, NNODE_1D> QSpectralElement<1, NNODE_1D>::integral;
45
46 //=======================================================================
47 /// The output function for general 1D QSpectralElements
48 //=======================================================================
49 template<unsigned NNODE_1D>
50 void QSpectralElement<1, NNODE_1D>::output(std::ostream& outfile)
51 {
52 // Tecplot header info
53 outfile << "ZONE I=" << NNODE_1D << std::endl;
54
55 // Find the dimension of the nodes
56 unsigned n_dim = this->nodal_dimension();
57
58 // Loop over element nodes
59 for (unsigned l = 0; l < NNODE_1D; l++)
60 {
61 // Loop over the dimensions and output the position
62 for (unsigned i = 0; i < n_dim; i++)
63 {
64 outfile << this->node_pt(l)->x(i) << " ";
65 }
66 // Find out how many data values at the node
67 unsigned initial_nvalue = this->node_pt(l)->nvalue();
68 // Lopp over the data and output whether pinned or not
69 for (unsigned i = 0; i < initial_nvalue; i++)
70 {
71 outfile << this->node_pt(l)->is_pinned(i) << " ";
72 }
73 outfile << std::endl;
74 }
75 outfile << std::endl;
76 }
77
78 //=======================================================================
79 /// The output function for nplot points in each coordinate direction
80 //=======================================================================
81 template<unsigned NNODE_1D>
82 void QSpectralElement<1, NNODE_1D>::output(std::ostream& outfile,
83 const unsigned& nplot)
84 {
85 // Local variables
86 Vector<double> s(1);
87 // Shape functions
88 Shape psi(NNODE_1D);
89
90 // Find the dimension of the nodes
91 unsigned n_dim = this->nodal_dimension();
92
93 // Tecplot header info
94 outfile << "ZONE I=" << nplot << std::endl;
95 // Loop over element nodes
96 for (unsigned l = 0; l < nplot; l++)
97 {
98 s[0] = -1.0 + l * 2.0 / (nplot - 1);
99 shape(s, psi);
100 for (unsigned i = 0; i < n_dim; i++)
101 {
102 // Output the x and y positions
103 outfile << this->interpolated_x(s, i) << " ";
104 }
105 for (unsigned i = 0; i < NNODE_1D; i++)
106 {
107 outfile << psi(i) << " ";
108 }
109 outfile << std::endl;
110 }
111 outfile << std::endl;
112 }
113
114 //===========================================================
115 /// Function to setup geometrical information for lower-dimensional
116 /// FaceElements (which are of type QSpectralElement<0,1>).
117 //===========================================================
118 template<unsigned NNODE_1D>
119 void QSpectralElement<1, NNODE_1D>::build_face_element(
120 const int& face_index, FaceElement* face_element_pt)
121 {
122 /*throw OomphLibError("Untested",
123 OOMPH_CURRENT_FUNCTION,
124 OOMPH_EXCEPTION_LOCATION);*/
125
126 // Overload the nodal dimension by reading out the value from the node
127 face_element_pt->set_nodal_dimension(this->node_pt(0)->ndim());
128
129 // Set the pointer to the "bulk" element
130 face_element_pt->bulk_element_pt() = this;
131
132#ifdef OOMPH_HAS_MPI
133 // If the bulk element is halo then the face element must be too
134 // if (this->is_halo())
135 {
136 face_element_pt->set_halo(Non_halo_proc_ID);
137 }
138#endif
139
140 // Resize the storage for the original number of values at the (one and
141 // only) node of the face element.
142 face_element_pt->nbulk_value_resize(1);
143
144 // Resize the storage for the bulk node number corresponding to the (one
145 // and only) node of the face element
146 face_element_pt->bulk_node_number_resize(1);
147
148 // Set the face index in the face element
149 face_element_pt->face_index() = face_index;
150
151 // Now set up the node pointer
152 // The convention is that the "normal", should always point
153 // out of the element
154 switch (face_index)
155 {
156 // Bottom, normal sign is negative (coordinate points into element)
157 case (-1):
158 face_element_pt->node_pt(0) = this->node_pt(0);
159 face_element_pt->bulk_node_number(0) = 0;
160 face_element_pt->normal_sign() = -1;
161
162 // Set the pointer to the function that determines the bulk coordinates
163 // in the face element
164 face_element_pt->face_to_bulk_coordinate_fct_pt() =
165 &QElement1FaceToBulkCoordinates::face0;
166
167 // Set the pointer to the function that determines the mapping of
168 // derivatives
169 face_element_pt->bulk_coordinate_derivatives_fct_pt() =
170 &QElement1BulkCoordinateDerivatives::faces0;
171
172 // Set the number of values stored when the node is part of the "bulk"
173 // element. The required_nvalue() must be used, rather than nvalue(),
174 // because otherwise nodes on boundaries will be resized multiple
175 // times. If you want any other behaviour, you MUST set nbulk_value()
176 // manually after construction of your specific face element.
177 face_element_pt->nbulk_value(0) = this->required_nvalue(0);
178 break;
179
180 // Top, normal sign is positive (coordinate points out of element)
181 case (1):
182 face_element_pt->node_pt(0) = this->node_pt(NNODE_1D - 1);
183 face_element_pt->bulk_node_number(0) = NNODE_1D - 1;
184 face_element_pt->normal_sign() = +1;
185
186
187 // Set the pointer to the function that determines the bulk coordinates
188 // in the face element
189 face_element_pt->face_to_bulk_coordinate_fct_pt() =
190 &QElement1FaceToBulkCoordinates::face1;
191
192 // Set the pointer to the function that determines the mapping of
193 // derivatives
194 face_element_pt->bulk_coordinate_derivatives_fct_pt() =
195 &QElement1BulkCoordinateDerivatives::faces0;
196
197
198 // Set the number of values stored when the node is part of the "bulk"
199 // element.
200 face_element_pt->nbulk_value(0) = this->required_nvalue(NNODE_1D - 1);
201 break;
202
203 // Other cases
204 default:
205 std::ostringstream error_message;
206 error_message << "Face_index should only take "
207 << "the values +/-1, not " << face_index << std::endl;
208
209 throw OomphLibError(error_message.str(),
210 OOMPH_CURRENT_FUNCTION,
211 OOMPH_EXCEPTION_LOCATION);
212 }
213 }
214
215
216 //////////////////////////////////////////////////////////////////////////
217 /// 2D QLegendreElements
218 //////////////////////////////////////////////////////////////////////////
219
220 //=======================================================================
221 /// Assign the static integral
222 //=======================================================================
223 template<unsigned NNODE_1D>
224 GaussLobattoLegendre<2, NNODE_1D> QSpectralElement<2, NNODE_1D>::integral;
225
226 //=======================================================================
227 /// The output function for general 1D QSpectralElements
228 //=======================================================================
229 template<unsigned NNODE_1D>
230 void QSpectralElement<2, NNODE_1D>::output(std::ostream& outfile)
231 {
232 // Tecplot header info
233 outfile << "ZONE I=" << NNODE_1D << ", J=" << NNODE_1D << std::endl;
234
235 // Find the dimension of the nodes
236 unsigned n_dim = this->nodal_dimension();
237
238 // Loop over element nodes
239 for (unsigned l2 = 0; l2 < NNODE_1D; l2++)
240 {
241 for (unsigned l1 = 0; l1 < NNODE_1D; l1++)
242 {
243 unsigned l = l2 * NNODE_1D + l1;
244
245 // Loop over the dimensions and output the position
246 for (unsigned i = 0; i < n_dim; i++)
247 {
248 outfile << this->node_pt(l)->x(i) << " ";
249 }
250 // Find out how many data values at the node
251 unsigned initial_nvalue = this->node_pt(l)->nvalue();
252 // Loop over the data and output whether pinned or not
253 for (unsigned i = 0; i < initial_nvalue; i++)
254 {
255 outfile << this->node_pt(l)->is_pinned(i) << " ";
256 }
257 outfile << std::endl;
258 }
259 }
260 outfile << std::endl;
261 }
262
263
264 //=======================================================================
265 /// The output function for n+plot points in each coordinate direction
266 //=======================================================================
267 template<unsigned NNODE_1D>
268 void QSpectralElement<2, NNODE_1D>::output(std::ostream& outfile,
269 const unsigned& n_plot)
270 {
271 // Local variables
272 Vector<double> s(2);
273
274 // Tecplot header info
275 outfile << "ZONE I=" << n_plot << ", J=" << n_plot << std::endl;
276
277 // Find the dimension of the first node
278 unsigned n_dim = this->nodal_dimension();
279
280 // Loop over plot points
281 for (unsigned l2 = 0; l2 < n_plot; l2++)
282 {
283 s[1] = -1.0 + l2 * 2.0 / (n_plot - 1);
284 for (unsigned l1 = 0; l1 < n_plot; l1++)
285 {
286 s[0] = -1.0 + l1 * 2.0 / (n_plot - 1);
287
288 // Output the coordinates
289 for (unsigned i = 0; i < n_dim; i++)
290 {
291 outfile << this->interpolated_x(s, i) << " ";
292 }
293 outfile << std::endl;
294 }
295 }
296 outfile << std::endl;
297 }
298
299
300 //===========================================================
301 /// Function to setup geometrical information for lower-dimensional
302 /// FaceElements (which are of type QSpectralElement<0,1>).
303 //===========================================================
304 template<unsigned NNODE_1D>
305 void QSpectralElement<2, NNODE_1D>::build_face_element(
306 const int& face_index, FaceElement* face_element_pt)
307 {
308 // Set the nodal dimension from the "bulk"
309 face_element_pt->set_nodal_dimension(this->node_pt(0)->ndim());
310
311 // Set the pointer to the "bulk" element
312 face_element_pt->bulk_element_pt() = this;
313
314#ifdef OOMPH_HAS_MPI
315 // If the bulk element is halo then the face element must be too
316 // if (this->is_halo())
317 {
318 face_element_pt->set_halo(Non_halo_proc_ID);
319 }
320#endif
321
322 // Resize the storage for the original number of values at
323 // NNODE_1D nodes of the FaceElement
324 face_element_pt->nbulk_value_resize(NNODE_1D);
325
326 // Resize the storage for the bulk node numbers corresponding
327 // to the NNODE_1D nodes of the FaceElement
328 face_element_pt->bulk_node_number_resize(NNODE_1D);
329
330 // Set the face index in the face element
331 // The faces are
332 // +1 East
333 // -1 West
334 // +2 North
335 // -2 South
336
337 // Set the face index in the face element
338 face_element_pt->face_index() = face_index;
339
340 // Now set up the node pointers
341 // =================================
342 // The convention here is that interior_tangent X tangent X tangent
343 // is the OUTWARD normal
344 switch (face_index)
345 {
346 unsigned bulk_number;
347 // West face, normal sign is positive
348 case (-1):
349 // Set the pointer to the bulk coordinate translation scheme
350 face_element_pt->face_to_bulk_coordinate_fct_pt() =
351 &QElement2FaceToBulkCoordinates::face0;
352
353 // Set the pointer to the derivative mappings
354 face_element_pt->bulk_coordinate_derivatives_fct_pt() =
355 &QElement2BulkCoordinateDerivatives::faces0;
356
357 for (unsigned i = 0; i < NNODE_1D; i++)
358 {
359 bulk_number = i * NNODE_1D;
360 face_element_pt->node_pt(i) = this->node_pt(bulk_number);
361 face_element_pt->bulk_node_number(i) = bulk_number;
362 face_element_pt->normal_sign() = 1;
363 // Set the number of values originally stored at this node
364 face_element_pt->nbulk_value(i) = this->required_nvalue(bulk_number);
365 }
366 break;
367 // South face, normal sign is positive
368 case (-2):
369 // Set the pointer to the bulk coordinate translation scheme
370 face_element_pt->face_to_bulk_coordinate_fct_pt() =
371 &QElement2FaceToBulkCoordinates::face1;
372
373 // Set the pointer to the derivative mappings
374 face_element_pt->bulk_coordinate_derivatives_fct_pt() =
375 &QElement2BulkCoordinateDerivatives::faces1;
376
377 for (unsigned i = 0; i < NNODE_1D; i++)
378 {
379 bulk_number = i;
380 face_element_pt->node_pt(i) = this->node_pt(bulk_number);
381 face_element_pt->bulk_node_number(i) = bulk_number;
382 face_element_pt->normal_sign() = 1;
383 // Set the number of values originally stored at this node
384 face_element_pt->nbulk_value(i) = this->required_nvalue(bulk_number);
385 }
386 break;
387 // East face, normal sign is negative
388 case (1):
389 // Set the pointer to the bulk coordinate translation scheme
390 face_element_pt->face_to_bulk_coordinate_fct_pt() =
391 &QElement2FaceToBulkCoordinates::face2;
392
393 // Set the pointer to the derivative mappings
394 face_element_pt->bulk_coordinate_derivatives_fct_pt() =
395 &QElement2BulkCoordinateDerivatives::faces0;
396
397 for (unsigned i = 0; i < NNODE_1D; i++)
398 {
399 bulk_number = NNODE_1D * i + NNODE_1D - 1;
400 face_element_pt->node_pt(i) = this->node_pt(bulk_number);
401 face_element_pt->bulk_node_number(i) = bulk_number;
402 face_element_pt->normal_sign() = -1;
403 // Set the number of values originally stored at this node
404 face_element_pt->nbulk_value(i) = this->required_nvalue(bulk_number);
405 }
406 break;
407 // North face, normal sign is negative
408 case (2):
409 // Set the pointer to the bulk coordinate translation scheme
410 face_element_pt->face_to_bulk_coordinate_fct_pt() =
411 &QElement2FaceToBulkCoordinates::face3;
412
413 // Set the pointer to the derivative mappings
414 face_element_pt->bulk_coordinate_derivatives_fct_pt() =
415 &QElement2BulkCoordinateDerivatives::faces1;
416
417 for (unsigned i = 0; i < NNODE_1D; i++)
418 {
419 bulk_number = NNODE_1D * (NNODE_1D - 1) + i;
420 face_element_pt->node_pt(i) = this->node_pt(bulk_number);
421 face_element_pt->bulk_node_number(i) = bulk_number;
422 face_element_pt->normal_sign() = -1;
423 // Set the number of values originally stored at this node
424 face_element_pt->nbulk_value(i) = this->required_nvalue(bulk_number);
425 }
426 break;
427
428 // Now cover the other cases
429 default:
430 std::ostringstream error_message;
431 error_message
432 << "Face index should only take the values +/- 1 or +/- 2,"
433 << " not " << face_index << std::endl;
434 throw OomphLibError(error_message.str(),
435 OOMPH_CURRENT_FUNCTION,
436 OOMPH_EXCEPTION_LOCATION);
437 }
438 }
439
440
441 //////////////////////////////////////////////////////////////////////////
442 /// 3D QLegendreElements
443 //////////////////////////////////////////////////////////////////////////
444
445 //=======================================================================
446 /// Assign the static integral
447 //=======================================================================
448 template<unsigned NNODE_1D>
449 GaussLobattoLegendre<3, NNODE_1D> QSpectralElement<3, NNODE_1D>::integral;
450
451 //=======================================================================
452 /// The output function for general 1D QSpectralElements
453 //=======================================================================
454 template<unsigned NNODE_1D>
455 void QSpectralElement<3, NNODE_1D>::output(std::ostream& outfile)
456 {
457 // Tecplot header info
458 outfile << "ZONE I=" << NNODE_1D << ", J=" << NNODE_1D << ", K=" << NNODE_1D
459 << std::endl;
460
461 // Find the dimension of the nodes
462 unsigned n_dim = this->nodal_dimension();
463
464 // Loop over element nodes
465 for (unsigned l3 = 0; l3 < NNODE_1D; l3++)
466 {
467 for (unsigned l2 = 0; l2 < NNODE_1D; l2++)
468 {
469 for (unsigned l1 = 0; l1 < NNODE_1D; l1++)
470 {
471 unsigned l = l3 * NNODE_1D * NNODE_1D + l2 * NNODE_1D + l1;
472
473 // Loop over the dimensions and output the position
474 for (unsigned i = 0; i < n_dim; i++)
475 {
476 outfile << this->node_pt(l)->x(i) << " ";
477 }
478 // Find out how many data values at the node
479 unsigned initial_nvalue = this->node_pt(l)->nvalue();
480 // Loop over the data and output whether pinned or not
481 for (unsigned i = 0; i < initial_nvalue; i++)
482 {
483 outfile << this->node_pt(l)->is_pinned(i) << " ";
484 }
485 outfile << std::endl;
486 }
487 }
488 }
489 outfile << std::endl;
490 }
491
492 //=======================================================================
493 /// The output function for n_plot points in each coordinate direction
494 //=======================================================================
495 template<unsigned NNODE_1D>
496 void QSpectralElement<3, NNODE_1D>::output(std::ostream& outfile,
497 const unsigned& n_plot)
498 {
499 // Local variables
500 Vector<double> s(3);
501
502 // Tecplot header info
503 outfile << "ZONE I=" << n_plot << ", J=" << n_plot << ", K=" << n_plot
504 << std::endl;
505
506 // Find the dimension of the first node
507 unsigned n_dim = this->nodal_dimension();
508
509 // Loop over element nodes
510 for (unsigned l3 = 0; l3 < n_plot; l3++)
511 {
512 s[2] = -1.0 + l3 * 2.0 / (n_plot - 1);
513 for (unsigned l2 = 0; l2 < n_plot; l2++)
514 {
515 s[1] = -1.0 + l2 * 2.0 / (n_plot - 1);
516 for (unsigned l1 = 0; l1 < n_plot; l1++)
517 {
518 s[0] = -1.0 + l1 * 2.0 / (n_plot - 1);
519
520 // Output the coordinates
521 for (unsigned i = 0; i < n_dim; i++)
522 {
523 outfile << this->interpolated_x(s, i) << " ";
524 }
525 outfile << std::endl;
526 }
527 }
528 }
529 outfile << std::endl;
530 }
531
532
533 //=======================================================================
534 /// Function to setup geometrical information for lower-dimensional
535 /// FaceElements (which are of type QElement<3,NNODE_1D>).
536 //=======================================================================
537 template<unsigned NNODE_1D>
538 void QSpectralElement<3, NNODE_1D>::build_face_element(
539 const int& face_index, FaceElement* face_element_pt)
540 {
541 oomph_info << " WARNING UNTESTED CODE" << std::endl;
542
543 // Set the nodal dimension from the "bulk"
544 face_element_pt->set_nodal_dimension(this->node_pt(0)->ndim());
545
546 // Set the pointer to the orginal "bulk" element
547 face_element_pt->bulk_element_pt() = this;
548
549#ifdef OOMPH_HAS_MPI
550 // If the bulk element is halo then the face element must be too
551 // if (this->is_halo())
552 {
553 face_element_pt->set_halo(Non_halo_proc_ID);
554 }
555#endif
556
557 // Resize storage for the number of values originally stored
558 // at the face element's NNODE_1D*NNODE_1D nodes.
559 face_element_pt->nbulk_value_resize(NNODE_1D * NNODE_1D);
560
561 // Set the face index in the element
562 // The faces are
563 // -3 : BACK (OLD: Bottom
564 // -2 : DOWN (OLD: Front
565 // -1 : LEFT (OLD: Left Side
566 // 1 : RIGHT (OLD: Right Side
567 // 2 : UP (OLD: Back
568 // 3 : FRONT (OLD: Top
569
570 face_element_pt->face_index() = face_index;
571
572 // Now set up the node pointers and the normal vectors
573 switch (face_index)
574 {
575 // BACK
576 //-----
577 case -3:
578 // Set the pointer to the bulk coordinate translation scheme
579 face_element_pt->face_to_bulk_coordinate_fct_pt() =
580 &QElement3FaceToBulkCoordinates::face2;
581
582 // Set the pointer to the derivative mappings
583 face_element_pt->bulk_coordinate_derivatives_fct_pt() =
584 &QElement3BulkCoordinateDerivatives::faces2;
585
586 // Copy nodes
587 for (unsigned i = 0; i < (NNODE_1D * NNODE_1D); i++)
588 {
589 face_element_pt->node_pt(i) = this->node_pt(i);
590 }
591 // Outer unit normal is negative of cross product of two in plane
592 // tangent vectors
593 face_element_pt->normal_sign() = -1;
594
595 break;
596
597 // FRONT
598 //------
599 case 3:
600
601 // Set the pointer to the bulk coordinate translation scheme
602 face_element_pt->face_to_bulk_coordinate_fct_pt() =
603 &QElement3FaceToBulkCoordinates::face5;
604
605 // Set the pointer to the derivative mappings
606 face_element_pt->bulk_coordinate_derivatives_fct_pt() =
607 &QElement3BulkCoordinateDerivatives::faces2;
608
609 // Copy nodes
610 for (unsigned i = 0; i < (NNODE_1D * NNODE_1D); i++)
611 {
612 face_element_pt->node_pt(i) =
613 this->node_pt(i + (NNODE_1D * NNODE_1D) * (NNODE_1D - 1));
614 }
615 // Outer unit normal is cross product of two in plane
616 // tangent vectors
617 face_element_pt->normal_sign() = 1;
618
619
620 break;
621
622 // DOWN:
623 //------
624 case -2:
625
626 {
627 // Set the pointer to the bulk coordinate translation scheme
628 face_element_pt->face_to_bulk_coordinate_fct_pt() =
629 &QElement3FaceToBulkCoordinates::face1;
630
631 // Set the pointer to the derivative mappings
632 face_element_pt->bulk_coordinate_derivatives_fct_pt() =
633 &QElement3BulkCoordinateDerivatives::faces1;
634
635 // Copy nodes
636 unsigned count = 0;
637 for (unsigned i = 0; i < NNODE_1D; i++)
638 {
639 for (unsigned j = 0; j < NNODE_1D; j++)
640 {
641 face_element_pt->node_pt(count) =
642 this->node_pt(j + i * (NNODE_1D * NNODE_1D));
643 count++;
644 }
645 }
646
647 // Outer unit normal is cross product of two in plane
648 // tangent vectors
649 face_element_pt->normal_sign() = 1;
650 }
651 break;
652
653
654 // UP:
655 //----
656 case 2:
657
658 {
659 // Set the pointer to the bulk coordinate translation scheme
660 face_element_pt->face_to_bulk_coordinate_fct_pt() =
661 &QElement3FaceToBulkCoordinates::face4;
662
663 // Set the pointer to the derivative mappings
664 face_element_pt->bulk_coordinate_derivatives_fct_pt() =
665 &QElement3BulkCoordinateDerivatives::faces1;
666
667 // Copy nodes
668 unsigned count = 0;
669 for (unsigned i = 0; i < NNODE_1D; i++)
670 {
671 for (unsigned j = 0; j < NNODE_1D; j++)
672 {
673 face_element_pt->node_pt(count) = this->node_pt(
674 j + i * (NNODE_1D * NNODE_1D) + (NNODE_1D * (NNODE_1D - 1)));
675 count++;
676 }
677 }
678
679 // Outer unit normal is negative of cross product of two in plane
680 // tangent vectors
681 face_element_pt->normal_sign() = -1;
682 }
683 break;
684
685 // LEFT:
686 //------
687 case -1:
688
689 {
690 // Set the pointer to the bulk coordinate translation scheme
691 face_element_pt->face_to_bulk_coordinate_fct_pt() =
692 &QElement3FaceToBulkCoordinates::face0;
693
694 // Set the pointer to the derivative mappings
695 face_element_pt->bulk_coordinate_derivatives_fct_pt() =
696 &QElement3BulkCoordinateDerivatives::faces0;
697
698 // Copy nodes
699 unsigned count = 0;
700 for (unsigned i = 0; i < NNODE_1D; i++)
701 {
702 for (unsigned j = 0; j < NNODE_1D; j++)
703 {
704 unsigned jj = j * NNODE_1D + i * (NNODE_1D * NNODE_1D);
705 face_element_pt->node_pt(count) = this->node_pt(jj);
706 count++;
707 }
708 }
709
710 // Outer unit normal is negative of cross product of two in plane
711 // tangent vectors
712 face_element_pt->normal_sign() = -1;
713 }
714 break;
715
716
717 // RIGHT:
718 //-------
719 case 1:
720
721 {
722 // Set the pointer to the bulk coordinate translation scheme
723 face_element_pt->face_to_bulk_coordinate_fct_pt() =
724 &QElement3FaceToBulkCoordinates::face3;
725
726 // Set the pointer to the derivative mappings
727 face_element_pt->bulk_coordinate_derivatives_fct_pt() =
728 &QElement3BulkCoordinateDerivatives::faces0;
729
730 // Copy nodes
731 unsigned count = 0;
732 for (unsigned i = 0; i < NNODE_1D; i++)
733 {
734 for (unsigned j = 0; j < NNODE_1D; j++)
735 {
736 unsigned jj =
737 j * NNODE_1D + i * (NNODE_1D * NNODE_1D) + (NNODE_1D - 1);
738 face_element_pt->node_pt(count) = this->node_pt(jj);
739 count++;
740 }
741 }
742
743 // Outer unit normal is cross product of two in plane
744 // tangent vectors
745 face_element_pt->normal_sign() = 1;
746 }
747 break;
748
749
750 // Cover all other cases
751 default:
752 std::ostringstream error_message;
753 error_message
754 << "Face index should only take the values +/- 1, +/- 2 or +/- 3,"
755 << " not " << face_index << std::endl;
756 throw OomphLibError(error_message.str(),
757 OOMPH_CURRENT_FUNCTION,
758 OOMPH_EXCEPTION_LOCATION);
759 } // end switch
760 }
761
762
763 template class QSpectralElement<1, 2>;
764 template class QSpectralElement<1, 3>;
765 template class QSpectralElement<1, 4>;
766 template class QSpectralElement<1, 5>;
767 template class QSpectralElement<1, 6>;
768 template class QSpectralElement<1, 7>;
769 template class QSpectralElement<1, 8>;
770 template class QSpectralElement<1, 9>;
771 template class QSpectralElement<1, 10>;
772 template class QSpectralElement<1, 11>;
773 template class QSpectralElement<1, 12>;
774 template class QSpectralElement<1, 13>;
775 template class QSpectralElement<1, 14>;
776
777 template class QSpectralElement<2, 2>;
778 template class QSpectralElement<2, 3>;
779 template class QSpectralElement<2, 4>;
780 template class QSpectralElement<2, 5>;
781 template class QSpectralElement<2, 6>;
782 template class QSpectralElement<2, 7>;
783 template class QSpectralElement<2, 8>;
784
785 template class QSpectralElement<3, 2>;
786 template class QSpectralElement<3, 3>;
787 template class QSpectralElement<3, 4>;
788 template class QSpectralElement<3, 5>;
789 template class QSpectralElement<3, 6>;
790 template class QSpectralElement<3, 7>;
791 template class QSpectralElement<3, 8>;
792
793} // namespace oomph
Qelement_face_coordinate_translation_schemes.h
Qspectral_elements.h
s
static char t char * s
Definition cfortran.h:568
i
cstr elem_len * i
Definition cfortran.h:603
oomph::Data::nvalue
unsigned nvalue() const
Return number of values stored in data object (incl pinned ones).
Definition nodes.h:483
oomph::Data::is_pinned
bool is_pinned(const unsigned &i) const
Test whether the i-th variable is pinned (1: true; 0: false).
Definition nodes.h:417
oomph::FaceElement
FaceElements are elements that coincide with the faces of higher-dimensional "bulk" elements....
Definition elements.h:4342
oomph::FaceElement::normal_sign
int & normal_sign()
Sign of outer unit normal (relative to cross-products of tangent vectors in the corresponding "bulk" ...
Definition elements.h:4616
oomph::FaceElement::face_index
int & face_index()
Index of the face (a number that uniquely identifies the face in the element)
Definition elements.h:4630
oomph::FaceElement::bulk_node_number_resize
void bulk_node_number_resize(const unsigned &i)
Resize the storage for the bulk node numbers.
Definition elements.h:4822
oomph::FaceElement::bulk_coordinate_derivatives_fct_pt
BulkCoordinateDerivativesFctPt & bulk_coordinate_derivatives_fct_pt()
Return the pointer to the function that returns the derivatives of the bulk coordinates wrt the face ...
Definition elements.h:4774
oomph::FaceElement::nbulk_value
unsigned & nbulk_value(const unsigned &n)
Return the number of values originally stored at local node n (before the FaceElement added additiona...
Definition elements.h:4849
oomph::FaceElement::bulk_node_number
unsigned & bulk_node_number(const unsigned &n)
Return the bulk node number that corresponds to the n-th local node number.
Definition elements.h:4829
oomph::FaceElement::bulk_element_pt
FiniteElement *& bulk_element_pt()
Pointer to higher-dimensional "bulk" element.
Definition elements.h:4739
oomph::FaceElement::face_to_bulk_coordinate_fct_pt
CoordinateMappingFctPt & face_to_bulk_coordinate_fct_pt()
Return the pointer to the function that maps the face coordinate to the bulk coordinate.
Definition elements.h:4759
oomph::FaceElement::nbulk_value_resize
void nbulk_value_resize(const unsigned &i)
Resize the storage for the number of values originally stored at the local nodes to i entries.
Definition elements.h:4864
oomph::FiniteElement::set_nodal_dimension
void set_nodal_dimension(const unsigned &nodal_dim)
Set the dimension of the nodes in the element. This will typically only be required when constructing...
Definition elements.h:1394
oomph::FiniteElement::shape
virtual void shape(const Vector< double > &s, Shape &psi) const =0
Calculate the geometric shape functions at local coordinate s. This function must be overloaded for e...
oomph::FiniteElement::interpolated_x
virtual double interpolated_x(const Vector< double > &s, const unsigned &i) const
Return FE interpolated coordinate x[i] at local coordinate s.
Definition elements.cc:3992
oomph::FiniteElement::node_pt
Node *& node_pt(const unsigned &n)
Return a pointer to the local node n.
Definition elements.h:2179
oomph::FiniteElement::nodal_dimension
unsigned nodal_dimension() const
Return the required Eulerian dimension of the nodes in this element.
Definition elements.h:2488
oomph::GeneralisedElement::Non_halo_proc_ID
int Non_halo_proc_ID
Non-halo processor ID for Data; -1 if it's not a halo.
Definition elements.h:128
oomph::GeneralisedElement::set_halo
void set_halo(const unsigned &non_halo_proc_ID)
Label the element as halo and specify processor that holds non-halo counterpart.
Definition elements.h:1138
oomph::GeomObject::ndim
unsigned ndim() const
Access function to # of Eulerian coordinates.
Definition geom_objects.h:177
oomph::Node::x
double & x(const unsigned &i)
Return the i-th nodal coordinate.
Definition nodes.h:1060
oomph::OneDLegendreShapeParam::z
static std::map< unsigned, Vector< double > > z
Definition Qspectral_elements.h:48
oomph::OomphLibError
An OomphLibError object which should be thrown when an run-time error is encountered....
Definition oomph_definitions.h:229
oomph::QSpectralElement
General QLegendreElement class.
Definition Qspectral_elements.h:363
oomph::Shape
A Class for shape functions. In simple cases, the shape functions have only one index that can be tho...
Definition shape.h:76
oomph::TAdvectionDiffusionReactionElement
TAdvectionDiffusionReactionElement<NREAGENT,DIM,NNODE_1D> elements are isoparametric triangular DIM-d...
Definition Tadvection_diffusion_reaction_elements.h:66
oomph::TAdvectionDiffusionReactionElement::required_nvalue
unsigned required_nvalue(const unsigned &n) const
Access function for Nvalue: # of ‘values’ (pinned or dofs) at node n (always returns the same value a...
Definition Tadvection_diffusion_reaction_elements.h:89
oomph::TAdvectionDiffusionReactionElement::TAdvectionDiffusionReactionElement
TAdvectionDiffusionReactionElement()
Constructor: Call constructors for TElement and AdvectionDiffusionReaction equations.
Definition Tadvection_diffusion_reaction_elements.h:70
oomph::QElement1BulkCoordinateDerivatives::faces0
void faces0(const Vector< double > &s, DenseMatrix< double > &dsbulk_dsface, unsigned &interior_direction)
Function for both faces – the bulk coordinate is fixed on both.
Definition Qelement_face_coordinate_translation_schemes.cc:59
oomph::QElement1FaceToBulkCoordinates::face1
void face1(const Vector< double > &s, Vector< double > &s_bulk)
The translation scheme for the face s0 = 1.0.
Definition Qelement_face_coordinate_translation_schemes.cc:44
oomph::QElement1FaceToBulkCoordinates::face0
void face0(const Vector< double > &s, Vector< double > &s_bulk)
The translation scheme for the face s0 = -1.0.
Definition Qelement_face_coordinate_translation_schemes.cc:38
oomph::QElement2BulkCoordinateDerivatives::faces1
void faces1(const Vector< double > &s, DenseMatrix< double > &dsbulk_dsface, unsigned &interior_direction)
Function for the north and south faces, along which s1 is fixed.
Definition Qelement_face_coordinate_translation_schemes.cc:131
oomph::QElement2BulkCoordinateDerivatives::faces0
void faces0(const Vector< double > &s, DenseMatrix< double > &dsbulk_dsface, unsigned &interior_direction)
Function for the east and west faces, along which s0 is fixed.
Definition Qelement_face_coordinate_translation_schemes.cc:116
oomph::QElement2FaceToBulkCoordinates::face2
void face2(const Vector< double > &s, Vector< double > &s_bulk)
The translation scheme for the east face (s0 = 1.0)
Definition Qelement_face_coordinate_translation_schemes.cc:93
oomph::QElement2FaceToBulkCoordinates::face1
void face1(const Vector< double > &s, Vector< double > &s_bulk)
The translation scheme for the south face (s1 = -1.0)
Definition Qelement_face_coordinate_translation_schemes.cc:86
oomph::QElement2FaceToBulkCoordinates::face0
void face0(const Vector< double > &s, Vector< double > &s_bulk)
The translation scheme for the west face (s0 = -1.0)
Definition Qelement_face_coordinate_translation_schemes.cc:79
oomph::QElement2FaceToBulkCoordinates::face3
void face3(const Vector< double > &s, Vector< double > &s_bulk)
The translation scheme for the north face (s1 = 1.0)
Definition Qelement_face_coordinate_translation_schemes.cc:100
oomph::QElement3BulkCoordinateDerivatives::faces2
void faces2(const Vector< double > &s, DenseMatrix< double > &dsbulk_dsface, unsigned &interior_direction)
Function for the left and right faces, along which s2 is fixed.
Definition Qelement_face_coordinate_translation_schemes.cc:249
oomph::QElement3BulkCoordinateDerivatives::faces0
void faces0(const Vector< double > &s, DenseMatrix< double > &dsbulk_dsface, unsigned &interior_direction)
Function for the back and front faces, along which s0 is fixed.
Definition Qelement_face_coordinate_translation_schemes.cc:210
oomph::QElement3BulkCoordinateDerivatives::faces1
void faces1(const Vector< double > &s, DenseMatrix< double > &dsbulk_dsface, unsigned &interior_direction)
Function for the up and down faces, along which s1 is fixed.
Definition Qelement_face_coordinate_translation_schemes.cc:230
oomph::QElement3FaceToBulkCoordinates::face4
void face4(const Vector< double > &s, Vector< double > &s_bulk)
The translation scheme for the up face (s1 = 1.0)
Definition Qelement_face_coordinate_translation_schemes.cc:185
oomph::QElement3FaceToBulkCoordinates::face1
void face1(const Vector< double > &s, Vector< double > &s_bulk)
The translation scheme for the down face (s1 = -1.0)
Definition Qelement_face_coordinate_translation_schemes.cc:161
oomph::QElement3FaceToBulkCoordinates::face2
void face2(const Vector< double > &s, Vector< double > &s_bulk)
The translation scheme for the back face (s2 = -1.0)
Definition Qelement_face_coordinate_translation_schemes.cc:169
oomph::QElement3FaceToBulkCoordinates::face0
void face0(const Vector< double > &s, Vector< double > &s_bulk)
The translation scheme for the left face (s0 = -1.0)
Definition Qelement_face_coordinate_translation_schemes.cc:153
oomph::QElement3FaceToBulkCoordinates::face5
void face5(const Vector< double > &s, Vector< double > &s_bulk)
The translation scheme for the front face (s2 = 1.0)
Definition Qelement_face_coordinate_translation_schemes.cc:193
oomph::QElement3FaceToBulkCoordinates::face3
void face3(const Vector< double > &s, Vector< double > &s_bulk)
The translation scheme for the right face (s0 = 1.0)
Definition Qelement_face_coordinate_translation_schemes.cc:177
oomph
DRAIG: Change all instances of (SPATIAL_DIM) to (DIM-1).
Definition advection_diffusion_elements.cc:30
oomph::oomph_info
OomphInfo oomph_info
Single (global) instantiation of the OomphInfo object – this is used throughout the library as a "rep...
Definition oomph_definitions.cc:326