time_harmonic_fourier_decomposed_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: 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 
31 //Include guards to prevent multiple inclusion of the header
32 #ifndef OOMPH_FOURIER_DECOMPOSED_TIME_HARMONIC_LINEAR_ELASTICITY_ELEMENTS_HEADER
33 #define OOMPH_FOURIER_DECOMPOSED_TIME_HARMONIC_LINEAR_ELASTICITY_ELEMENTS_HEADER
34 
35 // Config header generated by autoconfig
36 #ifdef HAVE_CONFIG_H
37  #include <oomph-lib-config.h>
38 #endif
39 
40 
41 #ifdef OOMPH_HAS_MPI
42 #include "mpi.h"
43 #endif
44 
45 #include<complex>
46 
47 
48 //OOMPH-LIB headers
49 #include "../generic/Qelements.h"
50 #include "../generic/Telements.h"
51 #include "../generic/projection.h"
52 #include "../generic/error_estimator.h"
53 
54 
55 namespace oomph
56 {
57 //=======================================================================
58 /// A base class for elements that solve the Fourier decomposed (in
59 /// cylindrical polars) equations of time-harmonic linear elasticity.
60 //=======================================================================
62  public virtual FiniteElement
63  {
64  public:
65 
66  /// \short Return the index at which the i-th (i=0: r, i=1: z; i=2: theta)
67  /// real or imag unknown displacement component is stored at the nodes.
68  /// The default assignment here (u_r_real, u_z_real, ..., u_theta_imag)
69  /// is appropriate for single-physics problems.
70  virtual inline std::complex<unsigned>
72  {
73  return std::complex<unsigned>(i,i+3);
74  }
75 
76  /// Compute vector of FE interpolated displacement u at local coordinate s
78  const Vector<double> &s,
79  Vector<std::complex<double> >& disp)
80  const
81  {
82  //Find number of nodes
83  unsigned n_node = nnode();
84 
85  //Local shape function
86  Shape psi(n_node);
87 
88  //Find values of shape function
89  shape(s,psi);
90 
91  for (unsigned i=0;i<3;i++)
92  {
93  //Index at which the nodal value is stored
94  std::complex<unsigned> u_nodal_index =
96 
97  //Initialise value of u
98  disp[i] = std::complex<double>(0.0,0.0);
99 
100  //Loop over the local nodes and sum
101  for(unsigned l=0;l<n_node;l++)
102  {
103  const std::complex<double> u_value(
104  nodal_value(l,u_nodal_index.real()),
105  nodal_value(l,u_nodal_index.imag()));
106 
107  disp[i] += u_value*psi[l];
108  }
109  }
110  }
111 
112  /// \short Return FE interpolated displacement u[i] (i=0: r, i=1: z; i=2:
113  /// theta) at local coordinate s
114  std::complex<double>
116  const Vector<double> &s,
117  const unsigned &i) const
118  {
119  //Find number of nodes
120  unsigned n_node = nnode();
121 
122  //Local shape function
123  Shape psi(n_node);
124 
125  //Find values of shape function
126  shape(s,psi);
127 
128  //Get nodal index at which i-th velocity is stored
129  std::complex<unsigned> u_nodal_index =
131 
132  //Initialise value of u
133  std::complex<double> interpolated_u(0.0,0.0);
134 
135  //Loop over the local nodes and sum
136  for(unsigned l=0;l<n_node;l++)
137  {
138  const std::complex<double> u_value(
139  nodal_value(l,u_nodal_index.real()),
140  nodal_value(l,u_nodal_index.imag()));
141 
142  interpolated_u += u_value*psi[l];
143  }
144 
145  return(interpolated_u);
146  }
147 
148 
149  /// \short Function pointer to function that specifies the body force
150  /// as a function of the Cartesian coordinates and time FCT(x,b) --
151  /// x and b are Vectors!
152  typedef void (*BodyForceFctPt)(const Vector<double>& x,
154 
155  /// \short Constructor: Set null pointers for constitutive law.
156  /// Set physical parameter values to
157  /// default values, and set body force to zero.
162 
163  ///Access function for square of non-dim frequency
164  const std::complex<double>& omega_sq() const {return *Omega_sq_pt;}
165 
166  /// Access function for square of non-dim frequency
167  std::complex<double>* &omega_sq_pt() {return Omega_sq_pt;}
168 
169  /// Return the pointer to Young's modulus
170  std::complex<double>*& youngs_modulus_pt() {return Youngs_modulus_pt;}
171 
172  /// Access function to Young's modulus
173  inline std::complex<double> youngs_modulus() const
174  {
175  return (*Youngs_modulus_pt);
176  }
177 
178  ///Access function for Poisson's ratio
179  std::complex<double>& nu() const
180  {
181 #ifdef PARANOID
182  if (Nu_pt==0)
183  {
184  std::ostringstream error_message;
185  error_message << "No pointer to Poisson's ratio set. Please set one!\n";
186  throw OomphLibError(
187  error_message.str(),
188  OOMPH_CURRENT_FUNCTION,
189  OOMPH_EXCEPTION_LOCATION);
190  }
191 #endif
192  return *Nu_pt;
193  }
194 
195  /// Access function for pointer to Poisson's ratio
196  std::complex<double>*& nu_pt() {return Nu_pt;}
197 
198  ///Access function for Fourier wavenumber
199  int& fourier_wavenumber() const
200  {
201 #ifdef PARANOID
202  if (Fourier_wavenumber_pt==0)
203  {
204  std::ostringstream error_message;
205  error_message
206  << "No pointer to Fourier wavenumber set. Please set one!\n";
207  throw OomphLibError(
208  error_message.str(),
209  OOMPH_CURRENT_FUNCTION,
210  OOMPH_EXCEPTION_LOCATION);
211  }
212 #endif
213  return *Fourier_wavenumber_pt;
214  }
215 
216  /// Access function for pointer to Fourier wavenumber
218 
219  /// Access function: Pointer to body force function
221 
222  /// Access function: Pointer to body force function (const version)
224 
225  /// \short Evaluate body force at Eulerian coordinate x at present time
226  /// (returns zero vector if no body force function pointer has been set)
227  inline void body_force(const Vector<double>& x,
228  Vector<std::complex<double> >& b) const
229  {
230  //If no function has been set, return zero vector
231  if(Body_force_fct_pt==0)
232  {
233  // Get spatial dimension of element
234  unsigned n=dim();
235  for (unsigned i=0;i<n;i++)
236  {
237  b[i] = std::complex<double>(0.0,0.0);
238  }
239  }
240  else
241  {
242  (*Body_force_fct_pt)(x,b);
243  }
244  }
245 
246  /// \short The number of "DOF types" that degrees of freedom in this element
247  /// are sub-divided into: for now lump them all into one DOF type.
248  /// Can be adjusted later
249  unsigned ndof_types() const
250  {
251  return 1;
252  }
253 
254  /// \short Create a list of pairs for all unknowns in this element,
255  /// so that the first entry in each pair contains the global equation
256  /// number of the unknown, while the second one contains the number
257  /// of the "DOF type" that this unknown is associated with.
258  /// (Function can obviously only be called if the equation numbering
259  /// scheme has been set up.)
261  std::list<std::pair<unsigned long,unsigned> >& dof_lookup_list) const
262  {
263 
264  // temporary pair (used to store dof lookup prior to being added
265  // to list)
266  std::pair<unsigned long,unsigned> dof_lookup;
267 
268  // number of nodes
269  const unsigned n_node = this->nnode();
270 
271  //Integer storage for local unknown
272  int local_unknown=0;
273 
274  //Loop over the nodes
275  for(unsigned n=0;n<n_node;n++)
276  {
277  //Loop over dimension (real and imag displacement components)
278  for(unsigned i=0;i<6;i++)
279  {
280  //If the variable is free
281  local_unknown = nodal_local_eqn(n,i);
282 
283  // ignore pinned values
284  if (local_unknown >= 0)
285  {
286  // store dof lookup in temporary pair: First entry in pair
287  // is global equation number; second entry is DOF type
288  dof_lookup.first = this->eqn_number(local_unknown);
289  dof_lookup.second = 0;
290 
291  // add to list
292  dof_lookup_list.push_front(dof_lookup);
293  }
294  }
295  }
296  }
297 
298 
299  protected:
300 
301  /// Square of nondim frequency
302  std::complex<double>* Omega_sq_pt;
303 
304  /// Pointer to the Young's modulus
305  std::complex<double>* Youngs_modulus_pt;
306 
307  /// Pointer to Poisson's ratio
308  std::complex<double>* Nu_pt;
309 
310  /// Pointer to Fourier wavenumber
312 
313  /// Pointer to body force function
315 
316  /// Static default value for squared frequency
317  static std::complex<double> Default_omega_sq_value;
318 
319  /// Static default value for Young's modulus (1.0 -- for natural
320  /// scaling, i.e. all stresses have been non-dimensionalised by
321  /// the same reference Young's modulus. Setting the "non-dimensional"
322  /// Young's modulus (obtained by de-referencing Youngs_modulus_pt)
323  /// to a number larger than one means that the material is stiffer
324  /// than assumed in the non-dimensionalisation.
325  static std::complex<double> Default_youngs_modulus_value;
326  };
327 
328 
329 ///////////////////////////////////////////////////////////////////////
330 ///////////////////////////////////////////////////////////////////////
331 ///////////////////////////////////////////////////////////////////////
332 
333 
334 //=======================================================================
335 /// A class for elements that solve the Fourier decomposed (in cylindrical
336 /// polars) equations of time-harmonic linear elasticity
337 //=======================================================================
340  {
341  public:
342 
343  /// \short Constructor
345 
346  /// Number of values required at node n.
347  unsigned required_nvalue(const unsigned &n) const {return 6;}
348 
349  /// \short Return the residuals for the equations (the discretised
350  /// principle of virtual displacements)
352  {
355  }
356 
357 
358  /// The jacobian is calculated by finite differences by default,
359  /// We need only to take finite differences w.r.t. positional variables
360  /// For this element
362  DenseMatrix<double> &jacobian)
363  {
364  //Add the contribution to the residuals
365  this->
367  residuals,jacobian,1);
368  }
369 
370 
371  /// Get strain tensor
372  void get_strain(const Vector<double>& s,
373  DenseMatrix<std::complex<double> >& strain);
374 
375  /// \short Compute norm of solution: square of the L2 norm
376  void compute_norm(double& norm);
377 
378  ///Output exact solution: r,z, u_r_real, u_z_real, ..., u_theta_imag
379  void output_fct(std::ostream &outfile,
380  const unsigned &nplot,
382 
383  /// Output: r,z, u_r_real, u_z_real, ..., u_theta_imag
384  void output(std::ostream &outfile)
385  {
386  unsigned n_plot=5;
387  output(outfile,n_plot);
388  }
389 
390  /// Output: r,z, u_r_real, u_z_real, ..., u_theta_imag
391  void output(std::ostream &outfile, const unsigned &n_plot);
392 
393  /// C-style output: r,z, u_r_real, u_z_real, ..., u_theta_imag
394  void output(FILE* file_pt)
395  {
396  unsigned n_plot=5;
397  output(file_pt,n_plot);
398  }
399 
400  /// Output: r,z, u_r_real, u_z_real, ..., u_theta_imag
401  void output(FILE* file_pt, const unsigned &n_plot);
402 
403  /// Validate against exact solution.
404  /// Solution is provided via function pointer.
405  /// Plot at a given number of plot points and compute L2 error
406  /// and L2 norm of displacement solution over element
407  void compute_error(std::ostream &outfile,
409  double& error, double& norm);
410 
411 
412  private:
413 
414  /// \short Private helper function to compute residuals and (if requested
415  /// via flag) also the Jacobian matrix.
417  Vector<double> &residuals,DenseMatrix<double> &jacobian,unsigned flag);
418 
419  };
420 
421 
422 ////////////////////////////////////////////////////////////////////////
423 ////////////////////////////////////////////////////////////////////////
424 ////////////////////////////////////////////////////////////////////////
425 
426 
427 //===========================================================================
428 /// An Element that solves the equations of Fourier decomposed (in cylindrical
429 /// polars) time-harmonic linear elasticity, using QElements for the geometry.
430 //============================================================================
431  template<unsigned NNODE_1D>
433  public virtual QElement<2,NNODE_1D>,
435  {
436  public:
437 
438  /// Constructor
440  QElement<2,NNODE_1D>(),
442 
443  /// Output function
444  void output(std::ostream &outfile)
446 
447  /// Output function
448  void output(std::ostream &outfile, const unsigned &n_plot)
450  output(outfile,n_plot);}
451 
452 
453  /// C-style output function
454  void output(FILE* file_pt)
456 
457  /// C-style output function
458  void output(FILE* file_pt, const unsigned &n_plot)
460  output(file_pt,n_plot);}
461 
462  };
463 
464 
465 //============================================================================
466 /// FaceGeometry of a linear
467 /// QTimeHarmonicFourierDecomposedLinearElasticityElement element
468 //============================================================================
469  template<unsigned NNODE_1D>
471  public virtual QElement<1,NNODE_1D>
472  {
473  public:
474  /// Constructor must call the constructor of the underlying element
475  FaceGeometry() : QElement<1,NNODE_1D>() {}
476  };
477 
478 
479 
480 
481 ////////////////////////////////////////////////////////////////////////
482 ////////////////////////////////////////////////////////////////////////
483 ////////////////////////////////////////////////////////////////////////
484 
485 
486 //===========================================================================
487 /// An Element that solves the equations of Fourier decomposed (in cylindrical
488 /// polars) time-harmonic linear elasticity, using TElements for the geometry.
489 //============================================================================
490  template<unsigned NNODE_1D>
492  public virtual TElement<2,NNODE_1D>,
494  public virtual ElementWithZ2ErrorEstimator
495  {
496  public:
497 
498  /// Constructor
500  TElement<2,NNODE_1D>(),
502 
503  /// Output function
504  void output(std::ostream &outfile)
506 
507  /// Output function
508  void output(std::ostream &outfile, const unsigned &n_plot)
510  output(outfile,n_plot);}
511 
512  /// C-style output function
513  void output(FILE* file_pt)
515 
516  /// C-style output function
517  void output(FILE* file_pt, const unsigned &n_plot)
519  output(file_pt,n_plot);}
520 
521 
522  /// \short Number of vertex nodes in the element
523  unsigned nvertex_node() const
525 
526  /// \short Pointer to the j-th vertex node in the element
527  Node* vertex_node_pt(const unsigned& j) const
528  {
530  }
531 
532  /// \short Order of recovery shape functions for Z2 error estimation:
533  /// Same order as shape functions.
534  unsigned nrecovery_order() {return NNODE_1D-1;}
535 
536  /// Number of 'flux' terms for Z2 error estimation
537  unsigned num_Z2_flux_terms()
538  {
539  // 3 Diagonal strain rates and 3 off diagonal terms for real and imag part
540  return 12;
541  }
542 
543  /// \short Get 'flux' for Z2 error recovery: Upper triangular entries
544  /// in strain tensor.
546  {
547 #ifdef PARANOID
548  unsigned num_entries=12;
549  if (flux.size()!=num_entries)
550  {
551  std::ostringstream error_message;
552  error_message << "The flux vector has the wrong number of entries, "
553  << flux.size() << ", whereas it should be "
554  << num_entries << std::endl;
555  throw OomphLibError(
556  error_message.str(),
557  OOMPH_CURRENT_FUNCTION,
558  OOMPH_EXCEPTION_LOCATION);
559  }
560 #endif
561 
562  // Get strain matrix
564  this->get_strain(s,strain);
565 
566  // Pack into flux Vector
567  unsigned icount=0;
568 
569  // Start with diagonal terms
570  for(unsigned i=0;i<3;i++)
571  {
572  flux[icount]=strain(i,i).real();
573  icount++;
574  flux[icount]=strain(i,i).imag();
575  icount++;
576  }
577 
578  //Off diagonals row by row
579  for(unsigned i=0;i<3;i++)
580  {
581  for(unsigned j=i+1;j<3;j++)
582  {
583  flux[icount]=strain(i,j).real();
584  icount++;
585  flux[icount]=strain(i,j).imag();
586  icount++;
587  }
588  }
589  }
590 
591 
592 
593  };
594 
595 
596 //============================================================================
597 /// FaceGeometry of a linear
598 /// TTimeHarmonicFourierDecomposedLinearElasticityElement element
599 //============================================================================
600  template<unsigned NNODE_1D>
602  public virtual TElement<1,NNODE_1D>
603  {
604  public:
605  /// Constructor must call the constructor of the underlying element
606  FaceGeometry() : TElement<1,NNODE_1D>() {}
607  };
608 
609 
610 ////////////////////////////////////////////////////////////////////
611 ////////////////////////////////////////////////////////////////////
612 ////////////////////////////////////////////////////////////////////
613 
614 
615 //==========================================================
616 /// Fourier-decomposed time-harmonic linear elasticity
617 /// upgraded to become projectable
618 //==========================================================
619  template<class TIME_HARMONIC_LINEAR_ELAST_ELEMENT>
621  public virtual ProjectableElement<TIME_HARMONIC_LINEAR_ELAST_ELEMENT>
622  {
623 
624  public:
625 
626  /// \short Constructor [this was only required explicitly
627  /// from gcc 4.5.2 onwards...]
629 
630  /// \short Specify the values associated with field fld.
631  /// The information is returned in a vector of pairs which comprise
632  /// the Data object and the value within it, that correspond to field fld.
633  /// In the underlying time-harmonic linear elasticity elemements the
634  /// real and complex parts of the displacements are stored
635  /// at the nodal values
637  {
638  // Create the vector
639  Vector<std::pair<Data*,unsigned> > data_values;
640 
641  // Loop over all nodes and extract the fld-th nodal value
642  unsigned nnod=this->nnode();
643  for (unsigned j=0;j<nnod;j++)
644  {
645  // Add the data value associated with the velocity components
646  data_values.push_back(std::make_pair(this->node_pt(j),fld));
647  }
648 
649  // Return the vector
650  return data_values;
651 
652  }
653 
654  /// \short Number of fields to be projected: 3*dim, corresponding to
655  /// real and imag parts of the displacement components
657  {
658  return 3*this->dim();
659  }
660 
661  /// \short Number of history values to be stored for fld-th field.
662  /// (includes present value!)
663  unsigned nhistory_values_for_projection(const unsigned &fld)
664  {
665 #ifdef PARANOID
666  if (fld>5)
667  {
668  std::stringstream error_stream;
669  error_stream
670  << "Elements only store six fields so fld can't be"
671  << " " << fld << std::endl;
672  throw OomphLibError(
673  error_stream.str(),
674  OOMPH_CURRENT_FUNCTION,
675  OOMPH_EXCEPTION_LOCATION);
676  }
677 #endif
678  return this->node_pt(0)->ntstorage();
679  }
680 
681  ///\short Number of positional history values: Read out from
682  /// positional timestepper
683  /// (Note: count includes current value!)
685  {
686  return this->node_pt(0)->position_time_stepper_pt()->ntstorage();
687  }
688 
689  /// \short Return Jacobian of mapping and shape functions of field fld
690  /// at local coordinate s
691  double jacobian_and_shape_of_field(const unsigned &fld,
692  const Vector<double> &s,
693  Shape &psi)
694  {
695  unsigned n_dim=this->dim();
696  unsigned n_node=this->nnode();
697  DShape dpsidx(n_node,n_dim);
698 
699  // Call the derivatives of the shape functions and return
700  // the Jacobian
701  return this->dshape_eulerian(s,psi,dpsidx);
702  }
703 
704 
705 
706  /// \short Return interpolated field fld at local coordinate s, at time level
707  /// t (t=0: present; t>0: history values)
708  double get_field(const unsigned &t,
709  const unsigned &fld,
710  const Vector<double>& s)
711  {
712  unsigned n_node=this->nnode();
713 
714 #ifdef PARANOID
715  unsigned n_dim=this->node_pt(0)->ndim();
716 #endif
717 
718  //Local shape function
719  Shape psi(n_node);
720 
721  //Find values of shape function
722  this->shape(s,psi);
723 
724  //Initialise value of u
725  double interpolated_u = 0.0;
726 
727  //Sum over the local nodes at that time
728  for(unsigned l=0;l<n_node;l++)
729  {
730 #ifdef PARANOID
731  unsigned nvalue=this->node_pt(l)->nvalue();
732  if (nvalue!=3*n_dim)
733  {
734  std::stringstream error_stream;
735  error_stream
736  << "Current implementation only works for non-resized nodes\n"
737  << "but nvalue= " << nvalue << "!= 3 dim = " << 3*n_dim << std::endl;
738  throw OomphLibError(
739  error_stream.str(),
740  OOMPH_CURRENT_FUNCTION,
741  OOMPH_EXCEPTION_LOCATION);
742  }
743 #endif
744  interpolated_u += this->nodal_value(t,l,fld)*psi[l];
745  }
746  return interpolated_u;
747  }
748 
749 
750  ///Return number of values in field fld
751  unsigned nvalue_of_field(const unsigned &fld)
752  {
753  return this->nnode();
754  }
755 
756 
757  ///Return local equation number of value j in field fld.
758  int local_equation(const unsigned &fld,
759  const unsigned &j)
760  {
761 #ifdef PARANOID
762  unsigned n_dim=this->node_pt(0)->ndim();
763  unsigned nvalue=this->node_pt(j)->nvalue();
764  if (nvalue!=3*n_dim)
765  {
766  std::stringstream error_stream;
767  error_stream
768  << "Current implementation only works for non-resized nodes\n"
769  << "but nvalue= " << nvalue << "!= 3 dim = " << 3*n_dim << std::endl;
770  throw OomphLibError(
771  error_stream.str(),
772  OOMPH_CURRENT_FUNCTION,
773  OOMPH_EXCEPTION_LOCATION);
774  }
775 #endif
776  return this->nodal_local_eqn(j,fld);
777  }
778 
779 
780  };
781 
782 
783 //=======================================================================
784 /// Face geometry for element is the same as that for the underlying
785 /// wrapped element
786 //=======================================================================
787  template<class ELEMENT>
789  : public virtual FaceGeometry<ELEMENT>
790  {
791  public:
792  FaceGeometry() : FaceGeometry<ELEMENT>() {}
793  };
794 
795 
796 //=======================================================================
797 /// Face geometry of the Face Geometry for element is the same as
798 /// that for the underlying wrapped element
799 //=======================================================================
800  template<class ELEMENT>
802  : public virtual FaceGeometry<FaceGeometry<ELEMENT> >
803  {
804  public:
806  };
807 
808 
809 }
810 
811 
812 #endif
813 
814 
815 
816 
ProjectableTimeHarmonicFourierDecomposedLinearElasticityElement()
Constructor [this was only required explicitly from gcc 4.5.2 onwards...].
void output_fct(std::ostream &outfile, const unsigned &nplot, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt)
Output exact solution: r,z, u_r_real, u_z_real, ..., u_theta_imag.
Base class for finite elements that can compute the quantities that are required for the Z2 error est...
unsigned nrecovery_order()
Order of recovery shape functions for Z2 error estimation: Same order as shape functions.
void get_dof_numbers_for_unknowns(std::list< std::pair< unsigned long, unsigned > > &dof_lookup_list) const
Create a list of pairs for all unknowns in this element, so that the first entry in each pair contain...
unsigned nfields_for_projection()
Number of fields to be projected: 3*dim, corresponding to real and imag parts of the displacement com...
TimeHarmonicFourierDecomposedLinearElasticityEquationsBase()
Constructor: Set null pointers for constitutive law. Set physical parameter values to default values...
unsigned ntstorage() const
Return total number of doubles stored per value to record time history of each value (one for steady ...
Definition: nodes.cc:847
const std::complex< double > & omega_sq() const
Access function for square of non-dim frequency.
double dshape_eulerian(const Vector< double > &s, Shape &psi, DShape &dpsidx) const
Compute the geometric shape functions and also first derivatives w.r.t. global coordinates at local c...
Definition: elements.cc:3225
int local_equation(const unsigned &fld, const unsigned &j)
Return local equation number of value j in field fld.
cstr elem_len * i
Definition: cfortran.h:607
unsigned long eqn_number(const unsigned &ieqn_local) const
Return the global equation number corresponding to the ieqn_local-th local equation number...
Definition: elements.h:709
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
A general Finite Element class.
Definition: elements.h:1271
void get_Z2_flux(const Vector< double > &s, Vector< double > &flux)
Get 'flux' for Z2 error recovery: Upper triangular entries in strain tensor.
char t
Definition: cfortran.h:572
unsigned required_nvalue(const unsigned &n) const
Number of values required at node n.
TimeStepper *& position_time_stepper_pt()
Return a pointer to the position timestepper.
Definition: nodes.h:969
Nodes are derived from Data, but, in addition, have a definite (Eulerian) position in a space of a gi...
Definition: nodes.h:852
unsigned ntstorage() const
Return the number of doubles required to represent history (one for steady)
Definition: timesteppers.h:562
unsigned nvalue() const
Return number of values stored in data object (incl pinned ones).
Definition: nodes.h:448
BodyForceFctPt body_force_fct_pt() const
Access function: Pointer to body force function (const version)
unsigned ndof_types() const
The number of "DOF types" that degrees of freedom in this element are sub-divided into: for now lump ...
virtual std::complex< unsigned > u_index_time_harmonic_fourier_decomposed_linear_elasticity(const unsigned i) const
Return the index at which the i-th (i=0: r, i=1: z; i=2: theta) real or imag unknown displacement com...
unsigned nhistory_values_for_coordinate_projection()
Number of positional history values: Read out from positional timestepper (Note: count includes curre...
std::complex< double > *& nu_pt()
Access function for pointer to Poisson's ratio.
double jacobian_and_shape_of_field(const unsigned &fld, const Vector< double > &s, Shape &psi)
Return Jacobian of mapping and shape functions of field fld at local coordinate s.
virtual void fill_in_generic_contribution_to_residuals_fourier_decomp_time_harmonic_linear_elasticity(Vector< double > &residuals, DenseMatrix< double > &jacobian, unsigned flag)
Private helper function to compute residuals and (if requested via flag) also the Jacobian matrix...
void(* BodyForceFctPt)(const Vector< double > &x, Vector< std::complex< double > > &b)
Function pointer to function that specifies the body force as a function of the Cartesian coordinates...
static char t char * s
Definition: cfortran.h:572
Vector< std::pair< Data *, unsigned > > data_values_of_field(const unsigned &fld)
Specify the values associated with field fld. The information is returned in a vector of pairs which ...
Node * vertex_node_pt(const unsigned &j) const
Pointer to the j-th vertex node in the element.
void compute_error(std::ostream &outfile, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt, double &error, double &norm)
void output(FILE *file_pt)
C-style output: r,z, u_r_real, u_z_real, ..., u_theta_imag.
std::complex< double > interpolated_u_time_harmonic_fourier_decomposed_linear_elasticity(const Vector< double > &s, const unsigned &i) const
Return FE interpolated displacement u[i] (i=0: r, i=1: z; i=2: theta) at local coordinate s...
std::complex< double > *& omega_sq_pt()
Access function for square of non-dim frequency.
unsigned nhistory_values_for_projection(const unsigned &fld)
Number of history values to be stored for fld-th field. (includes present value!) ...
void(* SteadyExactSolutionFctPt)(const Vector< double > &, Vector< double > &)
Function pointer for function that computes vector-valued steady "exact solution" as ...
Definition: elements.h:1720
double get_field(const unsigned &t, const unsigned &fld, const Vector< double > &s)
Return interpolated field fld at local coordinate s, at time level t (t=0: present; t>0: history valu...
void fill_in_contribution_to_residuals(Vector< double > &residuals)
Return the residuals for the equations (the discretised principle of virtual displacements) ...
Node *& node_pt(const unsigned &n)
Return a pointer to the local node n.
Definition: elements.h:2097
void output(std::ostream &outfile)
Output: r,z, u_r_real, u_z_real, ..., u_theta_imag.
int nodal_local_eqn(const unsigned &n, const unsigned &i) const
Return the local equation number corresponding to the i-th value at the n-th local node...
Definition: elements.h:1383
static std::complex< double > Default_omega_sq_value
Static default value for squared frequency.
unsigned dim() const
Return the spatial dimension of the element, i.e. the number of local coordinates required to paramet...
Definition: elements.h:2470
unsigned nnode() const
Return the number of nodes.
Definition: elements.h:2134
void interpolated_u_time_harmonic_fourier_decomposed_linear_elasticity(const Vector< double > &s, Vector< std::complex< double > > &disp) const
Compute vector of FE interpolated displacement u at local coordinate s.
unsigned ndim() const
Return (Eulerian) spatial dimension of the node.
Definition: nodes.h:992
void fill_in_contribution_to_jacobian(Vector< double > &residuals, DenseMatrix< double > &jacobian)
void body_force(const Vector< double > &x, Vector< std::complex< double > > &b) const
Evaluate body force at Eulerian coordinate x at present time (returns zero vector if no body force fu...
void get_strain(const Vector< double > &s, DenseMatrix< std::complex< double > > &strain)
Get strain tensor.
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...
static DenseMatrix< double > Dummy_matrix
Empty dense matrix used as a dummy argument to combined residual and jacobian functions in the case w...
Definition: elements.h:228