Public Types | Public Member Functions | Protected Member Functions | Static Protected Attributes | List of all members
oomph::RefineableQElement< 2 > Class Referenceabstract

#include <refineable_quad_element.h>

+ Inheritance diagram for oomph::RefineableQElement< 2 >:

Public Types

typedef void(RefineableQElement< 2 >::* VoidMemberFctPt) ()
 
- Public Types inherited from oomph::FiniteElement
typedef void(* SteadyExactSolutionFctPt) (const Vector< double > &, Vector< double > &)
 
typedef void(* UnsteadyExactSolutionFctPt) (const double &, const Vector< double > &, Vector< double > &)
 

Public Member Functions

 RefineableQElement ()
 Constructor: Pass refinement level (default 0 = root) More...
 
 RefineableQElement (const RefineableQElement< 2 > &dummy)=delete
 Broken copy constructor. More...
 
virtual ~RefineableQElement ()
 Broken assignment operator. More...
 
unsigned required_nsons () const
 A refineable quad element has four sons. More...
 
virtual Node * node_created_by_neighbour (const Vector< double > &s_fraction, bool &is_periodic)
 
virtual Node * node_created_by_son_of_neighbour (const Vector< double > &s_fraction, bool &is_periodic)
 
virtual void build (Mesh *&mesh_pt, Vector< Node * > &new_node_pt, bool &was_already_built, std::ofstream &new_nodes_file)
 
void check_integrity (double &max_error)
 
void output_corners (std::ostream &outfile, const std::string &colour) const
 Print corner nodes, use colour. More...
 
QuadTree * quadtree_pt ()
 Pointer to quadtree representation of this element. More...
 
QuadTree * quadtree_pt () const
 Pointer to quadtree representation of this element. More...
 
void setup_hanging_nodes (Vector< std::ofstream * > &output_stream)
 
virtual void further_setup_hanging_nodes ()=0
 
void get_boundaries (const int &edge, std::set< unsigned > &boundaries) const
 
void get_bcs (int bound, Vector< int > &bound_cons) const
 
void interpolated_zeta_on_edge (const unsigned &boundary, const int &edge, const Vector< double > &s, Vector< double > &zeta)
 
- Public Member Functions inherited from oomph::RefineableElement
 RefineableElement ()
 
virtual ~RefineableElement ()
 Destructor, delete the allocated storage for the hanging equations. More...
 
 RefineableElement (const RefineableElement &)=delete
 Broken copy constructor. More...
 
void operator= (const RefineableElement &)=delete
 Broken assignment operator. More...
 
Tree * tree_pt ()
 Access function: Pointer to quadtree representation of this element. More...
 
void set_tree_pt (Tree *my_tree_pt)
 Set pointer to quadtree representation of this element. More...
 
bool refinement_is_enabled ()
 Flag to indicate suppression of any refinement. More...
 
void disable_refinement ()
 Suppress of any refinement for this element. More...
 
void enable_refinement ()
 Emnable refinement for this element. More...
 
template<class ELEMENT >
void split (Vector< ELEMENT * > &son_pt) const
 
int local_hang_eqn (Node *const &node_pt, const unsigned &i)
 
void set_refinement_level (const int &refine_level)
 Set the refinement level. More...
 
unsigned refinement_level () const
 Return the Refinement level. More...
 
void select_for_refinement ()
 Select the element for refinement. More...
 
void deselect_for_refinement ()
 Deselect the element for refinement. More...
 
void select_sons_for_unrefinement ()
 Unrefinement will be performed by merging the four sons of this element. More...
 
void deselect_sons_for_unrefinement ()
 
bool to_be_refined ()
 Has the element been selected for refinement? More...
 
bool sons_to_be_unrefined ()
 Has the element been selected for unrefinement? More...
 
virtual void rebuild_from_sons (Mesh *&mesh_pt)=0
 
virtual void unbuild ()
 
virtual void deactivate_element ()
 
virtual bool nodes_built ()
 Return true if all the nodes have been built, false if not. More...
 
long number () const
 Element number (for debugging/plotting) More...
 
void set_number (const long &mynumber)
 Set element number (for debugging/plotting) More...
 
virtual unsigned ncont_interpolated_values () const =0
 
virtual void get_interpolated_values (const Vector< double > &s, Vector< double > &values)
 
virtual void get_interpolated_values (const unsigned &t, const Vector< double > &s, Vector< double > &values)=0
 
virtual Node * interpolating_node_pt (const unsigned &n, const int &value_id)
 
virtual double local_one_d_fraction_of_interpolating_node (const unsigned &n1d, const unsigned &i, const int &value_id)
 
virtual Node * get_interpolating_node_at_local_coordinate (const Vector< double > &s, const int &value_id)
 
virtual unsigned ninterpolating_node (const int &value_id)
 
virtual unsigned ninterpolating_node_1d (const int &value_id)
 
virtual void interpolating_basis (const Vector< double > &s, Shape &psi, const int &value_id) const
 
void identify_field_data_for_interactions (std::set< std::pair< Data *, unsigned >> &paired_field_data)
 
void assign_nodal_local_eqn_numbers (const bool &store_local_dof_pt)
 
virtual RefineableElement * root_element_pt ()
 
virtual RefineableElement * father_element_pt () const
 Return a pointer to the father element. More...
 
void get_father_at_refinement_level (unsigned &refinement_level, RefineableElement *&father_at_reflevel_pt)
 
virtual void initial_setup (Tree *const &adopted_father_pt=0, const unsigned &initial_p_order=0)
 
virtual void pre_build (Mesh *&mesh_pt, Vector< Node * > &new_node_pt)
 Pre-build the element. More...
 
virtual void further_build ()
 Further build: e.g. deal with interpolation of internal values. More...
 
void get_dresidual_dnodal_coordinates (RankThreeTensor< double > &dresidual_dnodal_coordinates)
 
unsigned nshape_controlling_nodes ()
 
std::map< Node *, unsigned > shape_controlling_node_lookup ()
 
- Public Member Functions inherited from oomph::FiniteElement
void set_dimension (const unsigned &dim)
 
void set_nodal_dimension (const unsigned &nodal_dim)
 
void set_nnodal_position_type (const unsigned &nposition_type)
 Set the number of types required to interpolate the coordinate. More...
 
void set_n_node (const unsigned &n)
 
int nodal_local_eqn (const unsigned &n, const unsigned &i) const
 
double dJ_eulerian_at_knot (const unsigned &ipt, Shape &psi, DenseMatrix< double > &djacobian_dX) const
 
 FiniteElement ()
 Constructor. More...
 
virtual ~FiniteElement ()
 
 FiniteElement (const FiniteElement &)=delete
 Broken copy constructor. More...
 
void get_centre_of_gravity_and_max_radius_in_terms_of_zeta (Vector< double > &cog, double &max_radius) const
 
virtual void local_coordinate_of_node (const unsigned &j, Vector< double > &s) const
 
virtual void local_fraction_of_node (const unsigned &j, Vector< double > &s_fraction)
 
virtual double local_one_d_fraction_of_node (const unsigned &n1d, const unsigned &i)
 
MacroElement * macro_elem_pt ()
 Access function to pointer to macro element. More...
 
void get_x (const Vector< double > &s, Vector< double > &x) const
 
void get_x (const unsigned &t, const Vector< double > &s, Vector< double > &x)
 
virtual void set_integration_scheme (Integral *const &integral_pt)
 Set the spatial integration scheme. More...
 
Integral *const & integral_pt () const
 Return the pointer to the integration scheme (const version) More...
 
virtual void shape (const Vector< double > &s, Shape &psi) const =0
 
virtual void shape_at_knot (const unsigned &ipt, Shape &psi) const
 
virtual void dshape_local (const Vector< double > &s, Shape &psi, DShape &dpsids) const
 
virtual void dshape_local_at_knot (const unsigned &ipt, Shape &psi, DShape &dpsids) const
 
virtual void d2shape_local (const Vector< double > &s, Shape &psi, DShape &dpsids, DShape &d2psids) const
 
virtual void d2shape_local_at_knot (const unsigned &ipt, Shape &psi, DShape &dpsids, DShape &d2psids) const
 
virtual double J_eulerian (const Vector< double > &s) const
 
virtual double J_eulerian_at_knot (const unsigned &ipt) const
 
void check_J_eulerian_at_knots (bool &passed) const
 
void check_jacobian (const double &jacobian) const
 
double dshape_eulerian (const Vector< double > &s, Shape &psi, DShape &dpsidx) const
 
virtual double dshape_eulerian_at_knot (const unsigned &ipt, Shape &psi, DShape &dpsidx) const
 
virtual double dshape_eulerian_at_knot (const unsigned &ipt, Shape &psi, DShape &dpsi, DenseMatrix< double > &djacobian_dX, RankFourTensor< double > &d_dpsidx_dX) const
 
double d2shape_eulerian (const Vector< double > &s, Shape &psi, DShape &dpsidx, DShape &d2psidx) const
 
virtual double d2shape_eulerian_at_knot (const unsigned &ipt, Shape &psi, DShape &dpsidx, DShape &d2psidx) const
 
virtual void describe_local_dofs (std::ostream &out, const std::string &current_string) const
 
virtual void describe_nodal_local_dofs (std::ostream &out, const std::string &current_string) const
 
virtual void assign_all_generic_local_eqn_numbers (const bool &store_local_dof_pt)
 
Node *& node_pt (const unsigned &n)
 Return a pointer to the local node n. More...
 
Node *const & node_pt (const unsigned &n) const
 Return a pointer to the local node n (const version) More...
 
unsigned nnode () const
 Return the number of nodes. More...
 
virtual unsigned nnode_1d () const
 
double raw_nodal_position (const unsigned &n, const unsigned &i) const
 
double raw_nodal_position (const unsigned &t, const unsigned &n, const unsigned &i) const
 
double raw_dnodal_position_dt (const unsigned &n, const unsigned &i) const
 
double raw_dnodal_position_dt (const unsigned &n, const unsigned &j, const unsigned &i) const
 
double raw_nodal_position_gen (const unsigned &n, const unsigned &k, const unsigned &i) const
 
double raw_nodal_position_gen (const unsigned &t, const unsigned &n, const unsigned &k, const unsigned &i) const
 
double raw_dnodal_position_gen_dt (const unsigned &n, const unsigned &k, const unsigned &i) const
 
double raw_dnodal_position_gen_dt (const unsigned &j, const unsigned &n, const unsigned &k, const unsigned &i) const
 
double nodal_position (const unsigned &n, const unsigned &i) const
 
double nodal_position (const unsigned &t, const unsigned &n, const unsigned &i) const
 
double dnodal_position_dt (const unsigned &n, const unsigned &i) const
 Return the i-th component of nodal velocity: dx/dt at local node n. More...
 
double dnodal_position_dt (const unsigned &n, const unsigned &j, const unsigned &i) const
 
double nodal_position_gen (const unsigned &n, const unsigned &k, const unsigned &i) const
 
double nodal_position_gen (const unsigned &t, const unsigned &n, const unsigned &k, const unsigned &i) const
 
double dnodal_position_gen_dt (const unsigned &n, const unsigned &k, const unsigned &i) const
 
double dnodal_position_gen_dt (const unsigned &j, const unsigned &n, const unsigned &k, const unsigned &i) const
 
virtual void disable_ALE ()
 
virtual void enable_ALE ()
 
virtual unsigned required_nvalue (const unsigned &n) const
 
unsigned nnodal_position_type () const
 
bool has_hanging_nodes () const
 
unsigned nodal_dimension () const
 Return the required Eulerian dimension of the nodes in this element. More...
 
virtual Node * construct_node (const unsigned &n)
 
virtual Node * construct_node (const unsigned &n, TimeStepper *const &time_stepper_pt)
 
virtual Node * construct_boundary_node (const unsigned &n)
 
virtual Node * construct_boundary_node (const unsigned &n, TimeStepper *const &time_stepper_pt)
 
int get_node_number (Node *const &node_pt) const
 
virtual Node * get_node_at_local_coordinate (const Vector< double > &s) const
 
double raw_nodal_value (const unsigned &n, const unsigned &i) const
 
double raw_nodal_value (const unsigned &t, const unsigned &n, const unsigned &i) const
 
double nodal_value (const unsigned &n, const unsigned &i) const
 
double nodal_value (const unsigned &t, const unsigned &n, const unsigned &i) const
 
unsigned dim () const
 
virtual double interpolated_x (const Vector< double > &s, const unsigned &i) const
 Return FE interpolated coordinate x[i] at local coordinate s. More...
 
virtual double interpolated_x (const unsigned &t, const Vector< double > &s, const unsigned &i) const
 
virtual void interpolated_x (const Vector< double > &s, Vector< double > &x) const
 Return FE interpolated position x[] at local coordinate s as Vector. More...
 
virtual void interpolated_x (const unsigned &t, const Vector< double > &s, Vector< double > &x) const
 
virtual double interpolated_dxdt (const Vector< double > &s, const unsigned &i, const unsigned &t)
 
virtual void interpolated_dxdt (const Vector< double > &s, const unsigned &t, Vector< double > &dxdt)
 
unsigned ngeom_data () const
 
Data * geom_data_pt (const unsigned &j)
 
void position (const Vector< double > &zeta, Vector< double > &r) const
 
void position (const unsigned &t, const Vector< double > &zeta, Vector< double > &r) const
 
void dposition_dt (const Vector< double > &zeta, const unsigned &t, Vector< double > &drdt)
 
virtual double zeta_nodal (const unsigned &n, const unsigned &k, const unsigned &i) const
 
void interpolated_zeta (const Vector< double > &s, Vector< double > &zeta) const
 
void locate_zeta (const Vector< double > &zeta, GeomObject *&geom_object_pt, Vector< double > &s, const bool &use_coordinate_as_initial_guess=false)
 
virtual void node_update ()
 
virtual void identify_geometric_data (std::set< Data * > &geometric_data_pt)
 
virtual double s_min () const
 Min value of local coordinate. More...
 
virtual double s_max () const
 Max. value of local coordinate. More...
 
double size () const
 
virtual double compute_physical_size () const
 
virtual void point_output_data (const Vector< double > &s, Vector< double > &data)
 
void point_output (std::ostream &outfile, const Vector< double > &s)
 
virtual unsigned nplot_points_paraview (const unsigned &nplot) const
 
virtual unsigned nsub_elements_paraview (const unsigned &nplot) const
 
void output_paraview (std::ofstream &file_out, const unsigned &nplot) const
 
virtual void write_paraview_output_offset_information (std::ofstream &file_out, const unsigned &nplot, unsigned &counter) const
 
virtual void write_paraview_type (std::ofstream &file_out, const unsigned &nplot) const
 
virtual void write_paraview_offsets (std::ofstream &file_out, const unsigned &nplot, unsigned &offset_sum) const
 
virtual unsigned nscalar_paraview () const
 
virtual void scalar_value_paraview (std::ofstream &file_out, const unsigned &i, const unsigned &nplot) const
 
virtual void scalar_value_fct_paraview (std::ofstream &file_out, const unsigned &i, const unsigned &nplot, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt) const
 
virtual void scalar_value_fct_paraview (std::ofstream &file_out, const unsigned &i, const unsigned &nplot, const double &time, FiniteElement::UnsteadyExactSolutionFctPt exact_soln_pt) const
 
virtual std::string scalar_name_paraview (const unsigned &i) const
 
virtual void output (std::ostream &outfile)
 
virtual void output (std::ostream &outfile, const unsigned &n_plot)
 
virtual void output (const unsigned &t, std::ostream &outfile, const unsigned &n_plot) const
 
virtual void output (FILE *file_pt)
 
virtual void output (FILE *file_pt, const unsigned &n_plot)
 
virtual void output_fct (std::ostream &outfile, const unsigned &n_plot, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt)
 Output an exact solution over the element. More...
 
virtual void output_fct (std::ostream &outfile, const unsigned &n_plot, const double &time, FiniteElement::UnsteadyExactSolutionFctPt exact_soln_pt)
 Output a time-dependent exact solution over the element. More...
 
virtual void output_fct (std::ostream &outfile, const unsigned &n_plot, const double &time, const SolutionFunctorBase &exact_soln) const
 Output a time-dependent exact solution over the element. More...
 
virtual void get_s_plot (const unsigned &i, const unsigned &nplot, Vector< double > &s, const bool &shifted_to_interior=false) const
 
virtual std::string tecplot_zone_string (const unsigned &nplot) const
 
virtual void write_tecplot_zone_footer (std::ostream &outfile, const unsigned &nplot) const
 
virtual void write_tecplot_zone_footer (FILE *file_pt, const unsigned &nplot) const
 
virtual unsigned nplot_points (const unsigned &nplot) const
 
virtual void compute_error (FiniteElement::SteadyExactSolutionFctPt exact_soln_pt, double &error, double &norm)
 Calculate the norm of the error and that of the exact solution. More...
 
virtual void compute_error (FiniteElement::UnsteadyExactSolutionFctPt exact_soln_pt, const double &time, double &error, double &norm)
 Calculate the norm of the error and that of the exact solution. More...
 
