SimpleSpineMesh< ELEMENT > Class Template Reference
+ Inheritance diagram for SimpleSpineMesh< ELEMENT >:

Public Types

typedef double(* HeightFctPt) (const double &x)
 
- 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)
 

Public Member Functions

 SimpleSpineMesh (const unsigned &nx, const unsigned &ny, const double &lx, const double &h, GeomObject *substrate_pt, TimeStepper *time_stepper_pt=&Mesh::Default_TimeStepper)
 
HeightFctPtheight_fct_pt ()
 Access function: Pointer to height function. More...
 
double height_fct (const double &x)
 
void update_spine_heights ()
 
virtual void spine_node_update (SpineNode *spine_node_pt)
 
- 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 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...
 
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)
 
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...
 
- Public Member Functions inherited from oomph::SpineMesh
virtual ~SpineMesh ()
 Destructor to clean up the memory allocated to the spines. More...
 
Spine *& spine_pt (const unsigned long &i)
 Return the i-th spine in the mesh. More...
 
const Spinespine_pt (const unsigned long &i) const
 Return the i-th spine in the mesh (const version) More...
 
unsigned long nspine () const
 Return the number of spines in the mesh. More...
 
void add_spine_pt (Spine *const &spine_pt)
 Add a spine to the mesh. More...
 
SpineNodenode_pt (const unsigned long &n)
 Return a pointer to the n-th global SpineNode. More...
 
SpineNodeelement_node_pt (const unsigned long &e, const unsigned &n)
 
unsigned long assign_global_spine_eqn_numbers (Vector< double * > &Dof_pt)
 Assign spines to Spine_pt vector of element. More...
 
void describe_spine_dofs (std::ostream &out, const std::string &current_string) const
 
void set_mesh_level_time_stepper (TimeStepper *const &time_stepper_pt, const bool &preserve_existing_data)
 
void set_spine_time_stepper (TimeStepper *const &time_stepper_pt, const bool &preserve_existing_data)
 Assign time stepper to spines data. More...
 
void set_consistent_pinned_spine_values_for_continuation (ContinuationStorageScheme *const &continuation_stepper_pt)
 
bool does_pointer_correspond_to_spine_data (double *const &parameter_pt)
 
void node_update (const bool &update_all_solid_nodes=false)
 
void dump (std::ofstream &dump_file) const
 Overload the dump function so that the spine data is dumped. More...
 
void read (std::ifstream &restart_file)
 Overload the read function so that the spine data is also read. More...
 

Private Attributes

double Default_height
 Default height. More...
 
HeightFctPt Height_fct_pt
 Pointer to height function. More...
 
GeomObjectSubstrate_pt
 Pointer to GeomObject that specifies the "substrate" (the lower wall) More...
 

Additional Inherited Members

- 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)
 
- 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
 
- Protected Attributes inherited from oomph::SpineMesh
Vector< Spine * > Spine_pt
 A Spine mesh contains a Vector of pointers to spines. More...
 

Detailed Description

template<class ELEMENT>
class SimpleSpineMesh< ELEMENT >

//////////////////////////////////////////////////////////////////// //////////////////////////////////////////////////////////////////// Simple spine mesh class derived from standard 2D mesh. Vertical lines of nodes are located on spines.

Member Typedef Documentation

◆ HeightFctPt

template<class ELEMENT >
typedef double(* SimpleSpineMesh< ELEMENT >::HeightFctPt) (const double &x)

Function pointer to function, h(x), that may be used to specify the "height" of the domain, by assigning the function values to the spine heights.

Constructor & Destructor Documentation

◆ SimpleSpineMesh()

template<class ELEMENT >
SimpleSpineMesh< ELEMENT >::SimpleSpineMesh ( const unsigned nx,
const unsigned ny,
const double lx,
const double h,
GeomObject substrate_pt,
TimeStepper time_stepper_pt = &Mesh::Default_TimeStepper 
)

Constructor: Pass number of elements in x-direction, number of elements in y-direction, length in x direction, initial height of mesh, and pointer to timestepper (defaults to Steady timestepper)

Constructor for spine 2D mesh: Pass number of elements in x-direction, number of elements in y-direction, axial length and height of layer, pointer to geometric object that specifies the substrate (the lower wall) and pointer to timestepper (defaults to Static timestepper).

