oomph::AxisymmetricLinearElasticityEquationsBase Class Reference

#include <axisym_linear_elasticity_elements.h>

+ Inheritance diagram for oomph::AxisymmetricLinearElasticityEquationsBase:

Public Types

typedef void(* BodyForceFctPt) (const double &time, const Vector< double > &x, Vector< double > &b)
 
- Public Types inherited from oomph::FiniteElement
typedef void(* SteadyExactSolutionFctPt) (const Vector< double > &, Vector< double > &)
 
typedef void(* UnsteadyExactSolutionFctPt) (const double &, const Vector< double > &, Vector< double > &)
 

Public Member Functions

virtual unsigned u_index_axisymmetric_linear_elasticity (const unsigned &i) const
 
double d2u_dt2_axisymmetric_linear_elasticity (const unsigned &n, const unsigned &i) const
 d^2u/dt^2 at local node n More...
 
double du_dt_axisymmetric_linear_elasticity (const unsigned &n, const unsigned &i) const
 du/dt at local node n More...
 
void interpolated_u_axisymmetric_linear_elasticity (const Vector< double > &s, Vector< double > &disp) const
 Compute vector of FE interpolated displacement u at local coordinate s. More...
 
double interpolated_u_axisymmetric_linear_elasticity (const Vector< double > &s, const unsigned &i) const
 
void interpolated_du_dt_axisymmetric_linear_elasticity (const Vector< double > &s, Vector< double > &du_dt) const
 Compute vector of FE interpolated velocity du/dt at local coordinate s. More...
 
void interpolated_d2u_dt2_axisymmetric_linear_elasticity (const Vector< double > &s, Vector< double > &d2u_dt2) const
 Compute vector of FE interpolated accel d2u/dt2 at local coordinate s. More...
 
 AxisymmetricLinearElasticityEquationsBase ()
 
double *& youngs_modulus_pt ()
 Return the pointer to Young's modulus. More...
 
double youngs_modulus () const
 Access function to Young's modulus. More...
 
doublenu () const
 Access function for Poisson's ratio. More...
 
double *& nu_pt ()
 Access function for pointer to Poisson's ratio. More...
 
double *& lambda_sq_pt ()
 Access function for pointer to timescale ratio (nondim density) More...
 
const doublelambda_sq () const
 Access function for timescale ratio (nondim density) More...
 
BodyForceFctPtbody_force_fct_pt ()
 Access function: Pointer to body force function. More...
 
BodyForceFctPt body_force_fct_pt () const
 Access function: Pointer to body force function (const version) More...
 
void body_force (const double &time, const Vector< double > &x, Vector< double > &b) const
 
unsigned ndof_types () const
 
void get_dof_numbers_for_unknowns (std::list< std::pair< unsigned long, unsigned >> &dof_lookup_list) const
 
- Public Member Functions inherited from oomph::FiniteElement
void set_dimension (const unsigned &dim)
 
void set_nodal_dimension (const unsigned &nodal_dim)
 
void set_nnodal_position_type (const unsigned &nposition_type)
 Set the number of types required to interpolate the coordinate. More...
 
void set_n_node (const unsigned &n)
 
int nodal_local_eqn (const unsigned &n, const unsigned &i) const
 
double dJ_eulerian_at_knot (const unsigned &ipt, Shape &psi, DenseMatrix< double > &djacobian_dX) const
 
 FiniteElement ()
 Constructor. More...
 
virtual ~FiniteElement ()
 
 FiniteElement (const FiniteElement &)=delete
 Broken copy constructor. More...
 
virtual bool local_coord_is_valid (const Vector< double > &s)
 Broken assignment operator. More...
 
virtual void move_local_coord_back_into_element (Vector< double > &s) const
 
void get_centre_of_gravity_and_max_radius_in_terms_of_zeta (Vector< double > &cog, double &max_radius) const
 
virtual void local_coordinate_of_node (const unsigned &j, Vector< double > &s) const
 
virtual void local_fraction_of_node (const unsigned &j, Vector< double > &s_fraction)
 
virtual double local_one_d_fraction_of_node (const unsigned &n1d, const unsigned &i)
 
virtual void set_macro_elem_pt (MacroElement *macro_elem_pt)
 
MacroElementmacro_elem_pt ()
 Access function to pointer to macro element. More...
 
void get_x (const Vector< double > &s, Vector< double > &x) const
 
void get_x (const unsigned &t, const Vector< double > &s, Vector< double > &x)
 
virtual void get_x_from_macro_element (const Vector< double > &s, Vector< double > &x) const
 
virtual void get_x_from_macro_element (const unsigned &t, const Vector< double > &s, Vector< double > &x)
 
virtual void set_integration_scheme (Integral *const &integral_pt)
 Set the spatial integration scheme. More...
 
Integral *const & integral_pt () const
 Return the pointer to the integration scheme (const version) More...
 
virtual void shape (const Vector< double > &s, Shape &psi) const =0
 
virtual void shape_at_knot (const unsigned &ipt, Shape &psi) const
 
virtual void dshape_local (const Vector< double > &s, Shape &psi, DShape &dpsids) const
 
virtual void dshape_local_at_knot (const unsigned &ipt, Shape &psi, DShape &dpsids) const
 
virtual void d2shape_local (const Vector< double > &s, Shape &psi, DShape &dpsids, DShape &d2psids) const
 
virtual void d2shape_local_at_knot (const unsigned &ipt, Shape &psi, DShape &dpsids, DShape &d2psids) const
 
virtual double J_eulerian (const Vector< double > &s) const
 
virtual double J_eulerian_at_knot (const unsigned &ipt) const
 
void check_J_eulerian_at_knots (bool &passed) const
 
void check_jacobian (const double &jacobian) const
 
double dshape_eulerian (const Vector< double > &s, Shape &psi, DShape &dpsidx) const
 
virtual double dshape_eulerian_at_knot (const unsigned &ipt, Shape &psi, DShape &dpsidx) const
 
virtual double dshape_eulerian_at_knot (const unsigned &ipt, Shape &psi, DShape &dpsi, DenseMatrix< double > &djacobian_dX, RankFourTensor< double > &d_dpsidx_dX) const
 
double d2shape_eulerian (const Vector< double > &s, Shape &psi, DShape &dpsidx, DShape &d2psidx) const
 
virtual double d2shape_eulerian_at_knot (const unsigned &ipt, Shape &psi, DShape &dpsidx, DShape &d2psidx) const
 
virtual void assign_nodal_local_eqn_numbers (const bool &store_local_dof_pt)
 
virtual void describe_local_dofs (std::ostream &out, const std::string &current_string) const
 
virtual void describe_nodal_local_dofs (std::ostream &out, const std::string &current_string) const
 
virtual void assign_all_generic_local_eqn_numbers (const bool &store_local_dof_pt)
 
Node *& node_pt (const unsigned &n)
 Return a pointer to the local node n. More...
 
Node *const & node_pt (const unsigned &n) const
 Return a pointer to the local node n (const version) More...
 