virtual void compute_error (FiniteElement::SteadyExactSolutionFctPt exact_soln_pt, Vector< double > &error, Vector< double > &norm)
 
virtual void compute_error (FiniteElement::UnsteadyExactSolutionFctPt exact_soln_pt, const double &time, Vector< double > &error, Vector< double > &norm)
 
virtual void compute_error (std::ostream &outfile, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt, double &error, double &norm)
 
virtual void compute_error (std::ostream &outfile, FiniteElement::UnsteadyExactSolutionFctPt exact_soln_pt, const double &time, double &error, double &norm)
 
virtual void compute_error (std::ostream &outfile, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt, Vector< double > &error, Vector< double > &norm)
 
virtual void compute_error (std::ostream &outfile, FiniteElement::UnsteadyExactSolutionFctPt exact_soln_pt, const double &time, Vector< double > &error, Vector< double > &norm)
 
virtual void compute_abs_error (std::ostream &outfile, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt, double &error)
 
void integrate_fct (FiniteElement::SteadyExactSolutionFctPt integrand_fct_pt, Vector< double > &integral)
 Integrate Vector-valued function over element. More...
 
void integrate_fct (FiniteElement::UnsteadyExactSolutionFctPt integrand_fct_pt, const double &time, Vector< double > &integral)
 Integrate Vector-valued time-dep function over element. More...
 
virtual void build_face_element (const int &face_index, FaceElement *face_element_pt)
 
virtual unsigned self_test ()
 
virtual unsigned get_bulk_node_number (const int &face_index, const unsigned &i) const
 
virtual int face_outer_unit_normal_sign (const int &face_index) const
 Get the sign of the outer unit normal on the face given by face_index. More...
 
void face_node_number_error_check (const unsigned &i) const
 Range check for face node numbers. More...
 
virtual CoordinateMappingFctPt face_to_bulk_coordinate_fct_pt (const int &face_index) const
 
virtual BulkCoordinateDerivativesFctPt bulk_coordinate_derivatives_fct_pt (const int &face_index) const
 
- Public Member Functions inherited from oomph::GeneralisedElement
 GeneralisedElement ()
 Constructor: Initialise all pointers and all values to zero. More...
 
virtual ~GeneralisedElement ()
 Virtual destructor to clean up any memory allocated by the object. More...
 
 GeneralisedElement (const GeneralisedElement &)=delete
 Broken copy constructor. More...
 
void operator= (const GeneralisedElement &)=delete
 Broken assignment operator. More...
 
Data *& internal_data_pt (const unsigned &i)
 Return a pointer to i-th internal data object. More...
 
Data *const & internal_data_pt (const unsigned &i) const
 Return a pointer to i-th internal data object (const version) More...
 
Data *& external_data_pt (const unsigned &i)
 Return a pointer to i-th external data object. More...
 
Data *const & external_data_pt (const unsigned &i) const
 Return a pointer to i-th external data object (const version) More...
 
unsigned long eqn_number (const unsigned &ieqn_local) const
 
int local_eqn_number (const unsigned long &ieqn_global) const
 
unsigned add_external_data (Data *const &data_pt, const bool &fd=true)
 
bool external_data_fd (const unsigned &i) const
 
void exclude_external_data_fd (const unsigned &i)
 
void include_external_data_fd (const unsigned &i)
 
void flush_external_data ()
 Flush all external data. More...
 
void flush_external_data (Data *const &data_pt)
 Flush the object addressed by data_pt from the external data array. More...
 
unsigned ninternal_data () const
 Return the number of internal data objects. More...
 
unsigned nexternal_data () const
 Return the number of external data objects. More...
 
unsigned ndof () const
 Return the number of equations/dofs in the element. More...
 
void dof_vector (const unsigned &t, Vector< double > &dof)
 Return the vector of dof values at time level t. More...
 
void dof_pt_vector (Vector< double * > &dof_pt)
 Return the vector of pointers to dof values. More...
 
void set_internal_data_time_stepper (const unsigned &i, TimeStepper *const &time_stepper_pt, const bool &preserve_existing_data)
 
void assign_internal_eqn_numbers (unsigned long &global_number, Vector< double * > &Dof_pt)
 
void describe_dofs (std::ostream &out, const std::string &current_string) const
 
void add_internal_value_pt_to_map (std::map< unsigned, double * > &map_of_value_pt)
 
virtual void assign_local_eqn_numbers (const bool &store_local_dof_pt)
 
virtual void complete_setup_of_dependencies ()
 
virtual void get_residuals (Vector< double > &residuals)
 
virtual void get_jacobian (Vector< double > &residuals, DenseMatrix< double > &jacobian)
 
virtual void get_mass_matrix (Vector< double > &residuals, DenseMatrix< double > &mass_matrix)
 
virtual void get_jacobian_and_mass_matrix (Vector< double > &residuals, DenseMatrix< double > &jacobian, DenseMatrix< double > &mass_matrix)
 
virtual void get_dresiduals_dparameter (double *const &parameter_pt, Vector< double > &dres_dparam)
 
virtual void get_djacobian_dparameter (double *const &parameter_pt, Vector< double > &dres_dparam, DenseMatrix< double > &djac_dparam)
 
virtual void get_djacobian_and_dmass_matrix_dparameter (double *const &parameter_pt, Vector< double > &dres_dparam, DenseMatrix< double > &djac_dparam, DenseMatrix< double > &dmass_matrix_dparam)
 
virtual void get_hessian_vector_products (Vector< double > const &Y, DenseMatrix< double > const &C, DenseMatrix< double > &product)
 
virtual void get_inner_products (Vector< std::pair< unsigned, unsigned >> const &history_index, Vector< double > &inner_product)
 
virtual void get_inner_product_vectors (Vector< unsigned > const &history_index, Vector< Vector< double >> &inner_product_vector)
 
virtual void compute_norm (Vector< double > &norm)
 
virtual void compute_norm (double &norm)
 
virtual unsigned ndof_types () const
 
virtual void get_dof_numbers_for_unknowns (std::list< std::pair< unsigned long, unsigned >> &dof_lookup_list) const
 
- Public Member Functions inherited from oomph::GeomObject
 GeomObject ()
 Default constructor. More...
 
 GeomObject (const unsigned &ndim)
 
 GeomObject (const unsigned &nlagrangian, const unsigned &ndim)
 
 GeomObject (const unsigned &nlagrangian, const unsigned &ndim, TimeStepper *time_stepper_pt)
 
 GeomObject (const GeomObject &dummy)=delete
 Broken copy constructor. More...
 
void operator= (const GeomObject &)=delete
 Broken assignment operator. More...
 
virtual ~GeomObject ()
 (Empty) destructor More...
 
unsigned nlagrangian () const
 Access function to # of Lagrangian coordinates. More...
 
unsigned ndim () const
 Access function to # of Eulerian coordinates. More...
 
void set_nlagrangian_and_ndim (const unsigned &n_lagrangian, const unsigned &n_dim)
 Set # of Lagrangian and Eulerian coordinates. More...
 
TimeStepper *& time_stepper_pt ()
 
TimeStepper * time_stepper_pt () const
 
virtual void position (const double &t, const Vector< double > &zeta, Vector< double > &r) const
 
virtual void dposition (const Vector< double > &zeta, DenseMatrix< double > &drdzeta) const
 
virtual void d2position (const Vector< double > &zeta, RankThreeTensor< double > &ddrdzeta) const
 
virtual void d2position (const Vector< double > &zeta, Vector< double > &r, DenseMatrix< double > &drdzeta, RankThreeTensor< double > &ddrdzeta) const
 
- Public Member Functions inherited from oomph::QuadElementBase
 QuadElementBase ()
 Constructor. Empty. More...
 
virtual unsigned nvertex_node () const =0
 Number of vertex nodes in the element. More...
 
virtual Node * vertex_node_pt (const unsigned &j) const =0
 Pointer to the j-th vertex node in the element. More...
 
- Public Member Functions inherited from oomph::QElementBase
 QElementBase ()
 Constructor: Initialise pointers to macro element reference coords. More...
 
 QElementBase (const QElementBase &)=delete
 Broken copy constructor. More...
 
virtual ~QElementBase ()
 Broken assignment operator. More...
 
bool local_coord_is_valid (const Vector< double > &s)
 Check whether the local coordinate are valid or not. More...
 
void move_local_coord_back_into_element (Vector< double > &s) const
 
virtual void set_macro_elem_pt (MacroElement *macro_elem_pt)
 
double & s_macro_ll (const unsigned &i)
 
double & s_macro_ur (const unsigned &i)
 
double s_macro_ll (const unsigned &i) const
 
double s_macro_ur (const unsigned &i) const
 
void get_x_from_macro_element (const Vector< double > &s, Vector< double > &x) const
 
void get_x_from_macro_element (const unsigned &t, const Vector< double > &s, Vector< double > &x)
 
unsigned nnode_on_face () const
 
ElementGeometry::ElementGeometry element_geometry () const
 It's a Q element! More...
 
- Public Member Functions inherited from oomph::QElementGeometricBase
 QElementGeometricBase ()
 Empty default constructor. More...
 
 QElementGeometricBase (const QElementGeometricBase &)=delete
 Broken copy constructor. More...
 

Protected Member Functions

void setup_father_bounds ()
 
void get_edge_bcs (const int &edge, Vector< int > &bound_cons) const
 Determine Vector of boundary conditions along edge (N/S/W/E) More...
 
void setup_hang_for_value (const int &value_id)
 
virtual void quad_hang_helper (const int &value_id, const int &my_edge, std::ofstream &output_hangfile)
 
- Protected Member Functions inherited from oomph::RefineableElement
void assemble_local_to_eulerian_jacobian (const DShape &dpsids, DenseMatrix< double > &jacobian) const
 
void assemble_local_to_eulerian_jacobian2 (const DShape &d2psids, DenseMatrix< double > &jacobian2) const
 
void assemble_eulerian_base_vectors (const DShape &dpsids, DenseMatrix< double > &interpolated_G) const
 
double local_to_eulerian_mapping_diagonal (const DShape &dpsids, DenseMatrix< double > &jacobian, DenseMatrix< double > &inverse_jacobian) const
 
void assign_hanging_local_eqn_numbers (const bool &store_local_dof_pt)
 Assign the local equation numbers for hanging node variables. More...
 
virtual void fill_in_jacobian_from_nodal_by_fd (Vector< double > &residuals, DenseMatrix< double > &jacobian)
 
- Protected Member Functions inherited from oomph::FiniteElement
template<unsigned DIM>
double invert_jacobian (const DenseMatrix< double > &jacobian, DenseMatrix< double > &inverse_jacobian) const
 
virtual double invert_jacobian_mapping (const DenseMatrix< double > &jacobian, DenseMatrix< double > &inverse_jacobian) const
 
virtual double local_to_eulerian_mapping (const DShape &dpsids, DenseMatrix< double > &jacobian, DenseMatrix< double > &inverse_jacobian) const
 
double local_to_eulerian_mapping (const DShape &dpsids, DenseMatrix< double > &inverse_jacobian) const
 
virtual void dJ_eulerian_dnodal_coordinates (const DenseMatrix< double > &jacobian, const DShape &dpsids, DenseMatrix< double > &djacobian_dX) const
 
template<unsigned DIM>
void dJ_eulerian_dnodal_coordinates_templated_helper (const DenseMatrix< double > &jacobian, const DShape &dpsids, DenseMatrix< double > &djacobian_dX) const
 
virtual void d_dshape_eulerian_dnodal_coordinates (const double &det_jacobian, const DenseMatrix< double > &jacobian, const DenseMatrix< double > &djacobian_dX, const DenseMatrix< double > &inverse_jacobian, const DShape &dpsids, RankFourTensor< double > &d_dpsidx_dX) const
 
template<unsigned DIM>
void d_dshape_eulerian_dnodal_coordinates_templated_helper (const double &det_jacobian, const DenseMatrix< double > &jacobian, const DenseMatrix< double > &djacobian_dX, const DenseMatrix< double > &inverse_jacobian, const DShape &dpsids, RankFourTensor< double > &d_dpsidx_dX) const
 
virtual void transform_derivatives (const DenseMatrix< double > &inverse_jacobian, DShape &dbasis) const
 
void transform_derivatives_diagonal (const DenseMatrix< double > &inverse_jacobian, DShape &dbasis) const
 
virtual void transform_second_derivatives (const DenseMatrix< double > &jacobian, const DenseMatrix< double > &inverse_jacobian, const DenseMatrix< double > &jacobian2, DShape &dbasis, DShape &d2basis) const
 
template<unsigned DIM>
void transform_second_derivatives_template (const DenseMatrix< double > &jacobian, const DenseMatrix< double > &inverse_jacobian, const DenseMatrix< double > &jacobian2, DShape &dbasis, DShape &d2basis) const
 
template<unsigned DIM>
void transform_second_derivatives_diagonal (const DenseMatrix< double > &jacobian, const DenseMatrix< double > &inverse_jacobian, const DenseMatrix< double > &jacobian2, DShape &dbasis, DShape &d2basis) const
 
void fill_in_jacobian_from_nodal_by_fd (DenseMatrix< double > &jacobian)
 
virtual void update_before_nodal_fd ()
 
virtual void reset_after_nodal_fd ()
 
virtual void update_in_nodal_fd (const unsigned &i)
 
virtual void reset_in_nodal_fd (const unsigned &i)
 
void fill_in_contribution_to_jacobian (Vector< double > &residuals, DenseMatrix< double > &jacobian)
 
template<>
double invert_jacobian (const DenseMatrix< double > &jacobian, DenseMatrix< double > &inverse_jacobian) const
 Zero-d specialisation of function to calculate inverse of jacobian mapping. More...
 
template<>
double invert_jacobian (const DenseMatrix< double > &jacobian, DenseMatrix< double > &inverse_jacobian) const
 One-d specialisation of function to calculate inverse of jacobian mapping. More...
 
template<>
double invert_jacobian (const DenseMatrix< double > &jacobian, DenseMatrix< double > &inverse_jacobian) const
 Two-d specialisation of function to calculate inverse of jacobian mapping. More...
 
template<>
double invert_jacobian (const DenseMatrix< double > &jacobian, DenseMatrix< double > &inverse_jacobian) const
 
template<>
void dJ_eulerian_dnodal_coordinates_templated_helper (const DenseMatrix< double > &jacobian, const DShape &dpsids, DenseMatrix< double > &djacobian_dX) const
 
template<>
void dJ_eulerian_dnodal_coordinates_templated_helper (const DenseMatrix< double > &jacobian, const DShape &dpsids, DenseMatrix< double > &djacobian_dX) const
 
template<>
void dJ_eulerian_dnodal_coordinates_templated_helper (const DenseMatrix< double > &jacobian, const DShape &dpsids, DenseMatrix< double > &djacobian_dX) const
 
template<>
void dJ_eulerian_dnodal_coordinates_templated_helper (const DenseMatrix< double > &jacobian, const DShape &dpsids, DenseMatrix< double > &djacobian_dX) const
 
template<>
void d_dshape_eulerian_dnodal_coordinates_templated_helper (const double &det_jacobian, const DenseMatrix< double > &jacobian, const DenseMatrix< double > &djacobian_dX, const DenseMatrix< double > &inverse_jacobian, const DShape &dpsids, RankFourTensor< double > &d_dpsidx_dX) const
 
template<>
void d_dshape_eulerian_dnodal_coordinates_templated_helper (const double &det_jacobian, const DenseMatrix< double > &jacobian, const DenseMatrix< double > &djacobian_dX, const DenseMatrix< double > &inverse_jacobian, const DShape &dpsids, RankFourTensor< double > &d_dpsidx_dX) const
 
template<>
void d_dshape_eulerian_dnodal_coordinates_templated_helper (const double &det_jacobian, const DenseMatrix< double > &jacobian, const DenseMatrix< double > &djacobian_dX, const DenseMatrix< double > &inverse_jacobian, const DShape &dpsids, RankFourTensor< double > &d_dpsidx_dX) const
 
template<>
void d_dshape_eulerian_dnodal_coordinates_templated_helper (const double &det_jacobian, const DenseMatrix< double > &jacobian, const DenseMatrix< double > &djacobian_dX, const DenseMatrix< double > &inverse_jacobian, const DShape &dpsids, RankFourTensor< double > &d_dpsidx_dX) const
 
template<>
void transform_second_derivatives_template (const DenseMatrix< double > &jacobian, const DenseMatrix< double > &inverse_jacobian, const DenseMatrix< double > &jacobian2, DShape &dbasis, DShape &d2basis) const
 
template<>
void transform_second_derivatives_template (const DenseMatrix< double > &jacobian, const DenseMatrix< double > &inverse_jacobian, const DenseMatrix< double > &jacobian2, DShape &dbasis, DShape &d2basis) const
 
