31 #ifndef OOMPH_FULL_CIRCLE_MESH_TEMPLATE_CC 32 #define OOMPH_FULL_CIRCLE_MESH_TEMPLATE_CC 44 template <
class ELEMENT>
54 if(radius_box.size() != 4 || theta_positions.size() != 4)
56 std::string err =
"This mesh is hard coded to only work for the case when there are 5 elements: the central square and 4 surrounding elements, but you gave vectors inconsistent with this.";
58 OOMPH_EXCEPTION_LOCATION);
63 MeshChecker::assert_geometric_element<QElementGeometricBase,ELEMENT>(2);
75 const unsigned nelem=5;
79 ELEMENT* dummy_el_pt=
new ELEMENT;
82 unsigned n_p = dummy_el_pt->nnode_1d();
88 unsigned nnodes_total=
89 (n_p*n_p + 4*(n_p-1)*(n_p -1));
99 for (
unsigned ielem=0;ielem<nelem;ielem++)
105 for (
unsigned i1=0;i1<n_p;i1++)
108 for (
unsigned i0=0;i0<n_p;i0++)
111 unsigned jnod_local=i0+i1*n_p;
115 construct_node(jnod_local,time_stepper_pt);
118 s[0]=-1.0+2.0*double(i0)/double(n_p-1);
119 s[1]=-1.0+2.0*double(i1)/double(n_p-1);
120 Domain_pt->macro_element_pt(ielem)->macro_map(s,r);
122 node_pt->
x(0) = r[0];
123 node_pt->
x(1) = r[1];
129 unsigned node_count=0;
132 double node_kill_tol=1.0e-12;
141 for (
unsigned ielem=0;ielem<nelem;ielem++)
151 for (
unsigned i1=0;i1<n_p;i1++)
155 for (
unsigned i0=0;i0<n_p;i0++)
159 unsigned jnod_local=i0+i1*n_p;
176 for (
unsigned i1=0;i1<n_p;i1++)
180 for (
unsigned i0=0;i0<n_p;i0++)
184 unsigned jnod_local=i0+i1*n_p;
193 unsigned ielem_neigh=ielem-1;
196 unsigned i0_neigh=i0;
198 unsigned jnod_local_neigh = i0_neigh + i1_neigh*n_p;
201 for (
unsigned i=0;
i<2;
i++)
203 double error=std::fabs(
207 if (error>node_kill_tol)
210 <<
i <<
" " << error << std::endl;
239 zeta[0]= theta_positions[0] +
240 double(i1)/double(n_p-1)*0.5*
241 (theta_positions[1]-theta_positions[0]);
243 Node_pt[node_count]->set_coordinates_on_boundary(0,zeta);
261 for (
unsigned i1=0;i1<n_p;i1++)
265 for (
unsigned i0=0;i0<n_p;i0++)
269 unsigned jnod_local=i0+i1*n_p;
279 unsigned ielem_neigh=ielem-1;
282 unsigned i0_neigh=n_p-1;
283 unsigned i1_neigh= n_p-1 - i0;
285 unsigned jnod_local_neigh =i0_neigh+i1_neigh*n_p;
288 for (
unsigned i=0;
i<2;
i++)
290 double error=std::fabs(
294 if (error>node_kill_tol)
297 <<
i <<
" " << error << std::endl;
314 if ((i0==0)&&(i1!=0))
318 unsigned ielem_neigh=ielem-2;
321 unsigned i0_neigh=n_p-1;
322 unsigned i1_neigh=i1;
323 unsigned jnod_local_neigh = i0_neigh+i1_neigh*n_p;
326 for (
unsigned i=0;
i<2;
i++)
328 double error=std::fabs(
332 if (error>node_kill_tol)
335 <<
i <<
" " << error << std::endl;
365 zeta[0]=theta_positions[1] +
366 double(i1)/double(n_p-1)*0.5*
367 (theta_positions[2]-theta_positions[1]);
369 Node_pt[node_count]->set_coordinates_on_boundary(0,zeta);
388 for (
unsigned i1=0;i1<n_p;i1++)
392 for (
unsigned i0=0;i0<n_p;i0++)
396 unsigned jnod_local=i0+i1*n_p;
407 unsigned ielem_neigh=ielem-1;
410 unsigned i0_neigh=i1;
411 unsigned i1_neigh= n_p-1;
412 unsigned jnod_local_neigh=i0_neigh+i1_neigh*n_p;
415 for (
unsigned i=0;
i<2;
i++)
417 double error=std::fabs(
421 if (error>node_kill_tol)
424 <<
i <<
" " << error << std::endl;
440 if ((i1==0)&&(i0!=n_p-1))
444 unsigned ielem_neigh=ielem-3;
447 unsigned i0_neigh=i0;
448 unsigned i1_neigh=n_p-1;
449 unsigned jnod_local_neigh= i0_neigh+i1_neigh*n_p;
452 for (
unsigned i=0;
i<2;
i++)
454 double error=std::fabs(
458 if (error>node_kill_tol)
461 <<
i <<
" " << error << std::endl;
491 zeta[0]=theta_positions[3] +
492 double(i0)/double(n_p-1)*0.5*
493 (theta_positions[2]-theta_positions[3]);
495 Node_pt[node_count]->set_coordinates_on_boundary(0,zeta);
513 for (
unsigned i1=0;i1<n_p;i1++)
517 for (
unsigned i0=0;i0<n_p;i0++)
521 unsigned jnod_local=i0+i1*n_p;
531 unsigned ielem_neigh=ielem-1;
535 unsigned i1_neigh= n_p - 1 -i0;
536 unsigned jnod_local_neigh=i0_neigh+i1_neigh*n_p;
539 for (
unsigned i=0;
i<2;
i++)
541 double error=std::fabs(
545 if (error>node_kill_tol)
548 <<
i <<
" " << error << std::endl;
564 if ((i0==n_p-1)&&(i1!=n_p-1))
568 unsigned ielem_neigh=ielem-4;
572 unsigned i1_neigh=i1;
573 unsigned jnod_local_neigh = i0_neigh+i1_neigh*n_p;
576 for (
unsigned i=0;
i<2;
i++)
578 double error=std::fabs(
582 if (error>node_kill_tol)
585 <<
i <<
" " << error << std::endl;
601 if ((i1==0)&&(i0!=n_p-1))
605 unsigned ielem_neigh=ielem-3;
609 unsigned i1_neigh=i0;
610 unsigned jnod_local_neigh = i0_neigh+i1_neigh*n_p;
613 for (
unsigned i=0;
i<2;
i++)
615 double error=std::fabs(
619 if (error>node_kill_tol)
622 <<
i <<
" " << error << std::endl;
653 zeta[0]=theta_positions[0] + 2.0*pi +
654 double(i1)/double(n_p-1)*0.5*
655 (theta_positions[3]-theta_positions[0] + 2.0*pi);
657 Node_pt[node_count]->set_coordinates_on_boundary(0,zeta);
674 std::ostringstream error_stream;
675 error_stream <<
"Error in killing nodes\n" 676 <<
"The most probable cause is that the domain is not\n" 677 <<
"compatible with the mesh.\n" 678 <<
"For the FullCircleMesh, the domain must be\n" 679 <<
"topologically consistent with a quarter tube with a\n" 680 <<
"non-curved centreline.\n";
682 OOMPH_CURRENT_FUNCTION,
683 OOMPH_EXCEPTION_LOCATION);
Topologically circular domain, e.g. a tube cross section. The entire domain must be defined by a Geom...
Vector< Node * > Node_pt
Vector of pointers to nodes.
void convert_to_boundary_node(Node *&node_pt, const Vector< FiniteElement *> &finite_element_pt)
A function that upgrades an ordinary node to a boundary node We shouldn't ever really use this...
void add_boundary_node(const unsigned &b, Node *const &node_pt)
Add a (pointer to) a node to the b-th boundary.
std::vector< bool > Boundary_coordinate_exists
Vector of boolean data that indicates whether the boundary coordinates have been set for the boundary...
const double Pi
50 digits from maple
Nodes are derived from Data, but, in addition, have a definite (Eulerian) position in a space of a gi...
FullCircleDomain * Domain_pt
Pointer to domain.
double & x(const unsigned &i)
Return the i-th nodal coordinate.
void set_nboundary(const unsigned &nbound)
Set the number of boundaries in the mesh.
void setup_boundary_element_info()
Vector< GeneralisedElement * > Element_pt
Vector of pointers to generalised elements.
Node *& node_pt(const unsigned &n)
Return a pointer to the local node n.
Node *& node_pt(const unsigned long &n)
Return pointer to global node n.
FiniteElement * finite_element_pt(const unsigned &e) const
Upcast (downcast?) to FiniteElement (needed to access FiniteElement member functions).
std::string string(const unsigned &i)
Return the i-th string or "" if the relevant string hasn't been defined.
Base class for time-stepping schemes. Timestepper provides an approximation of the temporal derivativ...
FullCircleMesh(GeomObject *wall_pt, const Vector< double > &theta_positions, const Vector< double > &radius_box, TimeStepper *time_stepper_pt=&Mesh::Default_TimeStepper)
Constructor: Pass pointer to geometric object that specifies the area; values of theta at which divid...