unsigned nnode () const
 Return the number of nodes. More...
 
virtual unsigned nnode_1d () const
 
double raw_nodal_position (const unsigned &n, const unsigned &i) const
 
double raw_nodal_position (const unsigned &t, const unsigned &n, const unsigned &i) const
 
double raw_dnodal_position_dt (const unsigned &n, const unsigned &i) const
 
double raw_dnodal_position_dt (const unsigned &n, const unsigned &j, const unsigned &i) const
 
double raw_nodal_position_gen (const unsigned &n, const unsigned &k, const unsigned &i) const
 
double raw_nodal_position_gen (const unsigned &t, const unsigned &n, const unsigned &k, const unsigned &i) const
 
double raw_dnodal_position_gen_dt (const unsigned &n, const unsigned &k, const unsigned &i) const
 
double raw_dnodal_position_gen_dt (const unsigned &j, const unsigned &n, const unsigned &k, const unsigned &i) const
 
double nodal_position (const unsigned &n, const unsigned &i) const
 
double nodal_position (const unsigned &t, const unsigned &n, const unsigned &i) const
 
double dnodal_position_dt (const unsigned &n, const unsigned &i) const
 Return the i-th component of nodal velocity: dx/dt at local node n. More...
 
double dnodal_position_dt (const unsigned &n, const unsigned &j, const unsigned &i) const
 
double nodal_position_gen (const unsigned &n, const unsigned &k, const unsigned &i) const
 
double nodal_position_gen (const unsigned &t, const unsigned &n, const unsigned &k, const unsigned &i) const
 
double dnodal_position_gen_dt (const unsigned &n, const unsigned &k, const unsigned &i) const
 
double dnodal_position_gen_dt (const unsigned &j, const unsigned &n, const unsigned &k, const unsigned &i) const
 
virtual void get_dresidual_dnodal_coordinates (RankThreeTensor< double > &dresidual_dnodal_coordinates)
 
virtual void disable_ALE ()
 
virtual void enable_ALE ()
 
virtual unsigned required_nvalue (const unsigned &n) const
 
unsigned nnodal_position_type () const
 
bool has_hanging_nodes () const
 
unsigned nodal_dimension () const
 Return the required Eulerian dimension of the nodes in this element. More...
 
virtual unsigned nvertex_node () const
 
virtual Nodevertex_node_pt (const unsigned &j) const
 
virtual Nodeconstruct_node (const unsigned &n)
 
virtual Nodeconstruct_node (const unsigned &n, TimeStepper *const &time_stepper_pt)
 
virtual Nodeconstruct_boundary_node (const unsigned &n)
 
virtual Nodeconstruct_boundary_node (const unsigned &n, TimeStepper *const &time_stepper_pt)
 
int get_node_number (Node *const &node_pt) const
 
virtual Nodeget_node_at_local_coordinate (const Vector< double > &s) const
 
double raw_nodal_value (const unsigned &n, const unsigned &i) const
 
double raw_nodal_value (const unsigned &t, const unsigned &n, const unsigned &i) const
 
double nodal_value (const unsigned &n, const unsigned &i) const
 
double nodal_value (const unsigned &t, const unsigned &n, const unsigned &i) const
 
unsigned dim () const
 
virtual ElementGeometry::ElementGeometry element_geometry () const
 Return the geometry type of the element (either Q or T usually). More...
 
virtual double interpolated_x (const Vector< double > &s, const unsigned &i) const
 Return FE interpolated coordinate x[i] at local coordinate s. More...
 
virtual double interpolated_x (const unsigned &t, const Vector< double > &s, const unsigned &i) const
 
virtual void interpolated_x (const Vector< double > &s, Vector< double > &x) const
 Return FE interpolated position x[] at local coordinate s as Vector. More...
 
virtual void interpolated_x (const unsigned &t, const Vector< double > &s, Vector< double > &x) const
 
virtual double interpolated_dxdt (const Vector< double > &s, const unsigned &i, const unsigned &t)
 
virtual void interpolated_dxdt (const Vector< double > &s, const unsigned &t, Vector< double > &dxdt)
 
unsigned ngeom_data () const
 
Datageom_data_pt (const unsigned &j)
 
void position (const Vector< double > &zeta, Vector< double > &r) const
 
void position (const unsigned &t, const Vector< double > &zeta, Vector< double > &r) const
 
void dposition_dt (const Vector< double > &zeta, const unsigned &t, Vector< double > &drdt)
 
virtual double zeta_nodal (const unsigned &n, const unsigned &k, const unsigned &i) const
 
void interpolated_zeta (const Vector< double > &s, Vector< double > &zeta) const
 
void locate_zeta (const Vector< double > &zeta, GeomObject *&geom_object_pt, Vector< double > &s, const bool &use_coordinate_as_initial_guess=false)
 
virtual void node_update ()
 
virtual void identify_field_data_for_interactions (std::set< std::pair< Data *, unsigned >> &paired_field_data)
 
virtual void identify_geometric_data (std::set< Data * > &geometric_data_pt)
 
virtual double s_min () const
 Min value of local coordinate. More...
 
virtual double s_max () const
 Max. value of local coordinate. More...
 
double size () const
 
virtual double compute_physical_size () const
 
virtual void point_output_data (const Vector< double > &s, Vector< double > &data)
 
void point_output (std::ostream &outfile, const Vector< double > &s)
 
virtual unsigned nplot_points_paraview (const unsigned &nplot) const
 
virtual unsigned nsub_elements_paraview (const unsigned &nplot) const
 
void output_paraview (std::ofstream &file_out, const unsigned &nplot) const
 
virtual void write_paraview_output_offset_information (std::ofstream &file_out, const unsigned &nplot, unsigned &counter) const
 
virtual void write_paraview_type (std::ofstream &file_out, const unsigned &nplot) const
 
virtual void write_paraview_offsets (std::ofstream &file_out, const unsigned &nplot, unsigned &offset_sum) const
 
virtual unsigned nscalar_paraview () const
 
virtual void scalar_value_paraview (std::ofstream &file_out, const unsigned &i, const unsigned &nplot) const
 
virtual void scalar_value_fct_paraview (std::ofstream &file_out, const unsigned &i, const unsigned &nplot, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt) const
 
virtual void scalar_value_fct_paraview (std::ofstream &file_out, const unsigned &i, const unsigned &nplot, const double &time, FiniteElement::UnsteadyExactSolutionFctPt exact_soln_pt) const
 
virtual std::string scalar_name_paraview (const unsigned &i) const
 
virtual void output (std::ostream &outfile)
 
virtual void output (std::ostream &outfile, const unsigned &n_plot)
 
virtual void output (const unsigned &t, std::ostream &outfile, const unsigned &n_plot) const
 
virtual void output (FILE *file_pt)
 
virtual void output (FILE *file_pt, const unsigned &n_plot)
 
virtual void output_fct (std::ostream &outfile, const unsigned &n_plot, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt)
 Output an exact solution over the element. More...
 
