oomph::CollapsibleChannelMesh< ELEMENT > Class Template Reference

#include <collapsible_channel_mesh.template.h>

+ Inheritance diagram for oomph::CollapsibleChannelMesh< ELEMENT >:

Public Member Functions

 CollapsibleChannelMesh (const unsigned &nup, const unsigned &ncollapsible, const unsigned &ndown, const unsigned &ny, const double &lup, const double &lcollapsible, const double &ldown, const double &ly, GeomObject *wall_pt, TimeStepper *time_stepper_pt=&Mesh::Default_TimeStepper)
 
 ~CollapsibleChannelMesh ()
 destructor More...
 
GeomObject *& wall_pt ()
 Access function to GeomObject representing wall. More...
 
CollapsibleChannelDomaindomain_pt ()
 Access function to domain. More...
 
virtual CollapsibleChannelDomain::BLSquashFctPtbl_squash_fct_pt ()
 
CollapsibleChannelDomain::BLSquashFctPt bl_squash_fct_pt () const
 
virtual CollapsibleChannelDomain::AxialSpacingFctPtaxial_spacing_fct_pt ()
 
virtual CollapsibleChannelDomain::AxialSpacingFctPtaxial_spacing_fct_pt () const
 
- Public Member Functions inherited from oomph::SimpleRectangularQuadMesh< ELEMENT >
 SimpleRectangularQuadMesh (const unsigned &Nx, const unsigned &Ny, const double &Lx, const double &Ly, TimeStepper *time_stepper_pt=&Mesh::Default_TimeStepper)
 
const unsignednx () const
 Access function for number of elements in x directions. More...
 
const unsignedny () const
 Access function for number of elements in y directions. More...
 
- Public Member Functions inherited from oomph::QuadMeshBase
 QuadMeshBase ()
 Constructor (empty) More...
 
 QuadMeshBase (const QuadMeshBase &node)=delete
 Broken copy constructor. More...
 
void operator= (const QuadMeshBase &)=delete
 Broken assignment operator. More...
 
virtual ~QuadMeshBase ()
 Destructor (empty) More...
 
void setup_boundary_element_info ()
 
void setup_boundary_element_info (std::ostream &outfile)
 
- Public Member Functions inherited from oomph::Mesh
 Mesh ()
 Default constructor. More...
 
 Mesh (const Vector< Mesh * > &sub_mesh_pt)
 
void merge_meshes (const Vector< Mesh * > &sub_mesh_pt)
 
virtual void reset_boundary_element_info (Vector< unsigned > &ntmp_boundary_elements, Vector< Vector< unsigned >> &ntmp_boundary_elements_in_region, Vector< FiniteElement * > &deleted_elements)
 Virtual function to perform the reset boundary elements info rutines. More...
 
template<class BULK_ELEMENT >
void doc_boundary_coordinates (const unsigned &b, std::ofstream &the_file)
 
virtual void scale_mesh (const double &factor)
 
 Mesh (const Mesh &dummy)=delete
 Broken copy constructor. More...
 
void operator= (const Mesh &)=delete
 Broken assignment operator. More...
 
virtual ~Mesh ()
 Virtual Destructor to clean up all memory. More...
 
void flush_element_and_node_storage ()
 
void flush_element_storage ()
 
void flush_node_storage ()
 
Node *& node_pt (const unsigned long &n)
 Return pointer to global node n. More...
 
Nodenode_pt (const unsigned long &n) const
 Return pointer to global node n (const version) More...
 
GeneralisedElement *& element_pt (const unsigned long &e)
 Return pointer to element e. More...
 
GeneralisedElementelement_pt (const unsigned long &e) const
 Return pointer to element e (const version) More...
 
const Vector< GeneralisedElement * > & element_pt () const
 Return reference to the Vector of elements. More...
 
Vector< GeneralisedElement * > & element_pt ()
 Return reference to the Vector of elements. More...
 
FiniteElementfinite_element_pt (const unsigned &e) const
 
Node *& boundary_node_pt (const unsigned &b, const unsigned &n)
 Return pointer to node n on boundary b. More...
 
Nodeboundary_node_pt (const unsigned &b, const unsigned &n) const
 Return pointer to node n on boundary b. More...
 
void set_nboundary (const unsigned &nbound)
 Set the number of boundaries in the mesh. More...
 
