32 #ifndef OOMPH_FOURIER_DECOMPOSED_HELMHOLTZ_FLUX_ELEMENTS_HEADER 33 #define OOMPH_FOURIER_DECOMPOSED_HELMHOLTZ_FLUX_ELEMENTS_HEADER 38 #include <oomph-lib-config.h> 45 #include "../generic/Qelements.h" 57 template <
class ELEMENT>
80 "Don't call empty constructor for FourierDecomposedHelmholtzFluxElement",
81 OOMPH_CURRENT_FUNCTION,
82 OOMPH_EXCEPTION_LOCATION);
125 (residuals,jacobian,1);
134 void output(std::ostream &outfile,
const unsigned &n_plot)
145 void output(FILE* file_pt,
const unsigned &n_plot)
151 virtual inline std::complex<unsigned>
154 return std::complex<unsigned>(
170 unsigned n_node =
nnode();
176 for(
unsigned i=0;
i<n_node;
i++) {test[
i] = psi[
i];}
191 unsigned n_node =
nnode();
197 for(
unsigned i=0;
i<n_node;
i++) {test[
i] = psi[
i];}
211 flux = std::complex<double>(0.0,0.0);
216 (*Flux_fct_pt)(x,flux);
232 const unsigned& flag);
253 template<
class ELEMENT>
278 "Bulk element must inherit from FourierDecomposedHelmholtzEquations.";
281 OOMPH_CURRENT_FUNCTION,
282 OOMPH_EXCEPTION_LOCATION);
288 eqn_pt->u_index_fourier_decomposed_helmholtz();
297 template<
class ELEMENT>
301 const unsigned& flag)
304 const unsigned n_node =
nnode();
307 Shape psif(n_node), testf(n_node);
316 int local_eqn_real=0 ,local_eqn_imag=0;
320 for(
unsigned ipt=0;ipt<n_intpt;ipt++)
340 for(
unsigned l=0;l<n_node;l++)
343 for(
unsigned i=0;
i<2;
i++)
350 double r = interpolated_x[0];
353 std::complex<double> flux(0.0,0.0);
358 for(
unsigned l=0;l<n_node;l++)
364 if(local_eqn_real >= 0)
367 residuals[local_eqn_real] -= flux.real()*testf[l]*r*
W;
377 if(local_eqn_imag >= 0)
380 residuals[local_eqn_imag] -= flux.imag()*testf[l]*r*
W;
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 broken_copy(const std::string &class_name)
Issue error message and terminate execution.
void output(FILE *file_pt)
FourierDecomposedHelmholtzFluxElement(const FourierDecomposedHelmholtzFluxElement &dummy)
Broken copy constructor.
FourierDecomposedHelmholtzPrescribedFluxFctPt & flux_fct_pt()
Broken assignment operator.
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 ...
void output(std::ostream &outfile)
int & face_index()
Index of the face (a number that uniquely identifies the face in the element)
A general Finite Element class.
virtual std::complex< unsigned > u_index_fourier_decomposed_helmholtz() const
Return the index at which the unknown value is stored. (Real/imag part gives real index of real/imag ...
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...
virtual double weight(const unsigned &i) const =0
Return weight of i-th integration point.
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.
double J_eulerian(const Vector< double > &s) const
Return the Jacobian of mapping from local to global coordinates at local position s...
void fill_in_contribution_to_residuals(Vector< double > &residuals)
Add the element's contribution to its residual vector.
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.
FourierDecomposedHelmholtzFluxElement()
Broken empty constructor.
std::complex< unsigned > U_index_fourier_decomposed_helmholtz
The index at which the real and imag part of the unknown is stored at the nodes.
FourierDecomposedHelmholtzPrescribedFluxFctPt Flux_fct_pt
Function pointer to the (global) prescribed-flux function.
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(* FourierDecomposedHelmholtzPrescribedFluxFctPt)(const Vector< double > &x, std::complex< double > &flux)
Function pointer to the prescribed-flux function fct(x,f(x)) – x is a Vector and the flux is a compl...
A class for elements that allow the imposition of an applied flux on the boundaries of Fourier decomp...
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 unsigned nweight() const =0
Return the number of integration points of the scheme.
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 ...
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...
void get_flux(const Vector< double > &x, std::complex< double > &flux)
unsigned nnode() const
Return the number of nodes.
virtual void fill_in_generic_residual_contribution_fourier_decomposed_helmholtz_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...
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...
void output(std::ostream &outfile, const unsigned &n_plot)
Output function – forward to broken version in FiniteElement until somebody decides what exactly the...
static DenseMatrix< double > Dummy_matrix
Empty dense matrix used as a dummy argument to combined residual and jacobian functions in the case w...
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...