rectangle_with_hole_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-2025 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
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"
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
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
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 // Now the elements have been created, but there will be nodes in
110 // common, need to loop over the common edges and sort it, by reassigning
111 // pointers and the deleting excess nodes
112
113 // Read out the number of linear points in the element
114 unsigned np = dynamic_cast<ELEMENT*>(finite_element_pt(0))->nnode_1d();
115
116 // Edge between Elements 0 and 1
117 for (unsigned n = 0; n < np; n++)
118 {
119 // Set the nodes in element 1 to be the same as in element 0
121 finite_element_pt(0)->node_pt((np - 1) * np + np - 1 - n);
122
123 // Remove the nodes in element 1 from the temporary node list
124 delete Tmp_node_pt[np * np + n * np];
125 Tmp_node_pt[np * np + n * np] = 0;
126 }
127
128 // Edge between Elements 0 and 3
129 for (unsigned n = 0; n < np; n++)
130 {
131 // Set the nodes in element 3 to be the same as in element 0
134
135 // Remove the nodes in element 3 from the temporary node list
136 delete Tmp_node_pt[3 * np * np + n * np];
137 Tmp_node_pt[3 * np * np + n * np] = 0;
138 }
139
140 // Edge between Element 1 and 2
141 for (unsigned n = 0; n < np; n++)
142 {
143 // Set the nodes in element 2 to be the same as in element 1
144 finite_element_pt(2)->node_pt(np * (np - 1) + n) =
145 finite_element_pt(1)->node_pt(np * n + np - 1);
146
147 // Remove the nodes in element 2 from the temporary node list
148 delete Tmp_node_pt[2 * np * np + np * (np - 1) + n];
149 Tmp_node_pt[2 * np * np + np * (np - 1) + n] = 0;
150 }
151
152 // Edge between Element 3 and 2
153 for (unsigned n = 0; n < np; n++)
154 {
155 // Set the nodes in element 2 to be the same as in element 3
157 finite_element_pt(3)->node_pt(np * (np - n - 1) + np - 1);
158
159 // Remove the nodes in element 2 from the temporary node list
160 delete Tmp_node_pt[2 * np * np + n];
161 Tmp_node_pt[2 * np * np + n] = 0;
162 }
163
164 // Now set the actual true nodes
165 for (unsigned long n = 0; n < node_count; n++)
166 {
167 if (Tmp_node_pt[n] != 0)
168 {
169 Node_pt.push_back(Tmp_node_pt[n]);
170 }
171 }
172
173 // Finally set the nodes on the boundaries
174 set_nboundary(5);
175
176 for (unsigned n = 0; n < np; n++)
177 {
178 // Left hand side
182
183 // Right hand side
184 nod_pt = finite_element_pt(2)->node_pt(n * np + np - 1);
187
188 // First part of lower boundary
192
193 // First part of upper boundary
194 nod_pt = finite_element_pt(1)->node_pt(np * (np - 1) + n);
197
198 // First part of hole boundary
199 nod_pt = finite_element_pt(3)->node_pt(np * (np - 1) + n);
202 }
203
204 for (unsigned n = 1; n < np; n++)
205 {
206 // Next part of hole
210 }
211
212 for (unsigned n = 1; n < np; n++)
213 {
214 // Next part of hole
215 Node* nod_pt = finite_element_pt(1)->node_pt(np - n - 1);
218 }
219
220 for (unsigned n = 1; n < np - 1; n++)
221 {
222 // Final part of hole
223 Node* nod_pt =
224 finite_element_pt(0)->node_pt(np * (np - n - 1) + np - 1);
227 }
228 }
229
230 /// Access function to the domain
232 {
233 return Domain_pt;
234 }
235
236 protected:
237 /// Pointer to the domain
239 };
240
241
242 //////////////////////////////////////////////////////////////////////////
243 //////////////////////////////////////////////////////////////////////////
244 //////////////////////////////////////////////////////////////////////////
245
246 //===================================================================
247 /// Refineable version of RectangleWithHoleMesh. For some reason
248 /// this needs on uniform refinement to work...
249 //===================================================================
250 template<class ELEMENT>
252 public RefineableQuadMesh<ELEMENT>
253 {
254 public:
255 /// Constructor. Pass pointer to geometric object that
256 /// represents the cylinder, the length and height of the domain.
257 /// The GeomObject must be parametrised such that
258 /// \f$\zeta \in [0,2\pi]\f$ sweeps around the circumference
259 /// in anticlockwise direction. Timestepper defaults to Steady
260 /// default timestepper.
262 GeomObject* cylinder_pt,
263 const double& length,
264 TimeStepper* time_stepper_pt = &Mesh::Default_TimeStepper)
265 : RectangleWithHoleMesh<ELEMENT>(cylinder_pt, length, time_stepper_pt)
266 {
267 // Nodal positions etc. were created in constructor for
268 // Cylinder...<...>. Need to setup adaptive information.
269
270 // Loop over all elements and set macro element pointer
271 for (unsigned e = 0; e < 4; e++)
272 {
273 dynamic_cast<ELEMENT*>(this->element_pt(e))
274 ->set_macro_elem_pt(this->Domain_pt->macro_element_pt(e));
275 }
276
277 // Setup boundary element lookup schemes
279
280 // Nodal positions etc. were created in constructor for
281 // RectangularMesh<...>. Only need to setup quadtree forest
282 this->setup_quadtree_forest();
283 }
284
285 /// Destructor: Empty
287 };
288
289} // namespace oomph
290
292#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:2513
Node *& node_pt(const unsigned &n)
Return a pointer to the local node n.
Definition elements.h:2179
A geometric object is an object that provides a parametrised description of its shape via the functio...
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:477
void set_nboundary(const unsigned &nbound)
Set the number of boundaries in the mesh.
Definition mesh.h:509
virtual void setup_boundary_element_info()
Interface for function that is used to setup the boundary information (Empty virtual function – imple...
Definition mesh.h:279
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:464
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...
Refineable version of RectangleWithHoleMesh. For some reason this needs on uniform refinement to work...
virtual ~RefineableRectangleWithHoleMesh()
Destructor: Empty.
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 ...
TAdvectionDiffusionReactionElement<NREAGENT,DIM,NNODE_1D> elements are isoparametric triangular DIM-d...
Base class for time-stepping schemes. Timestepper provides an approximation of the temporal derivativ...
DRAIG: Change all instances of (SPATIAL_DIM) to (DIM-1).