void remove_boundary_nodes ()
 Clear all pointers to boundary nodes. More...
 
void remove_boundary_nodes (const unsigned &b)
 
void remove_boundary_node (const unsigned &b, Node *const &node_pt)
 Remove a node from the boundary b. More...
 
void add_boundary_node (const unsigned &b, Node *const &node_pt)
 Add a (pointer to) a node to the b-th boundary. More...
 
void copy_boundary_node_data_from_nodes ()
 
bool boundary_coordinate_exists (const unsigned &i) const
 Indicate whether the i-th boundary has an intrinsic coordinate. More...
 
unsigned long nelement () const
 Return number of elements in the mesh. More...
 
unsigned long nnode () const
 Return number of nodes in the mesh. More...
 
unsigned ndof_types () const
 Return number of dof types in mesh. More...
 
unsigned elemental_dimension () const
 Return number of elemental dimension in mesh. More...
 
unsigned nodal_dimension () const
 Return number of nodal dimension in mesh. More...
 
void add_node_pt (Node *const &node_pt)
 Add a (pointer to a) node to the mesh. More...
 
void add_element_pt (GeneralisedElement *const &element_pt)
 Add a (pointer to) an element to the mesh. More...
 
virtual void node_update (const bool &update_all_solid_nodes=false)
 
virtual void reorder_nodes (const bool &use_old_ordering=true)
 
virtual void get_node_reordering (Vector< Node * > &reordering, const bool &use_old_ordering=true) const
 
template<class BULK_ELEMENT , template< class > class FACE_ELEMENT>
void build_face_mesh (const unsigned &b, Mesh *const &face_mesh_pt)
 
unsigned self_test ()
 Self-test: Check elements and nodes. Return 0 for OK. More...
 
void max_and_min_element_size (double &max_size, double &min_size)
 
double total_size ()
 
void check_inverted_elements (bool &mesh_has_inverted_elements, std::ofstream &inverted_element_file)
 
void check_inverted_elements (bool &mesh_has_inverted_elements)
 
unsigned check_for_repeated_nodes (const double &epsilon=1.0e-12)
 
Vector< Node * > prune_dead_nodes ()
 
unsigned nboundary () const
 Return number of boundaries. More...
 
unsigned long nboundary_node (const unsigned &ibound) const
 Return number of nodes on a particular boundary. More...
 
FiniteElementboundary_element_pt (const unsigned &b, const unsigned &e) const
 Return pointer to e-th finite element on boundary b. More...
 
Nodeget_some_non_boundary_node () const
 
unsigned nboundary_element (const unsigned &b) const
 Return number of finite elements that are adjacent to boundary b. More...
 
int face_index_at_boundary (const unsigned &b, const unsigned &e) const
 
virtual void dump (std::ofstream &dump_file, const bool &use_old_ordering=true) const
 Dump the data in the mesh into a file for restart. More...
 
void dump (const std::string &dump_file_name, const bool &use_old_ordering=true) const
 Dump the data in the mesh into a file for restart. More...
 
virtual void read (std::ifstream &restart_file)
 Read solution from restart file. More...
 
void output_paraview (std::ofstream &file_out, const unsigned &nplot) const
 
void output_fct_paraview (std::ofstream &file_out, const unsigned &nplot, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt) const
 
void output_fct_paraview (std::ofstream &file_out, const unsigned &nplot, const double &time, FiniteElement::UnsteadyExactSolutionFctPt exact_soln_pt) const
 
void output (std::ostream &outfile)
 Output for all elements. More...
 
void output (std::ostream &outfile, const unsigned &n_plot)
 Output at f(n_plot) points in each element. More...
 
void output (FILE *file_pt)
 Output for all elements (C-style output) More...
 
void output (FILE *file_pt, const unsigned &nplot)
 Output at f(n_plot) points in each element (C-style output) More...
 
void output (const std::string &output_filename)
 Output for all elements. More...
 
void output (const std::string &output_filename, const unsigned &n_plot)
 Output at f(n_plot) points in each element. More...
 
void output_fct (std::ostream &outfile, const unsigned &n_plot, FiniteElement::SteadyExactSolutionFctPt)
 Output a given Vector function at f(n_plot) points in each element. More...
 
void output_fct (std::ostream &outfile, const unsigned &n_plot, const double &time, FiniteElement::UnsteadyExactSolutionFctPt)
 
