32 #ifndef OOMPH_POISSON_FLUX_ELEMENTS_HEADER 33 #define OOMPH_POISSON_FLUX_ELEMENTS_HEADER 38 #include <oomph-lib-config.h> 42 #include "../generic/Qelements.h" 53 template <
class ELEMENT>
74 "Don't call empty constructor for PoissonFluxElement",
75 OOMPH_CURRENT_FUNCTION,
76 OOMPH_EXCEPTION_LOCATION);
97 const unsigned &
i)
const 126 const unsigned n_plot=5;
131 void output(std::ostream &outfile,
const unsigned &nplot)
135 unsigned el_dim=
dim();
145 for (
unsigned iplot=0;iplot<num_plot_points;iplot++)
152 for(
unsigned i=0;
i<el_dim+1;
i++)
155 outfile << x[
i] <<
" ";
159 outfile << flux << std::endl;
175 void output(FILE* file_pt,
const unsigned &n_plot)
188 unsigned n_node =
nnode();
194 for(
unsigned i=0;
i<n_node;
i++) {test[
i] = psi[
i];}
209 unsigned n_node =
nnode();
215 for(
unsigned i=0;
i<n_node;
i++) {test[
i] = psi[
i];}
234 (*Flux_fct_pt)(x,flux);
246 const unsigned& flag);
274 template<
class ELEMENT>
288 ELEMENT* elem_pt =
dynamic_cast<ELEMENT*
>(bulk_el_pt);
290 if(elem_pt->dim()==3)
299 "This flux element will not work correctly if nodes are hanging\n",
300 OOMPH_CURRENT_FUNCTION,
301 OOMPH_EXCEPTION_LOCATION);
335 "Bulk element must inherit from PoissonEquations.";
337 "Nodes are one dimensional, but cannot cast the bulk element to\n";
338 error_string +=
"PoissonEquations<1>\n.";
340 "If you desire this functionality, you must implement it yourself\n";
343 OOMPH_CURRENT_FUNCTION,
344 OOMPH_EXCEPTION_LOCATION);
364 "Bulk element must inherit from PoissonEquations.";
366 "Nodes are two dimensional, but cannot cast the bulk element to\n";
367 error_string +=
"PoissonEquations<2>\n.";
369 "If you desire this functionality, you must implement it yourself\n";
372 OOMPH_CURRENT_FUNCTION,
373 OOMPH_EXCEPTION_LOCATION);
392 "Bulk element must inherit from PoissonEquations.";
394 "Nodes are three dimensional, but cannot cast the bulk element to\n";
395 error_string +=
"PoissonEquations<3>\n.";
397 "If you desire this functionality, you must implement it yourself\n";
400 OOMPH_CURRENT_FUNCTION,
401 OOMPH_EXCEPTION_LOCATION);
414 std::ostringstream error_stream;
415 error_stream <<
"Dimension of node is " <<
Dim 416 <<
". It should be 1,2, or 3!" << std::endl;
419 OOMPH_CURRENT_FUNCTION,
420 OOMPH_EXCEPTION_LOCATION);
429 template<
class ELEMENT>
433 const unsigned& flag)
436 const unsigned n_node =
nnode();
439 Shape psif(n_node), testf(n_node);
455 for(
unsigned ipt=0;ipt<n_intpt;ipt++)
475 for(
unsigned l=0;l<n_node;l++)
478 for(
unsigned i=0;
i<
Dim;
i++)
491 for(
unsigned l=0;l<n_node;l++)
498 residuals[local_eqn] -= flux*testf[l]*
W;
virtual std::string tecplot_zone_string(const unsigned &nplot) const
Return string for tecplot zone header (when plotting nplot points in each "coordinate direction") ...
double J_eulerian_at_knot(const unsigned &ipt) const
Return the Jacobian of the mapping from local to global coordinates at the ipt-th integration point O...
void output(FILE *file_pt, const unsigned &n_plot)
C-style output function – forward to broken version in FiniteElement until somebody decides what exa...
void fill_in_generic_residual_contribution_poisson_flux(Vector< double > &residuals, DenseMatrix< double > &jacobian, const unsigned &flag)
Add the element's contribution to its residual vector. flag=1(or 0): do (or don't) compute the contri...
bool has_hanging_nodes() const
Return boolean to indicate if any of the element's nodes are geometrically hanging.
void broken_copy(const std::string &class_name)
Issue error message and terminate execution.
double shape_and_test_at_knot(const unsigned &ipt, Shape &psi, Shape &test) const
Function to compute the shape and test functions and to return the Jacobian of mapping between local ...
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 void output(std::ostream &outfile)
Output the element data — typically the values at the nodes in a format suitable for post-processing...
double nodal_position(const unsigned &n, const unsigned &i) const
Return the i-th coordinate at local node n. If the node is hanging, the appropriate interpolation is ...
A class for elements that allow the imposition of an applied flux on the boundaries of Poisson elemen...
int & face_index()
Index of the face (a number that uniquely identifies the face in the element)
A general Finite Element class.
double zeta_nodal(const unsigned &n, const unsigned &k, const unsigned &i) const
In a FaceElement, the "global" intrinsic coordinate of the element along the boundary, when viewed as part of a compound geometric object is specified using the boundary coordinate defined by the mesh. Note: Boundary coordinates will have been set up when creating the underlying mesh, and their values will have been stored at the nodes.
double interpolated_x(const Vector< double > &s, const unsigned &i) const
Return FE interpolated coordinate x[i] at local coordinate s. Overloaded to get information from bulk...
unsigned U_index_poisson
The index at which the unknown is stored at the nodes.
void get_flux(const Vector< double > &x, double &flux)
void fill_in_contribution_to_residuals(Vector< double > &residuals)
Add the element's contribution to its residual vector.
virtual double weight(const unsigned &i) const =0
Return weight of i-th integration point.
void output(std::ostream &outfile)
Output function.
void(* PoissonPrescribedFluxFctPt)(const Vector< double > &x, double &flux)
Function pointer to the prescribed-flux function fct(x,f(x)) – x is a Vector!
virtual unsigned nplot_points(const unsigned &nplot) const
Return total number of plot points (when plotting nplot points in each "coordinate direction") ...
double J_eulerian(const Vector< double > &s) const
Return the Jacobian of mapping from local to global coordinates at local position s...
unsigned ndim() const
Return (Eulerian) spatial dimension of the node.
virtual unsigned u_index_poisson() const
Return the index at which the unknown value is stored. The default value, 0, is appropriate for singl...
void output(FILE *file_pt)
Integral *const & integral_pt() const
Return the pointer to the integration scheme (const version)
virtual void shape_at_knot(const unsigned &ipt, Shape &psi) const
Return the geometric shape function at the ipt-th integration point.
virtual double knot(const unsigned &i, const unsigned &j) const =0
Return local coordinate s[j] of i-th integration point.
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"...
double shape_and_test(const Vector< double > &s, Shape &psi, Shape &test) const
Function to compute the shape and test functions and to return the Jacobian of mapping between local ...
PoissonFluxElement()
Broken empty constructor.
PoissonPrescribedFluxFctPt Flux_fct_pt
Function pointer to the (global) prescribed-flux function.
Node *& node_pt(const unsigned &n)
Return a pointer to the local node n.
unsigned Dim
The spatial dimension of the problem.
PoissonPrescribedFluxFctPt & flux_fct_pt()
Access function for the prescribed-flux function pointer.
void operator=(const PoissonFluxElement &)
Broken assignment operator.
unsigned dim() const
Return the spatial dimension of the element, i.e. the number of local coordinates required to paramet...
void broken_assign(const std::string &class_name)
Issue error message and terminate execution.
void output(std::ostream &outfile, const unsigned &nplot)
Output function.
virtual unsigned nweight() const =0
Return the number of integration points of the scheme.
double zeta_nodal(const unsigned &n, const unsigned &k, const unsigned &i) const
Specify the value of nodal zeta from the face geometry The "global" intrinsic coordinate of the eleme...
std::string string(const unsigned &i)
Return the i-th string or "" if the relevant string hasn't been defined.
virtual void build_face_element(const int &face_index, FaceElement *face_element_pt)
Function for building a lower dimensional FaceElement on the specified face of the FiniteElement...
unsigned nnode() const
Return the number of nodes.
PoissonFluxElement(const PoissonFluxElement &dummy)
Broken copy constructor.
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...
static DenseMatrix< double > Dummy_matrix
Empty dense matrix used as a dummy argument to combined residual and jacobian functions in the case w...
void fill_in_contribution_to_jacobian(Vector< double > &residuals, DenseMatrix< double > &jacobian)
Add the element's contribution to its residual vector and its Jacobian matrix.
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...