Tpml_helmholtz_elements.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// Header file for THelmholtz elements
27#ifndef OOMPH_TPML_HELMHOLTZ_ELEMENTS_HEADER
28#define OOMPH_TPML_HELMHOLTZ_ELEMENTS_HEADER
29
30
31// Config header generated by autoconfig
32#ifdef HAVE_CONFIG_H
33#include <oomph-lib-config.h>
34#endif
35
36
37// OOMPH-LIB headers
38#include "../generic/nodes.h"
39#include "../generic/oomph_utilities.h"
40#include "../generic/Telements.h"
41#include "../generic/error_estimator.h"
43
44namespace oomph
45{
46 /// //////////////////////////////////////////////////////////////////////
47 /// //////////////////////////////////////////////////////////////////////
48 // TPMLHelmholtzElement
49 /// /////////////////////////////////////////////////////////////////////
50 /// /////////////////////////////////////////////////////////////////////
51
52
53 //======================================================================
54 /// TPMLHelmholtzElement<DIM,NNODE_1D> elements are
55 /// isoparametric triangular DIM-dimensional PMLHelmholtz elements
56 /// with NNODE_1D nodal points along each element edge. Inherits from
57 /// TElement and PMLHelmholtzEquations
58 //======================================================================
59 template<unsigned DIM, unsigned NNODE_1D>
60 class TPMLHelmholtzElement : public TElement<DIM, NNODE_1D>,
61 public PMLHelmholtzEquations<DIM>,
62 public virtual ElementWithZ2ErrorEstimator
63 {
64 public:
65 /// Constructor: Call constructors for TElement and
66 /// PMLHelmholtz equations
68 : TElement<DIM, NNODE_1D>(), PMLHelmholtzEquations<DIM>()
69 {
70 }
71
72
73 /// Broken copy constructor
75 delete;
76
77 /// Broken assignment operator
78 // Commented out broken assignment operator because this can lead to a
79 // conflict warning when used in the virtual inheritence hierarchy.
80 // Essentially the compiler doesn't realise that two separate
81 // implementations of the broken function are the same and so, quite
82 // rightly, it shouts.
83 /*void operator=(const TPMLHelmholtzElement<DIM,NNODE_1D>&) = delete;*/
84
85 /// Access function for Nvalue: # of `values' (pinned or dofs)
86 /// at node n (always returns the same value at every node, 1)
87 inline unsigned required_nvalue(const unsigned& n) const
88 {
89 return Initial_Nvalue;
90 }
91
92 /// Output function:
93 /// x,y,u or x,y,z,u
94 void output(std::ostream& outfile)
95 {
97 }
98
99 /// Output function:
100 /// x,y,u or x,y,z,u at n_plot^DIM plot points
101 void output(std::ostream& outfile, const unsigned& n_plot)
102 {
103 PMLHelmholtzEquations<DIM>::output(outfile, n_plot);
104 }
105
106
107 /// C-style output function:
108 /// x,y,u or x,y,z,u
109 void output(FILE* file_pt)
110 {
112 }
113
114
115 /// C-style output function:
116 /// x,y,u or x,y,z,u at n_plot^DIM plot points
117 void output(FILE* file_pt, const unsigned& n_plot)
118 {
119 PMLHelmholtzEquations<DIM>::output(file_pt, n_plot);
120 }
121
122
123 /// Output function for an exact solution:
124 /// x,y,u_exact
125 void output_fct(std::ostream& outfile,
126 const unsigned& n_plot,
128 {
129 PMLHelmholtzEquations<DIM>::output_fct(outfile, n_plot, exact_soln_pt);
130 }
131
132
133 /// Output function for a time-dependent exact solution.
134 /// x,y,u_exact (calls the steady version)
135 void output_fct(std::ostream& outfile,
136 const unsigned& n_plot,
137 const double& time,
139 {
141 outfile, n_plot, time, exact_soln_pt);
142 }
143
144 protected:
145 /// Shape, test functions & derivs. w.r.t. to global coords. Return
146 /// Jacobian.
148 Shape& psi,
149 DShape& dpsidx,
150 Shape& test,
151 DShape& dtestdx) const;
152
153
154 /// Shape, test functions & derivs. w.r.t. to global coords. Return
155 /// Jacobian.
157 const unsigned& ipt,
158 Shape& psi,
159 DShape& dpsidx,
160 Shape& test,
161 DShape& dtestdx) const;
162
163
164 /// Order of recovery shape functions for Z2 error estimation:
165 /// Same order as shape functions.
167 {
168 return (NNODE_1D - 1);
169 }
170
171 /// Number of 'flux' terms for Z2 error estimation
173 {
174 return 2 * DIM;
175 }
176
177 /// Get 'flux' for Z2 error recovery: Standard flux from
178 /// UnsteadyHeat equations
180 {
181 Vector<std::complex<double>> complex_flux(DIM);
182 this->get_flux(s, complex_flux);
183 unsigned count = 0;
184 for (unsigned i = 0; i < DIM; i++)
185 {
186 flux[count++] = complex_flux[i].real();
187 flux[count++] = complex_flux[i].imag();
188 }
189 }
190
191 /// Number of vertex nodes in the element
192 unsigned nvertex_node() const
193 {
195 }
196
197 /// Pointer to the j-th vertex node in the element
198 Node* vertex_node_pt(const unsigned& j) const
199 {
201 }
202
203 private:
204 /// Static unsigned that holds the (same) number of variables at every node
205 static const unsigned Initial_Nvalue;
206 };
207
208
209 //!! Cleanup - this was not here before!
210 template<unsigned DIM, unsigned NNODE_1D>
212
213 // Inline functions:
214
215
216 //======================================================================
217 /// Define the shape functions and test functions and derivatives
218 /// w.r.t. global coordinates and return Jacobian of mapping.
219 ///
220 /// Galerkin: Test functions = shape functions
221 //======================================================================
222 template<unsigned DIM, unsigned NNODE_1D>
225 Shape& psi,
226 DShape& dpsidx,
227 Shape& test,
228 DShape& dtestdx) const
229 {
230 unsigned n_node = this->nnode();
231
232 // Call the geometrical shape functions and derivatives
233 double J = this->dshape_eulerian(s, psi, dpsidx);
234
235 // Loop over the test functions and derivatives and set them equal to the
236 // shape functions
237 for (unsigned i = 0; i < n_node; i++)
238 {
239 test[i] = psi[i];
240 dtestdx(i, 0) = dpsidx(i, 0);
241 dtestdx(i, 1) = dpsidx(i, 1);
242 }
243
244 // Return the jacobian
245 return J;
246 }
247
248
249 //======================================================================
250 /// Define the shape functions and test functions and derivatives
251 /// w.r.t. global coordinates and return Jacobian of mapping.
252 ///
253 /// Galerkin: Test functions = shape functions
254 //======================================================================
255 template<unsigned DIM, unsigned NNODE_1D>
258 Shape& psi,
259 DShape& dpsidx,
260 Shape& test,
261 DShape& dtestdx) const
262 {
263 // Call the geometrical shape functions and derivatives
264 double J = this->dshape_eulerian_at_knot(ipt, psi, dpsidx);
265
266 // Set the pointers of the test functions
267 test = psi;
268 dtestdx = dpsidx;
269
270 // Return the jacobian
271 return J;
272 }
273
274
275 //=======================================================================
276 /// Face geometry for the TPMLHelmholtzElement elements: The spatial
277 /// dimension of the face elements is one lower than that of the
278 /// bulk element but they have the same number of points
279 /// along their 1D edges.
280 //=======================================================================
281 template<unsigned DIM, unsigned NNODE_1D>
282 class FaceGeometry<TPMLHelmholtzElement<DIM, NNODE_1D>>
283 : public virtual TElement<DIM - 1, NNODE_1D>
284 {
285 public:
286 /// Constructor: Call the constructor for the
287 /// appropriate lower-dimensional TElement
288 FaceGeometry() : TElement<DIM - 1, NNODE_1D>() {}
289 };
290
291 //=======================================================================
292 /// Face geometry for the 1D TPMLHelmholtzElement elements:
293 /// Point elements
294 //=======================================================================
295 template<unsigned NNODE_1D>
297 : public virtual PointElement
298 {
299 public:
300 /// Constructor: Call the constructor for the
301 /// appropriate lower-dimensional TElement
303 };
304
305
306 /// /////////////////////////////////////////////////////////////////////
307 /// /////////////////////////////////////////////////////////////////////
308 /// /////////////////////////////////////////////////////////////////////
309
310 //=======================================================================
311 /// Policy class defining the elements to be used in the actual
312 /// PML layers. It's the corresponding quads.
313 //=======================================================================
314 template<unsigned NNODE_1D>
316 : public virtual QPMLHelmholtzElement<2, NNODE_1D>
317 {
318 public:
319 /// Constructor: Call the constructor for the
320 /// appropriate QElement
322 };
323
324 //=======================================================================
325 /// Face geometry for the TPMLHelmholtzElement elements: The spatial
326 /// dimension of the face elements is one lower than that of the
327 /// bulk element but they have the same number of points
328 /// along their 1D edges.
329 //=======================================================================
330 template<unsigned DIM, unsigned NNODE_1D>
332 : public virtual QElement<DIM - 1, NNODE_1D>
333 {
334 public:
335 /// Constructor: Call the constructor for the
336 /// appropriate lower-dimensional TElement
337 FaceGeometry() : QElement<DIM - 1, NNODE_1D>() {}
338 };
339
340 /// /////////////////////////////////////////////////////////////////////
341 /// /////////////////////////////////////////////////////////////////////
342 /// /////////////////////////////////////////////////////////////////////
343
344 //=======================================================================
345 /// Policy class defining the elements to be used in the actual
346 /// PML layers. It's the corresponding quads.
347 //=======================================================================
348 template<unsigned NNODE_1D>
351 : public virtual QPMLHelmholtzElement<2, NNODE_1D>
352 {
353 public:
354 /// Constructor: Call the constructor for the
355 /// appropriate QElement
357 };
358
359} // namespace oomph
360
361
362#endif
static char t char * s
Definition: cfortran.h:568
cstr elem_len * i
Definition: cfortran.h:603
A Class for the derivatives of shape functions The class design is essentially the same as Shape,...
Definition: shape.h:278
Base class for finite elements that can compute the quantities that are required for the Z2 error est...
FaceGeometry()
Constructor: Call the constructor for the appropriate lower-dimensional TElement.
FaceGeometry()
Constructor: Call the constructor for the appropriate lower-dimensional TElement.
FaceGeometry()
Constructor: Call the constructor for the appropriate lower-dimensional TElement.
//////////////////////////////////////////////////////////////////// ////////////////////////////////...
Definition: elements.h:4998
void(* SteadyExactSolutionFctPt)(const Vector< double > &, Vector< double > &)
Function pointer for function that computes vector-valued steady "exact solution" as .
Definition: elements.h:1759
void(* UnsteadyExactSolutionFctPt)(const double &, const Vector< double > &, Vector< double > &)
Function pointer for function that computes Vector-valued time-dependent function as .
Definition: elements.h:1765
Nodes are derived from Data, but, in addition, have a definite (Eulerian) position in a space of a gi...
Definition: nodes.h:906
A class for all isoparametric elements that solve the Helmholtz equations with pml capabilities....
void output(std::ostream &outfile)
Output with default number of plot points.
void output_fct(std::ostream &outfile, const unsigned &n_plot, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt)
Output exact soln: x,y,u_re_exact,u_im_exact or x,y,z,u_re_exact,u_im_exact at n_plot^DIM plot points...
void get_flux(const Vector< double > &s, Vector< std::complex< double > > &flux) const
Get flux: flux[i] = du/dx_i for real and imag part.
PMLLayerElement()
Constructor: Call the constructor for the appropriate QElement.
General definition of policy class defining the elements to be used in the actual PML layers....
Definition: pml_meshes.h:48
/////////////////////////////////////////////////////////////////////// /////////////////////////////...
Definition: elements.h:3439
///////////////////////////////////////////////////////////////////// ///////////////////////////////...
/////////////////////////////////////////////////////////////////////// /////////////////////////////...
Definition: Qelements.h:459
//////////////////////////////////////////////////////////////////////// ////////////////////////////...
A Class for shape functions. In simple cases, the shape functions have only one index that can be tho...
Definition: shape.h:76
General TElement class.
Definition: Telements.h:1208
//////////////////////////////////////////////////////////////////////
void output_fct(std::ostream &outfile, const unsigned &n_plot, const double &time, FiniteElement::UnsteadyExactSolutionFctPt exact_soln_pt)
Output function for a time-dependent exact solution. x,y,u_exact (calls the steady version)
unsigned nvertex_node() const
Number of vertex nodes in the element.
Node * vertex_node_pt(const unsigned &j) const
Pointer to the j-th vertex node in the element.
void output(FILE *file_pt)
C-style output function: x,y,u or x,y,z,u.
static const unsigned Initial_Nvalue
Static unsigned that holds the (same) number of variables at every node.
double dshape_and_dtest_eulerian_at_knot_helmholtz(const unsigned &ipt, Shape &psi, DShape &dpsidx, Shape &test, DShape &dtestdx) const
Shape, test functions & derivs. w.r.t. to global coords. Return Jacobian.
unsigned nrecovery_order()
Order of recovery shape functions for Z2 error estimation: Same order as shape functions.
unsigned required_nvalue(const unsigned &n) const
Broken assignment operator.
void output(std::ostream &outfile)
Output function: x,y,u or x,y,z,u.
double dshape_and_dtest_eulerian_helmholtz(const Vector< double > &s, Shape &psi, DShape &dpsidx, Shape &test, DShape &dtestdx) const
Shape, test functions & derivs. w.r.t. to global coords. Return Jacobian.
void output_fct(std::ostream &outfile, const unsigned &n_plot, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt)
Output function for an exact solution: x,y,u_exact.
void output(FILE *file_pt, const unsigned &n_plot)
C-style output function: x,y,u or x,y,z,u at n_plot^DIM plot points.
TPMLHelmholtzElement(const TPMLHelmholtzElement< DIM, NNODE_1D > &dummy)=delete
Broken copy constructor.
void output(std::ostream &outfile, const unsigned &n_plot)
Output function: x,y,u or x,y,z,u at n_plot^DIM plot points.
TPMLHelmholtzElement()
Constructor: Call constructors for TElement and PMLHelmholtz equations.
unsigned num_Z2_flux_terms()
Number of 'flux' terms for Z2 error estimation.
void get_Z2_flux(const Vector< double > &s, Vector< double > &flux)
Get 'flux' for Z2 error recovery: Standard flux from UnsteadyHeat equations.
//////////////////////////////////////////////////////////////////// ////////////////////////////////...