The mesh contains a layer of spine-ified fluid elements (of type ELEMENT; e.g SpineElement<QCrouzeixRaviartElement<2>)

211  :
212  SimpleRectangularQuadMesh<ELEMENT>(nx,ny,lx,h,time_stepper_pt),
213  Default_height(h), Height_fct_pt(0), Substrate_pt(substrate_pt)
214 
215 {
216 
217  // We've already called the constructor for the SimpleRectangularQuadMesh
218  // which sets up nodal positons/topology etc. Now attach the spine
219  // information
220 
221 
222  // Allocate memory for the spines and fractions along spines
223 
224  // Read out number of linear points in the element
225  unsigned np = dynamic_cast<ELEMENT*>(finite_element_pt(0))->nnode_1d();
226 
227  // Number of spines
228  unsigned nspine = (np-1)*nx+1;
229  Spine_pt.reserve(nspine);
230 
231  // Set up x-increments between spine origins
232  double x_lo=0.0;
233  double dx=lx/double(nx);
234 
235  //FIRST SPINE
236  //-----------
237 
238  //Element 0
239  //Node 0
240  //Assign the new spine with specified height
241  Spine* new_spine_pt=new Spine(h);
242  Spine_pt.push_back(new_spine_pt);
243 
244  // Get pointer to node
245  SpineNode* nod_pt=element_node_pt(0,0);
246 
247  //Pass the pointer to the spine to the node
248  nod_pt->spine_pt() = new_spine_pt;
249 
250  //Set the node's fraction along the spine
251  nod_pt->fraction() = 0.0;
252 
253  // Set the pointer to the mesh that implements the update fct for this node
254  nod_pt->spine_mesh_pt() = this;
255 
256  // Set update fct id (not really needed here...)
257  nod_pt->node_update_fct_id()=0;
258 
259 
260  // Create vector that stores additional geometric parameters that
261  // are required during the execution of the remesh function:
262  // For our remesh function, we need the x-coordinate of the
263  // spine origin.
264 
265  // Vector with enough storage for one geometic parameter
266  Vector<double> parameters(1);
267 
268  // Store the x-coordinate in it
269  parameters[0]=x_lo;
270 
271  // Pass the parameter to the spine
272  new_spine_pt->set_geom_parameter(parameters);
273 
274 
275 
276  // The remesh function involving this spine only involves
277  // a single geometric object: The substrate
278  Vector<GeomObject*> geom_object_pt(1);
279  geom_object_pt[0] = Substrate_pt;
280 
281  // Pass geom object(s) to spine
282  new_spine_pt->set_geom_object_pt(geom_object_pt);
283 
284 
285  //Loop vertically along the first spine
286 
287  //Loop over the elements
288  for(unsigned i=0;i<ny;i++)
289  {
290  //Loop over the vertical nodes, apart from the first
291  for(unsigned l1=1;l1<np;l1++)
292  {
293  // Get pointer to node
294  SpineNode* nod_pt=element_node_pt(i*nx,l1*np);
295 
296  //Pass the pointer to the spine to the node
297  nod_pt->spine_pt() = new_spine_pt;
298 
299  //Set the fraction
300  nod_pt->fraction()=(double(i)+double(l1)/double(np-1))/double(ny);
301 
302  // Pointer to the mesh that implements the update fct
303  nod_pt->spine_mesh_pt() = this;
304 
305  // Set update fct id (not really needed here)
306  nod_pt->node_update_fct_id()=0;
307 
308  }
309  } // end loop over elements
310 
311 
312  //LOOP OVER OTHER SPINES
313  //----------------------
314 
315  //Now loop over the elements horizontally
316  for(unsigned j=0;j<nx;j++)
317  {
318  //Loop over the nodes in the elements horizontally, ignoring
319  //the first column
320  unsigned npmax=np;
321  for(unsigned l2=1;l2<npmax;l2++)
322  {
323  //Assign the new spine with specified height
324  new_spine_pt=new Spine(h);
325 
326  // Add to collection of spines
327  Spine_pt.push_back(new_spine_pt);
328 
329  // Get the node
330  SpineNode* nod_pt=element_node_pt(j,l2);
331 
332  //Pass the pointer to the spine to the node
333  nod_pt->spine_pt() = new_spine_pt;
334 
335  //Set the fraction
336  nod_pt->fraction() = 0.0;
337 
338  // Pointer to the mesh that implements the update fct
339  nod_pt->spine_mesh_pt() = this;
340 
341  // Set update fct id (not really needed here)
342  nod_pt->node_update_fct_id()=0;
343 
344  // Create vector that stores additional geometric parameters that
345  // are required during the execution of the remesh function:
346  // For our remesh function, we need the x-coordinate of the
347  // spine origin.
348 
349  // Vector with enough storage for one geometic parameter
350  Vector<double> parameters(1);
351 
352  // Store the x-coordinate in it
353  parameters[0]=x_lo+double(j)*dx + double(l2)*dx/double(np-1);
354 
355  // Pass the parameter to the spine
356  new_spine_pt->set_geom_parameter(parameters);
357 
358 
359  // The remesh function involving this spine only involves
360  // a single geometric object: The substrate
361  Vector<GeomObject*> geom_object_pt(1);
362  geom_object_pt[0] = Substrate_pt;
363 
364  // Pass geom object(s) to spine
365  new_spine_pt->set_geom_object_pt(geom_object_pt);
366 
367  //Loop vertically along the spine
368  //Loop over the elements
369  for(unsigned i=0;i<ny;i++)
370  {
371  //Loop over the vertical nodes, apart from the first
372  for(unsigned l1=1;l1<np;l1++)
373  {
374  // Get the node
375  SpineNode* nod_pt=element_node_pt(i*nx+j,l1*np+l2);
376 
377  //Set the pointer to the spine
378  nod_pt->spine_pt() = new_spine_pt;
379 
380  //Set the fraction
381  nod_pt->fraction()=(double(i)+double(l1)/double(np-1))/double(ny);
382 
383  // Pointer to the mesh that implements the update fct
384  nod_pt->spine_mesh_pt() = this;
385 
386  // Set update fct id (not really needed here)
387  nod_pt->node_update_fct_id()=0;
388  }
389  }
390  }
391  }
392 
393 } // end of constructor
int i
Definition: BiCGSTAB_step_by_step.cpp:9
double Default_height
Default height.
Definition: simple_spine_channel.cc:184
GeomObject * Substrate_pt
Pointer to GeomObject that specifies the "substrate" (the lower wall)
Definition: simple_spine_channel.cc:190
HeightFctPt Height_fct_pt
Pointer to height function.
Definition: simple_spine_channel.cc:187
FiniteElement * finite_element_pt(const unsigned &e) const
Definition: mesh.h:473
Definition: simple_rectangular_quadmesh.template.h:58
const unsigned & ny() const
Access function for number of elements in y directions.
Definition: simple_rectangular_quadmesh.template.h:77
const unsigned & nx() const
Access function for number of elements in x directions.
Definition: simple_rectangular_quadmesh.template.h:71
Vector< Spine * > Spine_pt
A Spine mesh contains a Vector of pointers to spines.
Definition: spines.h:616
unsigned long nspine() const
Return the number of spines in the mesh.
Definition: spines.h:635
SpineNode * element_node_pt(const unsigned long &e, const unsigned &n)
Definition: spines.h:669
Definition: spines.h:328
SpineMesh *& spine_mesh_pt()
Definition: spines.h:391
Spine *& spine_pt()
Access function to spine.
Definition: spines.h:372
unsigned & node_update_fct_id()
Access function to ID of node update function (within specific mesh)
Definition: spines.h:384
double & fraction()
Set reference to fraction along spine.
Definition: spines.h:378
Definition: spines.h:64
void set_geom_object_pt(const Vector< GeomObject * > &geom_object_pt)
Definition: spines.h:225
void set_geom_parameter(const Vector< double > &geom_parameter)
Definition: spines.h:273
const double lx
Definition: ConstraintElementsUnitTest.cpp:33
std::ptrdiff_t j
Definition: tut_arithmetic_redux_minmax.cpp:2