virtual void output_fct (std::ostream &outfile, const unsigned &n_plot, const double &time, FiniteElement::UnsteadyExactSolutionFctPt exact_soln_pt)
 Output a time-dependent exact solution over the element. More...
 
virtual void output_fct (std::ostream &outfile, const unsigned &n_plot, const double &time, const SolutionFunctorBase &exact_soln) const
 Output a time-dependent exact solution over the element. More...
 
virtual void get_s_plot (const unsigned &i, const unsigned &nplot, Vector< double > &s, const bool &shifted_to_interior=false) const
 
virtual std::string tecplot_zone_string (const unsigned &nplot) const
 
virtual void write_tecplot_zone_footer (std::ostream &outfile, const unsigned &nplot) const
 
virtual void write_tecplot_zone_footer (FILE *file_pt, const unsigned &nplot) const
 
virtual unsigned nplot_points (const unsigned &nplot) const
 
virtual void compute_error (FiniteElement::SteadyExactSolutionFctPt exact_soln_pt, double &error, double &norm)
 Calculate the norm of the error and that of the exact solution. More...
 
virtual void compute_error (FiniteElement::UnsteadyExactSolutionFctPt exact_soln_pt, const double &time, double &error, double &norm)
 Calculate the norm of the error and that of the exact solution. More...
 
virtual void compute_error (FiniteElement::SteadyExactSolutionFctPt exact_soln_pt, Vector< double > &error, Vector< double > &norm)
 
virtual void compute_error (FiniteElement::UnsteadyExactSolutionFctPt exact_soln_pt, const double &time, Vector< double > &error, Vector< double > &norm)
 
virtual void compute_error (std::ostream &outfile, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt, double &error, double &norm)
 
virtual void compute_error (std::ostream &outfile, FiniteElement::UnsteadyExactSolutionFctPt exact_soln_pt, const double &time, double &error, double &norm)
 
virtual void compute_error (std::ostream &outfile, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt, Vector< double > &error, Vector< double > &norm)
 
virtual void compute_error (std::ostream &outfile, FiniteElement::UnsteadyExactSolutionFctPt exact_soln_pt, const double &time, Vector< double > &error, Vector< double > &norm)
 
virtual void compute_abs_error (std::ostream &outfile, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt, double &error)
 
void integrate_fct (FiniteElement::SteadyExactSolutionFctPt integrand_fct_pt, Vector< double > &integral)
 Integrate Vector-valued function over element. More...
 
void integrate_fct (FiniteElement::UnsteadyExactSolutionFctPt integrand_fct_pt, const double &time, Vector< double > &integral)
 Integrate Vector-valued time-dep function over element. More...
 
virtual void build_face_element (const int &face_index, FaceElement *face_element_pt)
 
virtual unsigned self_test ()
 
virtual unsigned get_bulk_node_number (const int &face_index, const unsigned &i) const
 
virtual int face_outer_unit_normal_sign (const int &face_index) const
 Get the sign of the outer unit normal on the face given by face_index. More...
 
virtual unsigned nnode_on_face () const
 
void face_node_number_error_check (const unsigned &i) const
 Range check for face node numbers. More...
 
virtual CoordinateMappingFctPt face_to_bulk_coordinate_fct_pt (const int &face_index) const
 
virtual BulkCoordinateDerivativesFctPt bulk_coordinate_derivatives_fct_pt (const int &face_index) const
 
- Public Member Functions inherited from oomph::GeneralisedElement
 GeneralisedElement ()
 Constructor: Initialise all pointers and all values to zero. More...
 
virtual ~GeneralisedElement ()
 Virtual destructor to clean up any memory allocated by the object. More...
 
 GeneralisedElement (const GeneralisedElement &)=delete
 Broken copy constructor. More...
 
void operator= (const GeneralisedElement &)=delete
 Broken assignment operator. More...
 
Data *& internal_data_pt (const unsigned &i)
 Return a pointer to i-th internal data object. More...
 
Data *const & internal_data_pt (const unsigned &i) const
 Return a pointer to i-th internal data object (const version) More...
 
Data *& external_data_pt (const unsigned &i)
 Return a pointer to i-th external data object. More...
 
Data *const & external_data_pt (const unsigned &i) const
 Return a pointer to i-th external data object (const version) More...
 
unsigned long eqn_number (const unsigned &ieqn_local) const
 
int local_eqn_number (const unsigned long &ieqn_global) const
 
unsigned add_external_data (Data *const &data_pt, const bool &fd=true)
 
bool external_data_fd (const unsigned &i) const
 
void exclude_external_data_fd (const unsigned &i)
 
void include_external_data_fd (const unsigned &i)
 
void flush_external_data ()
 Flush all external data. More...
 
void flush_external_data (Data *const &data_pt)
 Flush the object addressed by data_pt from the external data array. More...
 
unsigned ninternal_data () const
 Return the number of internal data objects. More...
 
unsigned nexternal_data () const
 Return the number of external data objects. More...
 
unsigned ndof () const
 Return the number of equations/dofs in the element. More...
 
void dof_vector (const unsigned &t, Vector< double > &dof)
 Return the vector of dof values at time level t. More...
 
void dof_pt_vector (Vector< double * > &dof_pt)
 Return the vector of pointers to dof values. More...
 
void set_internal_data_time_stepper (const unsigned &i, TimeStepper *const &time_stepper_pt, const bool &preserve_existing_data)
 
void assign_internal_eqn_numbers (unsigned long &global_number, Vector< double * > &Dof_pt)
 
void describe_dofs (std::ostream &out, const std::string &current_string) const
 
void add_internal_value_pt_to_map (std::map< unsigned, double * > &map_of_value_pt)
 
virtual void assign_local_eqn_numbers (const bool &store_local_dof_pt)
 
virtual void complete_setup_of_dependencies ()
 
virtual void get_residuals (Vector< double > &residuals)
 
virtual void get_jacobian (Vector< double > &residuals, DenseMatrix< double > &jacobian)
 
virtual void get_mass_matrix (Vector< double > &residuals, DenseMatrix< double > &mass_matrix)
 
virtual void get_jacobian_and_mass_matrix (Vector< double > &residuals, DenseMatrix< double > &jacobian, DenseMatrix< double > &mass_matrix)
 
virtual void get_dresiduals_dparameter (double *const &parameter_pt, Vector< double > &dres_dparam)
 
virtual void get_djacobian_dparameter (double *const &parameter_pt, Vector< double > &dres_dparam, DenseMatrix< double > &djac_dparam)
 
virtual void get_djacobian_and_dmass_matrix_dparameter (double *const &parameter_pt, Vector< double > &dres_dparam, DenseMatrix< double > &djac_dparam, DenseMatrix< double > &dmass_matrix_dparam)
 
virtual void get_hessian_vector_products (Vector< double > const &Y, DenseMatrix< double > const &C, DenseMatrix< double > &product)
 
virtual void get_inner_products (Vector< std::pair< unsigned, unsigned >> const &history_index, Vector< double > &inner_product)
 
virtual void get_inner_product_vectors (Vector< unsigned > const &history_index, Vector< Vector< double >> &inner_product_vector)
 
