Mesh traversal

On Thu, Dec 24, 2009 at 03:52:33PM -0000, dbeacham wrote:
> New question #95085 on DOLFIN:
> https://answers.launchpad.net/dolfin/+question/95085
>
> Is there an easy way to systematically traverse a given mesh in either python or cpp?
>
> ie I want to create a spline represenation of the boundary, but cannot find an easy way to calculate the index of the adjacent/connected vertices/edges/cells within the mesh. Currently I have to create my own connectivity structure using cells(mesh), but was wondering if getting adjacency/connectivity information between and within dimensions could be more easily achieved?
>
> Thanks, David.

You can get all adjacency information you need by iterators. For a
tetrahedral mesh, there are 4 different kinds of entities (vertices,
edges, faces, cells) and you can get all 16 combinations of
connectivities between the entities.

Say you want to iterate over all neighboring vertices of all vertices
of all edges on all cells. Then just do

  for cell in cells(mesh):
      for edge in edges(cell):
          for vertex0 in vertices(edge):
              for vertex1 in vertices(vertex0):
                  ...

You can do the same in C++ (which will be much faster) using
CellIterator, EdgeIterator, VertexIterator.

You can also get the array of indices for the neighboring entities by
using the num_entities and entities member functions. For example

  cell.num_entities(0)

returns the number of vertices of a cell and

  cell.entities(0)

returns the array of indices for those vertices.

Note that the num_entities and entities functions require explicit
initialization of entities (which the iterators do automatically), so
to get connectivity between edges and faces, you would need to do

  mesh.init(1, 2)

--
Anders

Thanks for this - I hadn't appreciated the need to call mesh.init.

One quick clarification:

Say that a vertex, v, has v.entities(1) = [300 400 500] (ie is adjacent to edges 300, 400 and 500), then in order to find the information on any of those edges, I'd have to cycle through the edge iterator until I found the edge, e, such that e.index() = 300, 400 or 500?

dbeacham wrote:
> Question #95085 on DOLFIN changed:
> https://answers.launchpad.net/dolfin/+question/95085
>
> Status: Answered => Solved
>
> dbeacham confirmed that the question is solved:
> Thanks for this - I hadn't appreciated the need to call mesh.init.
>
> One quick clarification:
>
> Say that a vertex, v, has v.entities(1) = [300 400 500] (ie is adjacent
> to edges 300, 400 and 500), then in order to find the information on any
> of those edges, I'd have to cycle through the edge iterator until I
> found the edge, e, such that e.index() = 300, 400 or 500?
>

If you know the index, you can create the edge

  Edge edge(const Mesh& mesh, uint index);

Garth

On Mon, Jan 04, 2010 at 11:58:39AM +0100, Garth N. Wells wrote:
>
>
> dbeacham wrote:
> > Question #95085 on DOLFIN changed:
> > https://answers.launchpad.net/dolfin/+question/95085
> >
> > Status: Answered => Solved
> >
> > dbeacham confirmed that the question is solved:
> > Thanks for this - I hadn't appreciated the need to call mesh.init.
> >
> > One quick clarification:
> >
> > Say that a vertex, v, has v.entities(1) = [300 400 500] (ie is adjacent
> > to edges 300, 400 and 500), then in order to find the information on any
> > of those edges, I'd have to cycle through the edge iterator until I
> > found the edge, e, such that e.index() = 300, 400 or 500?
> >
>
> If you know the index, you can create the edge
>
> Edge edge(const Mesh& mesh, uint index);
>
> Garth

Furthermore, if you should need to know the local number of that
entity, you can do

  Edge edge(mesh, 400);
  uint local_index = cell.index(edge);

Then local_index should be 1.

--
Anders