refineable_advection_diffusion_reaction_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: 1097 $
7 //LIC//
8 //LIC// $LastChangedDate: 2015-12-17 11:53:17 +0000 (Thu, 17 Dec 2015) $
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 elements that solve the advection diffusion equation
31 //and that can be refined.
32 
33 #ifndef OOMPH_REFINEABLE_ADVECTION_DIFFUSION_REACTION_ELEMENTS_HEADER
34 #define OOMPH_REFINEABLE_ADVECTION_DIFFUSION_REACTION_ELEMENTS_HEADER
35 
36 // Config header generated by autoconfig
37 #ifdef HAVE_CONFIG_H
38  #include <oomph-lib-config.h>
39 #endif
40 
41 //oomph-lib headers
42 #include "../generic/refineable_quad_element.h"
43 #include "../generic/refineable_brick_element.h"
44 #include "../generic/error_estimator.h"
46 
47 namespace oomph
48 {
49 
50 //======================================================================
51 /// \short A version of the Advection Diffusion Reaction equations that can be
52 /// used with non-uniform mesh refinement. In essence, the class overloads
53 /// the fill_in_generic_residual_contribution_adv_diff_react()
54 /// function so that contributions
55 /// from hanging nodes (or alternatively in-compatible function values)
56 /// are taken into account.
57 //======================================================================
58 template <unsigned NREAGENT, unsigned DIM>
60  public virtual AdvectionDiffusionReactionEquations<NREAGENT,DIM>,
61  public virtual RefineableElement,
62  public virtual ElementWithZ2ErrorEstimator
63 {
64  public:
65 
66  /// \short Empty Constructor
71  {}
72 
73 
74  /// Broken copy constructor
77  {
78  BrokenCopy::broken_copy("RefineableAdvectionDiffusionReactionEquations");
79  }
80 
81  /// Broken assignment operator
82  void operator=(
84  {
85  BrokenCopy::broken_assign("RefineableAdvectionDiffusionReactionEquations");
86  }
87 
88  /// Number of 'flux' terms for Z2 error estimation
89  unsigned num_Z2_flux_terms() {return NREAGENT*DIM;}
90 
91  /// \short Get 'flux' for Z2 error recovery:
92  /// Standard flux.from AdvectionDiffusionReaction equations
94  {this->get_flux(s,flux);}
95 
96 
97  /// \short Get the function values c in Vector.
98  /// Note: Given the generality of the interface (this function
99  /// is usually called from black-box documentation or interpolation routines),
100  /// the values Vector sets its own size in here.
102  {
103  // Set size of Vector: c
104  values.resize(NREAGENT);
105 
106  //Find number of nodes
107  unsigned n_node = nnode();
108 
109  //Local shape function
110  Shape psi(n_node);
111 
112  //Find values of shape function
113  shape(s,psi);
114 
115  //Loop over the unknowns
116  for(unsigned r=0;r<NREAGENT;r++)
117  {
118  unsigned c_nodal_index = this->c_index_adv_diff_react(r);
119 
120  //Initialise value of c
121  values[r] = 0.0;
122 
123  //Loop over the local nodes and sum
124  for(unsigned l=0;l<n_node;l++)
125  {
126  values[r] += this->nodal_value(l,c_nodal_index)*psi[l];
127  }
128  }
129  }
130 
131  /// \short Get the function values c in Vector.
132  /// Note: Given the generality of the interface (this function
133  /// is usually called from black-box documentation or interpolation routines),
134  /// the values Vector sets its own size in here.
135  void get_interpolated_values(const unsigned& t, const Vector<double>&s,
136  Vector<double>& values)
137  {
138  // Set size of Vector:
139  values.resize(NREAGENT);
140 
141  //Find out how many nodes there are
142  const unsigned n_node = nnode();
143 
144  // Shape functions
145  Shape psi(n_node);
146  shape(s,psi);
147 
148  //Loop over the reagents
149  for(unsigned r=0;r<NREAGENT;r++)
150  {
151  //Find the nodal index at which the unknown is stored
152  unsigned c_nodal_index = this->c_index_adv_diff_react(r);
153 
154  // Initialise
155  values[r]=0.0;
156 
157  //Calculate value
158  for(unsigned l=0;l<n_node;l++)
159  {
160  values[r] += this->nodal_value(t,l,c_nodal_index)*psi[l];
161  }
162  }
163  }
164 
165 
166  /// Further build: Copy all pointers from the father
167  /// element
169  {
171  cast_father_element_pt
172  =
174  this->father_element_pt());
175 
176  //Set the values of the pointers from the father
177  this->Source_fct_pt = cast_father_element_pt->source_fct_pt();
178  this->Wind_fct_pt = cast_father_element_pt->wind_fct_pt();
179  this->Reaction_fct_pt = cast_father_element_pt->reaction_fct_pt();
180  this->Reaction_deriv_fct_pt =
181  cast_father_element_pt->reaction_deriv_fct_pt();
182 
183  this->Diff_pt = cast_father_element_pt->diff_pt();
184  this->Tau_pt = cast_father_element_pt->tau_pt();
185 
186  //Set the value of the ALE flag
187  this->ALE_is_disabled = cast_father_element_pt->ALE_is_disabled;
188  }
189 
190  protected:
191 
192  /// \short Add the element's contribution to the elemental residual vector
193 /// and/or Jacobian matrix
194 /// flag=1: compute both
195 /// flag=0: compute only residual vector
197  Vector<double> &residuals, DenseMatrix<double> &jacobian,
198  DenseMatrix<double> &mass_matrix, unsigned flag);
199 
200 };
201 
202 
203 //======================================================================
204 /// \short Refineable version of QAdvectionDiffusionReactionElement.
205 /// Inherit from the standard QAdvectionDiffusionReactionElement and the
206 /// appropriate refineable geometric element and the refineable equations.
207 //======================================================================
208 template <unsigned NREAGENT, unsigned DIM, unsigned NNODE_1D>
210  public QAdvectionDiffusionReactionElement<NREAGENT,DIM,NNODE_1D>,
211  public virtual RefineableAdvectionDiffusionReactionEquations<NREAGENT,DIM>,
212  public virtual RefineableQElement<DIM>
213 {
214  public:
215 
216  /// \short Empty Constructor:
220  RefineableQElement<DIM>(),
221  QAdvectionDiffusionReactionElement<NREAGENT,DIM,NNODE_1D>()
222  {}
223 
224 
225  /// Broken copy constructor
228  dummy)
229  {
230  BrokenCopy::broken_copy("RefineableQuadAdvectionDiffusionReactionElement");
231  }
232 
233  /// Broken assignment operator
234  void operator=(
236  {
237  BrokenCopy::broken_assign("RefineableQuadAdvectionDiffusionReactionElement");
238  }
239 
240  /// Number of continuously interpolated values: NREAGENT
241  unsigned ncont_interpolated_values() const {return NREAGENT;}
242 
243  /// \short Number of vertex nodes in the element
244  unsigned nvertex_node() const
245  {return
247 
248  /// \short Pointer to the j-th vertex node in the element
249  Node* vertex_node_pt(const unsigned& j) const
250  {return
252 
253  /// Rebuild from sons: empty
254  void rebuild_from_sons(Mesh* &mesh_pt) {}
255 
256  /// \short Order of recovery shape functions for Z2 error estimation:
257  /// Same order as shape functions.
258  unsigned nrecovery_order() {return (NNODE_1D-1);}
259 
260  /// \short Perform additional hanging node procedures for variables
261  /// that are not interpolated by all nodes. Empty.
263 
264 };
265 
266 ////////////////////////////////////////////////////////////////////////
267 ////////////////////////////////////////////////////////////////////////
268 ////////////////////////////////////////////////////////////////////////
269 
270 
271 
272 //=======================================================================
273 /// Face geometry for the RefineableQuadAdvectionDiffusionReactionElement elements: The spatial
274 /// dimension of the face elements is one lower than that of the
275 /// bulk element but they have the same number of points
276 /// along their 1D edges.
277 //=======================================================================
278 template<unsigned NREAGENT,unsigned DIM, unsigned NNODE_1D>
280  RefineableQAdvectionDiffusionReactionElement<NREAGENT,DIM,NNODE_1D> >:
281  public virtual QElement<DIM-1,NNODE_1D>
282 {
283 
284  public:
285 
286  /// \short Constructor: Call the constructor for the
287  /// appropriate lower-dimensional QElement
288  FaceGeometry() : QElement<DIM-1,NNODE_1D>() {}
289 
290 };
291 
292 }
293 
294 #endif
295 
A version of the Advection Diffusion Reaction equations that can be used with non-uniform mesh refine...
void broken_copy(const std::string &class_name)
Issue error message and terminate execution.
Base class for finite elements that can compute the quantities that are required for the Z2 error est...
void fill_in_generic_residual_contribution_adv_diff_react(Vector< double > &residuals, DenseMatrix< double > &jacobian, DenseMatrix< double > &mass_matrix, unsigned flag)
Add the element's contribution to the elemental residual vector.
unsigned nvertex_node() const
Number of vertex nodes in the element.
unsigned ncont_interpolated_values() const
Number of continuously interpolated values: NREAGENT.
void operator=(const RefineableQAdvectionDiffusionReactionElement< NREAGENT, DIM, NNODE_1D > &)
Broken assignment operator.
Refineable version of QAdvectionDiffusionReactionElement. Inherit from the standard QAdvectionDiffusi...
void operator=(const RefineableAdvectionDiffusionReactionEquations< NREAGENT, DIM > &)
Broken assignment operator.
QAdvectionDiffusionReactionElement elements are linear/quadrilateral/brick-shaped Advection Diffusion...
AdvectionDiffusionReactionReactionDerivFctPt Reaction_deriv_fct_pt
Pointer to reaction derivatives.
bool ALE_is_disabled
Boolean flag to indicate if ALE formulation is disabled when time-derivatives are computed...
double nodal_value(const unsigned &n, const unsigned &i) const
Return the i-th value stored at local node n. Produces suitably interpolated values for hanging nodes...
Definition: elements.h:2458
virtual unsigned c_index_adv_diff_react(const unsigned &i) const
Return the index at which the unknown i-th reagent is stored. The default value, i, is appropriate for single-physics problems, when there are only i variables, the values that satisfies the advection-diffusion-reaction equation. In derived multi-physics elements, this function should be overloaded to reflect the chosen storage scheme. Note that these equations require that the unknown is always stored at the same index at each node.
char t
Definition: cfortran.h:572
AdvectionDiffusionReactionReactionFctPt & reaction_fct_pt()
Access function: Pointer to reaction function.
Nodes are derived from Data, but, in addition, have a definite (Eulerian) position in a space of a gi...
Definition: nodes.h:852
void get_flux(const Vector< double > &s, Vector< double > &flux) const
Get flux: .
RefineableQAdvectionDiffusionReactionElement(const RefineableQAdvectionDiffusionReactionElement< NREAGENT, DIM, NNODE_1D > &dummy)
Broken copy constructor.
void further_setup_hanging_nodes()
Perform additional hanging node procedures for variables that are not interpolated by all nodes...
AdvectionDiffusionReactionReactionDerivFctPt & reaction_deriv_fct_pt()
Access function: Pointer to reaction derivatives function.
AdvectionDiffusionReactionSourceFctPt & source_fct_pt()
Access function: Pointer to source function.
unsigned nrecovery_order()
Order of recovery shape functions for Z2 error estimation: Same order as shape functions.
A class for all elements that solve the Advection Diffusion Reaction equations using isoparametric el...
AdvectionDiffusionReactionWindFctPt & wind_fct_pt()
Access function: Pointer to wind function.
unsigned num_Z2_flux_terms()
Number of 'flux' terms for Z2 error estimation.
void get_interpolated_values(const unsigned &t, const Vector< double > &s, Vector< double > &values)
Get the function values c in Vector. Note: Given the generality of the interface (this function is us...
AdvectionDiffusionReactionWindFctPt Wind_fct_pt
Pointer to wind function:
static char t char * s
Definition: cfortran.h:572
Node * vertex_node_pt(const unsigned &j) const
Pointer to the j-th vertex node in the element.
FaceGeometry()
Constructor: Call the constructor for the appropriate lower-dimensional QElement. ...
Vector< double > *& tau_pt()
Pointer to vector of dimensionless timescales.
Vector< double > * Diff_pt
Pointer to global diffusion coefficients.
AdvectionDiffusionReactionSourceFctPt Source_fct_pt
Pointer to source function:
virtual RefineableElement * father_element_pt() const
Return a pointer to the father element.
virtual Node * vertex_node_pt(const unsigned &j) const
Pointer to the j-th vertex node in the element. Broken virtual function in "pure" finite elements...
Definition: elements.h:2371
void broken_assign(const std::string &class_name)
Issue error message and terminate execution.
Vector< double > * Tau_pt
Pointer to global timescales.
void get_Z2_flux(const Vector< double > &s, Vector< double > &flux)
Get 'flux' for Z2 error recovery: Standard flux.from AdvectionDiffusionReaction equations.
Vector< double > *& diff_pt()
Pointer to vector of diffusion coefficients.
unsigned nnode() const
Return the number of nodes.
Definition: elements.h:2134
virtual unsigned nvertex_node() const
Definition: elements.h:2362
RefineableAdvectionDiffusionReactionEquations(const RefineableAdvectionDiffusionReactionEquations< NREAGENT, DIM > &dummy)
Broken copy constructor.
AdvectionDiffusionReactionReactionFctPt Reaction_fct_pt
Pointer to reaction function.
void get_interpolated_values(const Vector< double > &s, Vector< double > &values)
Get the function values c in Vector. Note: Given the generality of the interface (this function is us...
virtual void shape(const Vector< double > &s, Shape &psi) const =0
Calculate the geometric shape functions at local coordinate s. This function must be overloaded for e...
A general mesh class.
Definition: mesh.h:74