refineable_time_harmonic_linear_elasticity_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 refineable linear elasticity elements
31 
32 //Include guard to prevent multiple inclusions of this header
33 #ifndef OOMPH_REFINEABLE_TIME_HARMONIC_LINEAR_ELASTICITY_ELEMENTS_HEADER
34 #define OOMPH_REFINEABLE_TIME_HARMONIC_LINEAR_ELASTICITY_ELEMENTS_HEADER
35 
36 //oomph-lib headers
38 #include "../generic/refineable_quad_element.h"
39 #include "../generic/refineable_brick_element.h"
40 #include "../generic/error_estimator.h"
41 
42 namespace oomph
43 {
44 
45 //========================================================================
46 /// Class for Refineable TimeHarmonicLinearElasticity equations
47 //========================================================================
48  template<unsigned DIM>
50  public virtual TimeHarmonicLinearElasticityEquations<DIM>,
51  public virtual RefineableElement,
52  public virtual ElementWithZ2ErrorEstimator
53  {
54  public:
55 
56  /// \short Constructor
61  {}
62 
63 
64  /// \short Get the function value u in Vector.
65  /// Note: Given the generality of the interface (this function
66  /// is usually called from black-box documentation or interpolation
67  /// routines), the values Vector sets its own size in here.
68  void get_interpolated_values(const unsigned& t, const Vector<double>&s,
69  Vector<double>& values)
70  {
71  // Create enough initialised storage
72  values.resize(2*DIM,0.0);
73 
74  //Find out how many nodes there are
75  unsigned n_node = this->nnode();
76 
77  // Shape functions
78  Shape psi(n_node);
79  this->shape(s,psi);
80 
81  //Calculate displacements
82  for(unsigned i=0;i<DIM;i++)
83  {
84  //Get the index at which the i-th velocity is stored
85  std::complex<unsigned> u_nodal_index =
87  for(unsigned l=0;l<n_node;l++)
88  {
89  values[i]+=this->nodal_value(t,l,u_nodal_index.real())*psi(l);
90  values[i+DIM]+=this->nodal_value(t,l,u_nodal_index.imag())*psi(l);
91  }
92  }
93  }
94 
95  /// \short Get the current interpolated values (displacements).
96  /// Note: Given the generality of the interface (this function
97  /// is usually called from black-box documentation or interpolation
98  /// routines) ,the values Vector sets its own size in here.
100  Vector<double> & values)
101  {
102  this->get_interpolated_values(0,s,values);
103  }
104 
105  /// Number of 'flux' terms for Z2 error estimation
106  unsigned num_Z2_flux_terms()
107  {
108  // DIM Diagonal strain rates and DIM*(DIM-1)/2 off diagonal terms
109  return 2*(DIM + DIM*(DIM-1)/2);
110  }
111 
112  /// \short Get 'flux' for Z2 error recovery: Upper triangular entries
113  /// in strain tensor.
115  {
116 #ifdef PARANOID
117  unsigned num_entries=2*(DIM+((DIM*DIM)-DIM)/2);
118  if (flux.size()!=num_entries)
119  {
120  std::ostringstream error_message;
121  error_message << "The flux vector has the wrong number of entries, "
122  << flux.size() << ", whereas it should be "
123  << num_entries << std::endl;
124  throw OomphLibError(
125  error_message.str(),
126  OOMPH_CURRENT_FUNCTION,
127  OOMPH_EXCEPTION_LOCATION);
128  }
129 #endif
130 
131  // Get strain matrix
132  DenseMatrix<std::complex<double> > strain(DIM);
133  this->get_strain(s,strain);
134 
135  // Pack into flux Vector
136  unsigned icount=0;
137 
138  // Start with diagonal terms
139  for(unsigned i=0;i<DIM;i++)
140  {
141  flux[icount]=strain(i,i).real();
142  icount++;
143  flux[icount]=strain(i,i).imag();
144  icount++;
145  }
146 
147  //Off diagonals row by row
148  for(unsigned i=0;i<DIM;i++)
149  {
150  for(unsigned j=i+1;j<DIM;j++)
151  {
152  flux[icount]=strain(i,j).real();
153  icount++;
154  flux[icount]=strain(i,j).imag();
155  icount++;
156  }
157  }
158  }
159 
160  /// Number of continuously interpolated values: 2*DIM
161  unsigned ncont_interpolated_values() const {return 2*DIM;}
162 
163  /// Further build function, pass the pointers down to the sons
165  {
168  (this->father_element_pt());
169 
170  // Set pointer to body force function
171  this->Body_force_fct_pt = cast_father_element_pt->body_force_fct_pt();
172 
173  // Set pointer to the contitutive law
174  this->Elasticity_tensor_pt =
175  cast_father_element_pt->elasticity_tensor_pt();
176 
177  // Set the frequency
178  this->Omega_sq_pt = cast_father_element_pt->omega_sq_pt();
179  }
180 
181 
182  private:
183 
184  /// Overloaded helper function to take hanging nodes into account
186  Vector<double> &residuals,DenseMatrix<double> &jacobian,unsigned flag);
187 
188  };
189 
190 
191 //////////////////////////////////////////////////////////////////////
192 //////////////////////////////////////////////////////////////////////
193 //////////////////////////////////////////////////////////////////////
194 
195 //========================================================================
196 /// Class for refineable QTimeHarmonicLinearElasticityElement elements
197 //========================================================================
198  template<unsigned DIM, unsigned NNODE_1D>
200  public virtual QTimeHarmonicLinearElasticityElement<DIM,NNODE_1D>,
202  public virtual RefineableQElement<DIM>
203  {
204  public:
205 
206  /// Constructor:
208  QTimeHarmonicLinearElasticityElement<DIM,NNODE_1D>(),
211  RefineableQElement<DIM>() {}
212 
213  /// Empty rebuild from sons, no need to reconstruct anything here
214  void rebuild_from_sons(Mesh* &mesh_pt) {}
215 
216  /// \short Number of vertex nodes in the element
217  unsigned nvertex_node() const
219 
220  /// \short Pointer to the j-th vertex node in the element
221  Node* vertex_node_pt(const unsigned& j) const
222  {
223  return
225  }
226 
227  /// \short Order of recovery shape functions for Z2 error estimation:
228  /// Same order as shape functions.
229  unsigned nrecovery_order() {return NNODE_1D-1;}
230 
231  /// \short No additional hanging node procedures are required
233 
234  };
235 
236 
237 
238 //==============================================================
239 /// FaceGeometry of the 2D
240 /// RefineableQTimeHarmonicLinearElasticityElement elements
241 //==============================================================
242  template<unsigned NNODE_1D>
244  public virtual QElement<1,NNODE_1D>
245  {
246  public:
247  //Make sure that we call the constructor of the QElement
248  //Only the Intel compiler seems to need this!
249  FaceGeometry() : QElement<1,NNODE_1D>() {}
250  };
251 
252 //==============================================================
253 /// FaceGeometry of the FaceGeometry of the 2D
254 /// RefineableQTimeHarmonicLinearElasticityElement
255 //==============================================================
256  template<unsigned NNODE_1D>
259  public virtual PointElement
260  {
261  public:
262  //Make sure that we call the constructor of the QElement
263  //Only the Intel compiler seems to need this!
265  };
266 
267 
268 //==============================================================
269 /// FaceGeometry of the 3D RefineableQTimeHarmonicLinearElasticityElement
270 /// elements
271 //==============================================================
272  template<unsigned NNODE_1D>
274  public virtual QElement<2,NNODE_1D>
275  {
276  public:
277  //Make sure that we call the constructor of the QElement
278  //Only the Intel compiler seems to need this!
279  FaceGeometry() : QElement<2,NNODE_1D>() {}
280  };
281 
282 //==============================================================
283 /// FaceGeometry of the FaceGeometry of the 3D
284 /// RefineableQTimeHarmonicLinearElasticityElement
285 //==============================================================
286  template<unsigned NNODE_1D>
289  public virtual QElement<1,NNODE_1D>
290  {
291  public:
292  //Make sure that we call the constructor of the QElement
293  //Only the Intel compiler seems to need this!
294  FaceGeometry() : QElement<1,NNODE_1D>() {}
295  };
296 
297 
298 }
299 
300 #endif
void further_setup_hanging_nodes()
No additional hanging node procedures are required.
unsigned num_Z2_flux_terms()
Number of 'flux' terms for Z2 error estimation.
TimeHarmonicElasticityTensor *& elasticity_tensor_pt()
Return the pointer to the elasticity_tensor.
Base class for finite elements that can compute the quantities that are required for the Z2 error est...
BodyForceFctPt & body_force_fct_pt()
Access function: Pointer to body force function.
Class for refineable QTimeHarmonicLinearElasticityElement elements.
cstr elem_len * i
Definition: cfortran.h:607
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
char t
Definition: cfortran.h:572
Nodes are derived from Data, but, in addition, have a definite (Eulerian) position in a space of a gi...
Definition: nodes.h:852
Node * vertex_node_pt(const unsigned &j) const
Pointer to the j-th vertex node in the element.
void fill_in_generic_contribution_to_residuals_time_harmonic_linear_elasticity(Vector< double > &residuals, DenseMatrix< double > &jacobian, unsigned flag)
Overloaded helper function to take hanging nodes into account.
void rebuild_from_sons(Mesh *&mesh_pt)
Empty rebuild from sons, no need to reconstruct anything here.
void get_strain(const Vector< double > &s, DenseMatrix< std::complex< double > > &strain) const
Return the strain tensor.
double *& omega_sq_pt()
Access function for square of non-dim frequency.
static char t char * s
Definition: cfortran.h:572
void get_interpolated_values(const Vector< double > &s, Vector< double > &values)
Get the current interpolated values (displacements). Note: Given the generality of the interface (thi...
BodyForceFctPt Body_force_fct_pt
Pointer to body force function.
virtual RefineableElement * father_element_pt() const
Return a pointer to the father element.
TimeHarmonicElasticityTensor * Elasticity_tensor_pt
Pointer to the elasticity tensor.
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 get_interpolated_values(const unsigned &t, const Vector< double > &s, Vector< double > &values)
Get the function value u in Vector. Note: Given the generality of the interface (this function is usu...
unsigned nnode() const
Return the number of nodes.
Definition: elements.h:2134
virtual std::complex< unsigned > u_index_time_harmonic_linear_elasticity(const unsigned i) const
Return the index at which the i-th real or imag unknown displacement component is stored...
virtual unsigned nvertex_node() const
Definition: elements.h:2362
Class for dense matrices, storing all the values of the matrix as a pointer to a pointer with assorte...
Definition: communicator.h:50
void get_Z2_flux(const Vector< double > &s, Vector< double > &flux)
Get 'flux' for Z2 error recovery: Upper triangular entries in strain tensor.
unsigned nrecovery_order()
Order of recovery shape functions for Z2 error estimation: Same order as shape functions.
void further_build()
Further build function, pass the pointers down to the sons.
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
unsigned ncont_interpolated_values() const
Number of continuously interpolated values: 2*DIM.