void output_boundaries (std::ostream &outfile)
 Output the nodes on the boundaries (into separate tecplot zones) More...
 
void output_boundaries (const std::string &output_filename)
 
void assign_initial_values_impulsive ()
 Assign initial values for an impulsive start. More...
 
void shift_time_values ()
 
void calculate_predictions ()
 
void set_nodal_and_elemental_time_stepper (TimeStepper *const &time_stepper_pt, const bool &preserve_existing_data)
 
virtual void set_mesh_level_time_stepper (TimeStepper *const &time_stepper_pt, const bool &preserve_existing_data)
 
void set_consistent_pinned_values_for_continuation (ContinuationStorageScheme *const &continuation_stepper_pt)
 Set consistent values for pinned data in continuation. More...
 
bool does_pointer_correspond_to_mesh_data (double *const &parameter_pt)
 Does the double pointer correspond to any mesh data. More...
 
void set_nodal_time_stepper (TimeStepper *const &time_stepper_pt, const bool &preserve_existing_data)
 Set the timestepper associated with the nodal data in the mesh. More...
 
void set_elemental_internal_time_stepper (TimeStepper *const &time_stepper_pt, const bool &preserve_existing_data)
 
virtual void compute_norm (double &norm)
 
virtual void compute_norm (Vector< 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, double &error, double &norm)
 
virtual void compute_error (FiniteElement::SteadyExactSolutionFctPt exact_soln_pt, double &error, double &norm)
 
