single_layer_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-2023 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_SINGLE_LAYER_SPINE_MESH_TEMPLATE_CC
27 #define OOMPH_SINGLE_LAYER_SPINE_MESH_TEMPLATE_CC
28 
31 
32 
33 namespace 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 must be called with spinified elements and it
41  /// constructs the spines and contains the information on how to update
42  /// the nodal positions within the mesh as a function of the spine lengths.
43  /// Equations that determine the spine heights (even if they are pinned)
44  /// must be specified externally or else there will be problems.
45 
46  //===========================================================================
47  template<class ELEMENT>
49  const unsigned& nx,
50  const unsigned& ny,
51  const double& lx,
52  const double& h,
53  TimeStepper* time_stepper_pt)
54  : RectangularQuadMesh<ELEMENT>(
55  nx, ny, 0.0, lx, 0.0, h, false, false, time_stepper_pt)
56  {
57  // Mesh can only be built with 2D Qelements.
58  MeshChecker::assert_geometric_element<QElementGeometricBase, ELEMENT>(2);
59 
60  // Mesh can only be built with spine elements
61  MeshChecker::assert_geometric_element<SpineFiniteElement, ELEMENT>(2);
62 
63  // We've called the "generic" constructor for the RectangularQuadMesh
64  // which doesn't do much...
65 
66  // Now build the mesh:
67  build_single_layer_mesh(time_stepper_pt);
68  }
69 
70 
71  //===========================================================================
72  /// Constuctor for spine 2D mesh: Pass number of elements in x-direction,
73  /// number of elements in y-direction, axial length and height of layer,
74  /// a boolean flag to make the mesh periodic in the x-direction,
75  /// and pointer to timestepper (defaults to Static timestepper).
76  ///
77  /// The mesh must be called with spinified elements and it
78  /// constructs the spines and contains the information on how to update
79  /// the nodal positions within the mesh as a function of the spine lengths.
80  /// Equations that determine the spine heights (even if they are pinned)
81  /// must be specified externally or else there will be problems.
82  //===========================================================================
83  template<class ELEMENT>
85  const unsigned& nx,
86  const unsigned& ny,
87  const double& lx,
88  const double& h,
89  const bool& periodic_in_x,
90  TimeStepper* time_stepper_pt)
91  : RectangularQuadMesh<ELEMENT>(
92  nx, ny, 0.0, lx, 0.0, h, periodic_in_x, false, time_stepper_pt)
93  {
94  // Mesh can only be built with 2D Qelements.
95  MeshChecker::assert_geometric_element<QElementGeometricBase, ELEMENT>(2);
96 
97  // Mesh can only be built with spine elements
98  MeshChecker::assert_geometric_element<SpineFiniteElement, ELEMENT>(2);
99 
100 
101  // We've called the "generic" constructor for the RectangularQuadMesh
102  // which doesn't do much...
103 
104  // Now build the mesh:
105  build_single_layer_mesh(time_stepper_pt);
106  }
107 
108 
109  //===========================================================================
110  /// Helper function that actually builds the single-layer spine mesh
111  /// based on the parameters set in the various constructors
112  //===========================================================================
113  template<class ELEMENT>
115  TimeStepper* time_stepper_pt)
116  {
117  // Mesh can only be built with 2D Qelements.
118  MeshChecker::assert_geometric_element<QElementGeometricBase, ELEMENT>(2);
119 
120  // Build the underlying quad mesh:
122 
123  // Read out the number of elements in the x-direction
124  unsigned n_x = this->Nx;
125  unsigned n_y = this->Ny;
126 
127  // Allocate memory for the spines and fractions along spines
128  //---------------------------------------------------------
129 
130  // Read out number of linear points in the element
131  unsigned n_p = dynamic_cast<ELEMENT*>(finite_element_pt(0))->nnode_1d();
132 
133  // Allocate store for the spines:
134  if (this->Xperiodic)
135  {
136  Spine_pt.reserve((n_p - 1) * n_x);
137  }
138  else
139  {
140  Spine_pt.reserve((n_p - 1) * n_x + 1);
141  }
142 
143 
144  // FIRST SPINE
145  //-----------
146 
147  // Element 0
148  // Node 0
149  // Assign the new spine with unit length
150  Spine* new_spine_pt = new Spine(1.0);
151  Spine_pt.push_back(new_spine_pt);
152 
153 
154  // Get pointer to node
155  SpineNode* nod_pt = element_node_pt(0, 0);
156  // Set the pointer to the spine
157  nod_pt->spine_pt() = new_spine_pt;
158  // Set the fraction
159  nod_pt->fraction() = 0.0;
160  // Pointer to the mesh that implements the update fct
161  nod_pt->spine_mesh_pt() = this;
162 
163  // Loop vertically along the spine
164  // Loop over the elements
165  for (unsigned long i = 0; i < n_y; i++)
166  {
167  // Loop over the vertical nodes, apart from the first
168  for (unsigned l1 = 1; l1 < n_p; l1++)
169  {
170  // Get pointer to node
171  SpineNode* nod_pt = element_node_pt(i * n_x, l1 * n_p);
172  // Set the pointer to the spine
173  nod_pt->spine_pt() = new_spine_pt;
174  // Set the fraction
175  nod_pt->fraction() =
176  (double(i) + double(l1) / double(n_p - 1)) / double(n_y);
177  // Pointer to the mesh that implements the update fct
178  nod_pt->spine_mesh_pt() = this;
179  }
180  }
181 
182 
183  // LOOP OVER OTHER SPINES
184  //----------------------
185 
186  // Now loop over the elements horizontally
187  for (unsigned long j = 0; j < n_x; j++)
188  {
189  // Loop over the nodes in the elements horizontally, ignoring
190  // the first column
191 
192  // Last spine needs special treatment in x-periodic meshes:
193  unsigned n_pmax = n_p;
194  if ((this->Xperiodic) && (j == n_x - 1)) n_pmax = n_p - 1;
195 
196  for (unsigned l2 = 1; l2 < n_pmax; l2++)
197  {
198  // Assign the new spine with unit height
199  new_spine_pt = new Spine(1.0);
200  Spine_pt.push_back(new_spine_pt);
201 
202  // Get the node
203  SpineNode* nod_pt = element_node_pt(j, l2);
204  // Set the pointer to spine
205  nod_pt->spine_pt() = new_spine_pt;
206  // Set the fraction
207  nod_pt->fraction() = 0.0;
208  // Pointer to the mesh that implements the update fct
209  nod_pt->spine_mesh_pt() = this;
210 
211  // Loop vertically along the spine
212  // Loop over the elements
213  for (unsigned long i = 0; i < n_y; i++)
214  {
215  // Loop over the vertical nodes, apart from the first
216  for (unsigned l1 = 1; l1 < n_p; l1++)
217  {
218  // Get the node
219  SpineNode* nod_pt = element_node_pt(i * n_x + j, l1 * n_p + l2);
220  // Set the pointer to the spine
221  nod_pt->spine_pt() = new_spine_pt;
222  // Set the fraction
223  nod_pt->fraction() =
224  (double(i) + double(l1) / double(n_p - 1)) / double(n_y);
225  // Pointer to the mesh that implements the update fct
226  nod_pt->spine_mesh_pt() = this;
227  }
228  }
229  }
230  }
231 
232 
233  // Last spine needs special treatment for periodic meshes
234  // because it's the same as the first one...
235  if (this->Xperiodic)
236  {
237  // Last spine is the same as first one...
238  Spine* final_spine_pt = Spine_pt[0];
239 
240  // Get the node
241  SpineNode* nod_pt = element_node_pt((n_x - 1), (n_p - 1));
242 
243  // Set the pointer for the first node
244  nod_pt->spine_pt() = final_spine_pt;
245  // Set the fraction to be the same as for the nodes on the first row
246  nod_pt->fraction() = element_node_pt(0, 0)->fraction();
247  // Pointer to the mesh that implements the update fct
248  nod_pt->spine_mesh_pt() = element_node_pt(0, 0)->spine_mesh_pt();
249 
250  // Now loop vertically along the spine
251  for (unsigned i = 0; i < n_y; i++)
252  {
253  // Loop over the vertical nodes, apart from the first
254  for (unsigned l1 = 1; l1 < n_p; l1++)
255  {
256  // Get the node
257  SpineNode* nod_pt =
258  element_node_pt(i * n_x + (n_x - 1), l1 * n_p + (n_p - 1));
259 
260  // Set the pointer to the spine
261  nod_pt->spine_pt() = final_spine_pt;
262  // Set the fraction to be the same as in first row
263  nod_pt->fraction() = element_node_pt(i * n_x, l1 * n_p)->fraction();
264  // Pointer to the mesh that implements the update fct
265  nod_pt->spine_mesh_pt() =
266  element_node_pt(i * n_x, l1 * n_p)->spine_mesh_pt();
267  }
268  }
269  }
270  }
271 
272 } // namespace oomph
273 #endif
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.
SingleLayerSpineMesh(const unsigned &nx, const unsigned &ny, const double &lx, const double &h, TimeStepper *time_stepper_pt=&Mesh::Default_TimeStepper)
Constructor: Pass number of elements in x-direction, number of elements in y-direction,...
virtual void build_single_layer_mesh(TimeStepper *time_stepper_pt)
Helper function to actually build the single-layer spine mesh (called from various constructors)
////////////////////////////////////////////////////////////////////// //////////////////////////////...