On Mon, 11 Aug 2008, yunfei zhu wrote:
> Hi, could you guys tell me which function can be used to mesh a
> thick-walled cylinder? In the file of mesh_generation.h, there are
> functions: void build_cube(), void build_line (), void
> build_square(), void build_sphere(). Do you have some function like
> build_cylinder()? or some other method to mesh a cylinder? Many
> thanks.
I've got some code I wrote for another local libMesh user, but I
never tested it enough to want to put it in the library. I'll copy
the relevant parts here, but buyer beware. If you need child element
nodes to be snapped to the correct geometry after refinement, you'll
have to do that yourself. If you need documentation, please write
some and email us a copy. If this code breaks, you get to keep both
pieces.
---
else if (param.domaintype == "square" || param.domaintype == "cylinder")
{
if (param.dimension == 1)
MeshTools::Generation::build_line
(mesh, param.coarsegridx,
param.domain_xmin, param.domain_xmin + param.domain_edge_width,
EDGE2);
else if (param.dimension == 2)
{
if (param.elementtype == "unstructured")
{
MeshTools::Generation::build_delaunay_square
(mesh, param.coarsegridx, param.coarsegridy,
param.domain_xmin, param.domain_xmin +
param.domain_edge_width,
param.domain_ymin, param.domain_ymin +
param.domain_edge_length
,
TRI3);
}
else
MeshTools::Generation::build_square
(mesh, param.coarsegridx, param.coarsegridy,
param.domain_xmin, param.domain_xmin + param.domain_edge_width,
param.domain_ymin, param.domain_ymin + param.domain_edge_length,
elemtype);
}
else if (param.dimension == 3)
{
if (param.elementtype == "tri")
libmesh_error();
MeshTools::Generation::build_cube
(mesh, param.coarsegridx, param.coarsegridy, param.coarsegridz,
param.domain_xmin, param.domain_xmin + param.domain_edge_width,
param.domain_ymin, param.domain_ymin + param.domain_edge_length,
param.domain_zmin, param.domain_zmin + param.domain_edge_height,
HEX8);
}
else
libmesh_error();
}
...
if (param.domaintype == "cylinder")
{
if (param.domain_ymin != 0.0 ||
param.domain_edge_length != 1.0)
{
std::cerr << "Theta must range from 0 to 1"
<< std::endl;
libmesh_error();
}
PointLocatorTree point_locator(mesh);
MeshBase::element_iterator el = mesh.elements_begin();
const MeshBase::element_iterator end_el = mesh.elements_end();
for (; el != end_el; ++el)
{
Elem* elem = *el;
// The boundary we sew up causes problems with BoundaryInfo
// later, so for now remove all boundary ids
mesh.boundary_info->remove(elem);
for (unsigned int n = 0; n != elem->n_nodes(); ++n)
{
Node &node = *(elem->get_node(n));
if (std::abs(node(1) - param.domain_edge_length) < TOLERANCE)
{
Point wrapped_point = node;
wrapped_point(1) = 0;
Elem *wrapped_elem =
const_cast<Elem *>(point_locator(wrapped_point));
libmesh_assert(elem);
bool found_matching_point = false;
for (unsigned int wn = 0; wn != wrapped_elem->n_nodes();
++wn)
{
Node *new_node = wrapped_elem->get_node(wn);
if (new_node->absolute_fuzzy_equals(wrapped_point))
{
elem->set_node(n) = new_node;
found_matching_point = true;
break;
}
}
libmesh_assert(found_matching_point);
}
}
}
MeshBase::node_iterator nd = mesh.nodes_begin();
const MeshBase::node_iterator end_nd = mesh.nodes_end();
if (param.dimension == 3 || param.dimension == 2)
{
for (; nd != end_nd; ++nd)
{
Node &node = **nd;
Real r = node(0), theta = node(1); //, z = node(2);
node(0) = r * cos(2.0*M_PI*theta);
node(1) = r * sin(2.0*M_PI*theta);
// node(2) = z;
// The boundary we sew up causes problems with
// BoundaryInfo later, so for now remove all boundary
// ids
mesh.boundary_info->remove(&node);
}
}
else
{
libmesh_error();
}
}
-------------------------------------------------------------------------
This SF.Net email is sponsored by the Moblin Your Move Developer's challenge
Build the coolest Linux based applications with Moblin SDK & win great prizes
Grand prize is a trip for two to an Open Source event anywhere in the world
http://moblin-contest.org/redirect.php?banner_id=100&url=/
_______________________________________________
Libmesh-users mailing list
[email protected]
https://lists.sourceforge.net/lists/listinfo/libmesh-users
