28 #ifndef OOMPH_FOURIER_DECOMPOSED_TIME_HARMONIC_LINEAR_ELASTICITY_ELEMENTS_HEADER
29 #define OOMPH_FOURIER_DECOMPOSED_TIME_HARMONIC_LINEAR_ELASTICITY_ELEMENTS_HEADER
33 #include <oomph-lib-config.h>
45 #include "../generic/Qelements.h"
46 #include "../generic/Telements.h"
47 #include "../generic/projection.h"
48 #include "../generic/error_estimator.h"
66 const unsigned i)
const
68 return std::complex<unsigned>(
i,
i + 3);
76 unsigned n_node =
nnode();
84 for (
unsigned i = 0;
i < 3;
i++)
87 std::complex<unsigned> u_nodal_index =
91 disp[
i] = std::complex<double>(0.0, 0.0);
94 for (
unsigned l = 0; l < n_node; l++)
96 const std::complex<double> u_value(
100 disp[
i] += u_value * psi[l];
111 unsigned n_node =
nnode();
120 std::complex<unsigned> u_nodal_index =
124 std::complex<double> interpolated_u(0.0, 0.0);
127 for (
unsigned l = 0; l < n_node; l++)
129 const std::complex<double> u_value(
133 interpolated_u += u_value * psi[l];
136 return (interpolated_u);
183 std::complex<double>&
nu()
const
188 std::ostringstream error_message;
189 error_message <<
"No pointer to Poisson's ratio set. Please set one!\n";
191 OOMPH_CURRENT_FUNCTION,
192 OOMPH_EXCEPTION_LOCATION);
210 std::ostringstream error_message;
212 <<
"No pointer to Fourier wavenumber set. Please set one!\n";
214 OOMPH_CURRENT_FUNCTION,
215 OOMPH_EXCEPTION_LOCATION);
242 Vector<std::complex<double>>& b)
const
249 for (
unsigned i = 0;
i < n;
i++)
251 b[
i] = std::complex<double>(0.0, 0.0);
256 (*Body_force_fct_pt)(x, b);
275 std::list<std::pair<unsigned long, unsigned>>& dof_lookup_list)
const
279 std::pair<unsigned long, unsigned> dof_lookup;
282 const unsigned n_node = this->
nnode();
285 int local_unknown = 0;
288 for (
unsigned n = 0; n < n_node; n++)
291 for (
unsigned i = 0;
i < 6;
i++)
297 if (local_unknown >= 0)
301 dof_lookup.first = this->
eqn_number(local_unknown);
302 dof_lookup.second = 0;
305 dof_lookup_list.push_front(dof_lookup);
380 ->fill_in_generic_contribution_to_residuals_fourier_decomp_time_harmonic_linear_elasticity(
381 residuals, jacobian, 1);
394 const unsigned& nplot,
405 void output(std::ostream& outfile,
const unsigned& n_plot);
415 void output(FILE* file_pt,
const unsigned& n_plot);
444 template<
unsigned NNODE_1D>
446 :
public virtual QElement<2, NNODE_1D>,
464 void output(std::ostream& outfile,
const unsigned& n_plot)
478 void output(FILE* file_pt,
const unsigned& n_plot)
490 template<
unsigned NNODE_1D>
493 :
public virtual QElement<1, NNODE_1D>
510 template<
unsigned NNODE_1D>
512 :
public virtual TElement<2, NNODE_1D>,
531 void output(std::ostream& outfile,
const unsigned& n_plot)
544 void output(FILE* file_pt,
const unsigned& n_plot)
582 unsigned num_entries = 12;
583 if (flux.size() != num_entries)
585 std::ostringstream error_message;
586 error_message <<
"The flux vector has the wrong number of entries, "
587 << flux.size() <<
", whereas it should be " << num_entries
590 OOMPH_CURRENT_FUNCTION,
591 OOMPH_EXCEPTION_LOCATION);
603 for (
unsigned i = 0;
i < 3;
i++)
605 flux[icount] = strain(
i,
i).real();
607 flux[icount] = strain(
i,
i).imag();
612 for (
unsigned i = 0;
i < 3;
i++)
614 for (
unsigned j =
i + 1; j < 3; j++)
616 flux[icount] = strain(
i, j).real();
618 flux[icount] = strain(
i, j).imag();
630 template<
unsigned NNODE_1D>
633 :
public virtual TElement<1, NNODE_1D>
650 template<
class TIME_HARMONIC_LINEAR_ELAST_ELEMENT>
671 unsigned nnod = this->
nnode();
672 for (
unsigned j = 0; j < nnod; j++)
675 data_values.push_back(std::make_pair(this->
node_pt(j), fld));
686 return 3 * this->
dim();
696 std::stringstream error_stream;
697 error_stream <<
"Elements only store six fields so fld can't be"
698 <<
" " << fld << std::endl;
700 error_stream.str(), OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
720 unsigned n_dim = this->
dim();
721 unsigned n_node = this->
nnode();
722 DShape dpsidx(n_node, n_dim);
736 unsigned n_node = this->
nnode();
749 double interpolated_u = 0.0;
752 for (
unsigned l = 0; l < n_node; l++)
756 if (nvalue != 3 * n_dim)
758 std::stringstream error_stream;
760 <<
"Current implementation only works for non-resized nodes\n"
761 <<
"but nvalue= " << nvalue <<
"!= 3 dim = " << 3 * n_dim
764 OOMPH_CURRENT_FUNCTION,
765 OOMPH_EXCEPTION_LOCATION);
768 interpolated_u += this->
nodal_value(t, l, fld) * psi[l];
770 return interpolated_u;
777 return this->
nnode();
787 if (nvalue != 3 * n_dim)
789 std::stringstream error_stream;
791 <<
"Current implementation only works for non-resized nodes\n"
792 <<
"but nvalue= " << nvalue <<
"!= 3 dim = " << 3 * n_dim
795 error_stream.str(), OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
807 template<
class ELEMENT>
821 template<
class ELEMENT>
A Class for the derivatives of shape functions The class design is essentially the same as Shape,...
unsigned nvalue() const
Return number of values stored in data object (incl pinned ones).
unsigned ntstorage() const
Return total number of doubles stored per value to record time history of each value (one for steady ...
Base class for finite elements that can compute the quantities that are required for the Z2 error est...
FaceGeometry()
Constructor must call the constructor of the underlying element.
FaceGeometry()
Constructor must call the constructor of the underlying element.
//////////////////////////////////////////////////////////////////// ////////////////////////////////...
A general Finite Element class.
Node *& node_pt(const unsigned &n)
Return a pointer to the local node n.
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...
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...
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.
unsigned dim() const
Return the spatial dimension of the element, i.e. the number of local coordinates required to paramet...
unsigned nnode() const
Return the number of nodes.
void(* SteadyExactSolutionFctPt)(const Vector< double > &, Vector< double > &)
Function pointer for function that computes vector-valued steady "exact solution" as .
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...
unsigned long eqn_number(const unsigned &ieqn_local) const
Return the global equation number corresponding to the ieqn_local-th local equation number.
static DenseMatrix< double > Dummy_matrix
Empty dense matrix used as a dummy argument to combined residual and jacobian functions in the case w...
Nodes are derived from Data, but, in addition, have a definite (Eulerian) position in a space of a gi...
unsigned ndim() const
Return (Eulerian) spatial dimension of the node.
TimeStepper *& position_time_stepper_pt()
Return a pointer to the position timestepper.
An OomphLibError object which should be thrown when an run-time error is encountered....
Wrapper class for projectable elements. Adds "projectability" to the underlying ELEMENT.
///////////////////////////////////////////////////////////////// ///////////////////////////////////...
unsigned nvalue_of_field(const unsigned &fld)
Return number of values in field fld.
unsigned nhistory_values_for_projection(const unsigned &fld)
Number of history values to be stored for fld-th field. (includes present value!)
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.
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...
unsigned nhistory_values_for_coordinate_projection()
Number of positional history values: Read out from positional timestepper (Note: count includes curre...
int local_equation(const unsigned &fld, const unsigned &j)
Return local equation number of value j in field fld.
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 ...
unsigned nfields_for_projection()
Number of fields to be projected: 3*dim, corresponding to real and imag parts of the displacement com...
ProjectableTimeHarmonicFourierDecomposedLinearElasticityElement()
Constructor [this was only required explicitly from gcc 4.5.2 onwards...].
/////////////////////////////////////////////////////////////////////// /////////////////////////////...
///////////////////////////////////////////////////////////////////// ///////////////////////////////...
void output(std::ostream &outfile, const unsigned &n_plot)
Output function.
QTimeHarmonicFourierDecomposedLinearElasticityElement()
Constructor.
void output(std::ostream &outfile)
Output function.
void output(FILE *file_pt, const unsigned &n_plot)
C-style output function.
void output(FILE *file_pt)
C-style output function.
A Class for shape functions. In simple cases, the shape functions have only one index that can be tho...
///////////////////////////////////////////////////////////////////// ///////////////////////////////...
void output(FILE *file_pt, const unsigned &n_plot)
C-style output function.
unsigned nvertex_node() const
Number of vertex nodes in the element.
unsigned nrecovery_order()
Order of recovery shape functions for Z2 error estimation: Same order as shape functions.
void output(FILE *file_pt)
C-style output function.
void output(std::ostream &outfile, const unsigned &n_plot)
Output function.
unsigned num_Z2_flux_terms()
Number of 'flux' terms for Z2 error estimation.
Node * vertex_node_pt(const unsigned &j) const
Pointer to the j-th vertex node in the element.
void get_Z2_flux(const Vector< double > &s, Vector< double > &flux)
Get 'flux' for Z2 error recovery: Upper triangular entries in strain tensor.
TTimeHarmonicFourierDecomposedLinearElasticityElement()
Constructor.
void output(std::ostream &outfile)
Output function.
A base class for elements that solve the Fourier decomposed (in cylindrical polars) equations of time...
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.
BodyForceFctPt & body_force_fct_pt()
Access function: Pointer to body force function.
std::complex< double > *& youngs_modulus_pt()
Return the pointer to Young's modulus.
static std::complex< double > Default_youngs_modulus_value
Static default value for Young's modulus (1.0 – for natural scaling, i.e. all stresses have been non-...
int * Fourier_wavenumber_pt
Pointer to Fourier wavenumber.
static std::complex< double > Default_omega_sq_value
Static default value for squared frequency.
std::complex< double > * Nu_pt
Pointer to Poisson's ratio.
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...
int & fourier_wavenumber() const
Access function for Fourier wavenumber.
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.
std::complex< double > & nu() const
Access function for Poisson's ratio.
BodyForceFctPt Body_force_fct_pt
Pointer to body force function.
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...
std::complex< double > * Youngs_modulus_pt
Pointer to the Young's modulus.
unsigned ndof_types() const
The number of "DOF types" that degrees of freedom in this element are sub-divided into: for now lump ...
std::complex< double > youngs_modulus() const
Access function to Young's modulus.
std::complex< double > *& nu_pt()
Access function for pointer to Poisson's ratio.
TimeHarmonicFourierDecomposedLinearElasticityEquationsBase()
Constructor: Set null pointers for constitutive law. Set physical parameter values to default values,...
std::complex< double > *& omega_sq_pt()
Access function for square of non-dim frequency.
std::complex< double > * Omega_sq_pt
Square of nondim frequency.
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...
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...
BodyForceFctPt body_force_fct_pt() const
Access function: Pointer to body force function (const version)
int *& fourier_wavenumber_pt()
Access function for pointer to Fourier wavenumber.
const std::complex< double > & omega_sq() const
Access function for square of non-dim frequency.
//////////////////////////////////////////////////////////////////// ////////////////////////////////...
TimeHarmonicFourierDecomposedLinearElasticityEquations()
Constructor.
void output(std::ostream &outfile)
Output: r,z, u_r_real, u_z_real, ..., u_theta_imag.
void fill_in_contribution_to_residuals(Vector< double > &residuals)
Return the residuals for the equations (the discretised principle of virtual displacements)
void compute_error(std::ostream &outfile, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt, double &error, double &norm)
Validate against exact solution. Solution is provided via function pointer. Plot at a given number of...
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.
void compute_norm(double &norm)
Compute norm of solution: square of the L2 norm.
unsigned required_nvalue(const unsigned &n) const
Number of values required at node n.
void get_strain(const Vector< double > &s, DenseMatrix< std::complex< double >> &strain)
Get strain tensor.
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 output(FILE *file_pt)
C-style output: r,z, u_r_real, u_z_real, ..., u_theta_imag.
void fill_in_contribution_to_jacobian(Vector< double > &residuals, DenseMatrix< double > &jacobian)
The jacobian is calculated by finite differences by default, We need only to take finite differences ...
unsigned ntstorage() const
Return the number of doubles required to represent history (one for steady)
//////////////////////////////////////////////////////////////////// ////////////////////////////////...