//////////////////////////////////////////////////////////////////// //////////////////////////////////////////////////////////////////// More...
#include <specific_node_update_interface_elements.h>
Public Member Functions | |
ElasticUpdateFluidInterfaceElement (FiniteElement *const &element_pt, const int &face_index, const unsigned &id=0) | |
Constructor, pass a pointer to the bulk element and the face index of the bulk element to which the element is to be attached to. The optional identifier can be used to distinguish the additional nodal value (Lagr mult) created by this element from those created by other FaceElements. More... | |
double | zeta_nodal (const unsigned &n, const unsigned &k, const unsigned &i) const |
The "global" intrinsic coordinate of the element when viewed as part of a geometric object should be given by the FaceElement representation, by default. More... | |
double & | lagrange (const unsigned &n) |
Return the lagrange multiplier at local node n. More... | |
void | fill_in_contribution_to_jacobian (Vector< double > &residuals, DenseMatrix< double > &jacobian) |
Fill in contribution to residuals and Jacobian. More... | |
void | output (std::ostream &outfile) |
Overload the output function. More... | |
void | output (std::ostream &outfile, const unsigned &n_plot) |
Output the element. More... | |
void | output (FILE *file_pt) |
Overload the C-style output function. More... | |
void | output (FILE *file_pt, const unsigned &n_plot) |
C-style Output function. More... | |
void | add_additional_residual_contributions_interface (Vector< double > &residuals, DenseMatrix< double > &jacobian, const unsigned &flag, const Shape &psif, const DShape &dpsifds, const DShape &dpsifdS, const DShape &dpsifdS_div, const Vector< double > &s, const Vector< double > &interpolated_x, const Vector< double > &interpolated_n, const double &W, const double &J) |
Helper function to calculate the additional contributions to be added at each integration point. This deals with Lagrange multiplier contribution, as well as any additional contributions by the other equations. More... | |
virtual FluidInterfaceBoundingElement * | make_bounding_element (const int &face_index) |
Create an "bounding" element (here actually a 2D line element of type ElasticLineFluidInterfaceBoundingElement<ELEMENT> that allows the application of a contact angle boundary condition on the the specified face. More... | |
Public Member Functions inherited from oomph::Hijacked< FaceGeometry< ELEMENT > > | |
Hijacked () | |
Constructor, call the constructors of the base elements. More... | |
Hijacked (FiniteElement *const &element_pt, const int &face_index) | |
Constructor used for hijacking face elements. More... | |
Hijacked (FiniteElement *const &element_pt, const int &face_index, const unsigned &id=0) | |
Constructor used for hijacking face elements with specification of ID of additional variables. More... | |
Data * | hijack_internal_value (const unsigned &n, const unsigned &i, const bool &return_data=true) |
Hijack the i-th value stored at internal data n. Optionally return a custom-made (copied) data object that contains only the hijacked value. This can be used as the input to other elements. Note that the calling program assumes responsibility for this data object and must clean it up. The default is that the data object is returned. More... | |
Data * | hijack_external_value (const unsigned &n, const unsigned &i, const bool &return_data=true) |
Hijack the i-th value stored at external data n. Optionally return a custom-made (copied) data object that contains only the hijacked value. Note that the calling program assumes responsibility for this data object and must clean it up. The default is that the data object is returned. More... | |
Data * | hijack_nodal_value (const unsigned &n, const unsigned &i, const bool &return_data=true) |
Hijack the i-th value stored at node n. Optionally return a custom-made (copied) data object that contains only the hijacked value. Once again, the calling program must clean up the allocated Data object. The default is that the data object is returned. More... | |
Data * | hijack_nodal_position_value (const unsigned &n, const unsigned &i, const bool &return_data=true) |
Hijack the i-th positional value stored at node n. Optionaly return a custom-made (copied) data object that contains only the hijacked value. Again, responsibility for the memory allocated lies with the calling function. The default is that the data object is returned. More... | |
Data * | hijack_nodal_spine_value (const unsigned &n, const unsigned &i, const bool &return_data=true) |
Hijack the i-th value stored at the spine that affects local node n. Optionally return a custom-made (copied) data object that contains only the hijacked value. Deletion must be handled at the higher level The default is that the data object is returned. More... | |
void | assign_local_eqn_numbers (const bool &store_local_dof_pt) |
Set up the local equation numbers for the underlying element, then set up the local arrays to hold the hijacked variables. If the boolean argument is true then pointers to the associated degrees of freedom are stored in the array Dof_pt. More... | |
void | get_residuals (Vector< double > &residuals) |
Get the residuals from the underlying element, but then wipe the entries in the residual vector that correspond to hijacked values – they will be computed by other elements. More... | |
void | get_jacobian (Vector< double > &residuals, DenseMatrix< double > &jacobian) |
Get the residuals and Jacobian matrix from the underlying element, but then wipe the entries in the residual vector and the rows in the Jacobian matrix that correspond to hijacked values – they will be computed by other elements. More... | |
Public Member Functions inherited from oomph::HijackedElementBase | |
HijackedElementBase () | |
Constructor, initialise the pointer to the equation numbers for the storage to zero. More... | |
virtual | ~HijackedElementBase () |
Destructor, destroy the storage for the equation numbers. More... | |
void | unhijack_all_data () |
Reset the hijacked data pt, so that none of the equations in the element are hijacked. More... | |
const double & | residual_multiplier () const |
Return the value of the residual multiplier. More... | |
double *& | residual_multiplier_pt () |
Return the pointer to the residual multiplier. More... | |
Protected Member Functions | |
double | compute_surface_derivatives (const Shape &psi, const DShape &dpsids, const DenseMatrix< double > &interpolated_t, const Vector< double > &interpolated_x, DShape &surface_gradient, DShape &surface_divergence) |
Fill in the specific surface derivative calculations by calling the appropriate function from the derivative class. More... | |
Protected Member Functions inherited from oomph::HijackedElementBase | |
void | hijack_global_eqn (long *const &global_eqn_pt) |
Mark the global equation, addressed by global_eqn_pt, as hijacked by this element. More... | |
void | unhijack_global_eqn (long *const &global_eqn_pt) |
The global equation, addressed by global_eqn_pt, is no longer hijacked by this element. More... | |
Private Member Functions | |
unsigned | lagrange_index (const unsigned &n) |
Return the index at which the lagrange multiplier is stored at the n-th node. More... | |
int | kinematic_local_eqn (const unsigned &n) |
Equation number of the kinematic BC associated with node j. (This is the equation for the Lagrange multiplier) More... | |
void | hijack_kinematic_conditions (const Vector< unsigned > &bulk_node_number) |
Hijacking the kinematic condition corresponds to hijacking the variables associated with the Lagrange multipliers that are assigned on construction of this element. More... | |
Private Attributes | |
Vector< unsigned > | Lagrange_index |
Storage for the location of the Lagrange multiplier (If other additional values have been added we need to add the Lagrange multiplier at the end) More... | |
Additional Inherited Members | |
Protected Attributes inherited from oomph::HijackedElementBase | |
std::set< long * > * | Hijacked_global_eqn_number_pt |
Pointer to a Set of pointers to the equation numbers that will be hijacked by this element. Note that these MUST be pointers because hijacking information is set BEFORE global equation numbers have been assigned. More... | |
Vector< int > * | Hijacked_local_eqn_number_pt |
Pointer to a vector of integers containing the local equation numbers of any hijacked variables in the element. More... | |
double * | Residual_multiplier_pt |
Pointer to a double that multiplies the contribution to the residuals from the original element. This is usually used as a homotopy parameter to permit a smooth transition between different types of boundary conditions, rather than switching them on or off abruptly. More... | |
Static Protected Attributes inherited from oomph::HijackedElementBase | |
static double | Default_residual_multiplier = 0.0 |
Static default value for the double that multiplies the original residuals. More... | |
//////////////////////////////////////////////////////////////////// ////////////////////////////////////////////////////////////////////
Generic Elastic node update interface template class that can be combined with a given surface equations class and surface derivative class to provide a concrete implementation of any surface element that uses elastic node updates.
Definition at line 672 of file specific_node_update_interface_elements.h.
|
inline |
Constructor, pass a pointer to the bulk element and the face index of the bulk element to which the element is to be attached to. The optional identifier can be used to distinguish the additional nodal value (Lagr mult) created by this element from those created by other FaceElements.
Definition at line 739 of file specific_node_update_interface_elements.h.
References oomph::FiniteElement::build_face_element(), i, oomph::ElasticUpdateFluidInterfaceElement< EQUATION_CLASS, DERIVATIVE_CLASS, ELEMENT >::Lagrange_index, oomph::FluidInterfaceAdditionalValues< ELEMENT >::nadditional_values(), and oomph::FluidInterfaceAdditionalValues< ELEMENT >::setup_equation_indices().
|
inline |
Helper function to calculate the additional contributions to be added at each integration point. This deals with Lagrange multiplier contribution, as well as any additional contributions by the other equations.
Definition at line 877 of file specific_node_update_interface_elements.h.
References i, oomph::ElasticUpdateFluidInterfaceElement< EQUATION_CLASS, DERIVATIVE_CLASS, ELEMENT >::kinematic_local_eqn(), oomph::ElasticUpdateFluidInterfaceElement< EQUATION_CLASS, DERIVATIVE_CLASS, ELEMENT >::lagrange(), s, and oomph::QuadTreeNames::W.
|
inlineprotected |
Fill in the specific surface derivative calculations by calling the appropriate function from the derivative class.
Definition at line 716 of file specific_node_update_interface_elements.h.
|
inline |
Fill in contribution to residuals and Jacobian.
Definition at line 837 of file specific_node_update_interface_elements.h.
|
inlineprivate |
Hijacking the kinematic condition corresponds to hijacking the variables associated with the Lagrange multipliers that are assigned on construction of this element.
Definition at line 701 of file specific_node_update_interface_elements.h.
References oomph::Hijacked< FaceGeometry< ELEMENT > >::hijack_nodal_value(), and oomph::ElasticUpdateFluidInterfaceElement< EQUATION_CLASS, DERIVATIVE_CLASS, ELEMENT >::lagrange_index().
|
inlineprivate |
Equation number of the kinematic BC associated with node j. (This is the equation for the Lagrange multiplier)
Definition at line 692 of file specific_node_update_interface_elements.h.
References oomph::ElasticUpdateFluidInterfaceElement< EQUATION_CLASS, DERIVATIVE_CLASS, ELEMENT >::lagrange_index().
|
inline |
Return the lagrange multiplier at local node n.
Definition at line 830 of file specific_node_update_interface_elements.h.
References oomph::ElasticUpdateFluidInterfaceElement< EQUATION_CLASS, DERIVATIVE_CLASS, ELEMENT >::lagrange_index().
|
inlineprivate |
Return the index at which the lagrange multiplier is stored at the n-th node.
Definition at line 685 of file specific_node_update_interface_elements.h.
Referenced by oomph::ElasticUpdateFluidInterfaceElement< EQUATION_CLASS, DERIVATIVE_CLASS, ELEMENT >::hijack_kinematic_conditions(), oomph::ElasticUpdateFluidInterfaceElement< EQUATION_CLASS, DERIVATIVE_CLASS, ELEMENT >::kinematic_local_eqn(), and oomph::ElasticUpdateFluidInterfaceElement< EQUATION_CLASS, DERIVATIVE_CLASS, ELEMENT >::lagrange().
|
inlinevirtual |
Create an "bounding" element (here actually a 2D line element of type ElasticLineFluidInterfaceBoundingElement<ELEMENT> that allows the application of a contact angle boundary condition on the the specified face.
Definition at line 960 of file specific_node_update_interface_elements.h.
|
inline |
Overload the C-style output function.
Definition at line 861 of file specific_node_update_interface_elements.h.
References oomph::CommandLineArgs::output().
|
inline |
C-style Output function.
Definition at line 867 of file specific_node_update_interface_elements.h.
References oomph::CommandLineArgs::output().
|
inline |
Overload the output function.
Definition at line 849 of file specific_node_update_interface_elements.h.
References oomph::CommandLineArgs::output().
|
inline |
Output the element.
Definition at line 855 of file specific_node_update_interface_elements.h.
References oomph::CommandLineArgs::output().
|
inline |
The "global" intrinsic coordinate of the element when viewed as part of a geometric object should be given by the FaceElement representation, by default.
Definition at line 822 of file specific_node_update_interface_elements.h.
References i, and oomph::FaceElement::zeta_nodal().
|
private |
Storage for the location of the Lagrange multiplier (If other additional values have been added we need to add the Lagrange multiplier at the end)
Definition at line 681 of file specific_node_update_interface_elements.h.
Referenced by oomph::ElasticUpdateFluidInterfaceElement< EQUATION_CLASS, DERIVATIVE_CLASS, ELEMENT >::ElasticUpdateFluidInterfaceElement().