36 template<
unsigned NNODE_1D>
49 const unsigned n_node =
nnode();
52 Shape psi(n_node), test(n_node);
53 DShape dpsidr(n_node, 1), dtestdr(n_node, 1);
64 const double nu_local =
nu();
65 const double eta_local =
eta();
71 for (
unsigned ipt = 0; ipt < n_intpt; ipt++)
78 ipt, psi, dpsidr, test, dtestdr);
81 double interpolated_r = 0.0;
83 double interpolated_w = 0.0;
84 double interpolated_laplacian_w = 0.0;
85 double interpolated_u_r = 0.0;
87 double interpolated_dwdr = 0.0;
88 double interpolated_dlaplacian_wdr = 0.0;
89 double interpolated_du_rdr = 0.0;
96 for (
unsigned l = 0; l < n_node; l++)
105 interpolated_laplacian_w +=
nodal_value[1] * psi(l);
110 interpolated_dwdr +=
nodal_value[0] * dpsidr(l, 0);
111 interpolated_dlaplacian_wdr +=
nodal_value[1] * dpsidr(l, 0);
112 interpolated_du_rdr +=
nodal_value[2] * dpsidr(l, 0);
117 double W = w * interpolated_r * J;
121 double pressure = 0.0;
127 double sigma_r_r = 0.0;
128 double sigma_phi_phi = 0.0;
133 1.0 / (1.0 - nu_local * nu_local) *
134 (interpolated_du_rdr + 0.5 * interpolated_dwdr * interpolated_dwdr +
135 nu_local * 1.0 / interpolated_r * interpolated_u_r);
138 1.0 / (1.0 - nu_local * nu_local) *
139 (1.0 / interpolated_r * interpolated_u_r +
140 nu_local * (interpolated_du_rdr +
141 0.5 * interpolated_dwdr * interpolated_dwdr));
145 sigma_r_r = 1.0 / (1.0 - nu_local * nu_local) *
146 (interpolated_du_rdr +
147 nu_local * 1.0 / interpolated_r * interpolated_u_r);
149 sigma_phi_phi = 1.0 / (1.0 - nu_local * nu_local) *
150 (1.0 / interpolated_r * interpolated_u_r +
151 nu_local * interpolated_du_rdr);
158 for (
unsigned l = 0; l < n_node; l++)
166 residuals[local_eqn] +=
167 (pressure * test(l) +
168 (dtestdr(l, 0)) * interpolated_dlaplacian_wdr) *
172 residuals[local_eqn] -=
173 eta_local * sigma_r_r * (dtestdr(l, 0)) * interpolated_dwdr *
W;
183 residuals[local_eqn] += (test(l) * interpolated_laplacian_w +
184 (dtestdr(l, 0)) * interpolated_dwdr) *
194 residuals[local_eqn] +=
195 (sigma_r_r * dtestdr(l, 0) +
196 1.0 / interpolated_r * sigma_phi_phi * test(l)) *
235 const Vector<double>&
s,
double& sigma_r_r,
double& sigma_phi_phi)
const
246 unsigned n_dim = this->
dim();
247 unsigned n_node = this->
nnode();
249 DShape dpsidr(n_node, n_dim);
255 double interpolated_r = 0.0;
256 double interpolated_u_r = 0.0;
258 double interpolated_dwdr = 0.0;
259 double interpolated_du_rdr = 0.0;
261 double nu_local =
nu();
267 for (
unsigned l = 0; l < n_node; l++)
275 interpolated_du_rdr +=
289 if (interpolated_r == 0.0)
292 sigma_r_r = interpolated_du_rdr / (1.0 - nu_local);
293 sigma_phi_phi = sigma_r_r;
300 1.0 / (1.0 - nu_local * nu_local) *
301 (interpolated_du_rdr + 0.5 * interpolated_dwdr * interpolated_dwdr +
302 nu_local * 1.0 / interpolated_r * interpolated_u_r);
305 1.0 / (1.0 - nu_local * nu_local) *
306 (1.0 / interpolated_r * interpolated_u_r +
307 nu_local * (interpolated_du_rdr +
308 0.5 * interpolated_dwdr * interpolated_dwdr));
320 const unsigned& nplot)
330 for (
unsigned iplot = 0; iplot < num_plot_points; iplot++)
336 double sigma_r_r = 0.0;
337 double sigma_phi_phi = 0.0;
343 << sigma_phi_phi << std::endl;
353 const unsigned& nplot)
363 for (
unsigned iplot = 0; iplot < num_plot_points; iplot++)
387 std::ostream& outfile,
388 const unsigned& nplot,
406 for (
unsigned iplot = 0; iplot < num_plot_points; iplot++)
415 (*exact_soln_pt)(r, exact_soln);
418 outfile << r[0] <<
" ";
419 outfile << exact_soln[0] << std::endl;
435 std::ostream& outfile,
451 unsigned n_node =
nnode();
459 outfile <<
"ZONE" << std::endl;
466 for (
unsigned ipt = 0; ipt < n_intpt; ipt++)
487 (*exact_soln_pt)(r, exact_soln);
490 outfile << r[0] <<
" ";
491 outfile << exact_soln[0] <<
" " << exact_soln[0] - w_fe << std::endl;
494 norm += exact_soln[0] * exact_soln[0] *
W;
495 error += (exact_soln[0] - w_fe) * (exact_soln[0] - w_fe) *
W;
//////////////////////////////////////////////////////////////////////// ////////////////////////////...
static const unsigned Initial_Nvalue
Static int that holds the number of variables at nodes: always the same.
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.
unsigned self_test()
Self-test: Return 0 for OK.
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 interpolated_stress(const Vector< double > &s, double &sigma_r_r, double &sigma_phi_phi) const
Compute in-plane stresses.
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.
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...
const double & nu() const
Poisson's ratio.
A Class for the derivatives of shape functions The class design is essentially the same as Shape,...
virtual double J_eulerian(const Vector< double > &s) const
Return the Jacobian of mapping from local to global coordinates at local position s.
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 std::string tecplot_zone_string(const unsigned &nplot) const
Return string for tecplot zone header (when plotting nplot points in each "coordinate direction")
virtual double interpolated_x(const Vector< double > &s, const unsigned &i) const
Return FE interpolated coordinate x[i] at local coordinate s.
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 .
Integral *const & integral_pt() const
Return the pointer to the integration scheme (const version)
virtual void get_s_plot(const unsigned &i, const unsigned &nplot, Vector< double > &s, const bool &shifted_to_interior=false) const
Get cector of local coordinates of plot point i (when plotting nplot points in each "coordinate direc...
virtual unsigned nplot_points(const unsigned &nplot) const
Return total number of plot points (when plotting nplot points in each "coordinate direction")
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...
double raw_nodal_value(const unsigned &n, const unsigned &i) const
Return the i-th value stored at local node n but do NOT take hanging nodes into account.
double raw_nodal_position(const unsigned &n, const unsigned &i) const
Return the i-th coordinate at local node n. Do not use the hanging node representation....
virtual void write_tecplot_zone_footer(std::ostream &outfile, const unsigned &nplot) const
Add tecplot zone "footer" to output stream (when plotting nplot points in each "coordinate direction"...
virtual unsigned self_test()
Self-test: Check inversion of element & do self-test for GeneralisedElement. Return 0 if OK.
virtual double knot(const unsigned &i, const unsigned &j) const =0
Return local coordinate s[j] of i-th integration point.
virtual unsigned nweight() const =0
Return the number of integration points of the scheme.
virtual double weight(const unsigned &i) const =0
Return weight of i-th integration point.
A Class for shape functions. In simple cases, the shape functions have only one index that can be tho...
//////////////////////////////////////////////////////////////////// ////////////////////////////////...