channel_spine_mesh.template.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-2022 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#ifndef OOMPH_CHANNEL_SPINE_MESH_TEMPLATE_CC
27#define OOMPH_CHANNEL_SPINE_MESH_TEMPLATE_CC
28
31
32
33namespace oomph
34{
35 //===========================================================================
36 /// Constructor for spine 2D mesh: Pass number of elements in x-direction,
37 /// number of elements in y-direction, axial length and height of layer,
38 /// and pointer to timestepper (defaults to Static timestepper).
39 ///
40 /// The mesh contains a layer of spinified fluid elements (of type ELEMENT;
41 /// e.g SpineElement<QCrouzeixRaviartElement<2>)
42 /// and a surface layer of corresponding Spine interface elements
43 /// of type INTERFACE_ELEMENT, e.g.
44 /// SpineLineFluidInterfaceElement<ELEMENT> for 2D planar
45 /// problems.
46 //===========================================================================
47 template<class ELEMENT>
49 const unsigned& nx1,
50 const unsigned& nx2,
51 const unsigned& ny,
52 const double& lx0,
53 const double& lx1,
54 const double& lx2,
55 const double& h,
56 GeomObject* wall_pt,
57 TimeStepper* time_stepper_pt)
58 : RectangularQuadMesh<ELEMENT>(nx0 + nx1 + nx2,
59 ny,
60 0.0,
61 lx0 + lx1 + lx2,
62 0.0,
63 h,
64 false,
65 false,
66 time_stepper_pt),
67 Nx0(nx0),
68 Nx1(nx1),
69 Nx2(nx2),
70 Lx0(lx0),
71 Lx1(lx1),
72 Lx2(lx2),
73 Wall_pt(wall_pt)
74 {
75 // Mesh can only be built with 2D Qelements.
76 MeshChecker::assert_geometric_element<QElementGeometricBase, ELEMENT>(2);
77
78 // Mesh can only be built with spine elements
79 MeshChecker::assert_geometric_element<SpineFiniteElement, ELEMENT>(2);
80
81 // We've called the "generic" constructor for the RectangularQuadMesh
82 // which doesn't do much...
83
84 // Build the straight line object
86
87 // Now build the mesh:
88 build_channel_spine_mesh(time_stepper_pt);
89 }
90
91
92 //===========================================================================
93 /// Constuctor for spine 2D mesh: Pass number of elements in x-direction,
94 /// number of elements in y-direction, axial length and height of layer,
95 /// a boolean flag to make the mesh periodic in the x-direction,
96 /// and pointer to timestepper (defaults to Static timestepper).
97 ///
98 /// The mesh contains a layer of elements (of type ELEMENT)
99 /// and a surface layer of corresponding Spine interface elements (of type
100 /// SpineLineFluidInterfaceElement<ELEMENT>).
101 //===========================================================================
102 template<class ELEMENT>
104 const unsigned& nx1,
105 const unsigned& nx2,
106 const unsigned& ny,
107 const double& lx0,
108 const double& lx1,
109 const double& lx2,
110 const double& h,
111 GeomObject* wall_pt,
112 const bool& periodic_in_x,
113 TimeStepper* time_stepper_pt)
114 : RectangularQuadMesh<ELEMENT>(nx0 + nx1 + nx2,
115 ny,
116 0.0,
117 lx0 + lx1 + lx2,
118 0.0,
119 h,
120 periodic_in_x,
121 false,
122 time_stepper_pt),
123 Nx0(nx0),
124 Nx1(nx1),
125 Nx2(nx2),
126 Lx0(lx0),
127 Lx1(lx1),
128 Lx2(lx2),
129 Wall_pt(wall_pt)
130 {
131 // Mesh can only be built with 2D Qelements.
132 MeshChecker::assert_geometric_element<QElementGeometricBase, ELEMENT>(2);
133
134 // Mesh can only be built with spine elements
135 MeshChecker::assert_geometric_element<SpineFiniteElement, ELEMENT>(2);
136
137 // We've called the "generic" constructor for the RectangularQuadMesh
138 // which doesn't do much...
139
140 // Build the straight line object
142
143 // Now build the mesh:
144 build_channel_spine_mesh(time_stepper_pt);
145 }
146
147
148 //===========================================================================
149 /// Helper function that actually builds the channel-spine mesh
150 /// based on the parameters set in the various constructors
151 //===========================================================================
152 template<class ELEMENT>
154 TimeStepper* time_stepper_pt)
155 {
156 // Build the underlying quad mesh:
158
159 // Read out the number of elements in the x-direction and y-direction
160 // and in each of the left, centre and right regions
161 unsigned n_x = this->Nx;
162 unsigned n_y = this->Ny;
163 unsigned n_x0 = this->Nx0;
164 unsigned n_x1 = this->Nx1;
165 unsigned n_x2 = this->Nx2;
166
167 // Set up the pointers to elements in the left region
168 unsigned nleft = n_x0 * n_y;
169 ;
170 Left_element_pt.reserve(nleft);
171 unsigned ncentre = n_x1 * n_y;
172 ;
173 Centre_element_pt.reserve(ncentre);
174 unsigned nright = n_x2 * n_y;
175 ;
176 Right_element_pt.reserve(nright);
177 for (unsigned irow = 0; irow < n_y; irow++)
178 {
179 for (unsigned e = 0; e < n_x0; e++)
180 {
181 Left_element_pt.push_back(this->finite_element_pt(irow * n_x + e));
182 }
183 for (unsigned e = 0; e < n_x1; e++)
184 {
185 Centre_element_pt.push_back(
186 this->finite_element_pt(irow * n_x + (n_x0 + e)));
187 }
188 for (unsigned e = 0; e < n_x2; e++)
189 {
190 Right_element_pt.push_back(
191 this->finite_element_pt(irow * n_x + (n_x0 + n_x1 + e)));
192 }
193 }
194
195#ifdef PARANOID
196 // Check that we have the correct number of elements
197 if (nelement() != nleft + ncentre + nright)
198 {
199 throw OomphLibError("Incorrect number of element pointers!",
200 OOMPH_CURRENT_FUNCTION,
201 OOMPH_EXCEPTION_LOCATION);
202 }
203#endif
204
205 // Allocate memory for the spines and fractions along spines
206 //---------------------------------------------------------
207
208 // Read out number of linear points in the element
209 unsigned n_p = dynamic_cast<ELEMENT*>(finite_element_pt(0))->nnode_1d();
210
211 unsigned nspine;
212 // Allocate store for the spines:
213 if (this->Xperiodic)
214 {
215 nspine = (n_p - 1) * n_x;
216 Spine_pt.reserve(nspine);
217 // Number of spines in each region
218 // NOTE that boundary spines are in both regions
219 Nleft_spine = (n_p - 1) * n_x0 + 1;
220 Ncentre_spine = (n_p - 1) * n_x1 + 1;
221 Nright_spine = (n_p - 1) * n_x2;
222 }
223 else
224 {
225 nspine = (n_p - 1) * n_x + 1;
226 Spine_pt.reserve(nspine);
227 // Number of spines in each region
228 // NOTE that boundary spines are in both regions
229 Nleft_spine = (n_p - 1) * n_x0 + 1;
230 Ncentre_spine = (n_p - 1) * n_x1 + 1;
231 Nright_spine = (n_p - 1) * n_x2 + 1;
232 }
233
234 // end Allocating memory
235
236
237 // set up the vectors of geometric data & objects for building spines
238 Vector<double> r_wall(2), zeta(1), s_wall(1);
239 GeomObject* geometric_object_pt = 0;
240
241 // LEFT REGION
242 // ===========
243
244 // SPINES IN LEFT REGION
245 // ---------------------
246
247 // Set up zeta increments
248 double zeta_lo = 0.0;
249 double dzeta = Lx0 / n_x0;
250
251 // Initialise number of elements in previous regions:
252 unsigned n_prev_elements = 0;
253
254
255 // FIRST SPINE
256 //-----------
257
258 // Element 0
259 // Node 0
260 // Assign the new spine with unit length
261 Spine* new_spine_pt = new Spine(1.0);
262 new_spine_pt->spine_height_pt()->pin(0);
263 Spine_pt.push_back(new_spine_pt);
264
265 // Get pointer to node
266 SpineNode* nod_pt = element_node_pt(0, 0);
267 // Set the pointer to the spine
268 nod_pt->spine_pt() = new_spine_pt;
269 // Set the fraction
270 nod_pt->fraction() = 0.0;
271 // Pointer to the mesh that implements the update fct
272 nod_pt->spine_mesh_pt() = this;
273 // Set update fct id
274 nod_pt->node_update_fct_id() = 0;
275
276 // Provide spine with additional storage for wall coordinate
277 // and wall geom object:
278
279 {
280 // Get the Lagrangian coordinate in the Lower Wall
281 zeta[0] = 0.0;
282 // Get the geometric object and local coordinate
283 Straight_wall_pt->locate_zeta(zeta, geometric_object_pt, s_wall);
284
285 // The local coordinate is a geometric parameter
286 // This needs to be set (rather than added) because the
287 // same spine may be visited more than once
288 Vector<double> parameters(1, s_wall[0]);
289 nod_pt->spine_pt()->set_geom_parameter(parameters);
290
291 // Get position of wall
292 Straight_wall_pt->position(s_wall, r_wall);
293
294 // Adjust spine height
295 nod_pt->spine_pt()->height() = r_wall[1];
296
297 // The sub geom object is one (and only) geom object
298 // for spine:
299 Vector<GeomObject*> geom_object_pt(1);
300 geom_object_pt[0] = geometric_object_pt;
301
302 // Pass geom object(s) to spine
303 nod_pt->spine_pt()->set_geom_object_pt(geom_object_pt);
304 }
305
306 // Loop vertically along the spine
307 // Loop over the elements
308 for (unsigned long i = 0; i < n_y; i++)
309 {
310 // Loop over the vertical nodes, apart from the first
311 for (unsigned l1 = 1; l1 < n_p; l1++)
312 {
313 // Get pointer to node
314 SpineNode* nod_pt = element_node_pt(i * n_x, l1 * n_p);
315 // Set the pointer to the spine
316 nod_pt->spine_pt() = new_spine_pt;
317 // Set the fraction
318 nod_pt->fraction() =
319 (double(i) + double(l1) / double(n_p - 1)) / double(n_y);
320 // Pointer to the mesh that implements the update fct
321 nod_pt->spine_mesh_pt() = this;
322 // Set update fct id
323 nod_pt->node_update_fct_id() = 0;
324 }
325 } // end loop over elements
326
327
328 // LOOP OVER OTHER SPINES
329 //----------------------
330
331 // Now loop over the elements horizontally
332 for (unsigned long j = 0; j < n_x0; j++)
333 {
334 // Loop over the nodes in the elements horizontally, ignoring
335 // the first column
336 unsigned n_pmax = n_p;
337 for (unsigned l2 = 1; l2 < n_pmax; l2++)
338 {
339 // Assign the new spine with unit height
340 new_spine_pt = new Spine(1.0);
341 new_spine_pt->spine_height_pt()->pin(0);
342 Spine_pt.push_back(new_spine_pt);
343
344 // Get the node
345 SpineNode* nod_pt = element_node_pt(j, l2);
346 // Set the pointer to spine
347 nod_pt->spine_pt() = new_spine_pt;
348 // Set the fraction
349 nod_pt->fraction() = 0.0;
350 // Pointer to the mesh that implements the update fct
351 nod_pt->spine_mesh_pt() = this;
352 // Set update fct id
353 nod_pt->node_update_fct_id() = 0;
354
355 {
356 // Provide spine with additional storage for wall coordinate
357 // and wall geom object:
358
359 // Get the Lagrangian coordinate in the Lower Wall
360 zeta[0] =
361 zeta_lo + double(j) * dzeta + double(l2) * dzeta / (n_p - 1.0);
362 // Get the geometric object and local coordinate
363 Straight_wall_pt->locate_zeta(zeta, geometric_object_pt, s_wall);
364
365 // The local coordinate is a geometric parameter
366 // This needs to be set (rather than added) because the
367 // same spine may be visited more than once
368 Vector<double> parameters(1, s_wall[0]);
369 nod_pt->spine_pt()->set_geom_parameter(parameters);
370
371 // Get position of wall
372 Straight_wall_pt->position(s_wall, r_wall);
373
374 // Adjust spine height
375 nod_pt->spine_pt()->height() = r_wall[1];
376
377 // The sub geom object is one (and only) geom object
378 // for spine:
379 Vector<GeomObject*> geom_object_pt(1);
380 geom_object_pt[0] = geometric_object_pt;
381
382 // Pass geom object(s) to spine
383 nod_pt->spine_pt()->set_geom_object_pt(geom_object_pt);
384 }
385
386 // Loop vertically along the spine
387 // Loop over the elements
388 for (unsigned long i = 0; i < n_y; i++)
389 {
390 // Loop over the vertical nodes, apart from the first
391 for (unsigned l1 = 1; l1 < n_p; l1++)
392 {
393 // Get the node
394 SpineNode* nod_pt = element_node_pt(i * n_x + j, l1 * n_p + l2);
395 // Set the pointer to the spine
396 nod_pt->spine_pt() = new_spine_pt;
397 // Set the fraction
398 nod_pt->fraction() =
399 (double(i) + double(l1) / double(n_p - 1)) / double(n_y);
400 // Pointer to the mesh that implements the update fct
401 nod_pt->spine_mesh_pt() = this;
402 // Set update fct id
403 nod_pt->node_update_fct_id() = 0;
404 }
405 }
406 }
407 }
408
409 // Increment number of previous elements
410 n_prev_elements += n_x0 * n_y;
411
412
413 // CENTRE REGION
414 // ===========
415
416 zeta_lo = Lx0;
417 dzeta = Lx1 / n_x1;
418
419 // SPINES IN LEFT REGION
420 // ---------------------
421
422 // LOOP OVER OTHER SPINES
423 //----------------------
424
425 // Now loop over the elements horizontally
426 for (unsigned long j = n_x0; j < n_x0 + n_x1; j++)
427 {
428 // Loop over the nodes in the elements horizontally, ignoring
429 // the first column
430 unsigned n_pmax = n_p;
431 for (unsigned l2 = 1; l2 < n_pmax; l2++)
432 {
433 // Assign the new spine with unit height
434 new_spine_pt = new Spine(1.0);
435 new_spine_pt->spine_height_pt()->pin(0);
436 Spine_pt.push_back(new_spine_pt);
437
438 // Get the node
439 SpineNode* nod_pt = element_node_pt(j, l2);
440 // Set the pointer to spine
441 nod_pt->spine_pt() = new_spine_pt;
442 // Set the fraction
443 nod_pt->fraction() = 0.0;
444 // Pointer to the mesh that implements the update fct
445 nod_pt->spine_mesh_pt() = this;
446 // Set update fct id
447 nod_pt->node_update_fct_id() = 1;
448
449 {
450 // Provide spine with additional storage for wall coordinate
451 // and wall geom object:
452
453 // Get the Lagrangian coordinate in the Lower Wall
454 zeta[0] = zeta_lo + double(j - n_x0) * dzeta +
455 double(l2) * dzeta / (n_p - 1.0);
456 // Get the geometric object and local coordinate
457 Wall_pt->locate_zeta(zeta, geometric_object_pt, s_wall);
458
459 // The local coordinate is a geometric parameter
460 // This needs to be set (rather than added) because the
461 // same spine may be visited more than once
462 Vector<double> parameters(1, s_wall[0]);
463 nod_pt->spine_pt()->set_geom_parameter(parameters);
464
465 // Get position of wall
466 Wall_pt->position(s_wall, r_wall);
467
468 // Adjust spine height
469 nod_pt->spine_pt()->height() = r_wall[1];
470
471 // The sub geom object is one (and only) geom object
472 // for spine:
473 Vector<GeomObject*> geom_object_pt(1);
474 geom_object_pt[0] = geometric_object_pt;
475
476 // Pass geom object(s) to spine
477 nod_pt->spine_pt()->set_geom_object_pt(geom_object_pt);
478 }
479
480 // Loop vertically along the spine
481 // Loop over the elements
482 for (unsigned long i = 0; i < n_y; i++)
483 {
484 // Loop over the vertical nodes, apart from the first
485 for (unsigned l1 = 1; l1 < n_p; l1++)
486 {
487 // Get the node
488 SpineNode* nod_pt = element_node_pt(i * n_x + j, l1 * n_p + l2);
489 // Set the pointer to the spine
490 nod_pt->spine_pt() = new_spine_pt;
491 // Set the fraction
492 nod_pt->fraction() =
493 (double(i) + double(l1) / double(n_p - 1)) / double(n_y);
494 // Pointer to the mesh that implements the update fct
495 nod_pt->spine_mesh_pt() = this;
496 // Set update fct id
497 nod_pt->node_update_fct_id() = 1;
498 }
499 }
500 }
501 }
502
503 // Increment number of previous elements
504 n_prev_elements += n_x1 * n_y;
505
506
507 // RIGHT REGION
508 // ===========
509
510 // SPINES IN RIGHT REGION
511 // ----------------------
512
513 // Set up zeta increments
514 zeta_lo = Lx0 + Lx1;
515 dzeta = Lx2 / n_x2;
516
517 // LOOP OVER OTHER SPINES
518 //----------------------
519
520 // Now loop over the elements horizontally
521 for (unsigned long j = n_x0 + n_x1; j < n_x0 + n_x1 + n_x2; j++)
522 {
523 // Need to copy last spine in previous element over to first spine
524 // in next elements, for all elements other than the first
525 if (j > 0)
526 {
527 // For first spine, re-assign the geometric objects, since
528 // we treat it as part of the right region.
529 if (j == n_x0 + n_x1)
530 {
531 SpineNode* nod_pt = element_node_pt(j, 0);
532 // Set update fct id
533 nod_pt->node_update_fct_id() = 2;
534 {
535 // Provide spine with additional storage for wall coordinate
536 // and wall geom object:
537
538 // Get the Lagrangian coordinate in the Lower Wall
539 zeta[0] = zeta_lo + double(j - n_x0 - n_x1) * dzeta;
540 // Get the geometric object and local coordinate
541 Straight_wall_pt->locate_zeta(zeta, geometric_object_pt, s_wall);
542
543 // The local coordinate is a geometric parameter
544 // This needs to be set (rather than added) because the
545 // same spine may be visited more than once
546 Vector<double> parameters(1, s_wall[0]);
547 nod_pt->spine_pt()->set_geom_parameter(parameters);
548
549 // Get position of wall
550 Straight_wall_pt->position(s_wall, r_wall);
551
552 // Adjust spine height
553 nod_pt->spine_pt()->height() = r_wall[1];
554
555 // The sub geom object is one (and only) geom object
556 // for spine:
557 Vector<GeomObject*> geom_object_pt(1);
558 geom_object_pt[0] = geometric_object_pt;
559
560 // Pass geom object(s) to spine
561 nod_pt->spine_pt()->set_geom_object_pt(geom_object_pt);
562 }
563 }
564 }
565 // Loop over the nodes in the elements horizontally, ignoring
566 // the first column
567
568 // Last spine needs special treatment in x-periodic meshes:
569 unsigned n_pmax = n_p;
570 if ((this->Xperiodic) && (j == n_x - 1)) n_pmax = n_p - 1;
571
572 for (unsigned l2 = 1; l2 < n_pmax; l2++)
573 {
574 // Assign the new spine with unit height
575 new_spine_pt = new Spine(1.0);
576 new_spine_pt->spine_height_pt()->pin(0);
577 Spine_pt.push_back(new_spine_pt);
578
579 // Get the node
580 SpineNode* nod_pt = element_node_pt(j, l2);
581 // Set the pointer to spine
582 nod_pt->spine_pt() = new_spine_pt;
583 // Set the fraction
584 nod_pt->fraction() = 0.0;
585 // Pointer to the mesh that implements the update fct
586 nod_pt->spine_mesh_pt() = this;
587 // Set update fct id
588 nod_pt->node_update_fct_id() = 2;
589
590 {
591 // Provide spine with additional storage for wall coordinate
592 // and wall geom object:
593
594 // Get the Lagrangian coordinate in the Lower Wall
595 zeta[0] = zeta_lo + double(j - n_x0 - n_x1) * dzeta +
596 double(l2) * dzeta / (n_p - 1.0);
597 // Get the geometric object and local coordinate
598 Straight_wall_pt->locate_zeta(zeta, geometric_object_pt, s_wall);
599
600 // The local coordinate is a geometric parameter
601 // This needs to be set (rather than added) because the
602 // same spine may be visited more than once
603 Vector<double> parameters(1, s_wall[0]);
604 nod_pt->spine_pt()->set_geom_parameter(parameters);
605
606 // Get position of wall
607 Straight_wall_pt->position(s_wall, r_wall);
608
609 // Adjust spine height
610 nod_pt->spine_pt()->height() = r_wall[1];
611
612 // The sub geom object is one (and only) geom object
613 // for spine:
614 Vector<GeomObject*> geom_object_pt(1);
615 geom_object_pt[0] = geometric_object_pt;
616
617 // Pass geom object(s) to spine
618 nod_pt->spine_pt()->set_geom_object_pt(geom_object_pt);
619 }
620
621 // Loop vertically along the spine
622 // Loop over the elements
623 for (unsigned long i = 0; i < n_y; i++)
624 {
625 // Loop over the vertical nodes, apart from the first
626 for (unsigned l1 = 1; l1 < n_p; l1++)
627 {
628 // Get the node
629 SpineNode* nod_pt = element_node_pt(i * n_x + j, l1 * n_p + l2);
630 // Set the pointer to the spine
631 nod_pt->spine_pt() = new_spine_pt;
632 // Set the fraction
633 nod_pt->fraction() =
634 (double(i) + double(l1) / double(n_p - 1)) / double(n_y);
635 // Pointer to the mesh that implements the update fct
636 nod_pt->spine_mesh_pt() = this;
637 // Set update fct id
638 nod_pt->node_update_fct_id() = 2;
639 }
640 }
641 }
642 }
643
644 // Increment number of previous elements
645 n_prev_elements += n_x2 * n_y;
646
647 // Last spine needs special treatment for periodic meshes
648 // because it's the same as the first one...
649 if (this->Xperiodic)
650 {
651 // Last spine is the same as first one...
652 Spine* final_spine_pt = Spine_pt[0];
653
654 // Get the node
655 SpineNode* nod_pt = element_node_pt((n_x - 1), (n_p - 1));
656
657 // Set the pointer for the first node
658 nod_pt->spine_pt() = final_spine_pt;
659 // Set the fraction to be the same as for the nodes on the first row
660 nod_pt->fraction() = element_node_pt(0, 0)->fraction();
661 // Pointer to the mesh that implements the update fct
662 nod_pt->spine_mesh_pt() = element_node_pt(0, 0)->spine_mesh_pt();
663
664 // Now loop vertically along the spine
665 for (unsigned i = 0; i < n_y; i++)
666 {
667 // Loop over the vertical nodes, apart from the first
668 for (unsigned l1 = 1; l1 < n_p; l1++)
669 {
670 // Get the node
671 SpineNode* nod_pt =
672 element_node_pt(i * n_x + (n_x - 1), l1 * n_p + (n_p - 1));
673
674 // Set the pointer to the spine
675 nod_pt->spine_pt() = final_spine_pt;
676 // Set the fraction to be the same as in first row
677 nod_pt->fraction() = element_node_pt(i * n_x, l1 * n_p)->fraction();
678 // Pointer to the mesh that implements the update fct
679 nod_pt->spine_mesh_pt() =
680 element_node_pt(i * n_x, l1 * n_p)->spine_mesh_pt();
681 }
682 }
683 }
684
685
686 } // end of build_channel_spine_mesh
687
688
689 //======================================================================
690 /// Reorder the elements, so we loop over them vertically first
691 /// (advantageous in "wide" domains if a frontal solver is used).
692 //======================================================================
693 template<class ELEMENT>
695 {
696 unsigned n_x = this->Nx;
697 unsigned n_y = this->Ny;
698 // Find out how many elements there are
699 unsigned long Nelement = nelement();
700 // Find out how many fluid elements there are
701 unsigned long Nfluid = n_x * n_y;
702 // Create a dummy array of elements
704
705 // Loop over the elements in horizontal order
706 for (unsigned long j = 0; j < n_x; j++)
707 {
708 // Loop over the elements in lower layer vertically
709 for (unsigned long i = 0; i < n_y; i++)
710 {
711 // Push back onto the new stack
712 dummy.push_back(finite_element_pt(n_x * i + j));
713 }
714
715 // Push back the line element onto the stack
716 dummy.push_back(finite_element_pt(Nfluid + j));
717 }
718
719 // Now copy the reordered elements into the element_pt
720 for (unsigned long e = 0; e < Nelement; e++)
721 {
722 Element_pt[e] = dummy[e];
723 }
724
725 } // end of element_reorder
726
727
728} // namespace oomph
729
730#endif
e
Definition: cfortran.h:571
cstr elem_len * i
Definition: cfortran.h:603
virtual void build_channel_spine_mesh(TimeStepper *time_stepper_pt)
Helper function to actually build the channel-spine mesh (called from various constructors)
GeomObject * Straight_wall_pt
GeomObject for the straight upper wall.
ChannelSpineMesh(const unsigned &nx0, const unsigned &nx1, const unsigned &nx2, const unsigned &ny, const double &lx0, const double &lx1, const double &lx2, const double &h, GeomObject *wall_pt, TimeStepper *time_stepper_pt=&Mesh::Default_TimeStepper)
Constructor: Pass number of elements in x-direction in regions 0,1 and 2, number of elements in y-dir...
void element_reorder()
Reorder the elements so we loop over them vertically first (advantageous in "wide" domains if a front...
void pin(const unsigned &i)
Pin the i-th stored variable.
Definition: nodes.h:385
/////////////////////////////////////////////////////////////////////
Definition: geom_objects.h:101
An OomphLibError object which should be thrown when an run-time error is encountered....
RectangularQuadMesh is a two-dimensional mesh of Quad elements with Nx elements in the "x" (horizonal...
void build_mesh(TimeStepper *time_stepper_pt=&Mesh::Default_TimeStepper)
Generic mesh construction function: contains all the hard work.
Class for nodes that live on spines. The assumption is that each Node lies at a fixed fraction on a s...
Definition: spines.h:328
double & fraction()
Set reference to fraction along spine.
Definition: spines.h:378
SpineMesh *& spine_mesh_pt()
Access function to Pointer to SpineMesh that this node is a part of and which implements the node upd...
Definition: spines.h:391
Spine *& spine_pt()
Access function to spine.
Definition: spines.h:372
unsigned & node_update_fct_id()
Access function to ID of node update function (within specific mesh)
Definition: spines.h:384
Spines are used for algebraic node update operations in free-surface fluid problems: They form the ba...
Definition: spines.h:64
Data *& spine_height_pt()
Access function to Data object that stores the spine height.
Definition: spines.h:156
void set_geom_object_pt(const Vector< GeomObject * > &geom_object_pt)
Set vector of (pointers to) geometric objects that is involved in the node update operations for this...
Definition: spines.h:225
void set_geom_parameter(const Vector< double > &geom_parameter)
Set vector of geometric parameters that are involved in the node update operations for this Spine....
Definition: spines.h:273
double & height()
Access function to spine height.
Definition: spines.h:150
////////////////////////////////////////////////////////////////////
Definition: geom_objects.h:452
////////////////////////////////////////////////////////////////////// //////////////////////////////...
Definition: timesteppers.h:231
//////////////////////////////////////////////////////////////////// ////////////////////////////////...