spectral_poisson_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: 1182 $
7 //LIC//
8 //LIC// $LastChangedDate: 2016-05-20 16:50:20 +0100 (Fri, 20 May 2016) $
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 Spectral Poisson elements
31 #ifndef OOMPH_SPECTRAL_POISSON_ELEMENTS_HEADER
32 #define OOMPH_SPECTRAL_POISSON_ELEMENTS_HEADER
33 
34 // Config header generated by autoconfig
35 #ifdef HAVE_CONFIG_H
36  #include <oomph-lib-config.h>
37 #endif
38 
39 //OOMPH-LIB headers
40 #include "poisson_elements.h"
41 #include "../generic/Qspectral_elements.h"
42 
43 namespace oomph
44 {
45 
46 //======================================================================
47 /// QSpectralPoissonElement elements are linear/quadrilateral/brick-shaped
48 /// Poisson elements with isoparametric spectral interpolation for the
49 /// function. Note that the implementation is PoissonEquations<DIM> does
50 /// not use sum factorisation for the evaluation of the residuals and is,
51 /// therefore, not optimal for higher dimensions.
52 //======================================================================
53 template <unsigned DIM, unsigned NNODE_1D>
54  class QSpectralPoissonElement : public virtual QSpectralElement<DIM,NNODE_1D>,
55  public virtual PoissonEquations<DIM>
56 {
57 
58  private:
59 
60  /// \short Static array of ints to hold number of variables at
61  /// nodes: Initial_Nvalue[n]
62  static const unsigned Initial_Nvalue;
63 
64  public:
65 
66 
67  ///\short Constructor: Call constructors for QSpectralElement and
68  /// Poisson equations
70  PoissonEquations<DIM>()
71  {}
72 
73  /// Broken copy constructor
75  {
76  BrokenCopy::broken_copy("QSpectralPoissonElement");
77  }
78 
79  /// Broken assignment operator
80 //Commented out broken assignment operator because this can lead to a conflict warning
81 //when used in the virtual inheritence hierarchy. Essentially the compiler doesn't
82 //realise that two separate implementations of the broken function are the same and so,
83 //quite rightly, it shouts.
84  /*void operator=(const QSpectralPoissonElement<DIM,NNODE_1D>&)
85  {
86  BrokenCopy::broken_assign("QSpectralPoissonElement");
87  }*/
88 
89  /// \short Required # of `values' (pinned or dofs)
90  /// at node n
91  inline unsigned required_nvalue(const unsigned &n) const
92  {return Initial_Nvalue;}
93 
94  /// \short Output function:
95  /// x,y,u or x,y,z,u
96  void output(std::ostream &outfile)
98 
99  /// \short 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  {PoissonEquations<DIM>::output(outfile,n_plot);}
103 
104 
105  /// \short C-style output function:
106  /// x,y,u or x,y,z,u
107  void output(FILE* file_pt)
109 
110 
111  /// \short C-style output function:
112  /// x,y,u or x,y,z,u at n_plot^DIM plot points
113  void output(FILE* file_pt, const unsigned &n_plot)
114  {PoissonEquations<DIM>::output(file_pt,n_plot);}
115 
116 
117  /// \short Output function for an exact solution:
118  /// x,y,u_exact or x,y,z,u_exact at n_plot^DIM plot points
119  void output_fct(std::ostream &outfile, const unsigned &n_plot,
121  {PoissonEquations<DIM>::output_fct(outfile,n_plot,exact_soln_pt);}
122 
123 
124 
125  /// \short Output function for a time-dependent exact solution.
126  /// x,y,u_exact or x,y,z,u_exact at n_plot^DIM plot points
127  /// (Calls the steady version)
128  void output_fct(std::ostream &outfile, const unsigned &n_plot,
129  const double& time,
131  {PoissonEquations<DIM>::output_fct(outfile,n_plot,time,exact_soln_pt);}
132 
133 
134 protected:
135 
136 /// Shape, test functions & derivs. w.r.t. to global coords. Return Jacobian.
137  inline double dshape_and_dtest_eulerian_poisson(
138  const Vector<double> &s, Shape &psi, DShape &dpsidx,
139  Shape &test, DShape &dtestdx) const;
140 
141 
142  /// \short Shape, test functions & derivs. w.r.t. to global coords. at
143  /// integration point ipt. Return Jacobian.
145  const unsigned& ipt, Shape &psi, DShape &dpsidx,
146  Shape &test, DShape &dtestdx)
147  const;
148 
149  /// \short Shape/test functions and derivs w.r.t. to global coords at
150  /// integration point ipt; return Jacobian of mapping (J). Also compute
151  /// derivatives of dpsidx, dtestdx and J w.r.t. nodal coordinates.
153  const unsigned &ipt,
154  Shape &psi,
155  DShape &dpsidx,
156  RankFourTensor<double> &d_dpsidx_dX,
157  Shape &test,
158  DShape &dtestdx,
159  RankFourTensor<double> &d_dtestdx_dX,
160  DenseMatrix<double> &djacobian_dX) const;
161 
162 };
163 
164 
165 //Inline functions:
166 
167 
168 //======================================================================
169 /// Define the shape functions and test functions and derivatives
170 /// w.r.t. global coordinates and return Jacobian of mapping.
171 ///
172 /// Galerkin: Test functions = shape functions
173 //======================================================================
174 template<unsigned DIM, unsigned NNODE_1D>
177  const Vector<double> &s,
178  Shape &psi,
179  DShape &dpsidx,
180  Shape &test,
181  DShape &dtestdx) const
182 {
183  //Call the geometrical shape functions and derivatives
184  double J = this->dshape_eulerian(s,psi,dpsidx);
185 
186  //Loop over the test functions and derivatives and set them equal to the
187  //shape functions
188  unsigned nnod=this->nnode();
189  for(unsigned i=0;i<nnod;i++)
190  {
191  test[i] = psi[i];
192  for(unsigned j=0;j<DIM;j++)
193  {
194  dtestdx(i,j) = dpsidx(i,j);
195  }
196  }
197 
198  //Return the jacobian
199  return J;
200 }
201 
202 //======================================================================
203 /// Define the shape functions and test functions and derivatives
204 /// w.r.t. global coordinates and return Jacobian of mapping.
205 ///
206 /// Galerkin: Test functions = shape functions
207 //======================================================================
208 template<unsigned DIM, unsigned NNODE_1D>
211  const unsigned &ipt,
212  Shape &psi,
213  DShape &dpsidx,
214  Shape &test,
215  DShape &dtestdx) const
216 {
217  //Call the geometrical shape functions and derivatives
218  double J = this->dshape_eulerian_at_knot(ipt,psi,dpsidx);
219 
220  //Set the pointers of the test functions
221  test = psi;
222  dtestdx = dpsidx;
223 
224  //Return the jacobian
225  return J;
226 }
227 
228 //======================================================================
229 /// Define the shape functions (psi) and test functions (test) and
230 /// their derivatives w.r.t. global coordinates (dpsidx and dtestdx)
231 /// and return Jacobian of mapping (J). Additionally compute the
232 /// derivatives of dpsidx, dtestdx and J w.r.t. nodal coordinates.
233 ///
234 /// Galerkin: Test functions = shape functions
235 //======================================================================
236 template<unsigned DIM, unsigned NNODE_1D>
239  const unsigned &ipt,
240  Shape &psi,
241  DShape &dpsidx,
242  RankFourTensor<double> &d_dpsidx_dX,
243  Shape &test,
244  DShape &dtestdx,
245  RankFourTensor<double> &d_dtestdx_dX,
246  DenseMatrix<double> &djacobian_dX) const
247  {
248  // Call the geometrical shape functions and derivatives
249  const double J = this->dshape_eulerian_at_knot(ipt,psi,dpsidx,
250  djacobian_dX,d_dpsidx_dX);
251 
252  // Set the pointers of the test functions
253  test = psi;
254  dtestdx = dpsidx;
255  d_dtestdx_dX = d_dpsidx_dX;
256 
257  //Return the jacobian
258  return J;
259 }
260 
261 ////////////////////////////////////////////////////////////////////////
262 ////////////////////////////////////////////////////////////////////////
263 ////////////////////////////////////////////////////////////////////////
264 
265 
266 
267 //=======================================================================
268 /// Face geometry for the QSpectralPoissonElement elements: The spatial
269 /// dimension of the face elements is one lower than that of the
270 /// bulk element but they have the same number of points
271 /// along their 1D edges.
272 //=======================================================================
273 template<unsigned DIM, unsigned NNODE_1D>
274 class FaceGeometry<QSpectralPoissonElement<DIM,NNODE_1D> >:
275  public virtual QSpectralElement<DIM-1,NNODE_1D>
276 {
277 
278  public:
279 
280  /// \short Constructor: Call the constructor for the
281  /// appropriate lower-dimensional QElement
282  FaceGeometry() : QSpectralElement<DIM-1,NNODE_1D>() {}
283 
284 };
285 
286 
287 //=======================================================================
288 /// Face geometry for the 1D QPoissonElement elements: Point elements
289 //=======================================================================
290 template<unsigned NNODE_1D>
292  public virtual PointElement
293 {
294 
295  public:
296 
297  /// \short Constructor: Call the constructor for the
298  /// appropriate lower-dimensional QElement
300 
301 };
302 
303 }
304 
305 #endif
void output(FILE *file_pt)
C-style output function: x,y,u or x,y,z,u.
void broken_copy(const std::string &class_name)
Issue error message and terminate execution.
cstr elem_len * i
Definition: cfortran.h:607
FaceGeometry()
Constructor: Call the constructor for the appropriate lower-dimensional QElement. ...
double dshape_and_dtest_eulerian_poisson(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.
A Rank 4 Tensor class.
Definition: matrices.h:1625
void output(std::ostream &outfile)
Output function: x,y,u or x,y,z,u.
void output(std::ostream &outfile)
Output with default number of plot points.
QSpectralPoissonElement()
Constructor: Call constructors for QSpectralElement and Poisson equations.
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.
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.
static char t char * s
Definition: cfortran.h:572
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 or x,y,z,u_exact at n_plot^DIM plot points...
unsigned required_nvalue(const unsigned &n) const
Broken assignment operator.
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, const double &time, FiniteElement::UnsteadyExactSolutionFctPt exact_soln_pt)
Output function for a time-dependent exact solution. x,y,u_exact or x,y,z,u_exact at n_plot^DIM plot ...
FaceGeometry()
Constructor: Call the constructor for the appropriate lower-dimensional QElement. ...
QSpectralPoissonElement(const QSpectralPoissonElement< DIM, NNODE_1D > &dummy)
Broken copy constructor.
double dshape_and_dtest_eulerian_at_knot_poisson(const unsigned &ipt, Shape &psi, DShape &dpsidx, Shape &test, DShape &dtestdx) const
Shape, test functions & derivs. w.r.t. to global coords. at integration point ipt. Return Jacobian.
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
static const unsigned Initial_Nvalue
Static array of ints to hold number of variables at nodes: Initial_Nvalue[n].
void output_fct(std::ostream &outfile, const unsigned &n_plot, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt)
Output exact soln: x,y,u_exact or x,y,z,u_exact at n_plot^DIM plot points.