27 #ifndef OOMPH_AXISYM_DISPL_BASED_FOEPPLVONKARMAN_ELEMENTS_HEADER
28 #define OOMPH_AXISYM_DISPL_BASED_FOEPPLVONKARMAN_ELEMENTS_HEADER
33 #include <oomph-lib-config.h>
37 #include "../generic/nodes.h"
38 #include "../generic/Qelements.h"
39 #include "../generic/oomph_utilities.h"
70 const double&
nu()
const
75 std::stringstream error_stream;
76 error_stream <<
"Nu has not yet been set!" << std::endl;
78 error_stream.str(), OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
91 const double&
eta()
const
119 const unsigned n_plot = 5;
125 void output(std::ostream& outfile,
const unsigned& n_plot);
130 const unsigned n_plot = 5;
136 void output(FILE* file_pt,
const unsigned& n_plot);
140 const unsigned& n_plot,
147 std::ostream& outfile,
148 const unsigned& n_plot,
153 "There is no time-dependent output_fct() for Foeppl von Karman"
155 OOMPH_CURRENT_FUNCTION,
156 OOMPH_EXCEPTION_LOCATION);
174 "There is no time-dependent compute_error() for Foeppl von Karman"
176 OOMPH_CURRENT_FUNCTION,
177 OOMPH_EXCEPTION_LOCATION);
198 double& pressure)
const
208 (*Pressure_fct_pt)(r, pressure);
217 const unsigned n_node =
nnode();
233 for (
unsigned l = 0; l < n_node; l++)
235 gradient[0] += this->
nodal_value(l, w_nodal_index) * dpsidr(l, 0);
252 const unsigned n_node =
nnode();
264 double interpolated_w = 0.0;
267 for (
unsigned l = 0; l < n_node; l++)
269 interpolated_w += this->
nodal_value(l, w_nodal_index) * psi[l];
272 return (interpolated_w);
280 const unsigned n_node =
nnode();
292 double interpolated_u = 0.0;
295 for (
unsigned l = 0; l < n_node; l++)
297 interpolated_u += this->
nodal_value(l, u_nodal_index) * psi[l];
300 return (interpolated_u);
307 double& sigma_phi_phi)
const;
330 unsigned total_fvk_nodal_indices = 3;
333 unsigned n_node =
nnode();
336 for (
unsigned index = first_fvk_nodal_index + 2;
337 index < first_fvk_nodal_index + total_fvk_nodal_indices;
341 for (
unsigned inod = 0; inod < n_node; inod++)
358 DShape& dtestdr)
const = 0;
368 DShape& dtestdr)
const = 0;
399 template<
unsigned NNODE_1D>
401 :
public virtual QElement<1, NNODE_1D>,
441 void output(std::ostream& outfile,
const unsigned& n_plot)
455 void output(FILE* file_pt,
const unsigned& n_plot)
463 const unsigned& n_plot,
467 outfile, n_plot, exact_soln_pt);
474 const unsigned& n_plot,
479 outfile, n_plot, time, exact_soln_pt);
511 template<
unsigned NNODE_1D>
512 double AxisymFoepplvonKarmanElement<
521 const double J = this->dshape_eulerian(
s, psi, dpsidr);
538 template<
unsigned NNODE_1D>
547 const double J = this->dshape_eulerian_at_knot(ipt, psi, dpsidr);
//////////////////////////////////////////////////////////////////////// ////////////////////////////...
void output(FILE *file_pt, const unsigned &n_plot)
C-style output function: r,w at n_plot plot points.
double dshape_and_dtest_eulerian_axisym_fvk(const Vector< double > &s, Shape &psi, DShape &dpsidr, Shape &test, DShape &dtestdr) const
Shape, test functions & derivs. w.r.t. to global coords. Return Jacobian.
void output(FILE *file_pt)
C-style output function: r,w.
double dshape_and_dtest_eulerian_at_knot_axisym_fvk(const unsigned &ipt, Shape &psi, DShape &dpsidr, Shape &test, DShape &dtestdr) const
Shape, test functions & derivs. w.r.t. to global coords. at integration point ipt....
void output_fct(std::ostream &outfile, const unsigned &n_plot, const double &time, FiniteElement::UnsteadyExactSolutionFctPt exact_soln_pt)
Output function for a time-dependent exact solution. r,w_exact at n_plot plot points (Calls the stead...
void output(std::ostream &outfile, const unsigned &n_plot)
Output function: r,w,u,sigma_r_r,sigma_phi_phi at n_plot plot points.
AxisymFoepplvonKarmanElement(const AxisymFoepplvonKarmanElement< NNODE_1D > &dummy)=delete
Broken copy constructor.
void output(std::ostream &outfile)
Output function: r,w,u,sigma_r_r,sigma_phi_phi.
static const unsigned Initial_Nvalue
Static int that holds the number of variables at nodes: always the same.
unsigned required_nvalue(const unsigned &n) const
Required # of ‘values’ (pinned or dofs) at node n.
void operator=(const AxisymFoepplvonKarmanElement< NNODE_1D > &)=delete
Broken assignment operator.
AxisymFoepplvonKarmanElement()
Constructor: Call constructors for QElement and AxisymFoepplvonKarmanEquations.
void output_fct(std::ostream &outfile, const unsigned &n_plot, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt)
Output function for an exact solution: r,w_exact at n_plot plot points.
A class for all isoparametric elements that solve the axisYm Foeppl von Karman equations in a displac...
double * Nu_pt
Pointer to Poisson's ratio.
AxisymFoepplvonKarmanPressureFctPt pressure_fct_pt() const
Access function: Pointer to pressure function. Const version.
AxisymFoepplvonKarmanPressureFctPt Pressure_fct_pt
Pointer to pressure function:
static double Default_Physical_Constant_Value
Default value for physical constants.
void operator=(const AxisymFoepplvonKarmanEquations &)=delete
Broken assignment operator.
virtual void output_fct(std::ostream &outfile, const unsigned &n_plot, const double &time, FiniteElement::UnsteadyExactSolutionFctPt exact_soln_pt)
Output exact soln: r,w_exact at n_plot plot points (dummy time-dependent version to keep intel compil...
void output(std::ostream &outfile)
Output with default number of plot points.
void fill_in_contribution_to_residuals(Vector< double > &residuals)
Fill in the residuals with this element's contribution.
void output(FILE *file_pt)
C_style output with default number of plot points.
void compute_error(std::ostream &outfile, FiniteElement::UnsteadyExactSolutionFctPt exact_soln_pt, const double &time, double &error, double &norm)
Dummy, time dependent error checker.
unsigned self_test()
Self-test: Return 0 for OK.
double * Eta_pt
Pointer to FvK parameter.
virtual double dshape_and_dtest_eulerian_at_knot_axisym_fvk(const unsigned &ipt, Shape &psi, DShape &dpsidr, Shape &test, DShape &dtestdr) const =0
Shape/test functions and derivs w.r.t. to global coords at integration point ipt; return Jacobian of ...
virtual unsigned nodal_index_fvk(const unsigned &i=0) const
Return the index at which the i-th unknown value is stored. The default value, i, is appropriate for ...
void(* AxisymFoepplvonKarmanPressureFctPt)(const double &r, double &f)
Function pointer to pressure function fct(r,f(r)) – r is a Vector!
void get_gradient_of_deflection(const Vector< double > &s, Vector< double > &gradient) const
Get gradient of deflection: gradient[i] = dw/dr_i */.
void interpolated_stress(const Vector< double > &s, double &sigma_r_r, double &sigma_phi_phi) const
Compute in-plane stresses.
AxisymFoepplvonKarmanPressureFctPt & pressure_fct_pt()
Access function: Pointer to pressure function.
AxisymFoepplvonKarmanEquations()
Constructor (must initialise the Pressure_fct_pt). Also set physical parameters to their default valu...
virtual double dshape_and_dtest_eulerian_axisym_fvk(const Vector< double > &s, Shape &psi, DShape &dpsidr, Shape &test, DShape &dtestdr) const =0
Shape/test functions and derivs w.r.t. to global coords at local coord. s; return Jacobian of mapping...
bool Linear_bending_model
Flag which stores whether we are using a linear, pure bending model instead of the full non-linear Fo...
const double & eta() const
FvK parameter.
void output_fct(std::ostream &outfile, const unsigned &n_plot, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt)
Output exact soln: r,w_exact at n_plot plot points.
void compute_error(std::ostream &outfile, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt, double &error, double &norm)
Get error against and norm of exact solution.
void use_linear_bending_model()
Sets a flag to signify that we are solving the linear, pure bending equations, and pin all the nodal ...
double interpolated_u_fvk(const Vector< double > &s) const
Return FE representation of radial displacement.
double interpolated_w_fvk(const Vector< double > &s) const
Return FE representation of transverse displacement.
virtual void get_pressure_fvk(const unsigned &ipt, const double &r, double &pressure) const
Get pressure term at (Eulerian) position r. This function is virtual to allow overloading in multi-ph...
double *& eta_pt()
Pointer to FvK parameter.
AxisymFoepplvonKarmanEquations(const AxisymFoepplvonKarmanEquations &dummy)=delete
Broken copy constructor.
const double & nu() const
Poisson's ratio.
double *& nu_pt()
Pointer to Poisson's ratio.
A Class for the derivatives of shape functions The class design is essentially the same as Shape,...
void pin(const unsigned &i)
Pin the i-th stored variable.
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...
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...
void(* UnsteadyExactSolutionFctPt)(const double &, const Vector< double > &, Vector< double > &)
Function pointer for function that computes Vector-valued time-dependent function as .
An OomphLibError object which should be thrown when an run-time error is encountered....
/////////////////////////////////////////////////////////////////////// /////////////////////////////...
A Class for shape functions. In simple cases, the shape functions have only one index that can be tho...
//////////////////////////////////////////////////////////////////// ////////////////////////////////...