annular_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// Version 1.0; svn revision $LastChangedRevision: 1097 $
7 //LIC//
8 //LIC// $LastChangedDate: 2015-12-17 11:53:17 +0000 (Thu, 17 Dec 2015) $
9 //LIC//
10 //LIC// Copyright (C) 2006-2016 Matthias Heil and Andrew Hazel
11 //LIC//
12 //LIC// This library is free software; you can redistribute it and/or
13 //LIC// modify it under the terms of the GNU Lesser General Public
14 //LIC// License as published by the Free Software Foundation; either
15 //LIC// version 2.1 of the License, or (at your option) any later version.
16 //LIC//
17 //LIC// This library is distributed in the hope that it will be useful,
18 //LIC// but WITHOUT ANY WARRANTY; without even the implied warranty of
19 //LIC// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
20 //LIC// Lesser General Public License for more details.
21 //LIC//
22 //LIC// You should have received a copy of the GNU Lesser General Public
23 //LIC// License along with this library; if not, write to the Free Software
24 //LIC// Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
25 //LIC// 02110-1301 USA.
26 //LIC//
27 //LIC// The authors may be contacted at oomph-lib@maths.man.ac.uk.
28 //LIC//
29 //LIC//====================================================================
30 
31 #ifndef OOMPH_ANNULAR_MESH_HEADER
32 #define OOMPH_ANNULAR_MESH_HEADER
33 
34 // Headers
37 
38 
39 //Include the headers file for domain
40 #include "annular_domain.h"
41 
42 namespace oomph
43 {
44 
45 
46 //===================================================================
47 /// 2D annular mesh with a unit circle in the middle and a layer
48 /// of thickness h surrounding it.
49 //==================================================================
50 template<class ELEMENT>
51 class TwoDAnnularMesh : public virtual RectangularQuadMesh<ELEMENT>
52 {
53 
54 public:
55 
56 
57  ///Constructor
58  TwoDAnnularMesh(const bool& periodic, const double& azimuthal_fraction,
59  const unsigned &ntheta, const unsigned &nr,
60  const double& a, const double &h, TimeStepper* time_stepper_pt
62  RectangularQuadMesh<ELEMENT>(ntheta,nr,1.0,1.0,
63  periodic,time_stepper_pt)
64  {
65  // Mesh can only be built with 2D Qelements.
66  MeshChecker::assert_geometric_element<QElementGeometricBase,ELEMENT>(2);
67 
68  // Wrap mesh into annular shape
69  double phi=0.0;
70  wrap_into_annular_shape(a,h,azimuthal_fraction,phi);
71  }
72 
73  ///Constructor; rotate mesh by angle phi.
74  TwoDAnnularMesh(const bool& periodic, const double& azimuthal_fraction,
75  const unsigned &ntheta, const unsigned &nr,
76  const double& a, const double &h,
77  const double& phi, TimeStepper* time_stepper_pt
79  RectangularQuadMesh<ELEMENT>(ntheta,nr,1.0,1.0,
80  periodic,time_stepper_pt)
81  {
82  // Mesh can only be built with 2D Qelements.
83  MeshChecker::assert_geometric_element<QElementGeometricBase,ELEMENT>(2);
84 
85  // Wrap mesh into annular shape
86  wrap_into_annular_shape(a,h,azimuthal_fraction,phi);
87  }
88 
89 
90  private:
91 
92  /// Wrap mesh into annular shape
93  void wrap_into_annular_shape(const double& a, const double& h,
94  const double& azimuthal_fraction,
95  const double& phi);
96 
97 };
98 
99 
100 //////////////////////////////////////////////////////////////////////
101 //////////////////////////////////////////////////////////////////////
102 //////////////////////////////////////////////////////////////////////
103 
104 
105 //===================================================================
106 /// Refineable 2D annular mesh with a unit circle in the middle and a layer
107 /// of thickness h surrounding it.
108 //==================================================================
109 template<class ELEMENT>
110  class RefineableTwoDAnnularMesh : public virtual TwoDAnnularMesh<ELEMENT>,
111  public virtual RefineableQuadMesh<ELEMENT>
112 {
113 
114 public:
115 
116  ///Constructor
117  RefineableTwoDAnnularMesh(const bool& periodic,
118  const double& azimuthal_fraction,
119  const unsigned &ntheta,
120  const unsigned &nr,
121  const double& a,
122  const double &h,
123  TimeStepper* time_stepper_pt
125  RectangularQuadMesh<ELEMENT>(ntheta,nr,1.0,1.0,
126  periodic,time_stepper_pt),
127  TwoDAnnularMesh<ELEMENT>(periodic,azimuthal_fraction,ntheta,
128  nr,a,h,time_stepper_pt)
129  {
130  // Mesh can only be built with 2D Qelements.
131  MeshChecker::assert_geometric_element<QElementGeometricBase,ELEMENT>(2);
132 
133  // Build macro element-based domain
134  double phi=0.0;
135  Domain_pt = new AnnularDomain(azimuthal_fraction,
136  ntheta,nr,a,h,phi);
137 
138  // Loop over all elements and set macro element pointer
139  unsigned nel=this->nelement();
140  for (unsigned ielem=0;ielem<nel;ielem++)
141  {
142  dynamic_cast<RefineableQElement<2>*>(this->element_pt(ielem))->
143  set_macro_elem_pt(this->Domain_pt->macro_element_pt(ielem));
144  }
145 
146  // Update nodal positions based on macro-element representation
147  this->node_update();
148 
149  // Nodal positions etc. were created in constructor for
150  // RectangularMesh<...>. Only need to setup quadtree forest
151  this->setup_quadtree_forest();
152 
153  //Setup the periodic neighbour information of the TreeRoots
154  //Cast to specific elements
155  if (periodic)
156  {
157  Vector<TreeRoot*> left_root_pt(nr);
158  Vector<TreeRoot*> right_root_pt(nr);
159  for(unsigned i=0;i<nr;i++)
160  {
161  left_root_pt[i] =
162  dynamic_cast<ELEMENT*>(this->element_pt(i*ntheta))->
163  tree_pt()->root_pt();
164 
165  right_root_pt[i] =
166  dynamic_cast<ELEMENT*>(this->element_pt((i+1)*ntheta-1))->
167  tree_pt()->root_pt();
168  }
169 
170  //Set the neighbour and periodicity
171  using namespace QuadTreeNames;
172  for(unsigned i=0;i<nr;i++)
173  {
174  left_root_pt[i]->neighbour_pt(W) = right_root_pt[i];
175  left_root_pt[i]->set_neighbour_periodic(W);
176 
177  right_root_pt[i]->neighbour_pt(E) = left_root_pt[i];
178  right_root_pt[i]->set_neighbour_periodic(E);
179  }
180  }
181  }
182 
183 
184 
185  ///Constructor; rotate mesh by angle phi
186  RefineableTwoDAnnularMesh(const bool& periodic,
187  const double& azimuthal_fraction,
188  const unsigned &ntheta,
189  const unsigned &nr,
190  const double& a,
191  const double &h,
192  const double &phi,
193  TimeStepper* time_stepper_pt
195  RectangularQuadMesh<ELEMENT>(ntheta,nr,1.0,1.0,
196  periodic,time_stepper_pt),
197  TwoDAnnularMesh<ELEMENT>(periodic,azimuthal_fraction,ntheta,
198  nr,a,h,phi,time_stepper_pt)
199  {
200  // Mesh can only be built with 2D Qelements.
201  MeshChecker::assert_geometric_element<QElementGeometricBase,ELEMENT>(2);
202 
203  // Build macro element-based domain
204  Domain_pt = new AnnularDomain(azimuthal_fraction,
205  ntheta,nr,a,h,phi);
206 
207  // Loop over all elements and set macro element pointer
208  unsigned nel=this->nelement();
209  for (unsigned ielem=0;ielem<nel;ielem++)
210  {
211  dynamic_cast<RefineableQElement<2>*>(this->element_pt(ielem))->
212  set_macro_elem_pt(this->Domain_pt->macro_element_pt(ielem));
213  }
214 
215  // Update nodal positions based on macro-element representation
216  this->node_update();
217 
218  // Nodal positions etc. were created in constructor for
219  // RectangularMesh<...>. Only need to setup quadtree forest
220  this->setup_quadtree_forest();
221 
222  //Setup the periodic neighbour information of the TreeRoots
223  //Cast to specific elements
224  if (periodic)
225  {
226  Vector<TreeRoot*> left_root_pt(nr);
227  Vector<TreeRoot*> right_root_pt(nr);
228  for(unsigned i=0;i<nr;i++)
229  {
230  left_root_pt[i] =
231  dynamic_cast<ELEMENT*>(this->element_pt(i*ntheta))->
232  tree_pt()->root_pt();
233 
234  right_root_pt[i] =
235  dynamic_cast<ELEMENT*>(this->element_pt((i+1)*ntheta-1))->
236  tree_pt()->root_pt();
237  }
238 
239  //Set the neighbour and periodicity
240  using namespace QuadTreeNames;
241  for(unsigned i=0;i<nr;i++)
242  {
243  left_root_pt[i]->neighbour_pt(W) = right_root_pt[i];
244  left_root_pt[i]->set_neighbour_periodic(W);
245 
246  right_root_pt[i]->neighbour_pt(E) = left_root_pt[i];
247  right_root_pt[i]->set_neighbour_periodic(E);
248  }
249  }
250  }
251 
252  private:
253 
254 
255  /// Pointer to domain
257 
258 };
259 
260 }
261 
262 #endif
static Steady< 0 > Default_TimeStepper
Default Steady Timestepper, to be used in default arguments to Mesh constructors. ...
Definition: mesh.h:85
void wrap_into_annular_shape(const double &a, const double &h, const double &azimuthal_fraction, const double &phi)
Wrap mesh into annular shape.
cstr elem_len * i
Definition: cfortran.h:607
void setup_quadtree_forest()
Set up QuadTreeForest. Wipes any existing tree structure below the minimum refinement level and regar...
MacroElement * macro_element_pt(const unsigned &i)
Access to i-th macro element.
Definition: domain.h:100
AnnularDomain * Domain_pt
Pointer to domain.
unsigned long nelement() const
Return number of elements in the mesh.
Definition: mesh.h:587
TwoDAnnularMesh(const bool &periodic, const double &azimuthal_fraction, const unsigned &ntheta, const unsigned &nr, const double &a, const double &h, TimeStepper *time_stepper_pt=&Mesh::Default_TimeStepper)
Constructor.
RefineableTwoDAnnularMesh(const bool &periodic, const double &azimuthal_fraction, const unsigned &ntheta, const unsigned &nr, const double &a, const double &h, const double &phi, TimeStepper *time_stepper_pt=&Mesh::Default_TimeStepper)
Constructor; rotate mesh by angle phi.
const Vector< GeneralisedElement * > & element_pt() const
Return reference to the Vector of elements.
Definition: mesh.h:470
TwoDAnnularMesh(const bool &periodic, const double &azimuthal_fraction, const unsigned &ntheta, const unsigned &nr, const double &a, const double &h, const double &phi, TimeStepper *time_stepper_pt=&Mesh::Default_TimeStepper)
Constructor; rotate mesh by angle phi.
virtual void node_update(const bool &update_all_solid_nodes=false)
Update nodal positions in response to changes in the domain shape. Uses the FiniteElement::get_x(...) function for FiniteElements and doesn't do anything for other element types. If a MacroElement pointer has been set for a FiniteElement, the MacroElement representation is used to update the nodal positions; if not get_x(...) uses the FE interpolation and thus leaves the nodal positions unchanged. Virtual, so it can be overloaded by specific meshes, such as AlgebraicMeshes or SpineMeshes. Generally, this function updates the position of all nodes in response to changes in the boundary position. However, we ignore all SolidNodes since their position is computed as part of the solution – unless the bool flag is set to true. Such calls are typically made when the initial mesh is created and/or after a mesh has been refined repeatedly before the start of the computation.
Definition: mesh.cc:290
Annular domain.
Base class for time-stepping schemes. Timestepper provides an approximation of the temporal derivativ...
Definition: timesteppers.h:219
RefineableTwoDAnnularMesh(const bool &periodic, const double &azimuthal_fraction, const unsigned &ntheta, const unsigned &nr, const double &a, const double &h, TimeStepper *time_stepper_pt=&Mesh::Default_TimeStepper)
Constructor.