33 #ifndef OOMPH_REFINEABLE_SPHERICAL_ADVECTION_DIFFUSION_ELEMENTS_HEADER 34 #define OOMPH_REFINEABLE_SPHERICAL_ADVECTION_DIFFUSION_ELEMENTS_HEADER 38 #include <oomph-lib-config.h> 42 #include "../generic/refineable_quad_element.h" 43 #include "../generic/refineable_brick_element.h" 44 #include "../generic/error_estimator.h" 106 unsigned n_node =
nnode();
121 for(
unsigned l=0;l<n_node;l++)
123 values[0] += this->
nodal_value(l,u_nodal_index)*psi[l];
138 const unsigned n_node =
nnode();
153 for(
unsigned l=0;l<n_node;l++)
155 values[0] += this->
nodal_value(t,l,u_nodal_index)*psi[l];
162 "Time-dependent version of get_interpolated_values() ";
163 error_message +=
"not implemented for this element \n";
167 "RefineableSphericalAdvectionDiffusionEquations::get_interpolated_values()",
168 OOMPH_EXCEPTION_LOCATION);
181 {
return x[0]*x[0]*sin(x[1]);}
193 this->
Pe_pt = cast_father_element_pt->
pe_pt();
211 unsigned n_node = this->
nnode();
227 for(
unsigned l=0;l<n_node;l++)
229 unsigned n_master = 1;
238 n_master = hang_info_pt->
nmaster();
247 for(
unsigned m=0;m<n_master;m++)
263 if (global_eqn >= 0) {++n_u_dof;}
268 du_ddata.resize(n_u_dof,0.0);
269 global_eqn_number.resize(n_u_dof,0);
274 for(
unsigned l=0;l<n_node;l++)
276 unsigned n_master = 1;
277 double hang_weight = 1.0;
286 n_master = hang_info_pt->
nmaster();
295 for(
unsigned m=0;m<n_master;m++)
325 global_eqn_number[count] = global_eqn;
327 du_ddata[count] = psi[l]*hang_weight;
354 template <
unsigned NNODE_1D>
423 template<
unsigned NNODE_1D>
Node * vertex_node_pt(const unsigned &j) const
Pointer to the j-th vertex node in the element.
void broken_copy(const std::string &class_name)
Issue error message and terminate execution.
Base class for finite elements that can compute the quantities that are required for the Z2 error est...
RefineableSphericalAdvectionDiffusionEquations()
Empty Constructor.
virtual RefineableElement * father_element_pt() const
Return a pointer to the father element.
unsigned num_Z2_flux_terms()
Number of 'flux' terms for Z2 error estimation.
HangInfo *const & hanging_pt() const
Return pointer to hanging node data (this refers to the geometric hanging node status) (const version...
void further_setup_hanging_nodes()
Perform additional hanging node procedures for variables that are not interpolated by all nodes...
SphericalAdvectionDiffusionSourceFctPt & source_fct_pt()
Access function: Pointer to source function.
void get_interpolated_values(const unsigned &t, const Vector< double > &s, Vector< double > &values)
Get the function value u in Vector. Note: Given the generality of the interface (this function is usu...
void get_flux(const Vector< double > &s, Vector< double > &flux) const
Get flux: .
double * PeSt_pt
Pointer to global Peclet number multiplied by Strouhal number.
void get_interpolated_values(const Vector< double > &s, Vector< double > &values)
Get the function value u in Vector. Note: Given the generality of the interface (this function is usu...
RefineableQSphericalAdvectionDiffusionElement(const RefineableQSphericalAdvectionDiffusionElement< NNODE_1D > &dummy)
Broken copy constructor.
unsigned nvertex_node() const
Number of vertex nodes in the element.
double geometric_jacobian(const Vector< double > &x)
Fill in the geometric Jacobian, which in this case is r*r*sin(theta)
double const & master_weight(const unsigned &i) const
Return weight for dofs on i-th master node.
Nodes are derived from Data, but, in addition, have a definite (Eulerian) position in a space of a gi...
unsigned nmaster() const
Return the number of master nodes.
void dinterpolated_u_adv_diff_ddata(const Vector< double > &s, Vector< double > &du_ddata, Vector< unsigned > &global_eqn_number)
Compute the derivatives of the i-th component of velocity at point s with respect to all data that ca...
bool is_hanging() const
Test whether the node is geometrically hanging.
long & eqn_number(const unsigned &i)
Return the equation number of the i-th stored variable.
void operator=(const RefineableQSphericalAdvectionDiffusionElement< NNODE_1D > &)
Broken assignment operator.
A class for all elements that solve the Advection Diffusion equations in a spherical polar coordinate...
void get_Z2_flux(const Vector< double > &s, Vector< double > &flux)
Get 'flux' for Z2 error recovery: Standard flux.from AdvectionDiffusion equations.
unsigned nvertex_node() const
Number of vertex nodes in the element.
FaceGeometry()
Constructor: Call the constructor for the appropriate lower-dimensional QElement. ...
Node * vertex_node_pt(const unsigned &j) const
Pointer to the j-th vertex node in the element.
unsigned long eqn_number(const unsigned &ieqn_local) const
Return the global equation number corresponding to the ieqn_local-th local equation number...
double * Pe_pt
Pointer to global Peclet number.
void fill_in_generic_residual_contribution_spherical_adv_diff(Vector< double > &residuals, DenseMatrix< double > &jacobian, DenseMatrix< double > &mass_matrix, unsigned flag)
Add the element's contribution to the elemental residual vector and/or Jacobian matrix flag=1: comput...
Refineable version of QSphericalAdvectionDiffusionElement. Inherit from the standard QSphericalAdvect...
SphericalAdvectionDiffusionSourceFctPt Source_fct_pt
Pointer to source function:
Node *const & master_node_pt(const unsigned &i) const
Return a pointer to the i-th master node.
Node *& node_pt(const unsigned &n)
Return a pointer to the local node n.
double *& pe_st_pt()
Pointer to Peclet number multipled by Strouha number.
A version of the Advection Diffusion in spherical coordinates equations that can be used with non-uni...
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...
QSphericalAdvectionDiffusionElement elements are linear/quadrilateral/brick-shaped Axisymmetric Advec...
Class that contains data for hanging nodes.
unsigned ncont_interpolated_values() const
Number of continuously interpolated values: 1.
void broken_assign(const std::string &class_name)
Issue error message and terminate execution.
std::string string(const unsigned &i)
Return the i-th string or "" if the relevant string hasn't been defined.
void rebuild_from_sons(Mesh *&mesh_pt)
Rebuild from sons: empty.
RefineableQSphericalAdvectionDiffusionElement()
Empty Constructor:
unsigned nrecovery_order()
Order of recovery shape functions for Z2 error estimation: Same order as shape functions.
SphericalAdvectionDiffusionWindFctPt Wind_fct_pt
Pointer to wind function:
void operator=(const RefineableSphericalAdvectionDiffusionEquations &)
Broken assignment operator.
unsigned nnode() const
Return the number of nodes.
void further_build()
Further build: Copy source function pointer from father element.
SphericalAdvectionDiffusionWindFctPt & wind_fct_pt()
Access function: Pointer to wind function.
double *& pe_pt()
Pointer to Peclet number.
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...
RefineableSphericalAdvectionDiffusionEquations(const RefineableSphericalAdvectionDiffusionEquations &dummy)
Broken copy constructor.
virtual unsigned u_index_spherical_adv_diff() const
Return the index at which the unknown value is stored. The default value, 0, is appropriate for singl...