References oomph::SpineMesh::element_node_pt(), oomph::Mesh::finite_element_pt(), oomph::SpineNode::fraction(), i, j, Mesh_Parameters::lx, oomph::SpineNode::node_update_fct_id(), oomph::SpineMesh::nspine(), oomph::SimpleRectangularQuadMesh< ELEMENT >::nx(), oomph::SimpleRectangularQuadMesh< ELEMENT >::ny(), oomph::Spine::set_geom_object_pt(), oomph::Spine::set_geom_parameter(), oomph::SpineNode::spine_mesh_pt(), oomph::SpineNode::spine_pt(), oomph::SpineMesh::Spine_pt, and SimpleSpineMesh< ELEMENT >::Substrate_pt.

Member Function Documentation

◆ height_fct()

template<class ELEMENT >
double SimpleSpineMesh< ELEMENT >::height_fct ( const double x)
inline

Height function – this is called by update_spine_heights() when spine heights are assigned.

122  {
123  // Resolve function pointer if non-NULL
124  if (Height_fct_pt==0)
125  {
126  return Default_height;
127  }
128  else
129  {
130  return Height_fct_pt(x);
131  }
132  }
list x
Definition: plotDoE.py:28

References plotDoE::x.

◆ height_fct_pt()

template<class ELEMENT >
HeightFctPt& SimpleSpineMesh< ELEMENT >::height_fct_pt ( )
inline