template<>
void transform_second_derivatives_diagonal (const DenseMatrix< double > &jacobian, const DenseMatrix< double > &inverse_jacobian, const DenseMatrix< double > &jacobian2, DShape &dbasis, DShape &d2basis) const
 
template<>
void transform_second_derivatives_diagonal (const DenseMatrix< double > &jacobian, const DenseMatrix< double > &inverse_jacobian, const DenseMatrix< double > &jacobian2, DShape &dbasis, DShape &d2basis) const
 
- Protected Member Functions inherited from oomph::GeneralisedElement
unsigned add_internal_data (Data *const &data_pt, const bool &fd=true)
 
bool internal_data_fd (const unsigned &i) const
 
void exclude_internal_data_fd (const unsigned &i)
 
void include_internal_data_fd (const unsigned &i)
 
void clear_global_eqn_numbers ()
 
void add_global_eqn_numbers (std::deque< unsigned long > const &global_eqn_numbers, std::deque< double * > const &global_dof_pt)
 
virtual void assign_internal_and_external_local_eqn_numbers (const bool &store_local_dof_pt)
 
virtual void assign_additional_local_eqn_numbers ()
 
int internal_local_eqn (const unsigned &i, const unsigned &j) const
 
int external_local_eqn (const unsigned &i, const unsigned &j)
 
virtual void fill_in_contribution_to_residuals (Vector< double > &residuals)
 
void fill_in_jacobian_from_internal_by_fd (Vector< double > &residuals, DenseMatrix< double > &jacobian, const bool &fd_all_data=false)
 
void fill_in_jacobian_from_internal_by_fd (DenseMatrix< double > &jacobian, const bool &fd_all_data=false)
 
void fill_in_jacobian_from_external_by_fd (Vector< double > &residuals, DenseMatrix< double > &jacobian, const bool &fd_all_data=false)
 
void fill_in_jacobian_from_external_by_fd (DenseMatrix< double > &jacobian, const bool &fd_all_data=false)
 
virtual void update_before_internal_fd ()
 
virtual void reset_after_internal_fd ()
 
virtual void update_in_internal_fd (const unsigned &i)
 
virtual void reset_in_internal_fd (const unsigned &i)
 
virtual void update_before_external_fd ()
 
virtual void reset_after_external_fd ()
 
virtual void update_in_external_fd (const unsigned &i)
 
virtual void reset_in_external_fd (const unsigned &i)
 
virtual void fill_in_contribution_to_mass_matrix (Vector< double > &residuals, DenseMatrix< double > &mass_matrix)
 
virtual void fill_in_contribution_to_jacobian_and_mass_matrix (Vector< double > &residuals, DenseMatrix< double > &jacobian, DenseMatrix< double > &mass_matrix)
 
virtual void fill_in_contribution_to_dresiduals_dparameter (double *const &parameter_pt, Vector< double > &dres_dparam)
 
virtual void fill_in_contribution_to_djacobian_dparameter (double *const &parameter_pt, Vector< double > &dres_dparam, DenseMatrix< double > &djac_dparam)
 
virtual void fill_in_contribution_to_djacobian_and_dmass_matrix_dparameter (double *const &parameter_pt, Vector< double > &dres_dparam, DenseMatrix< double > &djac_dparam, DenseMatrix< double > &dmass_matrix_dparam)
 
virtual void fill_in_contribution_to_hessian_vector_products (Vector< double > const &Y, DenseMatrix< double > const &C, DenseMatrix< double > &product)
 
virtual void fill_in_contribution_to_inner_products (Vector< std::pair< unsigned, unsigned >> const &history_index, Vector< double > &inner_product)
 
virtual void fill_in_contribution_to_inner_product_vectors (Vector< unsigned > const &history_index, Vector< Vector< double >> &inner_product_vector)
 

Static Protected Attributes

static std::map< unsigned, DenseMatrix< int > > Father_bound
 
- Static Protected Attributes inherited from oomph::RefineableElement
static double Max_integrity_tolerance = 1.0e-8
 Max. allowed discrepancy in element integrity check. More...
 
- Static Protected Attributes inherited from oomph::FiniteElement
static const unsigned Default_Initial_Nvalue = 0
 Default value for the number of values at a node. More...
 
static const double Node_location_tolerance = 1.0e-14
 
static const unsigned N2deriv [] = {0, 1, 3, 6}
 
- Static Protected Attributes inherited from oomph::GeneralisedElement
static DenseMatrix< double > Dummy_matrix
 
static std::deque< double * > Dof_pt_deque
 

Additional Inherited Members

- Static Public Member Functions inherited from oomph::RefineableElement
static double & max_integrity_tolerance ()
 Max. allowed discrepancy in element integrity check. More...
 
- Static Public Attributes inherited from oomph::FiniteElement
static double Tolerance_for_singular_jacobian = 1.0e-16
 Tolerance below which the jacobian is considered singular. More...
 
static bool Accept_negative_jacobian = false
 
static bool Suppress_output_while_checking_for_inverted_elements
 
- Static Public Attributes inherited from oomph::GeneralisedElement
static bool Suppress_warning_about_repeated_internal_data
 
static bool Suppress_warning_about_repeated_external_data = true
 
static double Default_fd_jacobian_step = 1.0e-8
 
- Static Protected Member Functions inherited from oomph::RefineableElement
static void check_value_id (const int &n_continuously_interpolated_values, const int &value_id)
 
- Protected Attributes inherited from oomph::RefineableElement
Tree * Tree_pt
 A pointer to a general tree object. More...
 
unsigned Refine_level
 Refinement level. More...
 
bool To_be_refined
 Flag for refinement. More...
 
bool Refinement_is_enabled
 Flag to indicate suppression of any refinement. More...
 
bool Sons_to_be_unrefined
 Flag for unrefinement. More...
 
long Number
 Global element number – for plotting/validation purposes. More...
 
- Protected Attributes inherited from oomph::FiniteElement
MacroElement * Macro_elem_pt
 Pointer to the element's macro element (NULL by default) More...
 
- Protected Attributes inherited from oomph::GeomObject
unsigned NLagrangian
 Number of Lagrangian (intrinsic) coordinates. More...
 
unsigned Ndim
 Number of Eulerian coordinates. More...
 
TimeStepper * Geom_object_time_stepper_pt
 

Detailed Description

Refineable version of QElement<2,NNODE_1D>.

Refinement is performed by quadtree procedures. When the element is subdivided, the geometry of its sons is established by calls to their father's

get_x(...)
oomph::FiniteElement::get_x
void get_x(const Vector< double > &s, Vector< double > &x) const
Definition: elements.h:1885

function which refers to

  • the father element's geometric FE mapping (this is the default)

or

  • to a MacroElement 's MacroElement::macro_map (if the pointer to the macro element is non-NULL)

The class provides a generic RefineableQElement<2>::build() function which deals with generic isoparametric QElements in which all values are associated with nodes. The RefineableQElement<2>::further_build() function provides an interface for any element-specific non-generic build operations.

Member Typedef Documentation

◆ VoidMemberFctPt

typedef void(RefineableQElement<2>::* oomph::RefineableQElement< 2 >::VoidMemberFctPt) ()

Shorthand for pointer to an argument-free void member function of the refineable element

Constructor & Destructor Documentation

◆ RefineableQElement() [1/2]

oomph::RefineableQElement< 2 >::RefineableQElement ( )
inline

Constructor: Pass refinement level (default 0 = root)

78  : RefineableElement()
79  {
80 #ifdef LEAK_CHECK
81  LeakCheckNames::RefineableQElement<2> _build += 1;
82 #endif
83  }
oomph::RefineableElement::RefineableElement
RefineableElement()
Definition: refineable_elements.h:188

◆ RefineableQElement() [2/2]

oomph::RefineableQElement< 2 >::RefineableQElement ( const RefineableQElement< 2 > &  dummy)
delete

Broken copy constructor.

◆ ~RefineableQElement()

virtual oomph::RefineableQElement< 2 >::~RefineableQElement ( )
inlinevirtual

Broken assignment operator.

Destructor

99  {
100 #ifdef LEAK_CHECK
101  LeakCheckNames::RefineableQElement<2> _build -= 1;
102 #endif
103  }

Member Function Documentation

◆ build()

void oomph::RefineableQElement< 2 >::build ( Mesh *&  mesh_pt,
Vector< Node * > &  new_node_pt,
bool &  was_already_built,
std::ofstream &  new_nodes_file 
)
virtual

Build the element, i.e. give it nodal positions, apply BCs, etc. Pointers to any new nodes will be returned in new_node_pt. If it is open, the positions of the new nodes will be written to the file stream new_nodes_file

