Tfourier_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// Copyright (C) 2006-2023 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 TFourierDecomposedHelmholtz elements
27 #ifndef OOMPH_TFOURIER_DECOMPOSED_HELMHOLTZ_ELEMENTS_HEADER
28 #define OOMPH_TFOURIER_DECOMPOSED_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 
44 namespace oomph
45 {
46  /// //////////////////////////////////////////////////////////////////////
47  /// //////////////////////////////////////////////////////////////////////
48  // TFourierDecomposedHelmholtzElement
49  /// /////////////////////////////////////////////////////////////////////
50  /// /////////////////////////////////////////////////////////////////////
51 
52 
53  //======================================================================
54  /// TFourierDecomposedHelmholtzElement<NNODE_1D> elements are
55  /// isoparametric triangular
56  /// FourierDecomposedHelmholtz elements with NNODE_1D nodal points along each
57  /// element edge. Inherits from TElement and
58  /// FourierDecomposedHelmholtzEquations
59  //======================================================================
60  template<unsigned NNODE_1D>
62  : public TElement<2, NNODE_1D>,
64  public virtual ElementWithZ2ErrorEstimator
65  {
66  public:
67  /// Constructor: Call constructors for TElement and
68  /// FourierDecomposedHelmholtz equations
71  {
72  }
73 
74 
75  /// Broken copy constructor
77  const TFourierDecomposedHelmholtzElement<NNODE_1D>& dummy) = delete;
78 
79  /// Broken assignment operator
80  // Commented out broken assignment operator because this can lead to a
81  // conflict warning when used in the virtual inheritence hierarchy.
82  // Essentially the compiler doesn't realise that two separate
83  // implementations of the broken function are the same and so, quite
84  // rightly, it shouts.
85  /*void operator=(const TFourierDecomposedHelmholtzElement<NNODE_1D>&) =
86  * delete;*/
87 
88  /// Access function for Nvalue: # of `values' (pinned or dofs)
89  /// at node n (always returns the same value at every node, 2)
90  inline unsigned required_nvalue(const unsigned& n) const
91  {
92  return Initial_Nvalue;
93  }
94 
95  /// Output function:
96  /// r,z,u
97  void output(std::ostream& outfile)
98  {
100  }
101 
102  /// Output function:
103  /// r,z,u n_plot^2 plot points
104  void output(std::ostream& outfile, const unsigned& n_plot)
105  {
107  }
108 
109 
110  /// C-style output function:
111  /// r,z,u or x,y,z,u
112  void output(FILE* file_pt)
113  {
115  }
116 
117 
118  /// C-style output function:
119  /// r,z,u at n_plot^2 plot points
120  void output(FILE* file_pt, const unsigned& n_plot)
121  {
123  }
124 
125 
126  /// Output function for an exact solution:
127  /// r,z,u_exact
128  void output_fct(std::ostream& outfile,
129  const unsigned& n_plot,
131  {
133  outfile, n_plot, exact_soln_pt);
134  }
135 
136 
137  /// Output function for a time-dependent exact solution.
138  /// x,y,u_exact (calls the steady version)
139  void output_fct(std::ostream& outfile,
140  const unsigned& n_plot,
141  const double& time,
143  {
145  outfile, n_plot, time, exact_soln_pt);
146  }
147 
148  protected:
149  /// Shape, test functions & derivs. w.r.t. to global coords. Return
150  /// Jacobian.
152  const Vector<double>& s,
153  Shape& psi,
154  DShape& dpsidx,
155  Shape& test,
156  DShape& dtestdx) const;
157 
158 
159  /// Shape, test functions & derivs. w.r.t. to global coords. Return
160  /// Jacobian.
162  const unsigned& ipt,
163  Shape& psi,
164  DShape& dpsidx,
165  Shape& test,
166  DShape& dtestdx) const;
167 
168 
169  /// Order of recovery shape functions for Z2 error estimation:
170  /// Same order as shape functions.
171  unsigned nrecovery_order()
172  {
173  return (NNODE_1D - 1);
174  }
175 
176  /// Number of 'flux' terms for Z2 error estimation
177  unsigned num_Z2_flux_terms()
178  {
179  return 2 * 2;
180  }
181 
182  /// Get 'flux' for Z2 error recovery: Standard flux from
183  /// UnsteadyHeat equations
185  {
186  Vector<std::complex<double>> complex_flux(2);
187  this->get_flux(s, complex_flux);
188  unsigned count = 0;
189  for (unsigned i = 0; i < 2; i++)
190  {
191  flux[count++] = complex_flux[i].real();
192  flux[count++] = complex_flux[i].imag();
193  }
194  }
195 
196  /// Number of vertex nodes in the element
197  unsigned nvertex_node() const
198  {
200  }
201 
202  /// Pointer to the j-th vertex node in the element
203  Node* vertex_node_pt(const unsigned& j) const
204  {
206  }
207 
208  private:
209  /// Static unsigned that holds the (same) number of variables at every node
210  static const unsigned Initial_Nvalue;
211  };
212 
213 
214  // Inline functions:
215 
216 
217  //======================================================================
218  /// Define the shape functions and test functions and derivatives
219  /// w.r.t. global coordinates and return Jacobian of mapping.
220  ///
221  /// Galerkin: Test functions = shape functions
222  //======================================================================
223  template<unsigned NNODE_1D>
226  const Vector<double>& s,
227  Shape& psi,
228  DShape& dpsidx,
229  Shape& test,
230  DShape& dtestdx) const
231  {
232  unsigned n_node = this->nnode();
233 
234  // Call the geometrical shape functions and derivatives
235  double J = this->dshape_eulerian(s, psi, dpsidx);
236 
237  // Loop over the test functions and derivatives and set them equal to the
238  // shape functions
239  for (unsigned i = 0; i < n_node; i++)
240  {
241  test[i] = psi[i];
242  dtestdx(i, 0) = dpsidx(i, 0);
243  dtestdx(i, 1) = dpsidx(i, 1);
244  }
245 
246  // Return the jacobian
247  return J;
248  }
249 
250 
251  //======================================================================
252  /// Define the shape functions and test functions and derivatives
253  /// w.r.t. global coordinates and return Jacobian of mapping.
254  ///
255  /// Galerkin: Test functions = shape functions
256  //======================================================================
257  template<unsigned NNODE_1D>
260  const unsigned& ipt,
261  Shape& psi,
262  DShape& dpsidx,
263  Shape& test,
264  DShape& dtestdx) const
265  {
266  // Call the geometrical shape functions and derivatives
267  double J = this->dshape_eulerian_at_knot(ipt, psi, dpsidx);
268 
269  // Set the pointers of the test functions
270  test = psi;
271  dtestdx = dpsidx;
272 
273  // Return the jacobian
274  return J;
275  }
276 
277 
278  //=======================================================================
279  /// Face geometry for the TFourierDecomposedHelmholtzElement elements:
280  /// The spatial dimension of the face elements is one lower than that of the
281  /// bulk element but they have the same number of points
282  /// along their 1D edges.
283  //=======================================================================
284  template<unsigned NNODE_1D>
286  : public virtual TElement<1, NNODE_1D>
287  {
288  public:
289  /// Constructor: Call the constructor for the
290  /// appropriate lower-dimensional TElement
291  FaceGeometry() : TElement<1, NNODE_1D>() {}
292  };
293 
294 
295 } // namespace oomph
296 
297 #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.
//////////////////////////////////////////////////////////////////// ////////////////////////////////...
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
//////////////////////////////////////////////////////////////////// ////////////////////////////////...
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.
void output_fct(std::ostream &outfile, const unsigned &n_plot, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt)
Output exact soln: r,z,u_re_exact,u_im_exact at n_plot^2 plot points.
void output(std::ostream &outfile)
Output with default number of plot points.
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 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)
void output(std::ostream &outfile, const unsigned &n_plot)
Output function: r,z,u n_plot^2 plot points.
unsigned nvertex_node() const
Number of vertex nodes in the element.
TFourierDecomposedHelmholtzElement(const TFourierDecomposedHelmholtzElement< NNODE_1D > &dummy)=delete
Broken copy constructor.
double dshape_and_dtest_eulerian_at_knot_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.
void output(FILE *file_pt)
C-style output function: r,z,u or x,y,z,u.
unsigned nrecovery_order()
Order of recovery shape functions for Z2 error estimation: Same order as shape functions.
void output(FILE *file_pt, const unsigned &n_plot)
C-style output function: r,z,u at n_plot^2 plot points.
static const unsigned Initial_Nvalue
Static unsigned that holds the (same) number of variables at every node.
TFourierDecomposedHelmholtzElement()
Constructor: Call constructors for TElement and FourierDecomposedHelmholtz equations.
void get_Z2_flux(const Vector< double > &s, Vector< double > &flux)
Get 'flux' for Z2 error recovery: Standard flux from UnsteadyHeat equations.
unsigned num_Z2_flux_terms()
Number of 'flux' terms for Z2 error estimation.
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(std::ostream &outfile)
Output function: r,z,u.
unsigned required_nvalue(const unsigned &n) const
Broken assignment operator.
double dshape_and_dtest_eulerian_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.
Node * vertex_node_pt(const unsigned &j) const
Pointer to the j-th vertex node in the element.
//////////////////////////////////////////////////////////////////// ////////////////////////////////...