Access function: Pointer to height function.

114  {
115  return Height_fct_pt;
116  }

◆ spine_node_update()

template<class ELEMENT >
virtual void SimpleSpineMesh< ELEMENT >::spine_node_update ( SpineNode spine_node_pt)
inlinevirtual

General node update function implements pure virtual function defined in SpineMesh base class and performs specific node update actions: Nodes are located along vertical "spines" that emanate from the "substrate" (the lower wall) specified by a GeomObject.

Implements oomph::SpineMesh.

158  {
159 
160  //Get fraction along the spine
161  double w = spine_node_pt->fraction();
162 
163  // Get height of spine
164  double h = spine_node_pt->spine_pt()->height();
165 
166  // Get x-position of spine origin (as vector)
167  Vector<double> x_spine_origin(1);
168  x_spine_origin[0]=spine_node_pt->spine_pt()->geom_parameter(0);
169 
170  // Get position vector to substrate (lower wall) for this spine
171  Vector<double> spine_base(2);
172  spine_node_pt->spine_pt()->geom_object_pt(0)->
173  position(x_spine_origin,spine_base);
174 
175  //Update the nodal position
176  spine_node_pt->x(0) = spine_base[0];
177  spine_node_pt->x(1) = spine_base[1]+w*h;
178  }
RowVector3d w
Definition: Matrix_resize_int.cpp:3
double & x(const unsigned &i)
Return the i-th nodal coordinate.
Definition: nodes.h:1060
GeomObject *& geom_object_pt(const unsigned &i)
Definition: spines.h:244
double & height()
Access function to spine height.
Definition: spines.h:150
double & geom_parameter(const unsigned &i)
Definition: spines.h:287

References oomph::SpineNode::fraction(), oomph::Spine::geom_object_pt(), oomph::Spine::geom_parameter(), oomph::Spine::height(), oomph::SpineNode::spine_pt(), w, and oomph::Node::x().

◆ update_spine_heights()

template<class ELEMENT >
void SimpleSpineMesh< ELEMENT >::update_spine_heights ( )
inline

Update the spine heights according to the function specified by height_fct_pt().

138  {
139  unsigned n_spine=nspine();
140  for (unsigned i=0;i<n_spine;i++)
141  {
142  // Zero-th geometric parameter of the spine is the x-coordinate
143  // of its basis
144  double x=spine_pt(i)->geom_parameter(0);
145 
146  // Set spine height
148  }
149  }
double height_fct(const double &x)
Definition: simple_spine_channel.cc:121
Spine *& spine_pt(const unsigned long &i)
Return the i-th spine in the mesh.
Definition: spines.h:623

References i, and plotDoE::x.

Member Data Documentation

◆ Default_height

template<class ELEMENT >
double SimpleSpineMesh< ELEMENT >::Default_height
private

Default height.

◆ Height_fct_pt

template<class ELEMENT >
HeightFctPt SimpleSpineMesh< ELEMENT >::Height_fct_pt
private

Pointer to height function.

◆ Substrate_pt

template<class ELEMENT >
GeomObject* SimpleSpineMesh< ELEMENT >::Substrate_pt
private

Pointer to GeomObject that specifies the "substrate" (the lower wall)

Referenced by SimpleSpineMesh< ELEMENT >::SimpleSpineMesh().


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