Build the element by doing the following:

  • Give it nodal positions (by establishing the pointers to its nodes)
  • In the process create new nodes where required (i.e. if they don't exist in father element or have already been created while building new neighbour elements). Node building involves the following steps:
    • Get nodal position from father element.
    • Establish the time-history of the newly created nodal point (its coordinates and the previous values) consistent with the father's history.
    • Determine the boundary conditions of the nodes (newly created nodes can only lie on the interior of any edges of the father element – this makes it possible to to figure out what their bc should be...)
    • Add node to the mesh's stoarge scheme for the boundary nodes.
    • Add the new node to the mesh itself
    • Doc newly created nodes in "new_nodes.dat" stored in the directory of the DocInfo object (only if it's open!)
  • Finally, excute the element-specific further_build() (empty by default – must be overloaded for specific elements). This deals with any build operations that are not included in the generic process outlined above. For instance, in Crouzeix Raviart elements we need to initialise the internal pressure values in manner consistent with the pressure distribution in the father element.

Implements oomph::RefineableElement.

Reimplemented in oomph::RefineableSolidQElement< 2 >.

646  {
647  using namespace QuadTreeNames;
648 
649  // Get the number of 1d nodes
650  unsigned n_p = nnode_1d();
651 
652  // Check whether static father_bound needs to be created
653  if (Father_bound[n_p].nrow() == 0)
654  {
655  setup_father_bounds();
656  }
657 
658  // Pointer to my father (in quadtree impersonation)
659  QuadTree* father_pt = dynamic_cast<QuadTree*>(quadtree_pt()->father_pt());
660 
661  // What type of son am I? Ask my quadtree representation...
662  int son_type = Tree_pt->son_type();
663 
664  // Has somebody build me already? (If any nodes have not been built)
665  if (!nodes_built())
666  {
667 #ifdef PARANOID
668  if (father_pt == 0)
669  {
670  std::string error_message =
671  "Something fishy here: I have no father and yet \n";
672  error_message += "I have no nodes. Who has created me then?!\n";
673 
674  throw OomphLibError(
675  error_message, OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
676  }
677 #endif
678 
679  // Indicate status:
680  was_already_built = false;
681 
682  // Return pointer to father element
683  RefineableQElement<2>* father_el_pt =
684  dynamic_cast<RefineableQElement<2>*>(father_pt->object_pt());
685 
686  // Timestepper should be the same for all nodes in father
687  // element -- use it create timesteppers for new nodes
688  TimeStepper* time_stepper_pt =
689  father_el_pt->node_pt(0)->time_stepper_pt();
690 
691  // Number of history values (incl. present)
692  unsigned ntstorage = time_stepper_pt->ntstorage();
693 
694  // Currently we can't handle the case of generalised coordinates
695  // since we haven't established how they should be interpolated
696  // Buffer this case:
697  if (father_el_pt->node_pt(0)->nposition_type() != 1)
698  {
699  throw OomphLibError("Can't handle generalised nodal positions (yet).",
700  OOMPH_CURRENT_FUNCTION,
701  OOMPH_EXCEPTION_LOCATION);
702  }
703 
704  Vector<double> s_lo(2);
705  Vector<double> s_hi(2);
706  Vector<double> s(2);
707  Vector<double> x(2);
708 
709  // Setup vertex coordinates in father element:
710  //--------------------------------------------
711  switch (son_type)
712  {
713  case SW:
714  s_lo[0] = -1.0;
715  s_hi[0] = 0.0;
716  s_lo[1] = -1.0;
717  s_hi[1] = 0.0;
718  break;
719 
720  case SE:
721  s_lo[0] = 0.0;
722  s_hi[0] = 1.0;
723  s_lo[1] = -1.0;
724  s_hi[1] = 0.0;
725  break;
726 
727  case NE:
728  s_lo[0] = 0.0;
729  s_hi[0] = 1.0;
730  s_lo[1] = 0.0;
731  s_hi[1] = 1.0;
732  break;
733 
734  case NW:
735  s_lo[0] = -1.0;
736  s_hi[0] = 0.0;
737  s_lo[1] = 0.0;
738  s_hi[1] = 1.0;
739  break;
740  }
741 
742  // Pass macro element pointer on to sons and
743  // set coordinates in macro element
744  // hierher why can I see this?
745  if (father_el_pt->Macro_elem_pt != 0)
746  {
747  set_macro_elem_pt(father_el_pt->Macro_elem_pt);
748  for (unsigned i = 0; i < 2; i++)
749  {
750  s_macro_ll(i) =
751  father_el_pt->s_macro_ll(i) +
752  0.5 * (s_lo[i] + 1.0) *
753  (father_el_pt->s_macro_ur(i) - father_el_pt->s_macro_ll(i));
754  s_macro_ur(i) =
755  father_el_pt->s_macro_ll(i) +
756  0.5 * (s_hi[i] + 1.0) *
757  (father_el_pt->s_macro_ur(i) - father_el_pt->s_macro_ll(i));
758  }
759  }
760 
761 
762  // If the father element hasn't been generated yet, we're stuck...
763  if (father_el_pt->node_pt(0) == 0)
764  {
765  throw OomphLibError(
766  "Trouble: father_el_pt->node_pt(0)==0\n Can't build son element!\n",
767  OOMPH_CURRENT_FUNCTION,
768  OOMPH_EXCEPTION_LOCATION);
769  }
770  else
771  {
772  unsigned jnod = 0;
773  Vector<double> x_small(2);
774  Vector<double> x_large(2);
775 
776  Vector<double> s_fraction(2);
777  // Loop over nodes in element
778  for (unsigned i0 = 0; i0 < n_p; i0++)
779  {
780  // Get the fractional position of the node in the direction of s[0]
781  s_fraction[0] = local_one_d_fraction_of_node(i0, 0);
782  // Local coordinate in father element
783  s[0] = s_lo[0] + (s_hi[0] - s_lo[0]) * s_fraction[0];
784 
785  for (unsigned i1 = 0; i1 < n_p; i1++)
786  {
787  // Get the fractional position of the node in the direction of s[1]
788  s_fraction[1] = local_one_d_fraction_of_node(i1, 1);
789  // Local coordinate in father element
790  s[1] = s_lo[1] + (s_hi[1] - s_lo[1]) * s_fraction[1];
791 
792  // Local node number
793  jnod = i0 + n_p * i1;
794 
795  // Check whether the father's node is periodic if so, complain
796  /* {
797  Node* father_node_pt = father_el_pt->node_pt(jnod);
798  if((father_node_pt->is_a_copy()) ||
799  (father_node_pt->position_is_a_copy()))
800  {
801  throw OomphLibError(
802  "Can't handle periodic nodes (yet).",
803  OOMPH_CURRENT_FUNCTION,
804  OOMPH_EXCEPTION_LOCATION);
805  }
806  }*/
807 
808  // Initialise flag: So far, this node hasn't been built
809  // or copied yet
810  bool node_done = false;
811 
812  // Get the pointer to the node in the father, returns NULL
813  // if there is not node
814  Node* created_node_pt =
815  father_el_pt->get_node_at_local_coordinate(s);
816 
817  // Does this node already exist in father element?
818  //------------------------------------------------
819  if (created_node_pt != 0)
820  {
821  // Copy node across
822  node_pt(jnod) = created_node_pt;
823 
824  // Make sure that we update the values at the node so that
825  // they are consistent with the present representation.
826  // This is only need for mixed interpolation where the value
827  // at the father could now become active.
828 
829  // Loop over all history values
830  for (unsigned t = 0; t < ntstorage; t++)
831  {
832  // Get values from father element
833  // Note: get_interpolated_values() sets Vector size itself.
834  Vector<double> prev_values;
835  father_el_pt->get_interpolated_values(t, s, prev_values);
836  // Find the minimum number of values
837  //(either those stored at the node, or those returned by
838  // the function)
839  unsigned n_val_at_node = created_node_pt->nvalue();
840  unsigned n_val_from_function = prev_values.size();
841  // Use the ternary conditional operator here
842  unsigned n_var = n_val_at_node < n_val_from_function ?
843  n_val_at_node :
844  n_val_from_function;
845  // Assign the values that we can
846  for (unsigned k = 0; k < n_var; k++)
847  {
848  created_node_pt->set_value(t, k, prev_values[k]);
849  }
850  }
851 
852  // Node has been created by copy
853  node_done = true;
854  }
855  // Node does not exist in father element but might already
856  //--------------------------------------------------------
857  // have been created by neighbouring elements
858  //-------------------------------------------
859  else
860  {
861  // Was the node created by one of its neighbours
862  // Whether or not the node lies on an edge can be calculated
863  // by from the fractional position
864  bool is_periodic = false;
865  ;
866  created_node_pt =
867  node_created_by_neighbour(s_fraction, is_periodic);
868 
869  // If the node was so created, assign the pointers
870  if (created_node_pt != 0)
871  {
872  // If the node is periodic
873  if (is_periodic)
874  {
875  // Now the node must be on a boundary, but we don't know which
876  // one
877  // The returned created_node_pt is actually the neighbouring
878  // periodic node
879  Node* neighbour_node_pt = created_node_pt;
880 
881  // Determine the edge on which the new node will live
882  int father_bound = Father_bound[n_p](jnod, son_type);
883 
884  // Storage for the set of Mesh boundaries on which the
885  // appropriate father edge lives.
886  // [New nodes should always be mid-edge nodes in father
887  // and therefore only live on one boundary but just to
888  // play it safe...]
889  std::set<unsigned> boundaries;
890  // Only get the boundaries if we are at the edge of
891  // an element. Nodes in the centre of an element cannot be
892  // on Mesh boundaries
893  if (father_bound != Tree::OMEGA)
894  {
895  father_el_pt->get_boundaries(father_bound, boundaries);
896  }
897 
898 #ifdef PARANOID
899  // Case where a new node lives on more than one boundary
900  // seems fishy enough to flag
901  if (boundaries.size() > 1)
902  {
903  throw OomphLibError(
904  "boundaries.size()!=1 seems a bit strange..\n",
905  OOMPH_CURRENT_FUNCTION,
906  OOMPH_EXCEPTION_LOCATION);
907  }
908 
909  // Case when there are no boundaries, we are in big trouble
910  if (boundaries.size() == 0)
911  {
912  std::ostringstream error_stream;
913  error_stream << "Periodic node is not on a boundary...\n"
914  << "Coordinates: " << created_node_pt->x(0)
915  << " " << created_node_pt->x(1) << "\n";
916  throw OomphLibError(error_stream.str(),
917  OOMPH_CURRENT_FUNCTION,
918  OOMPH_EXCEPTION_LOCATION);
919  }
920 #endif
921 
922  // Create node and set the pointer to it from the element
923  created_node_pt =
924  construct_boundary_node(jnod, time_stepper_pt);
925  // Make the node periodic from the neighbour
926  created_node_pt->make_periodic(neighbour_node_pt);
927  // Add to vector of new nodes
928  new_node_pt.push_back(created_node_pt);
929 
930  // Loop over # of history values
931  for (unsigned t = 0; t < ntstorage; t++)
932  {
933  // Get position from father element -- this uses the macro
934  // element representation if appropriate. If the node
935  // turns out to be a hanging node later on, then
936  // its position gets adjusted in line with its
937  // hanging node interpolation.
938  Vector<double> x_prev(2);
939  father_el_pt->get_x(t, s, x_prev);
940  // Set previous positions of the new node
941  for (unsigned i = 0; i < 2; i++)
942  {
943  created_node_pt->x(t, i) = x_prev[i];
944  }
945  }
946 
947  // Next, we Update the boundary lookup schemes
948  // Loop over the boundaries stored in the set
949  for (std::set<unsigned>::iterator it = boundaries.begin();
950  it != boundaries.end();
951  ++it)
952  {
953  // Add the node to the boundary
954  mesh_pt->add_boundary_node(*it, created_node_pt);
955 
956  // If we have set an intrinsic coordinate on this
957  // mesh boundary then it must also be interpolated on
958  // the new node
959  // Now interpolate the intrinsic boundary coordinate
960  if (mesh_pt->boundary_coordinate_exists(*it) == true)
961  {
962  Vector<double> zeta(1);
963  father_el_pt->interpolated_zeta_on_edge(
964  *it, father_bound, s, zeta);
965 
966  created_node_pt->set_coordinates_on_boundary(*it, zeta);
967  }
968  }
969 
970  // Make sure that we add the node to the mesh
971  mesh_pt->add_node_pt(created_node_pt);
972  } // End of periodic case
973  // Otherwise the node is not periodic, so just set the
974  // pointer to the neighbours node
975  else
976  {
977  node_pt(jnod) = created_node_pt;
978  }
979  // Node has been created
980  node_done = true;
981  }
982  // Node does not exist in neighbour element but might already
983  //-----------------------------------------------------------
984  // have been created by a son of a neighbouring element
985  //-----------------------------------------------------
986  else
987  {
988  // Was the node created by one of its neighbours' sons
989  // Whether or not the node lies on an edge can be calculated
990  // by from the fractional position
991  bool is_periodic = false;
992  ;
993  created_node_pt =
994  node_created_by_son_of_neighbour(s_fraction, is_periodic);
995 
996  // If the node was so created, assign the pointers
997  if (created_node_pt != 0)
998  {
999  // If the node is periodic
1000  if (is_periodic)
1001  {
1002  // Now the node must be on a boundary, but we don't know
1003  // which one The returned created_node_pt is actually the
1004  // neighbouring periodic node
1005  Node* neighbour_node_pt = created_node_pt;
1006 
1007  // Determine the edge on which the new node will live
1008  int father_bound = Father_bound[n_p](jnod, son_type);
1009 
1010  // Storage for the set of Mesh boundaries on which the
1011  // appropriate father edge lives.
1012  // [New nodes should always be mid-edge nodes in father
1013  // and therefore only live on one boundary but just to
1014  // play it safe...]
1015  std::set<unsigned> boundaries;
1016  // Only get the boundaries if we are at the edge of
1017  // an element. Nodes in the centre of an element cannot be
1018  // on Mesh boundaries
1019  if (father_bound != Tree::OMEGA)
1020  {
1021  father_el_pt->get_boundaries(father_bound, boundaries);
1022  }
1023 
1024 #ifdef PARANOID
1025  // Case where a new node lives on more than one boundary
1026  // seems fishy enough to flag
1027  if (boundaries.size() > 1)
1028  {
1029  throw OomphLibError(
1030  "boundaries.size()!=1 seems a bit strange..\n",
1031  OOMPH_CURRENT_FUNCTION,
1032  OOMPH_EXCEPTION_LOCATION);
1033  }
1034 
1035  // Case when there are no boundaries, we are in big trouble
1036  if (boundaries.size() == 0)
1037  {
1038  std::ostringstream error_stream;
1039  error_stream << "Periodic node is not on a boundary...\n"
1040  << "Coordinates: " << created_node_pt->x(0)
1041  << " " << created_node_pt->x(1) << "\n";
1042  throw OomphLibError(error_stream.str(),
1043  OOMPH_CURRENT_FUNCTION,
1044  OOMPH_EXCEPTION_LOCATION);
1045  }
1046 #endif
1047 
1048  // Create node and set the pointer to it from the element
1049  created_node_pt =
1050  construct_boundary_node(jnod, time_stepper_pt);
1051  // Make the node periodic from the neighbour
1052  created_node_pt->make_periodic(neighbour_node_pt);
1053  // Add to vector of new nodes
1054  new_node_pt.push_back(created_node_pt);
1055 
1056  // Loop over # of history values
1057  for (unsigned t = 0; t < ntstorage; t++)
1058  {
1059  // Get position from father element -- this uses the macro
1060  // element representation if appropriate. If the node
1061  // turns out to be a hanging node later on, then
1062  // its position gets adjusted in line with its
1063  // hanging node interpolation.
1064  Vector<double> x_prev(2);
1065  father_el_pt->get_x(t, s, x_prev);
1066  // Set previous positions of the new node
1067  for (unsigned i = 0; i < 2; i++)
1068  {
1069  created_node_pt->x(t, i) = x_prev[i];
1070  }
1071  }
1072 
1073  // Next, we Update the boundary lookup schemes
1074  // Loop over the boundaries stored in the set
1075  for (std::set<unsigned>::iterator it = boundaries.begin();
1076  it != boundaries.end();
1077  ++it)
1078  {
1079  // Add the node to the boundary
1080  mesh_pt->add_boundary_node(*it, created_node_pt);
1081 
1082  // If we have set an intrinsic coordinate on this
1083  // mesh boundary then it must also be interpolated on
1084  // the new node
1085  // Now interpolate the intrinsic boundary coordinate
1086  if (mesh_pt->boundary_coordinate_exists(*it) == true)
1087  {
1088  Vector<double> zeta(1);
1089  father_el_pt->interpolated_zeta_on_edge(
1090  *it, father_bound, s, zeta);
1091 
1092  created_node_pt->set_coordinates_on_boundary(*it, zeta);
1093  }
1094  }
1095 
1096  // Make sure that we add the node to the mesh
1097  mesh_pt->add_node_pt(created_node_pt);
1098  } // End of periodic case
1099  // Otherwise the node is not periodic, so just set the
1100  // pointer to the neighbours node
1101  else
1102  {
1103  node_pt(jnod) = created_node_pt;
1104  }
1105  // Node has been created
1106  node_done = true;
1107  } // Node does not exist in son of neighbouring element
1108  } // Node does not exist in neighbouring element
1109  } // Node does not exist in father element
1110 
1111  // Node has not been built anywhere ---> build it here
1112  if (!node_done)
1113  {
1114  // Firstly, we need to determine whether or not a node lies
1115  // on the boundary before building it, because
1116  // we actually assign a different type of node on boundaries.
1117 
1118  // The node can only be on a Mesh boundary if it
1119  // lives on an edge that is shared with an edge of its
1120  // father element; i.e. it is not created inside the father
1121  // element Determine the edge on which the new node will live
1122  int father_bound = Father_bound[n_p](jnod, son_type);
1123 
1124  // Storage for the set of Mesh boundaries on which the
1125  // appropriate father edge lives.
1126  // [New nodes should always be mid-edge nodes in father
1127  // and therefore only live on one boundary but just to
1128  // play it safe...]
1129  std::set<unsigned> boundaries;
1130  // Only get the boundaries if we are at the edge of
1131  // an element. Nodes in the centre of an element cannot be
1132  // on Mesh boundaries
1133  if (father_bound != Tree::OMEGA)
1134  {
1135  father_el_pt->get_boundaries(father_bound, boundaries);
1136  }
1137 
1138 #ifdef PARANOID
1139  // Case where a new node lives on more than one boundary
1140  // seems fishy enough to flag
1141  if (boundaries.size() > 1)
1142  {
1143  throw OomphLibError(
1144  "boundaries.size()!=1 seems a bit strange..\n",
1145  OOMPH_CURRENT_FUNCTION,
1146  OOMPH_EXCEPTION_LOCATION);
1147  }
1148 #endif
1149 
1150  // If the node lives on a mesh boundary,
1151  // then we need to create a boundary node
1152  if (boundaries.size() > 0)
1153  {
1154  // Create node and set the pointer to it from the element
1155  created_node_pt =
1156  construct_boundary_node(jnod, time_stepper_pt);
1157  // Add to vector of new nodes
1158  new_node_pt.push_back(created_node_pt);
1159 
1160  // Now we need to work out whether to pin the values at
1161  // the new node based on the boundary conditions applied at
1162  // its Mesh boundary
1163 
1164  // Get the boundary conditions from the father
1165  Vector<int> bound_cons(ncont_interpolated_values());
1166  father_el_pt->get_bcs(father_bound, bound_cons);
1167 
1168  // Loop over the values and pin, if necessary
1169  unsigned n_value = created_node_pt->nvalue();
1170  for (unsigned k = 0; k < n_value; k++)
1171  {
1172  if (bound_cons[k])
1173  {
1174  created_node_pt->pin(k);
1175  }
1176  }
1177 
1178  // Solid node? If so, deal with the positional boundary
1179  // conditions:
1180  SolidNode* solid_node_pt =
1181  dynamic_cast<SolidNode*>(created_node_pt);
1182  if (solid_node_pt != 0)
1183  {
1184  // Get the positional boundary conditions from the father:
1185  unsigned n_dim = created_node_pt->ndim();
1186  Vector<int> solid_bound_cons(n_dim);
1187  RefineableSolidQElement<2>* father_solid_el_pt =
1188  dynamic_cast<RefineableSolidQElement<2>*>(father_el_pt);
1189 #ifdef PARANOID
1190  if (father_solid_el_pt == 0)
1191  {
1192  std::string error_message =
1193  "We have a SolidNode outside a refineable SolidElement\n";
1194  error_message +=
1195  "during mesh refinement -- this doesn't make sense";
1196 
1197  throw OomphLibError(error_message,
1198  OOMPH_CURRENT_FUNCTION,
1199  OOMPH_EXCEPTION_LOCATION);
1200  }
1201 #endif
1202  father_solid_el_pt->get_solid_bcs(father_bound,
1203  solid_bound_cons);
1204 
1205  // Loop over the positions and pin, if necessary
1206  for (unsigned k = 0; k < n_dim; k++)
1207  {
1208  if (solid_bound_cons[k])
1209  {
1210  solid_node_pt->pin_position(k);
1211  }
1212  }
1213  } // End of if solid_node_pt
1214 
1215 
1216  // Next, we Update the boundary lookup schemes
1217  // Loop over the boundaries stored in the set
1218  for (std::set<unsigned>::iterator it = boundaries.begin();
1219  it != boundaries.end();
1220  ++it)
1221  {
1222  // Add the node to the boundary
1223  mesh_pt->add_boundary_node(*it, created_node_pt);
1224 
1225  // If we have set an intrinsic coordinate on this
1226  // mesh boundary then it must also be interpolated on
1227  // the new node
1228  // Now interpolate the intrinsic boundary coordinate
1229  if (mesh_pt->boundary_coordinate_exists(*it) == true)
1230  {
1231  Vector<double> zeta(1);
1232  father_el_pt->interpolated_zeta_on_edge(
1233  *it, father_bound, s, zeta);
1234 
1235  created_node_pt->set_coordinates_on_boundary(*it, zeta);
1236  }
1237  }
1238  }
1239  // Otherwise the node is not on a Mesh boundary and
1240  // we create a normal "bulk" node
1241  else
1242  {
1243  // Create node and set the pointer to it from the element
1244  created_node_pt = construct_node(jnod, time_stepper_pt);
1245  // Add to vector of new nodes
1246  new_node_pt.push_back(created_node_pt);
1247  }
1248 
1249  // Now we set the position and values at the newly created node
1250 
1251  // In the first instance use macro element or FE representation
1252  // to create past and present nodal positions.
1253  // (THIS STEP SHOULD NOT BE SKIPPED FOR ALGEBRAIC
1254  // ELEMENTS AS NOT ALL OF THEM NECESSARILY IMPLEMENT
1255  // NONTRIVIAL NODE UPDATE FUNCTIONS. CALLING
1256  // THE NODE UPDATE FOR SUCH ELEMENTS/NODES WILL LEAVE
1257  // THEIR NODAL POSITIONS WHERE THEY WERE (THIS IS APPROPRIATE
1258  // ONCE THEY HAVE BEEN GIVEN POSITIONS) BUT WILL
1259  // NOT ASSIGN SENSIBLE INITIAL POSITONS!
1260 
1261  // Loop over # of history values
1262  for (unsigned t = 0; t < ntstorage; t++)
1263  {
1264  // Get position from father element -- this uses the macro
1265  // element representation if appropriate. If the node
1266  // turns out to be a hanging node later on, then
1267  // its position gets adjusted in line with its
1268  // hanging node interpolation.
1269  Vector<double> x_prev(2);
1270  father_el_pt->get_x(t, s, x_prev);
1271 
1272  // Set previous positions of the new node
1273  for (unsigned i = 0; i < 2; i++)
1274  {
1275  created_node_pt->x(t, i) = x_prev[i];
1276  }
1277  }
1278 
1279  // Loop over all history values
1280  for (unsigned t = 0; t < ntstorage; t++)
1281  {
1282  // Get values from father element
1283  // Note: get_interpolated_values() sets Vector size itself.
1284  Vector<double> prev_values;
1285  father_el_pt->get_interpolated_values(t, s, prev_values);
1286  // Initialise the values at the new node
1287  unsigned n_value = created_node_pt->nvalue();
1288  for (unsigned k = 0; k < n_value; k++)
1289  {
1290  created_node_pt->set_value(t, k, prev_values[k]);
1291  }
1292  }
1293 
1294  // Add new node to mesh
1295  mesh_pt->add_node_pt(created_node_pt);
1296 
1297  } // End of case when we build the node ourselves
1298 
1299  // Check if the element is an algebraic element
1300  AlgebraicElementBase* alg_el_pt =
1301  dynamic_cast<AlgebraicElementBase*>(this);
1302 
1303  // If the element is an algebraic element, setup
1304  // node position (past and present) from algebraic node update
1305  // function. This over-writes previous assingments that
1306  // were made based on the macro-element/FE representation.
1307  // NOTE: YES, THIS NEEDS TO BE CALLED REPEATEDLY IF THE
1308  // NODE IS MEMBER OF MULTIPLE ELEMENTS: THEY ALL ASSIGN
1309  // THE SAME NODAL POSITIONS BUT WE NEED TO ADD THE REMESH
1310  // INFO FOR *ALL* ROOT ELEMENTS!
1311  if (alg_el_pt != 0)
1312  {
1313  // Build algebraic node update info for new node
1314  // This sets up the node update data for all node update
1315  // functions that are shared by all nodes in the father
1316  // element
1317  alg_el_pt->setup_algebraic_node_update(
1318  node_pt(jnod), s, father_el_pt);
1319  }
1320 
1321  // If we have built the node and we are documenting our progress
1322  // write the (hopefully consistent position) to the outputfile
1323  if ((!node_done) && (new_nodes_file.is_open()))
1324  {
1325  new_nodes_file << node_pt(jnod)->x(0) << " "
1326  << node_pt(jnod)->x(1) << std::endl;
1327  }
1328 
1329  } // End of vertical loop over nodes in element
1330 
1331  } // End of horizontal loop over nodes in element
1332 
1333 
1334  // If the element is a MacroElementNodeUpdateElement, set
1335  // the update parameters for the current element's nodes --
1336  // all this needs is the vector of (pointers to the)
1337  // geometric objects that affect the MacroElement-based
1338  // node update -- this is the same as that in the father element
1339  MacroElementNodeUpdateElementBase* father_m_el_pt =
1340  dynamic_cast<MacroElementNodeUpdateElementBase*>(father_el_pt);
1341  if (father_m_el_pt != 0)
1342  {
1343  // Get vector of geometric objects from father (construct vector
1344  // via copy operation)
1345  Vector<GeomObject*> geom_object_pt(father_m_el_pt->geom_object_pt());
1346 
1347  // Cast current element to MacroElementNodeUpdateElement:
1348  MacroElementNodeUpdateElementBase* m_el_pt =
1349  dynamic_cast<MacroElementNodeUpdateElementBase*>(this);
1350 
1351 #ifdef PARANOID
1352  if (m_el_pt == 0)
1353  {
1354  std::string error_message =
1355  "Failed to cast to MacroElementNodeUpdateElementBase*\n";
1356  error_message +=
1357  "Strange -- if the father is a MacroElementNodeUpdateElement\n";
1358  error_message += "the son should be too....\n";
1359 
1360  throw OomphLibError(
1361  error_message, OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
1362  }
1363 #endif
1364  // Build update info by passing vector of geometric objects:
1365  // This sets the current element to be the update element
1366  // for all of the element's nodes -- this is reversed
1367  // if the element is ever un-refined in the father element's
1368  // rebuild_from_sons() function which overwrites this
1369  // assignment to avoid nasty segmentation faults that occur
1370  // when a node tries to update itself via an element that no
1371  // longer exists...
1372  m_el_pt->set_node_update_info(geom_object_pt);
1373  }
1374 
1375 #ifdef OOMPH_HAS_MPI
1376  // Pass on non-halo proc id
1377  Non_halo_proc_ID =
1378  tree_pt()->father_pt()->object_pt()->non_halo_proc_ID();
1379 #endif
1380 
1381  // Is it an ElementWithMovingNodes?
1382  ElementWithMovingNodes* aux_el_pt =
1383  dynamic_cast<ElementWithMovingNodes*>(this);
1384 
1385  // Pass down the information re the method for the evaluation
1386  // of the shape derivatives
1387  if (aux_el_pt != 0)
1388  {
1389  ElementWithMovingNodes* aux_father_el_pt =
1390  dynamic_cast<ElementWithMovingNodes*>(father_el_pt);
1391 
1392 #ifdef PARANOID
1393  if (aux_father_el_pt == 0)
1394  {
1395  std::string error_message =
1396  "Failed to cast to ElementWithMovingNodes*\n";
1397  error_message +=
1398  "Strange -- if the son is a ElementWithMovingNodes\n";
1399  error_message += "the father should be too....\n";
1400 
1401  throw OomphLibError(
1402  error_message, OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
1403  }
1404 #endif
1405 
1406  // If evaluating the residuals by finite differences in the father
1407  // continue to do so in the child
1408  if (aux_father_el_pt
1409  ->are_dresidual_dnodal_coordinates_always_evaluated_by_fd())
1410  {
1411  aux_el_pt
1412  ->enable_always_evaluate_dresidual_dnodal_coordinates_by_fd();
1413  }
1414 
1415  aux_el_pt->method_for_shape_derivs() =
1416  aux_father_el_pt->method_for_shape_derivs();
1417 
1418  // If bypassing the evaluation of fill_in_jacobian_from_geometric_data
1419  // continue to do so
1420  if (aux_father_el_pt
1421  ->is_fill_in_jacobian_from_geometric_data_bypassed())
1422  {
1423  aux_el_pt->enable_bypass_fill_in_jacobian_from_geometric_data();
1424  }
1425  }
1426 
1427 
1428  // Now do further build (if any)
1429  further_build();
1430 
1431  } // Sanity check: Father element has been generated
1432  }
1433  // Element has already been built
1434  else
1435  {
1436  was_already_built = true;
1437  }
1438  }
i
int i
Definition: BiCGSTAB_step_by_step.cpp:9
oomph::Data::nvalue
unsigned nvalue() const
Return number of values stored in data object (incl pinned ones).
Definition: nodes.h:483
oomph::FiniteElement::set_macro_elem_pt
virtual void set_macro_elem_pt(MacroElement *macro_elem_pt)
Definition: elements.h:1872
oomph::FiniteElement::construct_node
virtual Node * construct_node(const unsigned &n)
Definition: elements.h:2509
oomph::FiniteElement::node_pt
Node *& node_pt(const unsigned &n)
Return a pointer to the local node n.
Definition: elements.h:2175
oomph::FiniteElement::nnode_1d
virtual unsigned nnode_1d() const
Definition: elements.h:2218
oomph::FiniteElement::construct_boundary_node
virtual Node * construct_boundary_node(const unsigned &n)
Definition: elements.h:2538
oomph::FiniteElement::local_one_d_fraction_of_node
virtual double local_one_d_fraction_of_node(const unsigned &n1d, const unsigned &i)
Definition: elements.h:1858
oomph::GeomObject::time_stepper_pt
TimeStepper *& time_stepper_pt()
Definition: geom_objects.h:192
oomph::Node::x
double & x(const unsigned &i)
Return the i-th nodal coordinate.
Definition: nodes.h:1060
oomph::QElementBase::s_macro_ll
double & s_macro_ll(const unsigned &i)
Definition: Qelements.h:186
oomph::QElementBase::s_macro_ur
double & s_macro_ur(const unsigned &i)
Definition: Qelements.h:202
oomph::RefineableElement::tree_pt
Tree * tree_pt()
Access function: Pointer to quadtree representation of this element.
Definition: refineable_elements.h:211
oomph::RefineableElement::further_build
virtual void further_build()
Further build: e.g. deal with interpolation of internal values.
Definition: refineable_elements.h:599
oomph::RefineableElement::Tree_pt
Tree * Tree_pt
A pointer to a general tree object.
Definition: refineable_elements.h:100
oomph::RefineableElement::nodes_built
virtual bool nodes_built()
Return true if all the nodes have been built, false if not.
Definition: refineable_elements.h:389
oomph::RefineableElement::ncont_interpolated_values
virtual unsigned ncont_interpolated_values() const =0
oomph::RefineableQElement< 2 >::quadtree_pt
QuadTree * quadtree_pt()
Pointer to quadtree representation of this element.
Definition: refineable_quad_element.h:152
oomph::RefineableQElement< 2 >::Father_bound
static std::map< unsigned, DenseMatrix< int > > Father_bound
Definition: refineable_quad_element.h:176
oomph::RefineableQElement< 2 >::node_created_by_neighbour
virtual Node * node_created_by_neighbour(const Vector< double > &s_fraction, bool &is_periodic)
Definition: refineable_quad_element.cc:516
oomph::RefineableQElement< 2 >::setup_father_bounds
void setup_father_bounds()
Definition: refineable_quad_element.cc:46
oomph::RefineableQElement< 2 >::node_created_by_son_of_neighbour
virtual Node * node_created_by_son_of_neighbour(const Vector< double > &s_fraction, bool &is_periodic)
Definition: refineable_quad_element.h:126
oomph::TimeStepper::ntstorage
unsigned ntstorage() const
Definition: timesteppers.h:601
oomph::Tree::object_pt
RefineableElement * object_pt() const
Definition: tree.h:88
oomph::Tree::father_pt
Tree * father_pt() const
Return pointer to father: NULL if it's a root node.
Definition: tree.h:235
oomph::Tree::son_type
int son_type() const
Return son type.
Definition: tree.h:214
oomph::Tree::OMEGA
static const int OMEGA
Default value for an unassigned neighbour.
Definition: tree.h:262
s
RealScalar s
Definition: level1_cplx_impl.h:130
k
char char char int int * k
Definition: level2_impl.h:374
Eigen::zeta
EIGEN_STRONG_INLINE const Eigen::CwiseBinaryOp< Eigen::internal::scalar_zeta_op< typename DerivedX::Scalar >, const DerivedX, const DerivedQ > zeta(const Eigen::ArrayBase< DerivedX > &x, const Eigen::ArrayBase< DerivedQ > &q)
Definition: SpecialFunctionsArrayAPI.h:152
oomph::Global_string_for_annotation::string
std::string string(const unsigned &i)
Definition: oomph_definitions.cc:286
oomph::QuadTreeNames::SE
@ SE
Definition: quadtree.h:57
oomph::QuadTreeNames::NW
@ NW
Definition: quadtree.h:58
oomph::QuadTreeNames::NE
@ NE
Definition: quadtree.h:59
oomph::QuadTreeNames::SW
@ SW
Definition: quadtree.h:56
plotDoE.x
list x
Definition: plotDoE.py:28
plotPSD.t
t
Definition: plotPSD.py:36
OOMPH_EXCEPTION_LOCATION
#define OOMPH_EXCEPTION_LOCATION
Definition: oomph_definitions.h:61
OOMPH_CURRENT_FUNCTION
#define OOMPH_CURRENT_FUNCTION
Definition: oomph_definitions.h:86

References oomph::Mesh::add_boundary_node(), oomph::Mesh::add_node_pt(), oomph::Mesh::boundary_coordinate_exists(), oomph::ElementWithMovingNodes::enable_always_evaluate_dresidual_dnodal_coordinates_by_fd(), oomph::ElementWithMovingNodes::enable_bypass_fill_in_jacobian_from_geometric_data(), oomph::MacroElementNodeUpdateElementBase::geom_object_pt(), get_bcs(), get_boundaries(), oomph::RefineableElement::get_interpolated_values(), oomph::FiniteElement::get_node_at_local_coordinate(), oomph::RefineableSolidQElement< 2 >::get_solid_bcs(), oomph::FiniteElement::get_x(), i, interpolated_zeta_on_edge(), k, oomph::FiniteElement::Macro_elem_pt, oomph::Node::make_periodic(), oomph::ElementWithMovingNodes::method_for_shape_derivs(), oomph::Node::ndim(), oomph::QuadTreeNames::NE, oomph::FiniteElement::node_pt(), oomph::Node::nposition_type(), oomph::Data::nvalue(), oomph::QuadTreeNames::NW, oomph::Tree::object_pt(), oomph::Tree::OMEGA, OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION, oomph::Data::pin(), oomph::SolidNode::pin_position(), s, oomph::QElementBase::s_macro_ll(), oomph::QElementBase::s_macro_ur(), oomph::QuadTreeNames::SE, oomph::Node::set_coordinates_on_boundary(), oomph::MacroElementNodeUpdateElementBase::set_node_update_info(), oomph::Data::set_value(), oomph::AlgebraicElementBase::setup_algebraic_node_update(), oomph::Tree::son_type(), oomph::Global_string_for_annotation::string(), oomph::QuadTreeNames::SW, plotPSD::t, oomph::Data::time_stepper_pt(), plotDoE::x, oomph::Node::x(), and Eigen::zeta().

◆ check_integrity()

void oomph::RefineableQElement< 2 >::check_integrity ( double &  max_error)
virtual

Check the integrity of the element: ensure that the position and values are continuous across the element edges

Check inter-element continuity of

  • nodal positions
  • (nodally) interpolated function values

Implements oomph::RefineableElement.

1838  {
1839  using namespace QuadTreeNames;
1840 
1841  // Number of nodes along edge
1842  unsigned n_p = nnode_1d();
1843 
1844  // Number of timesteps (incl. present) for which continuity is
1845  // to be checked.
1846  unsigned n_time = 1;
1847 
1848  // Initialise errors
1849  max_error = 0.0;
1850  Vector<double> max_error_x(2, 0.0);
1851  double max_error_val = 0.0;
1852 
1853  Vector<int> edges(4);
1854  edges[0] = S;
1855  edges[1] = N;
1856  edges[2] = W;
1857  edges[3] = E;
1858 
1859  // Loop over the edges
1860  for (unsigned edge_counter = 0; edge_counter < 4; edge_counter++)
1861  {
1862  Vector<unsigned> translate_s(2);
1863  Vector<double> s(2), s_lo_neigh(2), s_hi_neigh(2), s_fraction(2);
1864  int neigh_edge, diff_level;
1865  bool in_neighbouring_tree;
1866 
1867  // Find pointer to neighbour in this direction
1868  QuadTree* neigh_pt;
1869  neigh_pt = quadtree_pt()->gteq_edge_neighbour(edges[edge_counter],
1870  translate_s,
1871  s_lo_neigh,
1872  s_hi_neigh,
1873  neigh_edge,
1874  diff_level,
1875  in_neighbouring_tree);
1876 
1877  // Neighbour exists and has existing nodes
1878  if ((neigh_pt != 0) && (neigh_pt->object_pt()->nodes_built()))
1879  {
1880  // Need to exclude periodic nodes from this check
1881  // There are only periodic nodes if we are in a neighbouring tree
1882  bool is_periodic = false;
1883  if (in_neighbouring_tree)
1884  {
1885  // Is it periodic
1886  is_periodic =
1887  tree_pt()->root_pt()->is_neighbour_periodic(edges[edge_counter]);
1888  }
1889 
1890  // Loop over nodes along the edge
1891  for (unsigned i0 = 0; i0 < n_p; i0++)
1892  {
1893  // Storage for pointer to the local node
1894  Node* local_node_pt = 0;
1895 
1896  switch (edge_counter)
1897  {
1898  case 0:
1899  // Local fraction of node
1900  s_fraction[0] = local_one_d_fraction_of_node(i0, 0);
1901  s_fraction[1] = 0.0;
1902  // Get pointer to local node
1903  local_node_pt = node_pt(i0);
1904  break;
1905 
1906  case 1:
1907  // Local fraction of node
1908  s_fraction[0] = local_one_d_fraction_of_node(i0, 0);
1909  s_fraction[1] = 1.0;
1910  // Get pointer to local node
1911  local_node_pt = node_pt(i0 + n_p * (n_p - 1));
1912  break;
1913 
1914  case 2:
1915  // Local fraction of node
1916  s_fraction[0] = 0.0;
1917  s_fraction[1] = local_one_d_fraction_of_node(i0, 1);
1918  // Get pointer to local node
1919  local_node_pt = node_pt(n_p * i0);
1920  break;
1921 
1922  case 3:
1923  // Local fraction of node
1924  s_fraction[0] = 1.0;
1925  s_fraction[1] = local_one_d_fraction_of_node(i0, 1);
1926  // Get pointer to local node
1927  local_node_pt = node_pt(n_p - 1 + n_p * i0);
1928  break;
1929  }
1930 
1931  // Calculate the local coordinate and the local coordinate as viewed
1932  // from the neighbour
1933  Vector<double> s_in_neighb(2);
1934  for (unsigned i = 0; i < 2; i++)
1935  {
1936  // Local coordinate in this element
1937  s[i] = -1.0 + 2.0 * s_fraction[i];
1938  // Local coordinate in the neighbour
1939  s_in_neighb[i] = s_lo_neigh[i] + s_fraction[translate_s[i]] *
1940  (s_hi_neigh[i] - s_lo_neigh[i]);
1941  }
1942 
1943  // Loop over timesteps
1944  for (unsigned t = 0; t < n_time; t++)
1945  {
1946  // Get the nodal position from neighbour element
1947  Vector<double> x_in_neighb(2);
1948  neigh_pt->object_pt()->interpolated_x(t, s_in_neighb, x_in_neighb);
1949 
1950  // Check error only if the node is NOT periodic
1951  if (is_periodic == false)
1952  {
1953  for (int i = 0; i < 2; i++)
1954  {
1955  // Find the spatial error
1956  double err = std::fabs(local_node_pt->x(t, i) - x_in_neighb[i]);
1957 
1958  // If it's bigger than our tolerance, say so
1959  if (err > 1e-9)
1960  {
1961  oomph_info << "errx " << err << " " << t << " "
1962  << local_node_pt->x(t, i) << " " << x_in_neighb[i]
1963  << std::endl;
1964 
1965  oomph_info << "at " << local_node_pt->x(0) << " "
1966  << local_node_pt->x(1) << std::endl;
1967  }
1968 
1969  // If it's bigger than the previous max error, it is the
1970  // new max error!
1971  if (err > max_error_x[i])
1972  {
1973  max_error_x[i] = err;
1974  }
1975  }
1976  }
1977 
1978  // Get the values from neighbour element. Note: # of values
1979  // gets set by routine (because in general we don't know
1980  // how many interpolated values a certain element has
1981  Vector<double> values_in_neighb;
1982  neigh_pt->object_pt()->get_interpolated_values(
1983  t, s_in_neighb, values_in_neighb);
1984 
1985  // Get the values in current element.
1986  Vector<double> values;
1987  get_interpolated_values(t, s, values);
1988 
1989  // Now figure out how many continuously interpolated values there
1990  // are
1991  unsigned num_val =
1992  neigh_pt->object_pt()->ncont_interpolated_values();
1993 
1994  // Check error
1995  for (unsigned ival = 0; ival < num_val; ival++)
1996  {
1997  double err = std::fabs(values[ival] - values_in_neighb[ival]);
1998 
1999  if (err > 1.0e-10)
2000  {
2001  oomph_info << local_node_pt->x(0) << " " << local_node_pt->x(1)
2002  << " \n# "
2003  << "erru (S)" << err << " " << ival << " "
2004  << get_node_number(local_node_pt) << " "
2005  << values[ival] << " " << values_in_neighb[ival]
2006  << std::endl;
2007  }
2008 
2009  if (err > max_error_val)
2010  {
2011  max_error_val = err;
2012  }
2013  }
2014  }
2015  }
2016  }
2017  }
2018 
2019  max_error = max_error_x[0];
2020  if (max_error_x[1] > max_error) max_error = max_error_x[1];
2021  if (max_error_val > max_error) max_error = max_error_val;
2022 
2023  if (max_error > 1e-9)
2024  {
2025  oomph_info << "\n#------------------------------------ \n#Max error ";
2026  oomph_info << max_error_x[0] << " " << max_error_x[1] << " "
2027  << max_error_val << std::endl;
2028  oomph_info << "#------------------------------------ \n " << std::endl;
2029  }
2030  }
e
Array< double, 1, 3 > e(1./3., 0.5, 2.)
oomph::FiniteElement::get_node_number
int get_node_number(Node *const &node_pt) const
Definition: elements.cc:3814
oomph::QuadTree::gteq_edge_neighbour
QuadTree * gteq_edge_neighbour(const int &direction, Vector< unsigned > &translate_s, Vector< double > &s_lo, Vector< double > &s_hi, int &edge, int &diff_level, bool &in_neighbouring_tree) const
Definition: quadtree.cc:413
oomph::RefineableElement::get_interpolated_values
virtual void get_interpolated_values(const Vector< double > &s, Vector< double > &values)
Definition: refineable_elements.h:417
oomph::TreeRoot::is_neighbour_periodic
bool is_neighbour_periodic(const int &direction)
Definition: tree.h:364
oomph::Tree::root_pt
TreeRoot *& root_pt()
Return pointer to root of the tree.
Definition: tree.h:141
N
@ N
Definition: constructor.cpp:22
Global_Physical_Variables::E
double E
Elastic modulus.
Definition: TwenteMeshGluing.cpp:68
MeshRefinement::max_error
double max_error
Definition: MortaringCantileverCompareToNonMortaring.cpp:188
boost::multiprecision::fabs
Real fabs(const Real &a)
Definition: boostmultiprec.cpp:117
oomph::QuadTreeNames::S
@ S
Definition: quadtree.h:62
oomph::QuadTreeNames::W
@ W
Definition: quadtree.h:63
oomph::oomph_info
OomphInfo oomph_info
Definition: oomph_definitions.cc:319

References Global_Physical_Variables::E, e(), boost::multiprecision::fabs(), oomph::RefineableElement::get_interpolated_values(), oomph::QuadTree::gteq_edge_neighbour(), i, oomph::FiniteElement::interpolated_x(), MeshRefinement::max_error, N, oomph::RefineableElement::ncont_interpolated_values(), oomph::RefineableElement::nodes_built(), oomph::Tree::object_pt(), oomph::oomph_info, s, oomph::QuadTreeNames::S, plotPSD::t, oomph::QuadTreeNames::W, and oomph::Node::x().

◆ further_setup_hanging_nodes()

virtual void oomph::RefineableQElement< 2 >::further_setup_hanging_nodes ( )
pure virtual

Perform additional hanging node procedures for variables that are not interpolated by all nodes (e.g. lower order interpolations as for the pressure in Taylor Hood).

Reimplemented from oomph::RefineableElement.

Implemented in oomph::RefineableQYoungLaplaceElement< NNODE_1D >, oomph::RefineableQSphericalCrouzeixRaviartElement, oomph::RefineableQSphericalTaylorHoodElement, oomph::RefineableQSphericalAdvectionDiffusionElement< NNODE_1D >, oomph::RefineableQSphericalAdvectionDiffusionElement< 3 >, oomph::RefineablePolarCrouzeixRaviartElement, oomph::RefineablePolarTaylorHoodElement, oomph::RefineableLinearisedQTaylorHoodElement, oomph::RefineableLinearisedQCrouzeixRaviartElement, oomph::RefineableLinearisedAxisymmetricQTaylorHoodElement, oomph::RefineableLinearisedAxisymmetricQCrouzeixRaviartElement, oomph::PRefineableQElement< 2, INITIAL_NNODE_1D >, oomph::RefineableGeneralisedNewtonianAxisymmetricQCrouzeixRaviartElement, oomph::RefineableGeneralisedNewtonianAxisymmetricQTaylorHoodElement, oomph::RefineableAxisymmetricQCrouzeixRaviartElement, oomph::RefineableAxisymmetricQTaylorHoodElement, oomph::RefineableQGeneralisedAxisymAdvectionDiffusionElement< NNODE_1D >, oomph::RefineableQAxisymAdvectionDiffusionElement< NNODE_1D >, oomph::RefineableQAxisymAdvectionDiffusionElement< 3 >, oomph::RefineablePolarStreamfunctionElement, oomph::RefineableBuoyantQSphericalCrouzeixRaviartElement, oomph::RefineableBuoyantQAxisymCrouzeixRaviartElement, oomph::RefineableBuoyantQAxisymCrouzeixRaviartElement, oomph::RefineableLinearisedAxisymmetricQCrouzeixRaviartElement, and oomph::RefineableLinearisedAxisymmetricQTaylorHoodElement.

◆ get_bcs()

void oomph::RefineableQElement< 2 >::get_bcs ( int  bound,
Vector< int > &  bound_cons 
) const

Determine Vector of boundary conditions along edge (or on vertex) bound (S/W/N/E/SW/SE/NW/NE): For value ival on this boundary, bound_cons[ival]=1 if pinned and 0 if free.

Determine Vector of boundary conditions along the element's boundary (or vertex) bound (S/W/N/E/SW/SE/NW/NE).

This function assumes that the same boundary condition is applied along the entire length of an element's edge (of course, the vertices combine the boundary conditions of their two adjacent edges in the most restrictive combination. Hence, if we're at a vertex, we apply the most restrictive boundary condition of the two adjacent edges. If we're on an edge (in its proper interior), we apply the least restrictive boundary condition of all nodes along the edge.

Usual convention:

  • bound_cons[ival]=0 if value ival on this boundary is free
  • bound_cons[ival]=1 if value ival on this boundary is pinned
146  {
147  using namespace QuadTreeNames;
148 
149  // Max. number of nodal data values in element
150  unsigned nvalue = ncont_interpolated_values();
151  // Set up temporary vectors to hold edge boundary conditions
152  Vector<int> bound_cons1(nvalue), bound_cons2(nvalue);
153 
154  // Which boundary are we on?
155  switch (bound)
156  {
157  // If on edge, just get the bcs
158  case N:
159  case S:
160  case W:
161  case E:
162  get_edge_bcs(bound, bound_cons);
163  break;
164 
165  // Most restrictive boundary at SE corner
166  case SE:
167  get_edge_bcs(S, bound_cons1);
168  get_edge_bcs(E, bound_cons2);
169 
170  for (unsigned k = 0; k < nvalue; k++)
171  {
172  bound_cons[k] = (bound_cons1[k] || bound_cons2[k]);
173  }
174  break;
175 
176  // Most restrictive boundary at SW corner
177  case SW:
178  get_edge_bcs(S, bound_cons1);
179  get_edge_bcs(W, bound_cons2);
180 
181  for (unsigned k = 0; k < nvalue; k++)
182  {
183  bound_cons[k] = (bound_cons1[k] || bound_cons2[k]);
184  }
185  break;
186 
187  // Most restrictive boundary at NW corner
188  case NW:
189  get_edge_bcs(N, bound_cons1);
190  get_edge_bcs(W, bound_cons2);
191 
192  for (unsigned k = 0; k < nvalue; k++)
193  {
194  bound_cons[k] = (bound_cons1[k] || bound_cons2[k]);
195  }
196  break;
197 
198  // Most restrictive boundary at NE corner
199  case NE:
200  get_edge_bcs(N, bound_cons1);
201  get_edge_bcs(E, bound_cons2);
202 
203  for (unsigned k = 0; k < nvalue; k++)
204  {
205  bound_cons[k] = (bound_cons1[k] || bound_cons2[k]);
206  }
207  break;
208 
209  default:
210  throw OomphLibError(
211  "Wrong boundary", OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
212  }
213  }
oomph::RefineableQElement< 2 >::get_edge_bcs
void get_edge_bcs(const int &edge, Vector< int > &bound_cons) const
Determine Vector of boundary conditions along edge (N/S/W/E)
Definition: refineable_quad_element.cc:224

References Global_Physical_Variables::E, k, N, oomph::QuadTreeNames::NE, oomph::QuadTreeNames::NW, OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION, oomph::QuadTreeNames::S, oomph::QuadTreeNames::SE, oomph::QuadTreeNames::SW, and oomph::QuadTreeNames::W.

Referenced by build(), oomph::PRefineableQElement< 2, INITIAL_NNODE_1D >::p_refine(), oomph::PRefineableQElement< 2, INITIAL_NNODE_1D >::rebuild_from_sons(), and oomph::RefineableQSpectralElement< 2 >::rebuild_from_sons().

◆ get_boundaries()

void oomph::RefineableQElement< 2 >::get_boundaries ( const int &  edge,
std::set< unsigned > &  boundary 
) const

Determine set of (mesh) boundaries that the element edge/vertex lives on

Given an element edge/vertex, return a set that contains all the (mesh-)boundary numbers that this element edge/vertex lives on.

For proper edges, the boundary is the one (if any) that is shared by both vertex nodes). For vertex nodes, we just return their boundaries.

290  {
291  using namespace QuadTreeNames;
292 
293  // Number of 1d nodes along an edge
294  unsigned n_p = nnode_1d();
295  // Left and right-hand nodes
296  int left_node = -1, right_node = -1;
297 
298  // Set the left (lower) and right (upper) nodes for the edge
299  switch (edge)
300  {
301  case N:
302  left_node = n_p * (n_p - 1);
303  right_node = n_p * n_p - 1;
304  break;
305 
306  case S:
307  left_node = 0;
308  right_node = n_p - 1;
309  break;
310 
311  case W:
312  left_node = 0;
313  right_node = n_p * (n_p - 1);
314  break;
315 
316  case E:
317  left_node = n_p - 1;
318  right_node = n_p * n_p - 1;
319  break;
320 
321  // Vertices do not have left nodes!
322  case SE:
323  right_node = n_p - 1;
324  break;
325 
326  case SW:
327  right_node = 0;
328  break;
329 
330  case NE:
331  right_node = n_p * n_p - 1;
332  break;
333 
334  case NW:
335  right_node = n_p * (n_p - 1);
336  break;
337 
338  default:
339  std::ostringstream error_stream;
340  error_stream << "Wrong edge " << edge << " passed" << std::endl;
341 
342  throw OomphLibError(
343  error_stream.str(), OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
344  }
345 
346  // Empty the boundary set: Edge does not live on any boundary
347  boundary.clear();
348 
349  // Storage for the nodes that the right node lives on
350  std::set<unsigned>* right_boundaries_pt = 0;
351  // Get the boundaries that the right node lives on
352  node_pt(right_node)->get_boundaries_pt(right_boundaries_pt);
353 
354  // If the right node lives on some boundaries
355  if (right_boundaries_pt != 0)
356  {
357  // If the node is a vertex then add the boundaries at the node
358  // into the vector boundary
359  if (left_node < 0)
360  {
361  copy(right_boundaries_pt->begin(),
362  right_boundaries_pt->end(),
363  inserter(boundary, boundary.begin()));
364  }
365  // Otherwise only add if the boundary also exists at the left hand node
366  else
367  {
368  std::set<unsigned>* left_boundaries_pt = 0;
369  node_pt(left_node)->get_boundaries_pt(left_boundaries_pt);
370  // If the left node is on some boundaries
371  if (left_boundaries_pt != 0)
372  {
373  // Use the standard algorithms to compute the boundaries in
374  // common between the left and right nodes
375  std::set_intersection(right_boundaries_pt->begin(),
376  right_boundaries_pt->end(),
377  left_boundaries_pt->begin(),
378  left_boundaries_pt->end(),
379  inserter(boundary, boundary.begin()));
380  }
381  }
382  }
383  }
oomph::Node::get_boundaries_pt
virtual void get_boundaries_pt(std::set< unsigned > *&boundaries_pt)
Definition: nodes.h:1365
copy
EIGEN_BLAS_FUNC() copy(int *n, RealScalar *px, int *incx, RealScalar *py, int *incy)
Definition: level1_impl.h:32

References copy(), Global_Physical_Variables::E, N, oomph::QuadTreeNames::NE, oomph::QuadTreeNames::NW, OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION, oomph::QuadTreeNames::S, oomph::QuadTreeNames::SE, oomph::QuadTreeNames::SW, and oomph::QuadTreeNames::W.

Referenced by build(), oomph::PRefineableQElement< 2, INITIAL_NNODE_1D >::p_refine(), oomph::PRefineableQElement< 2, INITIAL_NNODE_1D >::rebuild_from_sons(), and oomph::RefineableQSpectralElement< 2 >::rebuild_from_sons().

◆ get_edge_bcs()

void oomph::RefineableQElement< 2 >::get_edge_bcs ( const int &  edge,
Vector< int > &  bound_cons 
) const
protected

Determine Vector of boundary conditions along edge (N/S/W/E)

Determine Vector of boundary conditions along the element's edge (S/N/W/E) – BC is the least restrictive combination of all the nodes on this edge

Usual convention:

  • bound_cons[ival]=0 if value ival on this boundary is free
  • bound_cons[ival]=1 if value ival on this boundary is pinned
226  {
227  using namespace QuadTreeNames;
228 
229  // Number of nodes along 1D edge
230  unsigned n_p = nnode_1d();
231  // Left- and Right-hand nodes
232  unsigned left_node, right_node;
233 
234  // Set the left (lower) and right (upper) nodes for the edge
235  switch (edge)
236  {
237  case N:
238  left_node = n_p * (n_p - 1);
239  right_node = n_p * n_p - 1;
240  break;
241 
242  case S:
243  left_node = 0;
244  right_node = n_p - 1;
245  break;
246 
247  case W:
248  left_node = 0;
249  right_node = n_p * (n_p - 1);
250  break;
251 
252  case E:
253  left_node = n_p - 1;
254  right_node = n_p * n_p - 1;
255  break;
256 
257  default:
258  std::ostringstream error_stream;
259  error_stream << "Wrong edge " << edge << " passed to get_edge_bcs(..)"
260  << std::endl;
261 
262  throw OomphLibError(
263  error_stream.str(), OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
264  }
265 
266  // Max. number of nodal data values in element
267  unsigned maxnvalue = ncont_interpolated_values();
268 
269  // Loop over all values, the least restrictive value is
270  // the boundary condition at the left multiplied by that at the right
271  // Assuming that free is always zero and pinned is one
272  for (unsigned k = 0; k < maxnvalue; k++)
273  {
274  bound_cons[k] =
275  node_pt(left_node)->is_pinned(k) * node_pt(right_node)->is_pinned(k);
276  }
277  }
oomph::Data::is_pinned
bool is_pinned(const unsigned &i) const
Test whether the i-th variable is pinned (1: true; 0: false).
Definition: nodes.h:417

References Global_Physical_Variables::E, k, N, OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION, oomph::QuadTreeNames::S, and oomph::QuadTreeNames::W.

◆ interpolated_zeta_on_edge()

void oomph::RefineableQElement< 2 >::interpolated_zeta_on_edge ( const unsigned &  boundary,
const int &  edge,
const Vector< double > &  s,
Vector< double > &  zeta 
)

Return the value of the intrinsic boundary coordinate interpolated along the edge (S/W/N/E)

395  {
396  using namespace QuadTreeNames;
397 
398  // Number of 1D nodes along an edge
399  unsigned n_p = nnode_1d();
400 
401  // Storage for the shape functions
402  Shape psi(n_p * n_p);
403  // Get the shape functions at the passed position
404  this->shape(s, psi);
405 
406  // Unsigned data that give starts and multipliers for the loop
407  // over the nodes on the edges.
408  unsigned start = 0, multiplier = 1;
409 
410  // Which edge?
411  switch (edge)
412  {
413  case S:
414 #ifdef PARANOID
415  if (s[1] != -1.0)
416  {
417  std::ostringstream error_stream;
418  error_stream << "Coordinate " << s[0] << " " << s[1]
419  << " is not on South edge\n";
420 
421  throw OomphLibError(error_stream.str(),
422  OOMPH_CURRENT_FUNCTION,
423  OOMPH_EXCEPTION_LOCATION);
424  }
425 #endif
426  // Start is zero and multiplier is one
427  break;
428 
429  case N:
430 #ifdef PARANOID
431  if (s[1] != 1.0)
432  {
433  std::ostringstream error_stream;
434  error_stream << "Coordinate " << s[0] << " " << s[1]
435  << " is not on North edge\n";
436 
437  throw OomphLibError(error_stream.str(),
438  OOMPH_CURRENT_FUNCTION,
439  OOMPH_EXCEPTION_LOCATION);
440  }
441 #endif
442  // Start from the top left corner of the element, multiplier still one
443  start = n_p * (n_p - 1);
444  break;
445 
446  case W:
447 #ifdef PARANOID
448  if (s[0] != -1.0)
449  {
450  std::ostringstream error_stream;
451  error_stream << "Coordinate " << s[0] << " " << s[1]
452  << " is not on West edge\n";
453 
454  throw OomphLibError(error_stream.str(),
455  OOMPH_CURRENT_FUNCTION,
456  OOMPH_EXCEPTION_LOCATION);
457  }
458 #endif
459  // Loop over left-hand edge of element (start from zero)
460  multiplier = n_p;
461  break;
462 
463  case E:
464 #ifdef PARANOID
465  if (s[0] != 1.0)
466  {
467  std::ostringstream error_stream;
468  error_stream << "Coordinate " << s[0] << " " << s[1]
469  << " is not on East edge\n";
470 
471  throw OomphLibError(error_stream.str(),
472  OOMPH_CURRENT_FUNCTION,
473  OOMPH_EXCEPTION_LOCATION);
474  }
475 #endif
476  // Start from the bottom right-hand corner
477  start = n_p - 1;
478  // Loop over the right-hand edge of the element
479  multiplier = n_p;
480  break;
481 
482 
483  default:
484  std::ostringstream error_stream;
485  error_stream << "Wrong edge " << edge << " passed" << std::endl;
486 
487  throw OomphLibError(
488  error_stream.str(), OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
489  }
490 
491  // Initialise the intrinsic coordinate
492  double inter_zeta = 0.0;
493  // Loop over the nodes on the edge
494  for (unsigned n = 0; n < n_p; n++)
495  {
496  // Get the node number
497  unsigned node_number = start + multiplier * n;
498  // Now get the intrinsic coordinate
499  node_pt(node_number)->get_coordinates_on_boundary(boundary, zeta);
500  // Now multiply by the shape function
501  inter_zeta += zeta[0] * psi(node_number);
502  }
503 
504  // Set the value of the intrinsic coordinate
505  zeta[0] = inter_zeta;
506  }
n
const unsigned n
Definition: CG3DPackingUnitTest.cpp:11
oomph::FiniteElement::shape
virtual void shape(const Vector< double > &s, Shape &psi) const =0
oomph::Node::get_coordinates_on_boundary
virtual void get_coordinates_on_boundary(const unsigned &b, const unsigned &k, Vector< double > &boundary_zeta)
Definition: nodes.cc:2379
Global_Physical_Variables::multiplier
double multiplier(const Vector< double > &xi)
Definition: disk_oscillation.cc:63
oomph::CumulativeTimings::start
void start(const unsigned &i)
(Re-)start i-th timer
Definition: oomph_utilities.cc:243

References Global_Physical_Variables::E, Global_Physical_Variables::multiplier(), n, N, OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION, s, oomph::QuadTreeNames::S, oomph::OneDimLagrange::shape(), oomph::CumulativeTimings::start(), oomph::QuadTreeNames::W, and Eigen::zeta().

Referenced by build(), oomph::PRefineableQElement< 2, INITIAL_NNODE_1D >::p_refine(), oomph::PRefineableQElement< 2, INITIAL_NNODE_1D >::rebuild_from_sons(), and oomph::RefineableQSpectralElement< 2 >::rebuild_from_sons().

◆ node_created_by_neighbour()

Node * oomph::RefineableQElement< 2 >::node_created_by_neighbour ( const Vector< double > &  s_fraction,
bool &  is_periodic 
)
virtual

If a neighbouring element has already created a node at a position corresponding to the local fractional position within the present element, s_fraction, return a pointer to that node. If not, return NULL (0). If the node is on a periodic boundary the flag is_periodic is true, otherwise it will be false.

If a neighbouring element has already created a node at a position corresponding to the local fractional position within the present element, s_fraction, return a pointer to that node. If not, return NULL (0). If the node is periodic the flag is_periodic will be true

Reimplemented in oomph::PRefineableQElement< 2, INITIAL_NNODE_1D >.

518  {
519  using namespace QuadTreeNames;
520 
521  // Calculate the edges on which the node lies
522  Vector<int> edges;
523  if (s_fraction[0] == 0.0)
524  {
525  edges.push_back(W);
526  }
527  if (s_fraction[0] == 1.0)
528  {
529  edges.push_back(E);
530  }
531  if (s_fraction[1] == 0.0)
532  {
533  edges.push_back(S);
534  }
535  if (s_fraction[1] == 1.0)
536  {
537  edges.push_back(N);
538  }
539 
540  // Find the number of edges
541  unsigned n_size = edges.size();
542  // If there are no edges, then there is no neighbour, return 0
543  if (n_size == 0)
544  {
545  return 0;
546  }
547 
548  Vector<unsigned> translate_s(2);
549  Vector<double> s_lo_neigh(2);
550  Vector<double> s_hi_neigh(2);
551  Vector<double> s(2);
552 
553  int neigh_edge, diff_level;
554  bool in_neighbouring_tree;
555 
556  // Loop over the edges
557  for (unsigned j = 0; j < n_size; j++)
558  {
559  // Find pointer to neighbouring element along edge
560  QuadTree* neigh_pt;
561  neigh_pt = quadtree_pt()->gteq_edge_neighbour(edges[j],
562  translate_s,
563  s_lo_neigh,
564  s_hi_neigh,
565  neigh_edge,
566  diff_level,
567  in_neighbouring_tree);
568 
569  // Neighbour exists
570  if (neigh_pt != 0)
571  {
572  // Have its nodes been created yet?
573  // if(neigh_pt->object_pt()->node_pt(0)!=0)
574  if (neigh_pt->object_pt()->nodes_built())
575  {
576  // We now need to translate the nodal location
577  // as defined in terms of the fractional coordinates of the present
578  // element into those of its neighbour
579 
580  // Calculate the local coordinate in the neighbour
581  // Note that we need to use the translation scheme in case
582  // the local coordinates are swapped in the neighbour.
583  for (unsigned i = 0; i < 2; i++)
584  {
585  s[i] = s_lo_neigh[i] +
586  s_fraction[translate_s[i]] * (s_hi_neigh[i] - s_lo_neigh[i]);
587  }
588 
589  // Find the node in the neighbour
590  Node* neighbour_node_pt =
591  neigh_pt->object_pt()->get_node_at_local_coordinate(s);
592 
593  // If there is a node, return it
594  if (neighbour_node_pt != 0)
595  {
596  // Now work out whether it's a periodic boundary
597  // only possible if we have moved into a neighbouring tree
598  if (in_neighbouring_tree)
599  {
600  // Return whether the neighbour is periodic
601  is_periodic =
602  quadtree_pt()->root_pt()->is_neighbour_periodic(edges[j]);
603  }
604  // Return the pointer to the neighbouring node
605  return neighbour_node_pt;
606  }
607  }
608  }
609  }
610  // Node not found, return null
611  return 0;
612  }
j
std::ptrdiff_t j
Definition: tut_arithmetic_redux_minmax.cpp:2

References Global_Physical_Variables::E, oomph::FiniteElement::get_node_at_local_coordinate(), oomph::QuadTree::gteq_edge_neighbour(), i, j, N, oomph::RefineableElement::nodes_built(), oomph::Tree::object_pt(), s, oomph::QuadTreeNames::S, and oomph::QuadTreeNames::W.

◆ node_created_by_son_of_neighbour()

virtual Node* oomph::RefineableQElement< 2 >::node_created_by_son_of_neighbour ( const Vector< double > &  s_fraction,
bool &  is_periodic 
)
inlinevirtual

If a son of a neighbouring element has already created a node at a position corresponding to the local fractional position within the present element, s_fraction, return a pointer to that node. If not, return NULL (0). If the node is on a periodic boundary the flag is_periodic is true, otherwise it will be false.

Reimplemented in oomph::PRefineableQElement< 2, INITIAL_NNODE_1D >.

128  {
129  // It is impossible for this situation to arise in meshes
130  // containing elements of uniform p-order. This is here so
131  // that it can be overloaded for p-refineable elements.
132  return 0;
133  }

◆ output_corners()

void oomph::RefineableQElement< 2 >::output_corners ( std::ostream &  outfile,
const std::string &  colour 
) const

Print corner nodes, use colour.

Print corner nodes, use colour (default "BLACK")

1445  {
1446  Vector<double> s(2);
1447  Vector<double> corner(2);
1448 
1449  outfile << "ZONE I=2,J=2, C=" << colour << std::endl;
1450 
1451  s[0] = -1.0;
1452  s[1] = -1.0;
1453  get_x(s, corner);
1454  outfile << corner[0] << " " << corner[1] << " " << Number << std::endl;
1455 
1456  s[0] = 1.0;
1457  s[1] = -1.0;
1458  get_x(s, corner);
1459  outfile << corner[0] << " " << corner[1] << " " << Number << std::endl;
1460 
1461  s[0] = -1.0;
1462  s[1] = 1.0;
1463  get_x(s, corner);
1464  outfile << corner[0] << " " << corner[1] << " " << Number << std::endl;
1465 
1466  s[0] = 1.0;
1467  s[1] = 1.0;
1468  get_x(s, corner);
1469  outfile << corner[0] << " " << corner[1] << " " << Number << std::endl;
1470 
1471 
1472  outfile << "TEXT CS = GRID, X = " << corner[0] << ",Y = " << corner[1]
1473  << ", HU = GRID, H = 0.01, AN = MIDCENTER, T =\"" << Number << "\""
1474  << std::endl;
1475  }
oomph::RefineableElement::Number
long Number
Global element number – for plotting/validation purposes.
Definition: refineable_elements.h:115
oomph::TecplotNames::colour
Vector< std::string > colour
Tecplot colours.
Definition: oomph_utilities.cc:1159

References oomph::TecplotNames::colour, oomph::Global_unsigned::Number, and s.

◆ quad_hang_helper()

void oomph::RefineableQElement< 2 >::quad_hang_helper ( const int &  value_id,
const int &  my_edge,
std::ofstream &  output_hangfile 
)
protectedvirtual

Internal helper function that is used to construct the hanging node schemes for the positions.

Internal function to set up the hanging nodes on a particular edge of the element

Reimplemented in oomph::PRefineableQElement< 2, INITIAL_NNODE_1D >.

1525  {
1526  using namespace QuadTreeNames;
1527 
1528  Vector<unsigned> translate_s(2);
1529  Vector<double> s_lo_neigh(2);
1530  Vector<double> s_hi_neigh(2);
1531  int neigh_edge, diff_level;
1532  bool in_neighbouring_tree;
1533 
1534  // Find pointer to neighbour in this direction
1535  QuadTree* neigh_pt;
1536  neigh_pt = quadtree_pt()->gteq_edge_neighbour(my_edge,
1537  translate_s,
1538  s_lo_neigh,
1539  s_hi_neigh,
1540  neigh_edge,
1541  diff_level,
1542  in_neighbouring_tree);
1543 
1544  // Neighbour exists and all nodes have been created
1545  if (neigh_pt != 0)
1546  {
1547  // Different sized element?
1548  if (diff_level != 0)
1549  {
1550  // Test for the periodic node case
1551  // Are we crossing a periodic boundary
1552  bool is_periodic = false;
1553  if (in_neighbouring_tree)
1554  {
1555  is_periodic = tree_pt()->root_pt()->is_neighbour_periodic(my_edge);
1556  }
1557 
1558  // If it is periodic we actually need to get the node in
1559  // the neighbour of the neighbour (which will be a parent of
1560  // the present element) so that the "fixed" coordinate is
1561  // correctly calculated.
1562  // The idea is to replace the neigh_pt and associated data
1563  // with those of the neighbour of the neighbour
1564  if (is_periodic)
1565  {
1566  // Required data for the neighbour finding routine
1567  Vector<unsigned> translate_s_in_neigh(2);
1568  Vector<double> s_lo_neigh_of_neigh(2);
1569  Vector<double> s_hi_neigh_of_neigh(2);
1570  int neigh_edge_of_neigh, diff_level_of_neigh;
1571  bool in_neighbouring_tree_of_neigh;
1572 
1573  // Find pointer to neighbour of the neighbour on the edge
1574  // that we are currently considering
1575  QuadTree* neigh_of_neigh_pt;
1576  neigh_of_neigh_pt =
1577  neigh_pt->gteq_edge_neighbour(neigh_edge,
1578  translate_s_in_neigh,
1579  s_lo_neigh_of_neigh,
1580  s_hi_neigh_of_neigh,
1581  neigh_edge_of_neigh,
1582  diff_level_of_neigh,
1583  in_neighbouring_tree_of_neigh);
1584 
1585  // Set the value of the NEW neighbour and edge
1586  neigh_pt = neigh_of_neigh_pt;
1587  neigh_edge = neigh_edge_of_neigh;
1588 
1589  // Set the values of the translation schemes
1590  // Need to find the values of s_lo and s_hi
1591  // in the neighbour of the neighbour
1592 
1593  // Get the minimum and maximum values of the coordinate
1594  // in the neighbour (don't like this, but I think it's
1595  // necessary) Note that these values are hardcoded
1596  // in the quadtrees at some point!!
1597  double s_min = neigh_pt->object_pt()->s_min();
1598  double s_max = neigh_pt->object_pt()->s_max();
1599  Vector<double> s_lo_frac(2), s_hi_frac(2);
1600  // Work out the fractional position of the low and high points
1601  // of the original element
1602  for (unsigned i = 0; i < 2; i++)
1603  {
1604  s_lo_frac[i] = (s_lo_neigh[i] - s_min) / (s_max - s_min);
1605  s_hi_frac[i] = (s_hi_neigh[i] - s_min) / (s_max - s_min);
1606  }
1607 
1608  // We should now be able to construct the low and high points in
1609  // the neighbour of the neighbour
1610  for (unsigned i = 0; i < 2; i++)
1611  {
1612  s_lo_neigh[i] = s_lo_neigh_of_neigh[i] +
1613  s_lo_frac[translate_s_in_neigh[i]] *
1614  (s_hi_neigh_of_neigh[i] - s_lo_neigh_of_neigh[i]);
1615  s_hi_neigh[i] = s_lo_neigh_of_neigh[i] +
1616  s_hi_frac[translate_s_in_neigh[i]] *
1617  (s_hi_neigh_of_neigh[i] - s_lo_neigh_of_neigh[i]);
1618  }
1619 
1620  // Finally we must sort out the translation scheme
1621  Vector<unsigned> temp_translate(2);
1622  for (unsigned i = 0; i < 2; i++)
1623  {
1624  temp_translate[i] = translate_s_in_neigh[translate_s[i]];
1625  }
1626  for (unsigned i = 0; i < 2; i++)
1627  {
1628  translate_s[i] = temp_translate[i];
1629  }
1630  } // End of special treatment for periodic hanging nodes
1631 
1632  // Number of nodes in one dimension
1633  unsigned n_p = ninterpolating_node_1d(value_id);
1634  // Storage for the local nodes along the edge of the quadtree
1635  Node* local_node_pt = 0;
1636  // Loop over nodes along the edge
1637  for (unsigned i0 = 0; i0 < n_p; i0++)
1638  {
1639  // Storage for the fractional position of the node in the element
1640  Vector<double> s_fraction(2);
1641 
1642  // Find the local node and the fractional position of the node
1643  // which depends on the edge, of course
1644  switch (my_edge)
1645  {
1646  case N:
1647  s_fraction[0] =
1648  local_one_d_fraction_of_interpolating_node(i0, 0, value_id);
1649  s_fraction[1] = 1.0;
1650  local_node_pt =
1651  interpolating_node_pt(i0 + n_p * (n_p - 1), value_id);
1652  break;
1653 
1654  case S:
1655  s_fraction[0] =
1656  local_one_d_fraction_of_interpolating_node(i0, 0, value_id);
1657  s_fraction[1] = 0.0;
1658  local_node_pt = interpolating_node_pt(i0, value_id);
1659  break;
1660 
1661  case E:
1662  s_fraction[0] = 1.0;
1663  s_fraction[1] =
1664  local_one_d_fraction_of_interpolating_node(i0, 1, value_id);
1665  local_node_pt =
1666  interpolating_node_pt(n_p - 1 + n_p * i0, value_id);
1667  break;
1668 
1669  case W:
1670  s_fraction[0] = 0.0;
1671  s_fraction[1] =
1672  local_one_d_fraction_of_interpolating_node(i0, 1, value_id);
1673  local_node_pt = interpolating_node_pt(n_p * i0, value_id);
1674  break;
1675 
1676  default:
1677  throw OomphLibError("my_edge not N, S, W, E\n",
1678  OOMPH_CURRENT_FUNCTION,
1679  OOMPH_EXCEPTION_LOCATION);
1680  }
1681 
1682  // Calculate the local coordinates of the node in the neighbour
1683  Vector<double> s_in_neighb(2);
1684  for (unsigned i = 0; i < 2; i++)
1685  {
1686  s_in_neighb[i] = s_lo_neigh[i] + s_fraction[translate_s[i]] *
1687  (s_hi_neigh[i] - s_lo_neigh[i]);
1688  }
1689 
1690  // Find the Node in the neighbouring element
1691  Node* neighbouring_node_pt =
1692  neigh_pt->object_pt()->get_interpolating_node_at_local_coordinate(
1693  s_in_neighb, value_id);
1694 
1695  // If the neighbour does not have a node at this point
1696  if (0 == neighbouring_node_pt)
1697  {
1698  // Do we need to make a hanging node, we assume that we don't
1699  // initially
1700  bool make_hanging_node = false;
1701 
1702  // If the node is not hanging geometrically, then we must make
1703  // it hang
1704  if (!local_node_pt->is_hanging())
1705  {
1706  make_hanging_node = true;
1707  }
1708  // Otherwise, it could be hanging geometrically, but still
1709  // require a different hanging scheme for this data value
1710  else
1711  {
1712  if (local_node_pt->hanging_pt(value_id) ==
1713  local_node_pt->hanging_pt())
1714  {
1715  make_hanging_node = true;
1716  }
1717  }
1718 
1719  // If we do need to make the hanging node, then let's do it
1720  if (make_hanging_node == true)
1721  {
1722  // Cache refineable element used here
1723  RefineableElement* const obj_pt = neigh_pt->object_pt();
1724 
1725  // Get shape functions in neighbour element
1726  Shape psi(obj_pt->ninterpolating_node(value_id));
1727  obj_pt->interpolating_basis(s_in_neighb, psi, value_id);
1728 
1729  // Allocate the storage for the Hang pointer
1730  // which contains n_p nodes
1731  HangInfo* hang_pt = new HangInfo(n_p);
1732 
1733  // Loop over nodes on edge in neighbour and mark them as nodes
1734  // that this node depends on
1735  unsigned n_neighbour;
1736 
1737  // Number of nodes along edge in neighbour element
1738  for (unsigned n_edge = 0; n_edge < n_p; n_edge++)
1739  {
1740  switch (neigh_edge)
1741  {
1742  case N:
1743  n_neighbour = n_p * (n_p - 1) + n_edge;
1744  break;
1745 
1746  case S:
1747  n_neighbour = n_edge;
1748  break;
1749 
1750  case W:
1751  n_neighbour = n_p * n_edge;
1752  break;
1753 
1754  case E:
1755  n_neighbour = n_p * n_edge + (n_p - 1);
1756  break;
1757 
1758  default:
1759  throw OomphLibError("neigh_edge not N, S, W, E\n",
1760  OOMPH_CURRENT_FUNCTION,
1761  OOMPH_EXCEPTION_LOCATION);
1762  }
1763 
1764  // Push back neighbouring node and weight into
1765  // Vector of (pointers to)
1766  // master nodes and weights
1767  // The weight is merely the value of the shape function
1768  // corresponding to the node in the neighbour
1769  hang_pt->set_master_node_pt(
1770  n_edge,
1771  obj_pt->interpolating_node_pt(n_neighbour, value_id),
1772  psi[n_neighbour]);
1773  }
1774 
1775  // Now set the hanging data for the position
1776  // This also constrains the data values associated with the
1777  // value id
1778  local_node_pt->set_hanging_pt(hang_pt, value_id);
1779  }
1780 
1781  // Dump the output if the file has been openeed
1782  if (output_hangfile.is_open())
1783  {
1784  output_hangfile << local_node_pt->x(0) << " "
1785  << local_node_pt->x(1) << std::endl;
1786  }
1787  }
1788  // Otherwise check that the nodes are the same
1789  else
1790  {
1791 #ifdef PARANOID
1792  if (local_node_pt != neighbouring_node_pt)
1793  {
1794  std::ostringstream warning_stream;
1795  warning_stream << "SANITY CHECK in quad_hang_helper \n"
1796  << "Current node " << local_node_pt << " at "
1797  << "(" << local_node_pt->x(0) << ", "
1798  << local_node_pt->x(1) << ")"
1799  << " is not hanging and has " << std::endl
1800  << "Neighbour's node " << neighbouring_node_pt
1801  << " at "
1802  << "(" << neighbouring_node_pt->x(0) << ", "
1803  << neighbouring_node_pt->x(1) << ")" << std::endl
1804  << "even though the two should be "
1805  << "identical" << std::endl;
1806  OomphLibWarning(warning_stream.str(),
1807  "RefineableQElement<2>::quad_hang_helper()",
1808  OOMPH_EXCEPTION_LOCATION);
1809  }
1810 #endif
1811  }
1812 
1813  // If we are doing the position, then
1814  if (value_id == -1)
1815  {
1816  // Get the nodal position from neighbour element
1817  Vector<double> x_in_neighb(2);
1818  neigh_pt->object_pt()->interpolated_x(s_in_neighb, x_in_neighb);
1819 
1820  // Fine adjust the coordinates (macro map will pick up boundary
1821  // accurately but will lead to different element edges)
1822  local_node_pt->x(0) = x_in_neighb[0];
1823  local_node_pt->x(1) = x_in_neighb[1];
1824  }
1825  }
1826  }
1827  }
1828  }
oomph::FiniteElement::s_min
virtual double s_min() const
Min value of local coordinate.
Definition: elements.h:2793
oomph::FiniteElement::s_max
virtual double s_max() const
Max. value of local coordinate.
Definition: elements.h:2803
oomph::RefineableElement::ninterpolating_node_1d
virtual unsigned ninterpolating_node_1d(const int &value_id)
Definition: refineable_elements.h:474
oomph::RefineableElement::interpolating_node_pt
virtual Node * interpolating_node_pt(const unsigned &n, const int &value_id)
Definition: refineable_elements.h:437
oomph::RefineableElement::local_one_d_fraction_of_interpolating_node
virtual double local_one_d_fraction_of_interpolating_node(const unsigned &n1d, const unsigned &i, const int &value_id)
Definition: refineable_elements.h:447

References Global_Physical_Variables::E, oomph::RefineableElement::get_interpolating_node_at_local_coordinate(), oomph::QuadTree::gteq_edge_neighbour(), oomph::Node::hanging_pt(), i, oomph::FiniteElement::interpolated_x(), oomph::RefineableElement::interpolating_basis(), oomph::RefineableElement::interpolating_node_pt(), oomph::Node::is_hanging(), N, oomph::RefineableElement::ninterpolating_node(), oomph::Tree::object_pt(), OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION, oomph::QuadTreeNames::S, oomph::FiniteElement::s_max(), oomph::FiniteElement::s_min(), oomph::Node::set_hanging_pt(), oomph::HangInfo::set_master_node_pt(), oomph::QuadTreeNames::W, and oomph::Node::x().

◆ quadtree_pt() [1/2]

QuadTree* oomph::RefineableQElement< 2 >::quadtree_pt ( )
inline

Pointer to quadtree representation of this element.

153  {
154  return dynamic_cast<QuadTree*>(Tree_pt);
155  }

Referenced by oomph::RefineableLinearisedAxisymmetricQCrouzeixRaviartElement::further_build(), oomph::RefineableAxisymmetricQCrouzeixRaviartElement::further_build(), oomph::RefineableGeneralisedNewtonianAxisymmetricQCrouzeixRaviartElement::further_build(), oomph::RefineableLinearisedQCrouzeixRaviartElement::further_build(), oomph::RefineablePolarCrouzeixRaviartElement::further_build(), oomph::RefineableQSphericalCrouzeixRaviartElement::further_build(), oomph::RefineableLinearisedAxisymmetricQCrouzeixRaviartElement::rebuild_from_sons(), oomph::RefineableAxisymmetricQCrouzeixRaviartElement::rebuild_from_sons(), oomph::RefineableGeneralisedNewtonianAxisymmetricQCrouzeixRaviartElement::rebuild_from_sons(), oomph::RefineableLinearisedQCrouzeixRaviartElement::rebuild_from_sons(), oomph::RefineablePolarCrouzeixRaviartElement::rebuild_from_sons(), and oomph::RefineableQSphericalCrouzeixRaviartElement::rebuild_from_sons().

◆ quadtree_pt() [2/2]

QuadTree* oomph::RefineableQElement< 2 >::quadtree_pt ( ) const
inline

Pointer to quadtree representation of this element.

159  {
160  return dynamic_cast<QuadTree*>(Tree_pt);
161  }

◆ required_nsons()

unsigned oomph::RefineableQElement< 2 >::required_nsons ( ) const
inlinevirtual

A refineable quad element has four sons.

Reimplemented from oomph::RefineableElement.

107  {
108  return 4;
109  }

◆ setup_father_bounds()

void oomph::RefineableQElement< 2 >::setup_father_bounds ( )
protected

Setup static matrix for coincidence between son nodal points and father boundaries

Setup static matrix for coincidence between son nodal points and father boundaries:

Father_boundd[nnode_1d](nnode_son,son_type)={SW/SE/NW/NE/S/E/N/W/OMEGA}

so that node nnode_son in element of type son_type lies on boundary/vertex Father_boundd[nnode_1d](nnode_son,son_type) in its father element. If the node doesn't lie on a boundary the value is OMEGA.

47  {
48  using namespace QuadTreeNames;
49 
50  // Find the number of nodes along a 1D edge
51  unsigned n_p = nnode_1d();
52  // Allocate space for the boundary information
53  Father_bound[n_p].resize(n_p * n_p, 4);
54 
55  // Initialise: By default points are not on the boundary
56  for (unsigned n = 0; n < n_p * n_p; n++)
57  {
58  for (unsigned ison = 0; ison < 4; ison++)
59  {
60  Father_bound[n_p](n, ison) = Tree::OMEGA;
61  }
62  }
63 
64  // Southwest son
65  //--------------
66  // SW node (0) is the SW node of the parent
67  Father_bound[n_p](0, SW) = SW;
68  // Southern boundary is the southern boundary of the parent
69  for (unsigned n = 1; n < n_p; n++)
70  {
71  Father_bound[n_p](n, SW) = S;
72  }
73  // Western boundary is the western boundary of the parent
74  for (unsigned n = 1; n < n_p; n++)
75  {
76  Father_bound[n_p](n_p * n, SW) = W;
77  }
78  // Other boundaries are in the interior
79 
80  // Northwest son
81  //--------------
82  // NW node (n_p*(n_p-1))is the NW node of the parent
83  Father_bound[n_p](n_p * (n_p - 1), NW) = NW;
84  // Northern boundary is the northern boundary of the parent
85  for (unsigned n = 1; n < n_p; n++)
86  {
87  Father_bound[n_p](n_p * (n_p - 1) + n, NW) = N;
88  }
89  // Western boundary is the western boundary of the parent
90  for (unsigned n = 0; n < (n_p - 1); n++)
91  {
92  Father_bound[n_p](n_p * n, NW) = W;
93  }
94  // Other boundaries are in the interior
95 
96  // Northeast son
97  //--------------
98  // NE node (n_p*n_p-1) is the NE node of the parent
99  Father_bound[n_p](n_p * n_p - 1, NE) = NE;
100  // Northern boundary is the northern boundary of the parent
101  for (unsigned n = 0; n < (n_p - 1); n++)
102  {
103  Father_bound[n_p](n_p * (n_p - 1) + n, NE) = N;
104  }
105  // Eastern boundary is the eastern boundary of the parent
106  for (unsigned n = 0; n < (n_p - 1); n++)
107  {
108  Father_bound[n_p](n_p - 1 + n_p * n, NE) = E;
109  }
110  // Other boundaries are in the interior
111 
112  // Southeast son
113  //--------------
114  // SE node (n_p-1) is the SE node of the parent
115  Father_bound[n_p](n_p - 1, SE) = SE;
116  // Southern boundary is the southern boundary of the parent
117  for (unsigned n = 0; n < (n_p - 1); n++)
118  {
119  Father_bound[n_p](n, SE) = S;
120  }
121  // Eastern boundary is the eastern boundary of the parent
122  for (unsigned n = 1; n < n_p; n++)
123  {
124  Father_bound[n_p](n_p - 1 + n_p * n, SE) = E;
125  }
126  }

References Global_Physical_Variables::E, n, N, oomph::QuadTreeNames::NE, oomph::QuadTreeNames::NW, oomph::Tree::OMEGA, oomph::QuadTreeNames::S, oomph::QuadTreeNames::SE, oomph::QuadTreeNames::SW, and oomph::QuadTreeNames::W.

◆ setup_hang_for_value()

void oomph::RefineableQElement< 2 >::setup_hang_for_value ( const int &  value_id)
protected

Internal helper function that is used to construct the hanging node schemes for the value_id-th interpolated value

Internal function that sets up the hanging node scheme for a particular continuously interpolated value

1507  {
1508  using namespace QuadTreeNames;
1509 
1510  std::ofstream dummy_hangfile;
1511  quad_hang_helper(value_id, S, dummy_hangfile);
1512  quad_hang_helper(value_id, N, dummy_hangfile);
1513  quad_hang_helper(value_id, W, dummy_hangfile);
1514  quad_hang_helper(value_id, E, dummy_hangfile);
1515  }
oomph::RefineableQElement< 2 >::quad_hang_helper
virtual void quad_hang_helper(const int &value_id, const int &my_edge, std::ofstream &output_hangfile)
Definition: refineable_quad_element.cc:1522

References Global_Physical_Variables::E, N, oomph::QuadTreeNames::S, and oomph::QuadTreeNames::W.

Referenced by oomph::RefineableLinearisedAxisymmetricQTaylorHoodElement::further_setup_hanging_nodes(), oomph::RefineableAxisymmetricQTaylorHoodElement::further_setup_hanging_nodes(), oomph::RefineableGeneralisedNewtonianAxisymmetricQTaylorHoodElement::further_setup_hanging_nodes(), oomph::RefineableLinearisedQTaylorHoodElement::further_setup_hanging_nodes(), oomph::RefineablePolarTaylorHoodElement::further_setup_hanging_nodes(), and oomph::RefineableQSphericalTaylorHoodElement::further_setup_hanging_nodes().

◆ setup_hanging_nodes()

void oomph::RefineableQElement< 2 >::setup_hanging_nodes ( Vector< std::ofstream * > &  output_stream)
virtual

Markup all hanging nodes & document the results in the output streams contained in the vector output_stream, if they are open.

Set up all hanging nodes. If we are documenting the output then open the output files and pass the open files to the helper function

Reimplemented from oomph::RefineableElement.

1483  {
1484 #ifdef PARANOID
1485  if (output_stream.size() != 4)
1486  {
1487  throw OomphLibError("There must be four output streams",
1488  OOMPH_CURRENT_FUNCTION,
1489  OOMPH_EXCEPTION_LOCATION);
1490  }
1491 #endif
1492 
1493  using namespace QuadTreeNames;
1494 
1495  // Setup hanging nodes on each edge of the element
1496  quad_hang_helper(-1, S, *(output_stream[0]));
1497  quad_hang_helper(-1, N, *(output_stream[1]));
1498  quad_hang_helper(-1, W, *(output_stream[2]));
1499  quad_hang_helper(-1, E, *(output_stream[3]));
1500  }

References Global_Physical_Variables::E, N, OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION, oomph::QuadTreeNames::S, and oomph::QuadTreeNames::W.

Member Data Documentation

◆ Father_bound

std::map< unsigned, DenseMatrix< int > > oomph::RefineableQElement< 2 >::Father_bound
staticprotected

Coincidence between son nodal points and father boundaries: Father_bound[node_1d](jnod_son,son_type)={SW/SE/NW/NE/S/E/N/W/OMEGA}

Static matrix for coincidence between son nodal points and father boundaries

The documentation for this class was generated from the following files: