Tpml_fourier_decomposed_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// Version 1.0; svn revision $LastChangedRevision: 1299 $
7 //LIC//
8 //LIC// $LastChangedDate: 2017-10-06 07:40:18 +0100 (Fri, 06 Oct 2017) $
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 //Header file for TPMLFourierDecomposedHelmholtz elements
31 #ifndef OOMPH_TPML_FOURIER_DECOMPOSED_HELMHOLTZ_ELEMENTS_HEADER
32 #define OOMPH_TPML_FOURIER_DECOMPOSED_HELMHOLTZ_ELEMENTS_HEADER
33 
34 
35 // Config header generated by autoconfig
36 #ifdef HAVE_CONFIG_H
37  #include <oomph-lib-config.h>
38 #endif
39 
40 
41 //OOMPH-LIB headers
42 #include "../generic/nodes.h"
43 #include "../generic/oomph_utilities.h"
44 #include "../generic/Telements.h"
45 #include "../generic/error_estimator.h"
47 
48 namespace oomph
49 {
50 
51 /////////////////////////////////////////////////////////////////////////
52 /////////////////////////////////////////////////////////////////////////
53 // TPMLFourierDecomposedHelmholtzElement
54 ////////////////////////////////////////////////////////////////////////
55 ////////////////////////////////////////////////////////////////////////
56 
57 
58 
59 //======================================================================
60 /// TPMLFourierDecomposedHelmholtzElement<NNODE_1D> elements are
61 /// isoparametric triangular
62 /// PMLFourierDecomposedHelmholtz elements with NNODE_1D nodal
63 /// points along each element edge. Inherits from TElement and
64 /// PMLFourierDecomposedHelmholtzEquations
65 //======================================================================
66 template <unsigned NNODE_1D>
68  public TElement<2,NNODE_1D>,
70  public virtual ElementWithZ2ErrorEstimator
71  {
72 
73  public:
74 
75  ///\short Constructor: Call constructors for TElement and
76  /// PMLFourierDecomposedHelmholtz equations
79  { }
80 
81 
82  /// Broken copy constructor
85  {
86  BrokenCopy::broken_copy("TPMLFourierDecomposedHelmholtzElement");
87  }
88 
89  /// Broken assignment operator
90 //Commented out broken assignment operator because this can lead to a conflict warning
91 //when used in the virtual inheritence hierarchy. Essentially the compiler doesn't
92 //realise that two separate implementations of the broken function are the same and so,
93 //quite rightly, it shouts.
94  /*void operator=(const TPMLFourierDecomposedHelmholtzElement<NNODE_1D>&)
95  {
96  BrokenCopy::broken_assign("TPMLFourierDecomposedHelmholtzElement");
97  }*/
98 
99  /// \short Access function for Nvalue: # of `values' (pinned or dofs)
100  /// at node n (always returns the same value at every node, 2)
101  inline unsigned required_nvalue(const unsigned &n) const
102  {return Initial_Nvalue;}
103 
104  /// \short Output function:
105  /// r,z,u
106  void output(std::ostream &outfile)
107  {
109  }
110 
111  /// \short Output function:
112  /// r,z,u n_plot^2 plot points
113  void output(std::ostream &outfile, const unsigned &n_plot)
114  {
116  }
117 
118 
119  /// \short C-style output function:
120  /// r,z,u or x,y,z,u
121  void output(FILE* file_pt)
122  {
124  }
125 
126 
127  /// \short C-style output function:
128  /// r,z,u at n_plot^2 plot points
129  void output(FILE* file_pt, const unsigned &n_plot)
130  {
132  }
133 
134 
135  /// \short Output function for an exact solution:
136  /// r,z,u_exact
137  void output_fct(std::ostream &outfile, const unsigned &n_plot,
139  {
141  exact_soln_pt);
142  }
143 
144 
145  /// \short Output function for a time-dependent exact solution.
146  /// x,y,u_exact (calls the steady version)
147  void output_fct(std::ostream &outfile, const unsigned &n_plot,
148  const double& time,
150  {
152  outfile,n_plot,time,exact_soln_pt);
153  }
154 
155  protected:
156 
157 /// Shape, test functions & derivs. w.r.t. to global coords. Return Jacobian.
158  inline double
160  const Vector<double> &s,
161  Shape &psi,
162  DShape &dpsidx,
163  Shape &test,
164  DShape &dtestdx) const;
165 
166 
167 /// Shape, test functions & derivs. w.r.t. to global coords. Return Jacobian.
168  inline double
170  const unsigned &ipt,
171  Shape &psi,
172  DShape &dpsidx,
173  Shape &test,
174  DShape &dtestdx)
175  const;
176 
177 
178  /// \short Order of recovery shape functions for Z2 error estimation:
179  /// Same order as shape functions.
180  unsigned nrecovery_order() {return (NNODE_1D-1);}
181 
182  /// Number of 'flux' terms for Z2 error estimation
183  unsigned num_Z2_flux_terms() {return 2*2;}
184 
185  /// \short Get 'flux' for Z2 error recovery: Standard flux from
186  /// UnsteadyHeat equations
188  {
189  Vector<std::complex <double> > complex_flux(2);
190  this->get_flux(s,complex_flux);
191  unsigned count=0;
192  for (unsigned i=0;i<2;i++)
193  {
194  flux[count++]=complex_flux[i].real();
195  flux[count++]=complex_flux[i].imag();
196  }
197  }
198 
199  /// \short Number of vertex nodes in the element
200  unsigned nvertex_node() const
202 
203  /// \short Pointer to the j-th vertex node in the element
204  Node* vertex_node_pt(const unsigned& j) const
206 
207  private:
208 
209  /// Static unsigned that holds the (same) number of variables at every node
210  static const unsigned Initial_Nvalue;
211 
212 
213 };
214 
215 
216 
217 
218 //Inline functions:
219 
220 
221 //======================================================================
222 /// Define the shape functions and test functions and derivatives
223 /// w.r.t. global coordinates and return Jacobian of mapping.
224 ///
225 /// Galerkin: Test functions = shape functions
226 //======================================================================
227 template<unsigned NNODE_1D>
230  const Vector<double> &s,
231  Shape &psi,
232  DShape &dpsidx,
233  Shape &test,
234  DShape &dtestdx) const
235  {
236  unsigned n_node = this->nnode();
237 
238  //Call the geometrical shape functions and derivatives
239  double J = this->dshape_eulerian(s,psi,dpsidx);
240 
241  //Loop over the test functions and derivatives and set them equal to the
242  //shape functions
243  for(unsigned i=0;i<n_node;i++)
244  {
245  test[i] = psi[i];
246  dtestdx(i,0) = dpsidx(i,0);
247  dtestdx(i,1) = dpsidx(i,1);
248  }
249 
250  //Return the jacobian
251  return J;
252  }
253 
254 
255 
256 //======================================================================
257 /// Define the shape functions and test functions and derivatives
258 /// w.r.t. global coordinates and return Jacobian of mapping.
259 ///
260 /// Galerkin: Test functions = shape functions
261 //======================================================================
262 template<unsigned NNODE_1D>
265  const unsigned &ipt,
266  Shape &psi,
267  DShape &dpsidx,
268  Shape &test,
269  DShape &dtestdx) const
270  {
271 
272  //Call the geometrical shape functions and derivatives
273  double J = this->dshape_eulerian_at_knot(ipt,psi,dpsidx);
274 
275  //Set the pointers of the test functions
276  test = psi;
277  dtestdx = dpsidx;
278 
279  //Return the jacobian
280  return J;
281 
282  }
283 
284 
285 
286 //=======================================================================
287 /// Face geometry for the TPMLFourierDecomposedHelmholtzElement
288 /// elements:
289 /// The spatial dimension of the face elements is one lower than that of the
290 /// bulk element but they have the same number of points
291 /// along their 1D edges.
292 //=======================================================================
293 template<unsigned NNODE_1D>
295  public virtual TElement<1,NNODE_1D>
296  {
297 
298  public:
299 
300  /// \short Constructor: Call the constructor for the
301  /// appropriate lower-dimensional TElement
302  FaceGeometry() : TElement<1,NNODE_1D>() {}
303 
304 };
305 
306 
307 ////////////////////////////////////////////////////////////////////////
308 ////////////////////////////////////////////////////////////////////////
309 ////////////////////////////////////////////////////////////////////////
310 
311 //=======================================================================
312 /// Policy class defining the elements to be used in the actual
313 /// PML layers. It's the corresponding quads.
314 //=======================================================================
315  template<unsigned NNODE_1D>
317  <TPMLFourierDecomposedHelmholtzElement<NNODE_1D> > > :
318  public virtual QPMLFourierDecomposedHelmholtzElement<NNODE_1D>
319 {
320 
321  public:
322 
323  /// \short Constructor: Call the constructor for the
324  /// appropriate QElement
326  {}
327 
328 };
329 
330 ////////////////////////////////////////////////////////////////////////
331 ////////////////////////////////////////////////////////////////////////
332 ////////////////////////////////////////////////////////////////////////
333 
334 //=======================================================================
335 /// Policy class defining the elements to be used in the actual
336 /// PML layers. It's the corresponding quad.
337 //=======================================================================
338  template<unsigned NNODE_1D>
341  public virtual QPMLFourierDecomposedHelmholtzElement<NNODE_1D>
342 {
343 
344  public:
345 
346  /// \short Constructor: Call the constructor for the
347  /// appropriate QElement
350  {}
351 
352 };
353 
354 }
355 
356 #endif
void broken_copy(const std::string &class_name)
Issue error message and terminate execution.
unsigned required_nvalue(const unsigned &n) const
Broken assignment operator.
Base class for finite elements that can compute the quantities that are required for the Z2 error est...
Fourier decomposed Helmholtz upgraded to become projectable.
cstr elem_len * i
Definition: cfortran.h:607
void output(std::ostream &outfile)
Output function: r,z,u.
static const unsigned Initial_Nvalue
Static unsigned that holds the (same) number of variables at every node.
Node * vertex_node_pt(const unsigned &j) const
Pointer to the j-th vertex node in the element.
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) ...
Nodes are derived from Data, but, in addition, have a definite (Eulerian) position in a space of a gi...
Definition: nodes.h:852
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.
double dshape_and_dtest_eulerian_pml_fourier_decomposed_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.
static char t char * s
Definition: cfortran.h:572
TPMLFourierDecomposedHelmholtzElement()
Constructor: Call constructors for TElement and PMLFourierDecomposedHelmholtz equations.
void(* SteadyExactSolutionFctPt)(const Vector< double > &, Vector< double > &)
Function pointer for function that computes vector-valued steady "exact solution" as ...
Definition: elements.h:1720
void output_fct(std::ostream &outfile, const unsigned &n_plot, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt)
void output_fct(std::ostream &outfile, const unsigned &n_plot, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt)
Output function for an exact solution: r,z,u_exact.
void output(FILE *file_pt)
C-style output function: r,z,u or x,y,z,u.
void output(FILE *file_pt, const unsigned &n_plot)
C-style output function: r,z,u at n_plot^2 plot points.
unsigned nvertex_node() const
Number of vertex nodes in the element.
void(* UnsteadyExactSolutionFctPt)(const double &, const Vector< double > &, Vector< double > &)
Function pointer for function that computes Vector-valued time-dependent function as ...
Definition: elements.h:1726
double dshape_and_dtest_eulerian_at_knot_pml_fourier_decomposed_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.
TPMLFourierDecomposedHelmholtzElement(const TPMLFourierDecomposedHelmholtzElement< NNODE_1D > &dummy)
Broken copy constructor.
void output(std::ostream &outfile)
Output with default number of plot points.
FaceGeometry()
Constructor: Call the constructor for the appropriate lower-dimensional TElement. ...
void output(std::ostream &outfile, const unsigned &n_plot)
Output function: r,z,u n_plot^2 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.
unsigned nrecovery_order()
Order of recovery shape functions for Z2 error estimation: Same order as shape functions.