channel_spine_mesh.template.h
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_HEADER
27#define OOMPH_CHANNEL_SPINE_MESH_HEADER
28
29// oomph-lib includes
30#include "../generic/spines.h"
32
33namespace oomph
34{
35 //======================================================================
36 /// Spine mesh class derived from standard 2D mesh.
37 /// The mesh contains a StraightLine GeomObject which defines the height
38 /// of the left and right regions (0,2) and another GeomObject is passed
39 /// to the constructor to define the height in the central region.
40 //======================================================================
41 template<class ELEMENT>
42 class ChannelSpineMesh : public RectangularQuadMesh<ELEMENT>, public SpineMesh
43 {
44 public:
45 /// Constructor: Pass number of elements in x-direction in regions
46 /// 0,1 and 2, number of elements in y-direction, length in x direction in
47 /// regions 0,1 and 2, height mesh, pointer to the GeomObject defining the
48 /// heightof the central region and pointer to timestepper (defaults to
49 /// Steady timestepper)
50 ChannelSpineMesh(const unsigned& nx0,
51 const unsigned& nx1,
52 const unsigned& nx2,
53 const unsigned& ny,
54 const double& lx0,
55 const double& lx1,
56 const double& lx2,
57 const double& h,
59 TimeStepper* time_stepper_pt = &Mesh::Default_TimeStepper);
60
61 /// Constructor: Pass number of elements in x-direction in regions
62 /// 0,1 and 2, number of elements in y-direction, length in x direction in
63 /// regions 0,1 and 2, height mesh, pointer to the GeomObject defining the
64 /// heightof the central region, a boolean flag to indicate whether or not
65 /// the mesh is periodic and pointer to timestepper (defaults to Steady
66 /// timestepper)
67 ChannelSpineMesh(const unsigned& nx0,
68 const unsigned& nx1,
69 const unsigned& nx2,
70 const unsigned& ny,
71 const double& lx0,
72 const double& lx1,
73 const double& lx2,
74 const double& h,
76 const bool& periodic_in_x,
77 TimeStepper* time_stepper_pt = &Mesh::Default_TimeStepper);
78
79 /// Access functions for pointers to the \f$ i \f$ -th element in
80 /// the left region.
81 FiniteElement*& left_element_pt(const unsigned long& i)
82 {
83 return Left_element_pt[i];
84 }
85
86 /// Access functions for pointers to the \f$ i \f$ -th element in
87 /// the centre region.
88 FiniteElement*& centre_element_pt(const unsigned long& i)
89 {
90 return Centre_element_pt[i];
91 }
92
93 /// Access functions for pointers to the \f$ i \f$ -th element in
94 /// the right region.
95 FiniteElement*& right_element_pt(const unsigned long& i)
96 {
97 return Right_element_pt[i];
98 }
99
100 /// Number of elements in left region
101 unsigned long nleft() const
102 {
103 return Left_element_pt.size();
104 }
105
106 /// Number of elements in centre region
107 unsigned long ncentre() const
108 {
109 return Centre_element_pt.size();
110 }
111
112 /// Number of elements in right region
113 unsigned long nright() const
114 {
115 return Right_element_pt.size();
116 }
117
118 /// Number of elements in bulk
119 unsigned long nbulk() const
120 {
121 unsigned long Nbulk = Left_element_pt.size() + Centre_element_pt.size() +
122 Right_element_pt.size();
123 return Nbulk;
124 }
125
126 /// Reorder the elements so we loop over them vertically first
127 /// (advantageous in "wide" domains if a frontal solver is used).
128 void element_reorder();
129
130 /// General node update function implements pure virtual function
131 /// defined in SpineMesh base class and performs specific node update
132 /// actions: along vertical spines
133 virtual void spine_node_update(SpineNode* spine_node_pt)
134 {
135 // Get spine node's fraction along the spine
136 double W = spine_node_pt->fraction();
137
138 // Get local coordinates
139 Vector<double> s_wall(1);
140 s_wall[0] = spine_node_pt->spine_pt()->geom_parameter(0);
141
142 // Get position vector to wall
143 Vector<double> position(2);
144 spine_node_pt->spine_pt()->geom_object_pt(0)->position(s_wall, position);
145
146 // Set the value of y
147 spine_node_pt->x(1) = this->Ymin + W * position[1];
148 }
149
150 /// Return the value of the x-coordinate at the node given by the
151 /// local node number (xnode, ynode) in the element (xelement,yelement).
152 /// The description is in a "psudeo" two-dimensional coordinate system,
153 /// so the range of xelement is [0,Nx-1], yelement is [0,Ny-1], and
154 /// that of xnode and ynode is [0,Np-1]. The default is to return
155 /// nodes that are equally spaced in the x coodinate.
156 virtual double x_spacing_function(unsigned xelement,
157 unsigned xnode,
158 unsigned yelement,
159 unsigned ynode)
160 {
161 // Calculate the values of equal increments in nodal values in left region
162 double xstep1 = Lx0 / ((this->Np - 1) * Nx0);
163 // Calculate the values of equal increments in nodal values in centre
164 // region
165 double xstep2 = Lx1 / ((this->Np - 1) * Nx1);
166 // Calculate the values of equal increments in nodal values in right
167 // region
168 double xstep3 = Lx2 / ((this->Np - 1) * Nx2);
169
170 // left region
171 if (xelement < Nx0)
172 {
173 // Return the appropriate value
174 return (this->Xmin + xstep1 * ((this->Np - 1) * xelement + xnode));
175 }
176 // centre region
177 else if (xelement < Nx0 + Nx1)
178 {
179 // Return the appropriate value
180 return (Lx0 + xstep2 * ((this->Np - 1) * (xelement - Nx0) + xnode));
181 }
182 // right region
183 else if (xelement < Nx0 + Nx1 + Nx2)
184 {
185 // Return the appropriate value
186 return (Lx0 + Lx1 +
187 xstep3 * ((this->Np - 1) * (xelement - Nx0 - Nx1) + xnode));
188 }
189 else
190 {
191 throw OomphLibError("Should not have got here",
192 OOMPH_CURRENT_FUNCTION,
193 OOMPH_EXCEPTION_LOCATION);
194 }
195 // Dummy return to keep compiler from barking
196 return 0.0;
197 }
198
199 /// Access function for spines in left region
200 Spine*& left_spine_pt(const unsigned long& i)
201 {
202#ifdef PARNOID
203 if (i > Nleft_spine)
204 {
205 throw OomphLibError("Arguemnt out of range",
206 OOMPH_CURRENT_FUNCTION,
207 OOMPH_EXCEPTION_LOCATION);
208 }
209#endif
210 return Spine_pt[i];
211 }
212
213 /// Access function for spines in centre region
214 Spine*& centre_spine_pt(const unsigned long& i)
215 {
216 if (i > Ncentre_spine)
217 {
218 throw OomphLibError("Arguemnt out of range",
219 OOMPH_CURRENT_FUNCTION,
220 OOMPH_EXCEPTION_LOCATION);
221 }
222 else
223 {
224 return Spine_pt[Nleft_spine + i];
225 }
226 }
227
228 /// Access function for spines in right region
229 Spine*& right_spine_pt(const unsigned long& i)
230 {
231 if (i > Nright_spine)
232 {
233 throw OomphLibError("Arguemnt out of range",
234 OOMPH_CURRENT_FUNCTION,
235 OOMPH_EXCEPTION_LOCATION);
236 }
237 else
238 {
239 return Spine_pt[Nleft_spine + Ncentre_spine - 1 + i];
240 }
241 }
242
243 /// Access function for the number of spines in the left region
244 unsigned nleft_spine()
245 {
246 return Nleft_spine;
247 }
248
249 /// Access function for the number of spines in the centre region
250 unsigned ncentre_spine()
251 {
252 return Ncentre_spine;
253 }
254
255 /// Access function for the number of spines in the right region
256 unsigned nright_spine()
257 {
258 return Nright_spine;
259 }
260
261 /// Access function to the GeomObject for upper wall
263 {
264 return Wall_pt;
265 }
266
267 /// Access function to the GeomObject for the straight upper wall
269 {
270 return Straight_wall_pt;
271 }
272
273 protected:
274 /// Vector of pointers to element in the left region
276
277 /// Vector of pointers to element in the centre region
279
280 /// Vector of pointers to element in the right region
282
283 /// Helper function to actually build the channel-spine mesh
284 /// (called from various constructors)
285 virtual void build_channel_spine_mesh(TimeStepper* time_stepper_pt);
286
287 /// Number of elements in the left region
288 unsigned Nx0;
289
290 /// Number of elements in the centre region
291 unsigned Nx1;
292
293 /// Number of elements in the right region
294 unsigned Nx2;
295
296 /// Length of left region
297 double Lx0;
298
299 /// Length of centre region
300 double Lx1;
301
302 /// Length of right region
303 double Lx2;
304
305 /// Number of spines in left region
306 unsigned Nleft_spine;
307
308 /// Number of spines in centre region
310
311 /// Number of spines in right region
312 unsigned Nright_spine;
313
314 /// GeomObject for upper wall
316
317 /// GeomObject for the straight upper wall
319 };
320
321} // namespace oomph
322
323#endif
cstr elem_len * i
Definition: cfortran.h:603
Spine mesh class derived from standard 2D mesh. The mesh contains a StraightLine GeomObject which def...
GeomObject * wall_pt()
Access function to the GeomObject for upper wall.
Vector< FiniteElement * > Centre_element_pt
Vector of pointers to element in the centre region.
FiniteElement *& left_element_pt(const unsigned long &i)
Access functions for pointers to the -th element in the left region.
virtual void build_channel_spine_mesh(TimeStepper *time_stepper_pt)
Helper function to actually build the channel-spine mesh (called from various constructors)
virtual double x_spacing_function(unsigned xelement, unsigned xnode, unsigned yelement, unsigned ynode)
Return the value of the x-coordinate at the node given by the local node number (xnode,...
unsigned Nx0
Number of elements in the left region.
Spine *& left_spine_pt(const unsigned long &i)
Access function for spines in left region.
GeomObject * Straight_wall_pt
GeomObject for the straight upper wall.
double Lx1
Length of centre region.
double Lx0
Length of left region.
unsigned long ncentre() const
Number of elements in centre region.
unsigned Nx2
Number of elements in the right region.
unsigned nright_spine()
Access function for the number of spines in the right region.
double Lx2
Length of right region.
unsigned long nright() const
Number of elements in right region.
unsigned Nright_spine
Number of spines in right region.
unsigned Nleft_spine
Number of spines in left region.
unsigned Nx1
Number of elements in the centre region.
FiniteElement *& right_element_pt(const unsigned long &i)
Access functions for pointers to the -th element in the right region.
unsigned long nleft() const
Number of elements in left region.
unsigned ncentre_spine()
Access function for the number of spines in the centre region.
unsigned long nbulk() const
Number of elements in bulk.
Spine *& right_spine_pt(const unsigned long &i)
Access function for spines in right region.
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...
Vector< FiniteElement * > Right_element_pt
Vector of pointers to element in the right region.
Spine *& centre_spine_pt(const unsigned long &i)
Access function for spines in centre region.
Vector< FiniteElement * > Left_element_pt
Vector of pointers to element in the left region.
FiniteElement *& centre_element_pt(const unsigned long &i)
Access functions for pointers to the -th element in the centre region.
void element_reorder()
Reorder the elements so we loop over them vertically first (advantageous in "wide" domains if a front...
GeomObject * Wall_pt
GeomObject for upper wall.
virtual void spine_node_update(SpineNode *spine_node_pt)
General node update function implements pure virtual function defined in SpineMesh base class and per...
unsigned nleft_spine()
Access function for the number of spines in the left region.
unsigned Ncentre_spine
Number of spines in centre region.
GeomObject * straight_wall_pt()
Access function to the GeomObject for the straight upper wall.
A general Finite Element class.
Definition: elements.h:1313
/////////////////////////////////////////////////////////////////////
Definition: geom_objects.h:101
virtual void position(const Vector< double > &zeta, Vector< double > &r) const =0
Parametrised position on object at current time: r(zeta).
static Steady< 0 > Default_TimeStepper
Default Steady Timestepper, to be used in default arguments to Mesh constructors.
Definition: mesh.h:75
double & x(const unsigned &i)
Return the i-th nodal coordinate.
Definition: nodes.h:1060
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...
const unsigned & ny() const
Return number of elements in y direction.
unsigned Np
Np: number of (linear) points in the element.
double Ymin
Minimum value of y coordinate.
double Xmin
Minimum value of x coordinate.
/////////////////////////////////////////////////////////////////////// /////////////////////////////...
Definition: spines.h:613
Vector< Spine * > Spine_pt
A Spine mesh contains a Vector of pointers to spines.
Definition: spines.h:616
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
Spine *& spine_pt()
Access function to spine.
Definition: spines.h:372
Spines are used for algebraic node update operations in free-surface fluid problems: They form the ba...
Definition: spines.h:64
double & geom_parameter(const unsigned &i)
Return i-th geometric parameter that is involved in the node update operations for this Spine.
Definition: spines.h:287
GeomObject *& geom_object_pt(const unsigned &i)
Return i-th geometric object that is involved in the node update operations for this Spine.
Definition: spines.h:244
////////////////////////////////////////////////////////////////////// //////////////////////////////...
Definition: timesteppers.h:231
//////////////////////////////////////////////////////////////////// ////////////////////////////////...