virtual void compute_error (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_error (std::ostream &outfile, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt, Vector< double > &error, Vector< double > &norm)
 
virtual void compute_error (FiniteElement::UnsteadyExactSolutionFctPt exact_soln_pt, const double &time, double &error, double &norm)
 Returns the norm of the error and that of the exact solution. More...
 
virtual void compute_error (FiniteElement::UnsteadyExactSolutionFctPt exact_soln_pt, const double &time, Vector< double > &error, Vector< double > &norm)
 
bool is_mesh_distributed () const
 Boolean to indicate if Mesh has been distributed. More...
 
OomphCommunicatorcommunicator_pt () const
 
void delete_all_external_storage ()
 Wipe the storage for all externally-based elements. More...
 

Protected Attributes

CollapsibleChannelDomainDomain_pt
 Pointer to domain. More...
 
unsigned Nup
 Number of element columns in upstream part. More...
 
unsigned Ncollapsible
 Number of element columns in collapsible part. More...
 
unsigned Ndown
 Number of element columns in downstream part. More...
 
unsigned Ny
 Number of element rows across channel. More...
 
GeomObjectWall_pt
 Pointer to geometric object that represents the moving wall. More...
 
- Protected Attributes inherited from oomph::Mesh
Vector< Vector< Node * > > Boundary_node_pt
 
bool Lookup_for_elements_next_boundary_is_setup
 
Vector< Vector< FiniteElement * > > Boundary_element_pt
 
Vector< Vector< int > > Face_index_at_boundary
 
Vector< Node * > Node_pt
 Vector of pointers to nodes. More...
 
Vector< GeneralisedElement * > Element_pt
 Vector of pointers to generalised elements. More...
 
std::vector< boolBoundary_coordinate_exists
 

Additional Inherited Members

- Public Types inherited from oomph::Mesh
typedef void(FiniteElement::* SteadyExactSolutionFctPt) (const Vector< double > &x, Vector< double > &soln)
 
typedef void(FiniteElement::* UnsteadyExactSolutionFctPt) (const double &time, const Vector< double > &x, Vector< double > &soln)
 
- Static Public Attributes inherited from oomph::Mesh
static Steady< 0 > Default_TimeStepper
 The Steady Timestepper. More...
 
static bool Suppress_warning_about_empty_mesh_level_time_stepper_function
 Static boolean flag to control warning about mesh level timesteppers. More...
 
- Protected Member Functions inherited from oomph::Mesh
unsigned long assign_global_eqn_numbers (Vector< double * > &Dof_pt)
 Assign (global) equation numbers to the nodes. More...
 
void describe_dofs (std::ostream &out, const std::string &current_string) const
 
void describe_local_dofs (std::ostream &out, const std::string &current_string) const
 
void assign_local_eqn_numbers (const bool &store_local_dof_pt)
 Assign local equation numbers in all elements. More...
 
void convert_to_boundary_node (Node *&node_pt, const Vector< FiniteElement * > &finite_element_pt)
 
void convert_to_boundary_node (Node *&node_pt)
 

Detailed Description

template<class ELEMENT>
class oomph::CollapsibleChannelMesh< ELEMENT >

Basic collapsible channel mesh. The mesh is derived from the SimpleRectangularQuadMesh so it's node and element numbering scheme is the same as in that mesh. Only the boundaries are numbered differently to allow the easy identification of the "collapsible" segment. Boundary coordinates are set up for all nodes located on boundary 3 (the collapsible segment). The curvilinear ("collapsible") segment is defined by a GeomObject.

Constructor & Destructor Documentation

◆ CollapsibleChannelMesh()

template<class ELEMENT >
oomph::CollapsibleChannelMesh< ELEMENT >::CollapsibleChannelMesh ( const unsigned nup,
const unsigned ncollapsible,
const unsigned ndown,
const unsigned ny,
const double lup,
const double lcollapsible,
const double ldown,
const double ly,
GeomObject wall_pt,
TimeStepper time_stepper_pt = &Mesh::Default_TimeStepper 
)

Constructor: Pass number of elements in upstream/collapsible/ downstream segment and across the channel; lengths of upstream/ collapsible/downstream segments and width of channel, pointer to GeomObject that defines the collapsible segment and pointer to TimeStepper (defaults to the default timestepper, Steady).

Constructor: Pass number of elements in upstream/collapsible/downstream segment and across the channel; lengths of upstream/collapsible/downstream segments and width of channel, pointer to GeomObject that defines the collapsible segment and pointer to TimeStepper (defaults to the default timestepper, Steady).

54  : SimpleRectangularQuadMesh<ELEMENT>(nup + ncollapsible + ndown,
55  ny,
56  lup + lcollapsible + ldown,
57  ly,
58  time_stepper_pt),
59  Nup(nup),
60  Ncollapsible(ncollapsible),
61  Ndown(ndown),
62  Ny(ny),
64  {
65  // Mesh can only be built with 2D Qelements.
66  MeshChecker::assert_geometric_element<QElementGeometricBase, ELEMENT>(2);
67 
68  // Create the CollapsibleChannelDomain with the wall represented
69  // by the geometric object
70  Domain_pt = new CollapsibleChannelDomain(
71  nup, ncollapsible, ndown, ny, lup, lcollapsible, ldown, ly, wall_pt);
72 
73  // Total number of (macro/finite)elements
74  unsigned nmacro = (nup + ncollapsible + ndown) * ny;
75 
76  // Loop over all elements and set macro element pointer
77  for (unsigned e = 0; e < nmacro; e++)
78  {
79  // Get pointer to finite element
80  FiniteElement* el_pt = this->finite_element_pt(e);
81 
82  // Set pointer to macro element so the curvlinear boundaries
83  // of the mesh/domain can get picked up during adaptive
84  // mesh refinement
85  el_pt->set_macro_elem_pt(this->Domain_pt->macro_element_pt(e));
86  }
87 
88  // Update the nodal positions corresponding to their
89  // macro element representations.
90  this->node_update();
91 
92  // Update the boundary numbering scheme and set boundary coordinate
93  //-----------------------------------------------------------------
94 
95  // (Note: The original SimpleRectangularQuadMesh had four boundaries.
96  // We need to overwrite the boundary lookup scheme for the current
97  // mesh so that the collapsible segment becomes identifiable).
98  // While we're doing this, we're also setting up a boundary
99  // coordinate for the nodes located on the collapsible segment.
100  // The boundary coordinate can be used to setup FSI.
101 
102  // How many boundaries does the mesh have now?
103  unsigned nbound = this->nboundary();
104  for (unsigned b = 0; b < nbound; b++)
105  {
106  // Remove all nodes on this boundary from the mesh's lookup scheme
107  // and also delete the reverse lookup scheme held by the nodes
108  this->remove_boundary_nodes(b);
109  }
110 
111 #ifdef PARANOID
112  // Sanity check
113  unsigned nnod = this->nnode();
114  for (unsigned j = 0; j < nnod; j++)
115  {
116  if (this->node_pt(j)->is_on_boundary())
117  {
118  std::ostringstream error_message;
119  error_message << "Node " << j << "is still on boundary " << std::endl;
120 
121  throw OomphLibError(error_message.str(),
124  }
125  }
126 #endif
127 
128  // Change the numbers of boundaries
129  this->set_nboundary(6);
130 
131  // Get the number of nodes along the element edge from first element
132  unsigned nnode_1d = this->finite_element_pt(0)->nnode_1d();
133 
134  // Vector of Lagrangian coordinates used as boundary coordinate
135  Vector<double> zeta(1);
136 
137  // Index of first element underneath the collapsible bit
138  unsigned first_collapsible = (ny - 1) * (nup + ncollapsible + ndown) + nup;
139 
140  // Zeta increment over elements (used for assignment of
141  // boundary coordinate)
142  double dzeta = lcollapsible / double(ncollapsible);
143 
144  // Manually loop over the elements near the boundaries and
145  // assign nodes to boundaries. Also set up boundary coordinate
146  unsigned nelem = this->nelement();
147  for (unsigned e = 0; e < nelem; e++)
148  {
149  // Bottom row of elements
150  if (e < nup + ncollapsible + ndown)
151  {
152  for (unsigned i = 0; i < nnode_1d; i++)
153  {
154  this->add_boundary_node(0, this->finite_element_pt(e)->node_pt(i));
155  }
156  }
157  // Upstream upper rigid bit
158  if ((e > (ny - 1) * (nup + ncollapsible + ndown) - 1) &&
159  (e < (ny - 1) * (nup + ncollapsible + ndown) + nup))
160  {
161  for (unsigned i = 0; i < nnode_1d; i++)
162  {
163  this->add_boundary_node(
164  4,
165  this->finite_element_pt(e)->node_pt((nnode_1d - 1) * nnode_1d + i));
166  }
167  }
168  // Collapsible bit
169  if ((e > (ny - 1) * (nup + ncollapsible + ndown) + nup - 1) &&
170  (e < (ny - 1) * (nup + ncollapsible + ndown) + nup + ncollapsible))
171  {
172  for (unsigned i = 0; i < nnode_1d; i++)
173  {
174  this->add_boundary_node(
175  3,
176  this->finite_element_pt(e)->node_pt((nnode_1d - 1) * nnode_1d + i));
177 
178  // What column of elements are we in?
179  unsigned ix = e - first_collapsible;
180 
181  // Zeta coordinate
182  zeta[0] =
183  double(ix) * dzeta + double(i) * dzeta / double(nnode_1d - 1);
184 
185  // Set boundary coordinate
186  this->finite_element_pt(e)
187  ->node_pt((nnode_1d - 1) * nnode_1d + i)
189  }
190  }
191  // Downstream upper rigid bit
192  if ((e >
193  (ny - 1) * (nup + ncollapsible + ndown) + nup + ncollapsible - 1) &&
194  (e < ny * (nup + ncollapsible + ndown)))
195  {
196  for (unsigned i = 0; i < nnode_1d; i++)
197  {
198  this->add_boundary_node(
199  2,
200  this->finite_element_pt(e)->node_pt((nnode_1d - 1) * nnode_1d + i));
201  }
202  }
203  // Left end
204  if (e % (nup + ncollapsible + ndown) == 0)
205  {
206  for (unsigned i = 0; i < nnode_1d; i++)
207  {
208  this->add_boundary_node(
209  5, this->finite_element_pt(e)->node_pt(i * nnode_1d));
210  }
211  }
212  // Right end
213  if (e % (nup + ncollapsible + ndown) == (nup + ncollapsible + ndown) - 1)
214  {
215  for (unsigned i = 0; i < nnode_1d; i++)
216  {
217  this->add_boundary_node(
218  1, this->finite_element_pt(e)->node_pt((i + 1) * nnode_1d - 1));
219  }
220  }
221  }
222 
223  // Re-setup lookup scheme that establishes which elements are located
224  // on the mesh boundaries (doesn't need to be wiped)
226 
227  // We have only bothered to parametrise boundary 3
228  this->Boundary_coordinate_exists[3] = true;
229  }
int i
Definition: BiCGSTAB_step_by_step.cpp:9
Array< double, 1, 3 > e(1./3., 0.5, 2.)
Scalar * b
Definition: benchVecAdd.cpp:17
GeomObject *& wall_pt()
Access function to GeomObject representing wall.
Definition: collapsible_channel_mesh.template.h:91
GeomObject * Wall_pt
Pointer to geometric object that represents the moving wall.
Definition: collapsible_channel_mesh.template.h:162
unsigned Ncollapsible
Number of element columns in collapsible part.
Definition: collapsible_channel_mesh.template.h:153
unsigned Nup
Number of element columns in upstream part.
Definition: collapsible_channel_mesh.template.h:150
unsigned Ndown
Number of element columns in downstream part.
Definition: collapsible_channel_mesh.template.h:156
CollapsibleChannelDomain * Domain_pt
Pointer to domain.
Definition: collapsible_channel_mesh.template.h:147
unsigned Ny
Number of element rows across channel.
Definition: collapsible_channel_mesh.template.h:159
MacroElement * macro_element_pt(const unsigned &i)
Access to i-th macro element.
Definition: domain.h:116
Node *& node_pt(const unsigned &n)
Return a pointer to the local node n.
Definition: elements.h:2175
virtual unsigned nnode_1d() const
Definition: elements.h:2218
void add_boundary_node(const unsigned &b, Node *const &node_pt)
Add a (pointer to) a node to the b-th boundary.
Definition: mesh.cc:243
std::vector< bool > Boundary_coordinate_exists
Definition: mesh.h:190
FiniteElement * finite_element_pt(const unsigned &e) const
Definition: mesh.h:473
void set_nboundary(const unsigned &nbound)
Set the number of boundaries in the mesh.
Definition: mesh.h:505
unsigned nboundary() const
Return number of boundaries.
Definition: mesh.h:827
void remove_boundary_nodes()
Clear all pointers to boundary nodes.
Definition: mesh.cc:204
unsigned long nnode() const
Return number of nodes in the mesh.
Definition: mesh.h:596
virtual void node_update(const bool &update_all_solid_nodes=false)
Definition: mesh.cc:287
Node *& node_pt(const unsigned long &n)
Return pointer to global node n.
Definition: mesh.h:436
unsigned long nelement() const
Return number of elements in the mesh.
Definition: mesh.h:590
virtual void set_coordinates_on_boundary(const unsigned &b, const unsigned &k, const Vector< double > &boundary_zeta)
Definition: nodes.cc:2394
void setup_boundary_element_info()
Definition: quad_mesh.h:73
const unsigned & ny() const
Access function for number of elements in y directions.
Definition: simple_rectangular_quadmesh.template.h:77
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
const double ly
Definition: ConstraintElementsUnitTest.cpp:34
#define OOMPH_EXCEPTION_LOCATION
Definition: oomph_definitions.h:61
#define OOMPH_CURRENT_FUNCTION
Definition: oomph_definitions.h:86
std::ptrdiff_t j
Definition: tut_arithmetic_redux_minmax.cpp:2

References oomph::Mesh::add_boundary_node(), b, oomph::Mesh::Boundary_coordinate_exists, oomph::CollapsibleChannelMesh< ELEMENT >::Domain_pt, e(), oomph::Mesh::finite_element_pt(), i, j, Mesh_Parameters::ly, oomph::Domain::macro_element_pt(), oomph::Mesh::nboundary(), oomph::Mesh::nelement(), oomph::Mesh::nnode(), oomph::FiniteElement::nnode_1d(), oomph::FiniteElement::node_pt(), oomph::Mesh::node_pt(), oomph::Mesh::node_update(), oomph::SimpleRectangularQuadMesh< ELEMENT >::ny(), OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION, oomph::Mesh::remove_boundary_nodes(), oomph::Node::set_coordinates_on_boundary(), oomph::FiniteElement::set_macro_elem_pt(), oomph::Mesh::set_nboundary(), oomph::QuadMeshBase::setup_boundary_element_info(), oomph::CollapsibleChannelMesh< ELEMENT >::wall_pt(), and Eigen::zeta().

◆ ~CollapsibleChannelMesh()

template<class ELEMENT >
oomph::CollapsibleChannelMesh< ELEMENT >::~CollapsibleChannelMesh ( )
inline

destructor

86  {
87  delete Domain_pt;
88  }

References oomph::CollapsibleChannelMesh< ELEMENT >::Domain_pt.

Member Function Documentation

◆ axial_spacing_fct_pt() [1/2]

template<class ELEMENT >
virtual CollapsibleChannelDomain::AxialSpacingFctPt& oomph::CollapsibleChannelMesh< ELEMENT >::axial_spacing_fct_pt ( )
inlinevirtual

Function pointer for function that redistributes the elements in the axial direction. Virtual so we can break it in derived classes (e.g. the Algebraic versions of this mesh where it doesn't make any sense to provide the bl_squash_fct after the mesh has been built).

Reimplemented in oomph::AlgebraicCollapsibleChannelMesh< ELEMENT >.

131  {
133  }
AxialSpacingFctPt & axial_spacing_fct_pt()
Definition: collapsible_channel_domain.h:186

References oomph::CollapsibleChannelDomain::axial_spacing_fct_pt(), and oomph::CollapsibleChannelMesh< ELEMENT >::Domain_pt.

◆ axial_spacing_fct_pt() [2/2]

template<class ELEMENT >
virtual CollapsibleChannelDomain::AxialSpacingFctPt& oomph::CollapsibleChannelMesh< ELEMENT >::axial_spacing_fct_pt ( ) const
inlinevirtual

Function pointer for function that redistributes the elements in the axial direction. Const version

140  {
142  }

References oomph::CollapsibleChannelDomain::axial_spacing_fct_pt(), and oomph::CollapsibleChannelMesh< ELEMENT >::Domain_pt.

◆ bl_squash_fct_pt() [1/2]

template<class ELEMENT >
virtual CollapsibleChannelDomain::BLSquashFctPt& oomph::CollapsibleChannelMesh< ELEMENT >::bl_squash_fct_pt ( )
inlinevirtual

Function pointer for function that squashes the mesh near the walls. Default trivial mapping (the identity) leaves vertical nodal positions unchanged. Mapping is used in underlying CollapsibleChannelDomain. Virtual so we can break it in derived classes (e.g. the Algebraic versions of this mesh where it doesn't make any sense to provide the bl_squash_fct after the mesh has been built).

Reimplemented in oomph::AlgebraicCollapsibleChannelMesh< ELEMENT >, and oomph::MyAlgebraicCollapsibleChannelMesh< ELEMENT >.

110  {
111  return Domain_pt->bl_squash_fct_pt();
112  }
BLSquashFctPt & bl_squash_fct_pt()
Definition: collapsible_channel_domain.h:165

References oomph::CollapsibleChannelDomain::bl_squash_fct_pt(), and oomph::CollapsibleChannelMesh< ELEMENT >::Domain_pt.

◆ bl_squash_fct_pt() [2/2]

template<class ELEMENT >
CollapsibleChannelDomain::BLSquashFctPt oomph::CollapsibleChannelMesh< ELEMENT >::bl_squash_fct_pt ( ) const
inline

Function pointer for function that squashes the mesh near the walls. Default trivial mapping (the identity) leaves vertical nodal positions unchanged. Mapping is used in underlying CollapsibleChannelDomain. Const version.

120  {
121  return Domain_pt->bl_squash_fct_pt();
122  }

References oomph::CollapsibleChannelDomain::bl_squash_fct_pt(), and oomph::CollapsibleChannelMesh< ELEMENT >::Domain_pt.

◆ domain_pt()

◆ wall_pt()

template<class ELEMENT >
GeomObject*& oomph::CollapsibleChannelMesh< ELEMENT >::wall_pt ( )
inline

Member Data Documentation

◆ Domain_pt

◆ Ncollapsible

template<class ELEMENT >
unsigned oomph::CollapsibleChannelMesh< ELEMENT >::Ncollapsible
protected

Number of element columns in collapsible part.

◆ Ndown

template<class ELEMENT >
unsigned oomph::CollapsibleChannelMesh< ELEMENT >::Ndown
protected

Number of element columns in downstream part.

◆ Nup

template<class ELEMENT >
unsigned oomph::CollapsibleChannelMesh< ELEMENT >::Nup
protected

Number of element columns in upstream part.

◆ Ny

template<class ELEMENT >
unsigned oomph::CollapsibleChannelMesh< ELEMENT >::Ny
protected

Number of element rows across channel.

◆ Wall_pt

template<class ELEMENT >
GeomObject* oomph::CollapsibleChannelMesh< ELEMENT >::Wall_pt
protected

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