virtual void compute_norm (Vector< double > &norm)
 
virtual void compute_norm (double &norm)
 
- Public Member Functions inherited from oomph::GeomObject
 GeomObject ()
 Default constructor. More...
 
 GeomObject (const unsigned &ndim)
 
 GeomObject (const unsigned &nlagrangian, const unsigned &ndim)
 
 GeomObject (const unsigned &nlagrangian, const unsigned &ndim, TimeStepper *time_stepper_pt)
 
 GeomObject (const GeomObject &dummy)=delete
 Broken copy constructor. More...
 
void operator= (const GeomObject &)=delete
 Broken assignment operator. More...
 
virtual ~GeomObject ()
 (Empty) destructor More...
 
unsigned nlagrangian () const
 Access function to # of Lagrangian coordinates. More...
 
unsigned ndim () const
 Access function to # of Eulerian coordinates. More...
 
void set_nlagrangian_and_ndim (const unsigned &n_lagrangian, const unsigned &n_dim)
 Set # of Lagrangian and Eulerian coordinates. More...
 
TimeStepper *& time_stepper_pt ()
 
TimeSteppertime_stepper_pt () const
 
virtual void position (const double &t, const Vector< double > &zeta, Vector< double > &r) const
 
virtual void dposition (const Vector< double > &zeta, DenseMatrix< double > &drdzeta) const
 
virtual void d2position (const Vector< double > &zeta, RankThreeTensor< double > &ddrdzeta) const
 
virtual void d2position (const Vector< double > &zeta, Vector< double > &r, DenseMatrix< double > &drdzeta, RankThreeTensor< double > &ddrdzeta) const
 

Protected Attributes

doubleYoungs_modulus_pt
 Pointer to the Young's modulus. More...
 
doubleNu_pt
 Pointer to Poisson's ratio. More...
 
doubleLambda_sq_pt
 Timescale ratio (non-dim. density) More...
 
BodyForceFctPt Body_force_fct_pt
 Pointer to body force function. More...
 
- Protected Attributes inherited from oomph::FiniteElement
MacroElementMacro_elem_pt
 Pointer to the element's macro element (NULL by default) More...
 
- Protected Attributes inherited from oomph::GeomObject
unsigned NLagrangian
 Number of Lagrangian (intrinsic) coordinates. More...
 
unsigned Ndim
 Number of Eulerian coordinates. More...
 
TimeStepperGeom_object_time_stepper_pt
 

Static Protected Attributes

static double Default_youngs_modulus_value
 
static double Default_lambda_sq_value
 Static default value for timescale ratio (1.0 for natural scaling) More...
 
- Static Protected Attributes inherited from oomph::FiniteElement
static const unsigned Default_Initial_Nvalue = 0
 Default value for the number of values at a node. More...
 
static const double Node_location_tolerance = 1.0e-14
 
static const unsigned N2deriv [] = {0, 1, 3, 6}
 
- Static Protected Attributes inherited from oomph::GeneralisedElement
static DenseMatrix< doubleDummy_matrix
 
static std::deque< double * > Dof_pt_deque
 

Additional Inherited Members

- Static Public Attributes inherited from oomph::FiniteElement
static double Tolerance_for_singular_jacobian = 1.0e-16
 Tolerance below which the jacobian is considered singular. More...
 
static bool Accept_negative_jacobian = false
 
static bool Suppress_output_while_checking_for_inverted_elements
 
- Static Public Attributes inherited from oomph::GeneralisedElement
static bool Suppress_warning_about_repeated_internal_data
 
static bool Suppress_warning_about_repeated_external_data = true
 
static double Default_fd_jacobian_step = 1.0e-8
 
- Protected Member Functions inherited from oomph::FiniteElement
virtual void assemble_local_to_eulerian_jacobian (const DShape &dpsids, DenseMatrix< double > &jacobian) const
 
virtual void assemble_local_to_eulerian_jacobian2 (const DShape &d2psids, DenseMatrix< double > &jacobian2) const
 
virtual void assemble_eulerian_base_vectors (const DShape &dpsids, DenseMatrix< double > &interpolated_G) const
 
template<unsigned DIM>
double invert_jacobian (const DenseMatrix< double > &jacobian, DenseMatrix< double > &inverse_jacobian) const
 
virtual double invert_jacobian_mapping (const DenseMatrix< double > &jacobian, DenseMatrix< double > &inverse_jacobian) const
 
virtual double local_to_eulerian_mapping (const DShape &dpsids, DenseMatrix< double > &jacobian, DenseMatrix< double > &inverse_jacobian) const
 
double local_to_eulerian_mapping (const DShape &dpsids, DenseMatrix< double > &inverse_jacobian) const
 
virtual double local_to_eulerian_mapping_diagonal (const DShape &dpsids, DenseMatrix< double > &jacobian, DenseMatrix< double > &inverse_jacobian) const
 
virtual void dJ_eulerian_dnodal_coordinates (const DenseMatrix< double > &jacobian, const DShape &dpsids, DenseMatrix< double > &djacobian_dX) const
 
template<unsigned DIM>
void dJ_eulerian_dnodal_coordinates_templated_helper (const DenseMatrix< double > &jacobian, const DShape &dpsids, DenseMatrix< double > &djacobian_dX) const
 
virtual void d_dshape_eulerian_dnodal_coordinates (const double &det_jacobian, const DenseMatrix< double > &jacobian, const DenseMatrix< double > &djacobian_dX, const DenseMatrix< double > &inverse_jacobian, const DShape &dpsids, RankFourTensor< double > &d_dpsidx_dX) const
 
template<unsigned DIM>
void d_dshape_eulerian_dnodal_coordinates_templated_helper (const double &det_jacobian, const DenseMatrix< double > &jacobian, const DenseMatrix< double > &djacobian_dX, const DenseMatrix< double > &inverse_jacobian, const DShape &dpsids, RankFourTensor< double > &d_dpsidx_dX) const
 
virtual void transform_derivatives (const DenseMatrix< double > &inverse_jacobian, DShape &dbasis) const
 
void transform_derivatives_diagonal (const DenseMatrix< double > &inverse_jacobian, DShape &dbasis) const
 
virtual void transform_second_derivatives (const DenseMatrix< double > &jacobian, const DenseMatrix< double > &inverse_jacobian, const DenseMatrix< double > &jacobian2, DShape &dbasis, DShape &d2basis) const
 
template<unsigned DIM>
void transform_second_derivatives_template (const DenseMatrix< double > &jacobian, const DenseMatrix< double > &inverse_jacobian, const DenseMatrix< double > &jacobian2, DShape &dbasis, DShape &d2basis) const
 
template<unsigned DIM>
void transform_second_derivatives_diagonal (const DenseMatrix< double > &jacobian, const DenseMatrix< double > &inverse_jacobian, const DenseMatrix< double > &jacobian2, DShape &dbasis, DShape &d2basis) const
 
virtual void fill_in_jacobian_from_nodal_by_fd (Vector< double > &residuals, DenseMatrix< double > &jacobian)
 
void fill_in_jacobian_from_nodal_by_fd (DenseMatrix< double > &jacobian)
 
virtual void update_before_nodal_fd ()
 
virtual void reset_after_nodal_fd ()
 
virtual void update_in_nodal_fd (const unsigned &i)
 
virtual void reset_in_nodal_fd (const unsigned &i)
 
void fill_in_contribution_to_jacobian (Vector< double > &residuals, DenseMatrix< double > &jacobian)
 
template<>
double invert_jacobian (const DenseMatrix< double > &jacobian, DenseMatrix< double > &inverse_jacobian) const
 Zero-d specialisation of function to calculate inverse of jacobian mapping. More...
 
template<>
double invert_jacobian (const DenseMatrix< double > &jacobian, DenseMatrix< double > &inverse_jacobian) const
 One-d specialisation of function to calculate inverse of jacobian mapping. More...
 
template<>
double invert_jacobian (const DenseMatrix< double > &jacobian, DenseMatrix< double > &inverse_jacobian) const
 Two-d specialisation of function to calculate inverse of jacobian mapping. More...
 
template<>
double invert_jacobian (const DenseMatrix< double > &jacobian, DenseMatrix< double > &inverse_jacobian) const
 
template<>
void dJ_eulerian_dnodal_coordinates_templated_helper (const DenseMatrix< double > &jacobian, const DShape &dpsids, DenseMatrix< double > &djacobian_dX) const
 
template<>
void dJ_eulerian_dnodal_coordinates_templated_helper (const DenseMatrix< double > &jacobian, const DShape &dpsids, DenseMatrix< double > &djacobian_dX) const
 
template<>
void dJ_eulerian_dnodal_coordinates_templated_helper (const DenseMatrix< double > &jacobian, const DShape &dpsids, DenseMatrix< double > &djacobian_dX) const
 
template<>
void dJ_eulerian_dnodal_coordinates_templated_helper (const DenseMatrix< double > &jacobian, const DShape &dpsids, DenseMatrix< double > &djacobian_dX) const
 
template<>
void d_dshape_eulerian_dnodal_coordinates_templated_helper (const double &det_jacobian, const DenseMatrix< double > &jacobian, const DenseMatrix< double > &djacobian_dX, const DenseMatrix< double > &inverse_jacobian, const DShape &dpsids, RankFourTensor< double > &d_dpsidx_dX) const
 
template<>
void d_dshape_eulerian_dnodal_coordinates_templated_helper (const double &det_jacobian, const DenseMatrix< double > &jacobian, const DenseMatrix< double > &djacobian_dX, const DenseMatrix< double > &inverse_jacobian, const DShape &dpsids, RankFourTensor< double > &d_dpsidx_dX) const
 
template<>
void d_dshape_eulerian_dnodal_coordinates_templated_helper (const double &det_jacobian, const DenseMatrix< double > &jacobian, const DenseMatrix< double > &djacobian_dX, const DenseMatrix< double > &inverse_jacobian, const DShape &dpsids, RankFourTensor< double > &d_dpsidx_dX) const
 
template<>
void d_dshape_eulerian_dnodal_coordinates_templated_helper (const double &det_jacobian, const DenseMatrix< double > &jacobian, const DenseMatrix< double > &djacobian_dX, const DenseMatrix< double > &inverse_jacobian, const DShape &dpsids, RankFourTensor< double > &d_dpsidx_dX) const
 
template<>
void transform_second_derivatives_template (const DenseMatrix< double > &jacobian, const DenseMatrix< double > &inverse_jacobian, const DenseMatrix< double > &jacobian2, DShape &dbasis, DShape &d2basis) const
 
template<>
void transform_second_derivatives_template (const DenseMatrix< double > &jacobian, const DenseMatrix< double > &inverse_jacobian, const DenseMatrix< double > &jacobian2, DShape &dbasis, DShape &d2basis) const
 
template<>
void transform_second_derivatives_diagonal (const DenseMatrix< double > &jacobian, const DenseMatrix< double > &inverse_jacobian, const DenseMatrix< double > &jacobian2, DShape &dbasis, DShape &d2basis) const
 
template<>
void transform_second_derivatives_diagonal (const DenseMatrix< double > &jacobian, const DenseMatrix< double > &inverse_jacobian, const DenseMatrix< double > &jacobian2, DShape &dbasis, DShape &d2basis) const
 
- Protected Member Functions inherited from oomph::GeneralisedElement
unsigned add_internal_data (Data *const &data_pt, const bool &fd=true)
 
bool internal_data_fd (const unsigned &i) const
 
void exclude_internal_data_fd (const unsigned &i)
 
void include_internal_data_fd (const unsigned &i)
 
void clear_global_eqn_numbers ()
 
void add_global_eqn_numbers (std::deque< unsigned long > const &global_eqn_numbers, std::deque< double * > const &global_dof_pt)
 
virtual void assign_internal_and_external_local_eqn_numbers (const bool &store_local_dof_pt)
 
virtual void assign_additional_local_eqn_numbers ()
 
int internal_local_eqn (const unsigned &i, const unsigned &j) const
 
int external_local_eqn (const unsigned &i, const unsigned &j)
 
virtual void fill_in_contribution_to_residuals (Vector< double > &residuals)
 
void fill_in_jacobian_from_internal_by_fd (Vector< double > &residuals, DenseMatrix< double > &jacobian, const bool &fd_all_data=false)
 
void fill_in_jacobian_from_internal_by_fd (DenseMatrix< double > &jacobian, const bool &fd_all_data=false)
 
void fill_in_jacobian_from_external_by_fd (Vector< double > &residuals, DenseMatrix< double > &jacobian, const bool &fd_all_data=false)
 
void fill_in_jacobian_from_external_by_fd (DenseMatrix< double > &jacobian, const bool &fd_all_data=false)
 
virtual void update_before_internal_fd ()
 
virtual void reset_after_internal_fd ()
 
virtual void update_in_internal_fd (const unsigned &i)
 
virtual void reset_in_internal_fd (const unsigned &i)
 
virtual void update_before_external_fd ()
 
virtual void reset_after_external_fd ()
 
virtual void update_in_external_fd (const unsigned &i)
 
virtual void reset_in_external_fd (const unsigned &i)
 
virtual void fill_in_contribution_to_mass_matrix (Vector< double > &residuals, DenseMatrix< double > &mass_matrix)
 
virtual void fill_in_contribution_to_jacobian_and_mass_matrix (Vector< double > &residuals, DenseMatrix< double > &jacobian, DenseMatrix< double > &mass_matrix)
 
virtual void fill_in_contribution_to_dresiduals_dparameter (double *const &parameter_pt, Vector< double > &dres_dparam)
 
virtual void fill_in_contribution_to_djacobian_dparameter (double *const &parameter_pt, Vector< double > &dres_dparam, DenseMatrix< double > &djac_dparam)
 
virtual void fill_in_contribution_to_djacobian_and_dmass_matrix_dparameter (double *const &parameter_pt, Vector< double > &dres_dparam, DenseMatrix< double > &djac_dparam, DenseMatrix< double > &dmass_matrix_dparam)
 
virtual void fill_in_contribution_to_hessian_vector_products (Vector< double > const &Y, DenseMatrix< double > const &C, DenseMatrix< double > &product)
 
virtual void fill_in_contribution_to_inner_products (Vector< std::pair< unsigned, unsigned >> const &history_index, Vector< double > &inner_product)
 
virtual void fill_in_contribution_to_inner_product_vectors (Vector< unsigned > const &history_index, Vector< Vector< double >> &inner_product_vector)
 

Detailed Description

A base class for elements that solve the axisymmetric (in cylindrical polars) equations of linear elasticity.

Member Typedef Documentation

◆ BodyForceFctPt

typedef void(* oomph::AxisymmetricLinearElasticityEquationsBase::BodyForceFctPt) (const double &time, const Vector< double > &x, Vector< double > &b)

Function pointer to function that specifies the body force as a function of the Cartesian coordinates and time FCT(x,b) – x and b are Vectors!

Constructor & Destructor Documentation

◆ AxisymmetricLinearElasticityEquationsBase()

oomph::AxisymmetricLinearElasticityEquationsBase::AxisymmetricLinearElasticityEquationsBase ( )
inline

Constructor: Set null pointers for constitutive law. Set physical parameter values to default values, and set body force to zero.

259  Nu_pt(0),
262  {
263  }
double * Youngs_modulus_pt
Pointer to the Young's modulus.
Definition: axisym_linear_elasticity_elements.h:399
double * Nu_pt
Pointer to Poisson's ratio.
Definition: axisym_linear_elasticity_elements.h:402
static double Default_lambda_sq_value
Static default value for timescale ratio (1.0 for natural scaling)
Definition: axisym_linear_elasticity_elements.h:419
static double Default_youngs_modulus_value
Definition: axisym_linear_elasticity_elements.h:416
BodyForceFctPt Body_force_fct_pt
Pointer to body force function.
Definition: axisym_linear_elasticity_elements.h:408
double * Lambda_sq_pt
Timescale ratio (non-dim. density)
Definition: axisym_linear_elasticity_elements.h:405

Member Function Documentation

◆ body_force()

void oomph::AxisymmetricLinearElasticityEquationsBase::body_force ( const double time,
const Vector< double > &  x,
Vector< double > &  b 
) const
inline

Evaluate body force at Eulerian coordinate x at present time (returns zero vector if no body force function pointer has been set)

328  {
329  // If no function has been set, return zero vector
330  if (Body_force_fct_pt == 0)
331  {
332  // Get spatial dimension of element
333  unsigned n = dim();
334  for (unsigned i = 0; i < n; i++)
335  {
336  b[i] = 0.0;
337  }
338  }
339  else
340  {
341  (*Body_force_fct_pt)(time, x, b);
342  }
343  }
int i
Definition: BiCGSTAB_step_by_step.cpp:9
const unsigned n
Definition: CG3DPackingUnitTest.cpp:11
Scalar * b
Definition: benchVecAdd.cpp:17
unsigned dim() const
Definition: elements.h:2611
list x
Definition: plotDoE.py:28

References b, Body_force_fct_pt, oomph::FiniteElement::dim(), i, n, and plotDoE::x.

Referenced by oomph::AxisymmetricLinearElasticityEquations::fill_in_generic_contribution_to_residuals_axisymmetric_linear_elasticity().

◆ body_force_fct_pt() [1/2]

BodyForceFctPt& oomph::AxisymmetricLinearElasticityEquationsBase::body_force_fct_pt ( )
inline

Access function: Pointer to body force function.

313  {
314  return Body_force_fct_pt;
315  }

References Body_force_fct_pt.

◆ body_force_fct_pt() [2/2]

BodyForceFctPt oomph::AxisymmetricLinearElasticityEquationsBase::body_force_fct_pt ( ) const
inline

Access function: Pointer to body force function (const version)

319  {
320  return Body_force_fct_pt;
321  }

References Body_force_fct_pt.

◆ d2u_dt2_axisymmetric_linear_elasticity()

double oomph::AxisymmetricLinearElasticityEquationsBase::d2u_dt2_axisymmetric_linear_elasticity ( const unsigned n,
const unsigned i 
) const
inline

d^2u/dt^2 at local node n

69  {
70  // Get the timestepper
71  TimeStepper* time_stepper_pt = node_pt(n)->time_stepper_pt();
72 
73  // Storage for the derivative - initialise to 0
74  double d2u_dt2 = 0.0;
75 
76  // If we are doing an unsteady solve then calculate the derivative
77  if (!time_stepper_pt->is_steady())
78  {
79  // Get the nodal index
80  const unsigned u_nodal_index =
82 
83  // Get the number of values required to represent history
84  const unsigned n_time = time_stepper_pt->ntstorage();
85 
86  // Loop over history values
87  for (unsigned t = 0; t < n_time; t++)
88  {
89  // Add the contribution to the derivative
90  d2u_dt2 +=
91  time_stepper_pt->weight(2, t) * nodal_value(t, n, u_nodal_index);
92  }
93  }
94 
95  return d2u_dt2;
96  }
virtual unsigned u_index_axisymmetric_linear_elasticity(const unsigned &i) const
Definition: axisym_linear_elasticity_elements.h:60
TimeStepper *& time_stepper_pt()
Return the pointer to the timestepper.
Definition: nodes.h:238
Node *& node_pt(const unsigned &n)
Return a pointer to the local node n.
Definition: elements.h:2175
double nodal_value(const unsigned &n, const unsigned &i) const
Definition: elements.h:2593
TimeStepper *& time_stepper_pt()
Definition: geom_objects.h:192
unsigned ntstorage() const
Definition: timesteppers.h:601
virtual double weight(const unsigned &i, const unsigned &j) const
Access function for j-th weight for the i-th derivative.
Definition: timesteppers.h:594
bool is_steady() const
Definition: timesteppers.h:389
t
Definition: plotPSD.py:36

References i, oomph::TimeStepper::is_steady(), n, oomph::FiniteElement::nodal_value(), oomph::FiniteElement::node_pt(), oomph::TimeStepper::ntstorage(), plotPSD::t, oomph::GeomObject::time_stepper_pt(), oomph::Data::time_stepper_pt(), u_index_axisymmetric_linear_elasticity(), and oomph::TimeStepper::weight().

Referenced by oomph::AxisymmetricLinearElasticityEquations::fill_in_generic_contribution_to_residuals_axisymmetric_linear_elasticity(), and interpolated_d2u_dt2_axisymmetric_linear_elasticity().

◆ du_dt_axisymmetric_linear_elasticity()

double oomph::AxisymmetricLinearElasticityEquationsBase::du_dt_axisymmetric_linear_elasticity ( const unsigned n,
const unsigned i 
) const
inline

du/dt at local node n

102  {
103  // Get the timestepper
104  TimeStepper* time_stepper_pt = node_pt(n)->time_stepper_pt();
105 
106  // Storage for the derivative - initialise to 0
107  double du_dt = 0.0;
108 
109  // If we are doing an unsteady solve then calculate the derivative
110  if (!time_stepper_pt->is_steady())
111  {
112  // Get the nodal index
113  const unsigned u_nodal_index =
115 
116  // Get the number of values required to represent history
117  const unsigned n_time = time_stepper_pt->ntstorage();
118 
119  // Loop over history values
120  for (unsigned t = 0; t < n_time; t++)
121  {
122  // Add the contribution to the derivative
123  du_dt +=
124  time_stepper_pt->weight(1, t) * nodal_value(t, n, u_nodal_index);
125  }
126  }
127  return du_dt;
128  }

References i, oomph::TimeStepper::is_steady(), n, oomph::FiniteElement::nodal_value(), oomph::FiniteElement::node_pt(), oomph::TimeStepper::ntstorage(), plotPSD::t, oomph::GeomObject::time_stepper_pt(), oomph::Data::time_stepper_pt(), u_index_axisymmetric_linear_elasticity(), and oomph::TimeStepper::weight().

Referenced by interpolated_du_dt_axisymmetric_linear_elasticity().

◆ get_dof_numbers_for_unknowns()

void oomph::AxisymmetricLinearElasticityEquationsBase::get_dof_numbers_for_unknowns ( std::list< std::pair< unsigned long, unsigned >> &  dof_lookup_list) const
inlinevirtual

Create a list of pairs for all unknowns in this element, so that the first entry in each pair contains the global equation number of the unknown, while the second one contains the number of the "DOF type" that this unknown is associated with. (Function can obviously only be called if the equation numbering scheme has been set up.)

Reimplemented from oomph::GeneralisedElement.

361  {
362  // temporary pair (used to store DOF lookup prior to being added
363  // to list)
364  std::pair<unsigned long, unsigned> dof_lookup;
365 
366  // number of nodes
367  const unsigned n_node = this->nnode();
368 
369  // Integer storage for local unknown
370  int local_unknown = 0;
371 
372  // Loop over the nodes
373  for (unsigned n = 0; n < n_node; n++)
374  {
375  // Loop over dimension
376  for (unsigned i = 0; i < 3; i++)
377  {
378  // If the variable is free
379  local_unknown = nodal_local_eqn(n, i);
380 
381  // ignore pinned values
382  if (local_unknown >= 0)
383  {
384  // store DOF type lookup in temporary pair: First entry in pair
385  // is global equation number; second entry is DOF type
386  dof_lookup.first = this->eqn_number(local_unknown);
387  dof_lookup.second = 0;
388 
389  // add to list
390  dof_lookup_list.push_front(dof_lookup);
391  }
392  }
393  }
394  }
int nodal_local_eqn(const unsigned &n, const unsigned &i) const
Definition: elements.h:1432
unsigned nnode() const
Return the number of nodes.
Definition: elements.h:2210
unsigned long eqn_number(const unsigned &ieqn_local) const
Definition: elements.h:704

References oomph::GeneralisedElement::eqn_number(), i, n, oomph::FiniteElement::nnode(), and oomph::FiniteElement::nodal_local_eqn().

◆ interpolated_d2u_dt2_axisymmetric_linear_elasticity()

void oomph::AxisymmetricLinearElasticityEquationsBase::interpolated_d2u_dt2_axisymmetric_linear_elasticity ( const Vector< double > &  s,
Vector< double > &  d2u_dt2 
) const
inline

Compute vector of FE interpolated accel d2u/dt2 at local coordinate s.

223  {
224  // Find number of nodes
225  unsigned n_node = nnode();
226 
227  // Local shape function
228  Shape psi(n_node);
229 
230  // Find values of shape function
231  shape(s, psi);
232 
233  // Loop over directions
234  for (unsigned i = 0; i < 3; i++)
235  {
236  // Initialise value of u
237  d2u_dt2[i] = 0.0;
238 
239  // Loop over the local nodes and sum
240  for (unsigned l = 0; l < n_node; l++)
241  {
242  d2u_dt2[i] += d2u_dt2_axisymmetric_linear_elasticity(l, i) * psi[l];
243  }
244  }
245  }
double d2u_dt2_axisymmetric_linear_elasticity(const unsigned &n, const unsigned &i) const
d^2u/dt^2 at local node n
Definition: axisym_linear_elasticity_elements.h:67
virtual void shape(const Vector< double > &s, Shape &psi) const =0
RealScalar s
Definition: level1_cplx_impl.h:130

References d2u_dt2_axisymmetric_linear_elasticity(), i, oomph::FiniteElement::nnode(), s, and oomph::FiniteElement::shape().

Referenced by oomph::AxisymmetricLinearElasticityEquations::output().

◆ interpolated_du_dt_axisymmetric_linear_elasticity()

void oomph::AxisymmetricLinearElasticityEquationsBase::interpolated_du_dt_axisymmetric_linear_elasticity ( const Vector< double > &  s,
Vector< double > &  du_dt 
) const
inline

Compute vector of FE interpolated velocity du/dt at local coordinate s.

196  {
197  // Find number of nodes
198  unsigned n_node = nnode();
199 
200  // Local shape function
201  Shape psi(n_node);
202 
203  // Find values of shape function
204  shape(s, psi);
205 
206  // Loop over directions
207  for (unsigned i = 0; i < 3; i++)
208  {
209  // Initialise value of u
210  du_dt[i] = 0.0;
211 
212  // Loop over the local nodes and sum
213  for (unsigned l = 0; l < n_node; l++)
214  {
215  du_dt[i] += du_dt_axisymmetric_linear_elasticity(l, i) * psi[l];
216  }
217  }
218  }
double du_dt_axisymmetric_linear_elasticity(const unsigned &n, const unsigned &i) const
du/dt at local node n
Definition: axisym_linear_elasticity_elements.h:100

References du_dt_axisymmetric_linear_elasticity(), i, oomph::FiniteElement::nnode(), s, and oomph::FiniteElement::shape().

Referenced by oomph::AxisymmetricLinearElasticityEquations::output().

◆ interpolated_u_axisymmetric_linear_elasticity() [1/2]

double oomph::AxisymmetricLinearElasticityEquationsBase::interpolated_u_axisymmetric_linear_elasticity ( const Vector< double > &  s,
const unsigned i 
) const
inline

Return FE interpolated displacement u[i] (i=0: r, i=1: z; i=2: theta) at local coordinate s

165  {
166  // Find number of nodes
167  unsigned n_node = nnode();
168 
169  // Local shape function
170  Shape psi(n_node);
171 
172  // Find values of shape function
173  shape(s, psi);
174 
175  // Get nodal index at which i-th velocity is stored
176  unsigned u_nodal_index = u_index_axisymmetric_linear_elasticity(i);
177 
178  // Initialise value of u
179  double interpolated_u = 0.0;
180 
181  // Loop over the local nodes and sum
182  for (unsigned l = 0; l < n_node; l++)
183  {
184  const double u_value = nodal_value(l, u_nodal_index);
185 
186  interpolated_u += u_value * psi[l];
187  }
188 
189  return (interpolated_u);
190  }

References i, oomph::FiniteElement::nnode(), oomph::FiniteElement::nodal_value(), s, oomph::FiniteElement::shape(), and u_index_axisymmetric_linear_elasticity().

◆ interpolated_u_axisymmetric_linear_elasticity() [2/2]

void oomph::AxisymmetricLinearElasticityEquationsBase::interpolated_u_axisymmetric_linear_elasticity ( const Vector< double > &  s,
Vector< double > &  disp 
) const
inline

Compute vector of FE interpolated displacement u at local coordinate s.

133  {
134  // Find number of nodes
135  unsigned n_node = nnode();
136 
137  // Local shape function
138  Shape psi(n_node);
139 
140  // Find values of shape function
141  shape(s, psi);
142 
143  for (unsigned i = 0; i < 3; i++)
144  {
145  // Index at which the nodal value is stored
146  unsigned u_nodal_index = u_index_axisymmetric_linear_elasticity(i);
147 
148  // Initialise value of u
149  disp[i] = 0.0;
150 
151  // Loop over the local nodes and sum
152  for (unsigned l = 0; l < n_node; l++)
153  {
154  const double u_value = nodal_value(l, u_nodal_index);
155 
156  disp[i] += u_value * psi[l];
157  }
158  }
159  }

References i, oomph::FiniteElement::nnode(), oomph::FiniteElement::nodal_value(), s, oomph::FiniteElement::shape(), and u_index_axisymmetric_linear_elasticity().

Referenced by oomph::AxisymmetricLinearElasticityEquations::compute_error(), and oomph::AxisymmetricLinearElasticityEquations::output().

◆ lambda_sq()

const double& oomph::AxisymmetricLinearElasticityEquationsBase::lambda_sq ( ) const
inline

Access function for timescale ratio (nondim density)

307  {
308  return *Lambda_sq_pt;
309  }

References Lambda_sq_pt.

Referenced by oomph::AxisymmetricLinearElasticityEquations::fill_in_generic_contribution_to_residuals_axisymmetric_linear_elasticity().

◆ lambda_sq_pt()

double*& oomph::AxisymmetricLinearElasticityEquationsBase::lambda_sq_pt ( )
inline

Access function for pointer to timescale ratio (nondim density)

301  {
302  return Lambda_sq_pt;
303  }

References Lambda_sq_pt.

◆ ndof_types()

unsigned oomph::AxisymmetricLinearElasticityEquationsBase::ndof_types ( ) const
inlinevirtual

The number of "DOF types" that degrees of freedom in this element are sub-divided into: for now lump them all into one DOF type. Can be adjusted later

Reimplemented from oomph::GeneralisedElement.

349  {
350  return 1;
351  }

◆ nu()

double& oomph::AxisymmetricLinearElasticityEquationsBase::nu ( ) const
inline

Access function for Poisson's ratio.

279  {
280 #ifdef PARANOID
281  if (Nu_pt == 0)
282  {
283  std::ostringstream error_message;
284  error_message << "No pointer to Poisson's ratio set. Please set one!\n";
285  throw OomphLibError(error_message.str(),
288  }
289 #endif
290  return *Nu_pt;
291  }
#define OOMPH_EXCEPTION_LOCATION
Definition: oomph_definitions.h:61
#define OOMPH_CURRENT_FUNCTION
Definition: oomph_definitions.h:86

References Nu_pt, OOMPH_CURRENT_FUNCTION, and OOMPH_EXCEPTION_LOCATION.

Referenced by oomph::AxisymmetricLinearElasticityEquations::fill_in_generic_contribution_to_residuals_axisymmetric_linear_elasticity().

◆ nu_pt()

double*& oomph::AxisymmetricLinearElasticityEquationsBase::nu_pt ( )
inline

Access function for pointer to Poisson's ratio.

295  {
296  return Nu_pt;
297  }

References Nu_pt.

◆ u_index_axisymmetric_linear_elasticity()

virtual unsigned oomph::AxisymmetricLinearElasticityEquationsBase::u_index_axisymmetric_linear_elasticity ( const unsigned i) const
inlinevirtual

Return the index at which the i-th (i=0: r, i=1: z; i=2: theta) unknown displacement component is stored at the nodes. The default assignment here (u_r, u_z, u_theta) is appropriate for single-physics problems.

62  {
63  return i;
64  }

References i.

Referenced by d2u_dt2_axisymmetric_linear_elasticity(), du_dt_axisymmetric_linear_elasticity(), oomph::AxisymmetricLinearElasticityEquations::fill_in_generic_contribution_to_residuals_axisymmetric_linear_elasticity(), oomph::AxisymmetricLinearElasticityEquations::get_strain(), and interpolated_u_axisymmetric_linear_elasticity().

◆ youngs_modulus()

double oomph::AxisymmetricLinearElasticityEquationsBase::youngs_modulus ( ) const
inline

◆ youngs_modulus_pt()

double*& oomph::AxisymmetricLinearElasticityEquationsBase::youngs_modulus_pt ( )
inline

Return the pointer to Young's modulus.

267  {
268  return Youngs_modulus_pt;
269  }

References Youngs_modulus_pt.

Member Data Documentation

◆ Body_force_fct_pt

BodyForceFctPt oomph::AxisymmetricLinearElasticityEquationsBase::Body_force_fct_pt
protected

Pointer to body force function.

Referenced by body_force(), and body_force_fct_pt().

◆ Default_lambda_sq_value

double oomph::AxisymmetricLinearElasticityEquationsBase::Default_lambda_sq_value
staticprotected
Initial value:
=
1.0

Static default value for timescale ratio (1.0 for natural scaling)

Static default value for timescale ratio (1.0 – for natural scaling)

◆ Default_youngs_modulus_value

double oomph::AxisymmetricLinearElasticityEquationsBase::Default_youngs_modulus_value
staticprotected
Initial value:
=
1.0

Static default value for Young's modulus (1.0 – for natural scaling, i.e. all stresses have been non-dimensionalised by the same reference Young's modulus. Setting the "non-dimensional" Young's modulus (obtained by de-referencing Youngs_modulus_pt) to a number larger than one means that the material is stiffer than assumed in the non-dimensionalisation.

◆ Lambda_sq_pt

double* oomph::AxisymmetricLinearElasticityEquationsBase::Lambda_sq_pt
protected

Timescale ratio (non-dim. density)

Referenced by lambda_sq(), and lambda_sq_pt().

◆ Nu_pt

double* oomph::AxisymmetricLinearElasticityEquationsBase::Nu_pt
protected

Pointer to Poisson's ratio.

Referenced by nu(), and nu_pt().

◆ Youngs_modulus_pt

double* oomph::AxisymmetricLinearElasticityEquationsBase::Youngs_modulus_pt
protected

Pointer to the Young's modulus.

Referenced by youngs_modulus(), and youngs_modulus_pt().


The documentation for this class was generated from the following files: