rectangle_with_hole_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_RECTANGLE_WITH_HOLE_MESH_HEADER
27#define OOMPH_RECTANGLE_WITH_HOLE_MESH_HEADER
28
29// Config header generated by autoconfig
30#ifdef HAVE_CONFIG_H
31#include <oomph-lib-config.h>
32#endif
33
34// OOMPH-LIB headers
35#include "../generic/mesh.h"
36#include "../generic/quadtree.h"
37#include "../generic/quad_mesh.h"
38#include "../generic/refineable_quad_mesh.h"
40
41namespace oomph
42{
43 //=============================================================
44 /// Domain-based mesh for rectangular mesh with circular hole
45 //=============================================================
46 template<class ELEMENT>
47 class RectangleWithHoleMesh : public virtual Mesh
48 {
49 public:
50 /// Constructor: Pass pointer to geometric object that
51 /// represents the cylinder, the length and height of the domain.
52 /// The GeomObject must be parametrised such that
53 /// \f$\zeta \in [0,2\pi]\f$ sweeps around the circumference
54 /// in anticlockwise direction. Timestepper defaults to Steady
55 /// default timestepper.
57 GeomObject* cylinder_pt,
58 const double& length,
59 TimeStepper* time_stepper_pt = &Mesh::Default_TimeStepper)
60 {
61 // Create the domain
62 Domain_pt = new RectangleWithHoleDomain(cylinder_pt, length);
63
64 // Initialise the node counter
65 unsigned long node_count = 0;
66
67 // Vectors used to get data from domains
68 Vector<double> s(2), r(2);
69
70 // Setup temporary storage for the Node
71 Vector<Node*> Tmp_node_pt;
72
73 // Now blindly loop over the macro elements and associate and finite
74 // element with each
75 unsigned nmacro_element = Domain_pt->nmacro_element();
76 for (unsigned e = 0; e < nmacro_element; e++)
77 {
78 // Create the FiniteElement and add to the Element_pt Vector
79 Element_pt.push_back(new ELEMENT);
80
81 // Read out the number of linear points in the element
82 unsigned np = dynamic_cast<ELEMENT*>(finite_element_pt(e))->nnode_1d();
83
84 // Loop over nodes in the column
85 for (unsigned l1 = 0; l1 < np; l1++)
86 {
87 // Loop over the nodes in the row
88 for (unsigned l2 = 0; l2 < np; l2++)
89 {
90 // Allocate the memory for the node
91 Tmp_node_pt.push_back(finite_element_pt(e)->construct_node(
92 l1 * np + l2, time_stepper_pt));
93
94 // Read out the position of the node from the macro element
95 s[0] = -1.0 + 2.0 * (double)l2 / (double)(np - 1);
96 s[1] = -1.0 + 2.0 * (double)l1 / (double)(np - 1);
98
99 // Set the position of the node
100 Tmp_node_pt[node_count]->x(0) = r[0];
101 Tmp_node_pt[node_count]->x(1) = r[1];
102
103 // Increment the node number
104 node_count++;
105 }
106 }
107 } // End of loop over macro elements
108
109
110 // Now the elements have been created, but there will be nodes in
111 // common, need to loop over the common edges and sort it, by reassigning
112 // pointers and the deleting excess nodes
113
114 // Read out the number of linear points in the element
115 unsigned np = dynamic_cast<ELEMENT*>(finite_element_pt(0))->nnode_1d();
116
117 // Edge between Elements 0 and 1
118 for (unsigned n = 0; n < np; n++)
119 {
120 // Set the nodes in element 1 to be the same as in element 0
121 finite_element_pt(1)->node_pt(n * np) =
122 finite_element_pt(0)->node_pt((np - 1) * np + np - 1 - n);
123
124 // Remove the nodes in element 1 from the temporary node list
125 delete Tmp_node_pt[np * np + n * np];
126 Tmp_node_pt[np * np + n * np] = 0;
127 }
128
129 // Edge between Elements 0 and 3
130 for (unsigned n = 0; n < np; n++)
131 {
132 // Set the nodes in element 3 to be the same as in element 0
133 finite_element_pt(3)->node_pt(n * np) =
135
136 // Remove the nodes in element 3 from the temporary node list
137 delete Tmp_node_pt[3 * np * np + n * np];
138 Tmp_node_pt[3 * np * np + n * np] = 0;
139 }
140
141 // Edge between Element 1 and 2
142 for (unsigned n = 0; n < np; n++)
143 {
144 // Set the nodes in element 2 to be the same as in element 1
145 finite_element_pt(2)->node_pt(np * (np - 1) + n) =
146 finite_element_pt(1)->node_pt(np * n + np - 1);
147
148 // Remove the nodes in element 2 from the temporary node list
149 delete Tmp_node_pt[2 * np * np + np * (np - 1) + n];
150 Tmp_node_pt[2 * np * np + np * (np - 1) + n] = 0;
151 }
152
153
154 // Edge between Element 3 and 2
155 for (unsigned n = 0; n < np; n++)
156 {
157 // Set the nodes in element 2 to be the same as in element 3
159 finite_element_pt(3)->node_pt(np * (np - n - 1) + np - 1);
160
161 // Remove the nodes in element 2 from the temporary node list
162 delete Tmp_node_pt[2 * np * np + n];
163 Tmp_node_pt[2 * np * np + n] = 0;
164 }
165
166
167 // Now set the actual true nodes
168 for (unsigned long n = 0; n < node_count; n++)
169 {
170 if (Tmp_node_pt[n] != 0)
171 {
172 Node_pt.push_back(Tmp_node_pt[n]);
173 }
174 }
175
176 // Finally set the nodes on the boundaries
177 set_nboundary(5);
178
179 for (unsigned n = 0; n < np; n++)
180 {
181 // Left hand side
182 Node* nod_pt = finite_element_pt(0)->node_pt(n * np);
184 add_boundary_node(3, nod_pt);
185
186 // Right hand side
187 nod_pt = finite_element_pt(2)->node_pt(n * np + np - 1);
189 add_boundary_node(1, nod_pt);
190
191 // First part of lower boundary
192 nod_pt = finite_element_pt(3)->node_pt(n);
194 add_boundary_node(0, nod_pt);
195
196 // First part of upper boundary
197 nod_pt = finite_element_pt(1)->node_pt(np * (np - 1) + n);
199 add_boundary_node(2, nod_pt);
200
201 // First part of hole boundary
202 nod_pt = finite_element_pt(3)->node_pt(np * (np - 1) + n);
204 add_boundary_node(4, nod_pt);
205 }
206
207 for (unsigned n = 1; n < np; n++)
208 {
209 // Next part of hole
210 Node* nod_pt = finite_element_pt(2)->node_pt(n * np);
212 add_boundary_node(4, nod_pt);
213 }
214
215 for (unsigned n = 1; n < np; n++)
216 {
217 // Next part of hole
218 Node* nod_pt = finite_element_pt(1)->node_pt(np - n - 1);
220 add_boundary_node(4, nod_pt);
221 }
222
223 for (unsigned n = 1; n < np - 1; n++)
224 {
225 // Final part of hole
226 Node* nod_pt =
227 finite_element_pt(0)->node_pt(np * (np - n - 1) + np - 1);
229 add_boundary_node(4, nod_pt);
230 }
231 }
232
233 /// Access function to the domain
235 {
236 return Domain_pt;
237 }
238
239 protected:
240 /// Pointer to the domain
242 };
243
244
245 /// ///////////////////////////////////////////////////////////////////////
246 /// ///////////////////////////////////////////////////////////////////////
247 /// ///////////////////////////////////////////////////////////////////////
248
249
250 //===================================================================
251 /// Refineable version of RectangleWithHoleMesh. For some reason
252 /// this needs on uniform refinement to work...
253 //===================================================================
254 template<class ELEMENT>
256 public RefineableQuadMesh<ELEMENT>
257 {
258 public:
259 /// Constructor. Pass pointer to geometric object that
260 /// represents the cylinder, the length and height of the domain.
261 /// The GeomObject must be parametrised such that
262 /// \f$\zeta \in [0,2\pi]\f$ sweeps around the circumference
263 /// in anticlockwise direction. Timestepper defaults to Steady
264 /// default timestepper.
266 GeomObject* cylinder_pt,
267 const double& length,
268 TimeStepper* time_stepper_pt = &Mesh::Default_TimeStepper)
269 : RectangleWithHoleMesh<ELEMENT>(cylinder_pt, length, time_stepper_pt)
270 {
271 // Nodal positions etc. were created in constructor for
272 // Cylinder...<...>. Need to setup adaptive information.
273
274 // Loop over all elements and set macro element pointer
275 for (unsigned e = 0; e < 4; e++)
276 {
277 dynamic_cast<ELEMENT*>(this->element_pt(e))
278 ->set_macro_elem_pt(this->Domain_pt->macro_element_pt(e));
279 }
280
281 // Setup boundary element lookup schemes
283
284 // Nodal positions etc. were created in constructor for
285 // RectangularMesh<...>. Only need to setup quadtree forest
286 this->setup_quadtree_forest();
287 }
288
289 /// Destructor: Empty
291 };
292
293
294} // namespace oomph
295
296
297#endif
e
Definition: cfortran.h:571
static char t char * s
Definition: cfortran.h:568
unsigned nmacro_element()
Number of macro elements in domain.
Definition: domain.h:123
MacroElement * macro_element_pt(const unsigned &i)
Access to i-th macro element.
Definition: domain.h:116
virtual Node * construct_node(const unsigned &n)
Construct the local node n and return a pointer to the newly created node object.
Definition: elements.h:2509
Node *& node_pt(const unsigned &n)
Return a pointer to the local node n.
Definition: elements.h:2175
/////////////////////////////////////////////////////////////////////
Definition: geom_objects.h:101
void macro_map(const Vector< double > &s, Vector< double > &r)
The mapping from local to global coordinates at the current time : r(s)
A general mesh class.
Definition: mesh.h:67
void add_boundary_node(const unsigned &b, Node *const &node_pt)
Add a (pointer to) a node to the b-th boundary.
Definition: mesh.cc:243
Vector< Node * > Node_pt
Vector of pointers to nodes.
Definition: mesh.h:183
static Steady< 0 > Default_TimeStepper
Default Steady Timestepper, to be used in default arguments to Mesh constructors.
Definition: mesh.h:75
FiniteElement * finite_element_pt(const unsigned &e) const
Upcast (downcast?) to FiniteElement (needed to access FiniteElement member functions).
Definition: mesh.h:473
void set_nboundary(const unsigned &nbound)
Set the number of boundaries in the mesh.
Definition: mesh.h:505
virtual void setup_boundary_element_info()
Interface for function that is used to setup the boundary information (Empty virtual function – imple...
Definition: mesh.h:275
void convert_to_boundary_node(Node *&node_pt, const Vector< FiniteElement * > &finite_element_pt)
A function that upgrades an ordinary node to a boundary node We shouldn't ever really use this,...
Definition: mesh.cc:2590
const Vector< GeneralisedElement * > & element_pt() const
Return reference to the Vector of elements.
Definition: mesh.h:460
Vector< GeneralisedElement * > Element_pt
Vector of pointers to generalised elements.
Definition: mesh.h:186
Nodes are derived from Data, but, in addition, have a definite (Eulerian) position in a space of a gi...
Definition: nodes.h:906
Rectangular domain with circular whole.
Domain-based mesh for rectangular mesh with circular hole.
RectangleWithHoleDomain * Domain_pt
Pointer to the domain.
RectangleWithHoleMesh(GeomObject *cylinder_pt, const double &length, TimeStepper *time_stepper_pt=&Mesh::Default_TimeStepper)
Constructor: Pass pointer to geometric object that represents the cylinder, the length and height of ...
RectangleWithHoleDomain * domain_pt()
Access function to the domain.
Intermediate mesh class that implements the mesh adaptation functions specified in the TreeBasedRefin...
void setup_quadtree_forest()
Set up QuadTreeForest. Wipes any existing tree structure below the minimum refinement level and regar...
/////////////////////////////////////////////////////////////////////// /////////////////////////////...
RefineableRectangleWithHoleMesh(GeomObject *cylinder_pt, const double &length, TimeStepper *time_stepper_pt=&Mesh::Default_TimeStepper)
Constructor. Pass pointer to geometric object that represents the cylinder, the length and height of ...
////////////////////////////////////////////////////////////////////// //////////////////////////////...
Definition: timesteppers.h:231
//////////////////////////////////////////////////////////////////// ////////////////////////////////...