my_alg_channel_mesh.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-2024 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_MY_ALGEBRAIC_COLLAPSIBLE_CHANNEL_MESH
27 #define OOMPH_MY_ALGEBRAIC_COLLAPSIBLE_CHANNEL_MESH
28 
29 //Include the mesh
30 #include "meshes/collapsible_channel_mesh.h"
31 
32 
33 
34 
35 
36 /// /////////////////////////////////////////////////////////////////
37 /// /////////////////////////////////////////////////////////////////
38 /// /////////////////////////////////////////////////////////////////
39 
40 namespace oomph
41 {
42 
43 
44 //===========start_algebraic_mesh==================================
45 /// Collapsible channel mesh with algebraic node update
46 //=================================================================
47 template<class ELEMENT>
49  public AlgebraicMesh,
50  public virtual CollapsibleChannelMesh<ELEMENT>
51 {
52 
53 public:
54 
55  /// Constructor: Pass number of elements in upstream/collapsible/
56  /// downstream segment and across the channel; lengths of upstream/
57  /// collapsible/downstream segments and width of channel, pointer to
58  /// GeomObject that defines the collapsible segment, and pointer to
59  /// TimeStepper (defaults to the default timestepper, Steady).
60  MyAlgebraicCollapsibleChannelMesh(const unsigned& nup,
61  const unsigned& ncollapsible,
62  const unsigned& ndown,
63  const unsigned& ny,
64  const double& lup,
65  const double& lcollapsible,
66  const double& ldown,
67  const double& ly,
68  GeomObject* wall_pt,
69  TimeStepper* time_stepper_pt=
70  &Mesh::Default_TimeStepper) :
71  CollapsibleChannelMesh<ELEMENT>(nup, ncollapsible, ndown, ny,
72  lup, lcollapsible, ldown, ly,
73  wall_pt,
74  time_stepper_pt)
75  {
76  // Setup algebraic node update operations
78  }
79 
80 
81 
82  /// Constructor: Pass number of elements in upstream/collapsible/
83  /// downstream segment and across the channel; lengths of upstream/
84  /// collapsible/downstream segments and width of channel, pointer to
85  /// GeomObject that defines the collapsible segment, function pointer
86  /// to "boundary layer squash function", and pointer to
87  /// TimeStepper (defaults to the default timestepper, Steady).
89  const unsigned& nup,
90  const unsigned& ncollapsible,
91  const unsigned& ndown,
92  const unsigned& ny,
93  const double& lup,
94  const double& lcollapsible,
95  const double& ldown,
96  const double& ly,
97  GeomObject* wall_pt,
98  CollapsibleChannelDomain::BLSquashFctPt bl_squash_function_pt,
99  TimeStepper* time_stepper_pt=
100  &Mesh::Default_TimeStepper) :
101  CollapsibleChannelMesh<ELEMENT>(nup, ncollapsible, ndown, ny,
102  lup, lcollapsible, ldown, ly,
103  wall_pt,
104  time_stepper_pt)
105  {
106  // Set boundary layer squash function
107  this->Domain_pt->bl_squash_fct_pt()=bl_squash_function_pt;
108 
109  // Do MacroElement-based node update
110  CollapsibleChannelMesh<ELEMENT>::node_update();
111 
112  // Setup algebraic node update operations
114  }
115 
116 
117  /// Destructor: empty
119 
120  /// Function pointer for function that squashes
121  /// the mesh near the walls. Default trivial mapping (the identity)
122  /// leaves vertical nodal positions unchanged. Mapping is
123  /// used in underlying CollapsibleChannelDomain.
124  /// This (deliberately broken) function overloads the one
125  /// in the CollapsibleChannelMesh base class.
126  CollapsibleChannelDomain::BLSquashFctPt& bl_squash_fct_pt()
127  {
128  std::ostringstream error_message;
129  error_message
130  << "It does not make sense to set the bl_squash_fct_pt \n"
131  << "outside the constructor as it's only used to set up the \n"
132  << "algebraic remesh data when the algebraic mesh is first built. \n";
133  std::string function_name =
134  "MyAlgebraicCollapsibleChannelMesh::bl_squash_fct_pt()\n";
135 
136  throw OomphLibError(error_message.str(),
137  OOMPH_CURRENT_FUNCTION,
138  OOMPH_EXCEPTION_LOCATION);
139 
140  // Dummy return
141  return Dummy_fct_pt;
142  }
143 
144 
145  /// Update nodal position at time level t (t=0: present;
146  /// t>0: previous)
147  void algebraic_node_update(const unsigned& t, AlgebraicNode*& node_pt);
148 
149  /// Update the geometric references that are used
150  /// to update node after mesh adaptation.
151  /// Empty -- no update of node update required
152  void update_node_update(AlgebraicNode*& node_pt) {}
153 
154 protected:
155 
156  /// Function to setup the algebraic node update
158 
159 
160  /// Dummy function pointer
161  CollapsibleChannelDomain::BLSquashFctPt Dummy_fct_pt;
162 
163 };
164 
165 
166 
167 
168 
169 /// ////////////////////////////////////////////////////////////////////////
170 /// ////////////////////////////////////////////////////////////////////////
171 /// ////////////////////////////////////////////////////////////////////////
172 
173 
174 
175 //===========start_refineable_algebraic_collapsible_channel_mesh======
176 /// Refineable version of the CollapsibleChannel mesh with
177 /// algebraic node update.
178 //====================================================================
179 template<class ELEMENT>
181  public RefineableQuadMesh<ELEMENT>,
182  public virtual MyAlgebraicCollapsibleChannelMesh<ELEMENT>
183 {
184 
185 public:
186 
187 
188  /// Constructor: Pass number of elements in upstream/collapsible/
189  /// downstream segment and across the channel; lengths of upstream/
190  /// collapsible/downstream segments and width of channel, pointer to
191  /// GeomObject that defines the collapsible segment, function pointer
192  /// to "boundary layer squash function", and pointer to
193  /// TimeStepper (defaults to the default timestepper, Steady).
195  const unsigned& ncollapsible,
196  const unsigned& ndown,
197  const unsigned& ny,
198  const double& lup,
199  const double& lcollapsible,
200  const double& ldown,
201  const double& ly,
202  GeomObject* wall_pt,
203  TimeStepper* time_stepper_pt=
204  &Mesh::Default_TimeStepper) :
205  CollapsibleChannelMesh<ELEMENT>(nup, ncollapsible, ndown, ny,
206  lup, lcollapsible, ldown, ly,
207  wall_pt,
208  time_stepper_pt),
209  MyAlgebraicCollapsibleChannelMesh<ELEMENT>(nup, ncollapsible, ndown, ny,
210  lup, lcollapsible, ldown, ly,
211  wall_pt,
212  time_stepper_pt)
213  {
214  // Build quadtree forest
215  this->setup_quadtree_forest();
216  }
217 
218 
219 
220 
221  /// Constructor: Pass number of elements in upstream/collapsible/
222  /// downstream segment and across the channel; lengths of upstream/
223  /// collapsible/downstream segments and width of channel, pointer to
224  /// GeomObject that defines the collapsible segment, function pointer
225  /// to "boundary layer squash function", and pointer to
226  /// TimeStepper (defaults to the default timestepper, Steady).
228  const unsigned& nup,
229  const unsigned& ncollapsible,
230  const unsigned& ndown,
231  const unsigned& ny,
232  const double& lup,
233  const double& lcollapsible,
234  const double& ldown,
235  const double& ly,
236  GeomObject* wall_pt,
237  CollapsibleChannelDomain::BLSquashFctPt bl_squash_function_pt,
238  TimeStepper* time_stepper_pt=
239  &Mesh::Default_TimeStepper) :
240  CollapsibleChannelMesh<ELEMENT>(nup, ncollapsible, ndown, ny,
241  lup, lcollapsible, ldown, ly,
242  wall_pt,
243  time_stepper_pt),
244  MyAlgebraicCollapsibleChannelMesh<ELEMENT>(nup, ncollapsible, ndown, ny,
245  lup, lcollapsible, ldown, ly,
246  wall_pt,
247  bl_squash_function_pt,
248  time_stepper_pt)
249  {
250  // Build quadtree forest
251  this->setup_quadtree_forest();
252  }
253 
254 };
255 
256 
257 
258 
259 
260 
261 
262 
263 //============start_setup_algebraic_node_update====================
264 /// Setup algebraic mesh update -- assumes that mesh has
265 /// initially been set up with the wall in its undeformed position.
266 //=================================================================
267 template<class ELEMENT>
269 {
270 
271  // Extract some reference lengths from the CollapsibleChannelDomain.
272  double l_up=this->domain_pt()->l_up();
273  double l_collapsible=this->domain_pt()->l_collapsible();
274 
275  // Loop over all nodes in mesh
276  unsigned nnod=this->nnode();
277  for (unsigned j=0;j<nnod;j++)
278  {
279  // Get pointer to node
280  AlgebraicNode* nod_pt=node_pt(j);
281 
282  // Get coordinates
283  double x=nod_pt->x(0);
284  double y=nod_pt->x(1);
285 
286  // Check if the node is in the collapsible part:
287  if ( (x>=l_up) && (x<=(l_up+l_collapsible)) )
288  {
289 
290  // Get zeta coordinate on the undeformed wall
291  Vector<double> zeta(1);
292  zeta[0]=x-l_up;
293 
294  // Get position vector to wall:
295  Vector<double> r_wall(2);
296  this->Wall_pt->position(zeta,r_wall);
297 
298 
299  // Sanity check: Confirm that the wall is in its undeformed position
300 #ifdef PARANOID
301  if ((std::abs(r_wall[0]-x)>1.0e-15)&&(std::abs(r_wall[1]-y)>1.0e-15))
302  {
303  std::ostringstream error_stream;
304  error_stream
305  << "Wall must be in its undeformed position when\n"
306  << "algebraic node update information is set up!\n "
307  << "x-discrepancy: " << std::abs(r_wall[0]-x) << std::endl
308  << "y-discrepancy: " << std::abs(r_wall[1]-y) << std::endl;
309 
310  throw OomphLibError(
311  error_stream.str(),
312  OOMPH_CURRENT_FUNCTION,
313  OOMPH_EXCEPTION_LOCATION);
314  }
315 #endif
316 
317  // Only a single geometric object is involved in the node update operation
318  Vector<GeomObject*> geom_object_pt(1);
319 
320  // The wall geometric object
321  geom_object_pt[0]=this->Wall_pt;
322 
323  // The update function requires three parameters:
324  Vector<double> ref_value(3);
325 
326  // First reference value: Original x-position on the lower wall
327  ref_value[0]=r_wall[0];
328 
329  // Second reference value: Fractional position along
330  // straight line from the bottom (at the original x-position)
331  // to the point on the wall
332  ref_value[1]=y/r_wall[1];
333 
334  // Third reference value: Zeta coordinate on wall
335  ref_value[2]=zeta[0];
336 
337  // Setup algebraic update for node: Pass update information
338  // to AlgebraicNode:
339  nod_pt->add_node_update_info(
340  this, // mesh
341  geom_object_pt, // vector of geom objects
342  ref_value); // vector of ref. values
343  }
344 
345  }
346 
347 } //end of setup_algebraic_node_update
348 
349 
350 //=============start_of_algebraic_node_update======================
351 /// Perform algebraic mesh update at time level t (t=0: present;
352 /// t>0: previous)
353 //=================================================================
354 template<class ELEMENT>
356  const unsigned& t, AlgebraicNode*& node_pt)
357 {
358 
359  // Extract reference values for update by copy construction
360  Vector<double> ref_value(node_pt->vector_ref_value());
361 
362  // Extract geometric objects for update by copy construction
363  Vector<GeomObject*> geom_object_pt(node_pt->vector_geom_object_pt());
364 
365  // First reference value: Original x-position
366  double x_bottom=ref_value[0];
367 
368  // Second reference value: Fractional position along
369  // straight line from the bottom (at the original x-position)
370  // to the point on the wall
371  double fract=ref_value[1];
372 
373  // Third reference value: Zeta coordinate on wall
374  Vector<double> zeta(1);
375  zeta[0]=ref_value[2];
376 
377  // Pointer to wall geom object
378  GeomObject* wall_pt=geom_object_pt[0];
379 
380  // Get position vector to wall at previous timestep t
381  Vector<double> r_wall(2);
382  wall_pt->position(t,zeta,r_wall);
383 
384  // Assign new nodal coordinates
385  node_pt->x(t,0)=x_bottom+fract*(r_wall[0]-x_bottom);
386  node_pt->x(t,1)= fract* r_wall[1];
387 
388 }
389 
390 
391 
392 
393 }
394 
395 
396 
397 
398 
399 
400 
401 
402 #endif
Collapsible channel mesh with algebraic node update.
MyAlgebraicCollapsibleChannelMesh(const unsigned &nup, const unsigned &ncollapsible, const unsigned &ndown, const unsigned &ny, const double &lup, const double &lcollapsible, const double &ldown, const double &ly, GeomObject *wall_pt, TimeStepper *time_stepper_pt=&Mesh::Default_TimeStepper)
Constructor: Pass number of elements in upstream/collapsible/ downstream segment and across the chann...
void algebraic_node_update(const unsigned &t, AlgebraicNode *&node_pt)
Update nodal position at time level t (t=0: present; t>0: previous)
MyAlgebraicCollapsibleChannelMesh(const unsigned &nup, const unsigned &ncollapsible, const unsigned &ndown, const unsigned &ny, const double &lup, const double &lcollapsible, const double &ldown, const double &ly, GeomObject *wall_pt, CollapsibleChannelDomain::BLSquashFctPt bl_squash_function_pt, TimeStepper *time_stepper_pt=&Mesh::Default_TimeStepper)
Constructor: Pass number of elements in upstream/collapsible/ downstream segment and across the chann...
virtual ~MyAlgebraicCollapsibleChannelMesh()
Destructor: empty.
CollapsibleChannelDomain::BLSquashFctPt & bl_squash_fct_pt()
Function pointer for function that squashes the mesh near the walls. Default trivial mapping (the ide...
void setup_algebraic_node_update()
Function to setup the algebraic node update.
CollapsibleChannelDomain::BLSquashFctPt Dummy_fct_pt
Dummy function pointer.
void update_node_update(AlgebraicNode *&node_pt)
Update the geometric references that are used to update node after mesh adaptation....
//////////////////////////////////////////////////////////////////////// ////////////////////////////...
MyRefineableAlgebraicCollapsibleChannelMesh(const unsigned &nup, const unsigned &ncollapsible, const unsigned &ndown, const unsigned &ny, const double &lup, const double &lcollapsible, const double &ldown, const double &ly, GeomObject *wall_pt, CollapsibleChannelDomain::BLSquashFctPt bl_squash_function_pt, TimeStepper *time_stepper_pt=&Mesh::Default_TimeStepper)
Constructor: Pass number of elements in upstream/collapsible/ downstream segment and across the chann...
MyRefineableAlgebraicCollapsibleChannelMesh(const unsigned &nup, const unsigned &ncollapsible, const unsigned &ndown, const unsigned &ny, const double &lup, const double &lcollapsible, const double &ldown, const double &ly, GeomObject *wall_pt, TimeStepper *time_stepper_pt=&Mesh::Default_TimeStepper)
Constructor: Pass number of elements in upstream/collapsible/ downstream segment and across the chann...
///////////////////////////////////////////////////////////////// ///////////////////////////////////...