Poster of Linux kernelThe best gift for a Linux geek
TriangleInterface

TriangleInterface

Section: C Library Functions (3) Updated: Thu Apr 7 2011
Local index Up
 

NAME

TriangleInterface -  

SYNOPSIS


#include <mesh_triangle_support.h>  

Classes


class ArbitraryHole

class Hole

class PolygonHole
 

Public Types


enum TriangulationType { GENERATE_CONVEX_HULL = 0, PSLG = 1, INVALID_TRIANGULATION_TYPE }
 

Public Member Functions


TriangleInterface (UnstructuredMesh &mesh)

~TriangleInterface ()

void triangulate ()

ElemType & elem_type ()

Real & desired_area ()

TriangulationType & triangulation_type ()

bool & insert_extra_points ()

bool & smooth_after_generating ()

void attach_hole_list (const std::vector< Hole * > *holes)
 

Public Attributes


std::vector< std::pair< unsigned int, unsigned int > > segments
 

Private Attributes


UnstructuredMesh & _mesh

const std::vector< Hole * > * _holes

ElemType _elem_type

Real _desired_area

TriangulationType _triangulation_type

bool _insert_extra_points

bool _smooth_after_generating
 

Detailed Description

A C++ interface between LibMesh and the Triangle library written by J.R. Shewchuk.

Author:

John W. Peterson

Definition at line 108 of file mesh_triangle_support.h.  

Member Enumeration Documentation

 

enum TriangleInterface::TriangulationTypeThe TriangulationType is used with the general triangulate function defind below.

Enumerator:

GENERATE_CONVEX_HULL
PSLG
Uses the triangle library to first generate a convex hull from the set of points passed in, and then triangulate this set of points. This is probably the most common type of usage.
INVALID_TRIANGULATION_TYPE
Use the triangle library to triangulate a Planar Straight Line Graph which is defined implicitly by the order of the 'points' vector: a straight line is assumed to lie between each successive pair of points, with an additional line joining the final and first points. In case your triangulation is a little too 'structured' looking (which can happen when the initial PSLG is really simple) you can try to make the resulting triangulation a little more 'unstructured' looking by setting insert_points to true in the triangulate() function. Does nothing, used as a default value.

Definition at line 136 of file mesh_triangle_support.h.

    {
      GENERATE_CONVEX_HULL = 0,
      PSLG = 1,
      INVALID_TRIANGULATION_TYPE
    };
 

Constructor & Destructor Documentation

 

TriangleInterface::TriangleInterface (UnstructuredMesh &mesh) [inline]The constructor. A reference to the mesh containing the points which are to be triangulated must be provided. Unless otherwise specified, a convex hull will be computed for set of input points and the convex hull will be meshed.

Definition at line 117 of file mesh_triangle_support.h.

    : _mesh(mesh),
      _holes(NULL),
      _elem_type(TRI3),
      _desired_area(0.1),
      _triangulation_type(GENERATE_CONVEX_HULL),
      _insert_extra_points(false),
      _smooth_after_generating(true)
  {}
 

TriangleInterface::~TriangleInterface () [inline]Empty destructor.

Definition at line 130 of file mesh_triangle_support.h.

{}
 

Member Function Documentation

 

void TriangleInterface::attach_hole_list (const std::vector< Hole * > *holes) [inline]Attaches a vector of Hole* pointers which will be meshed around.

Definition at line 209 of file mesh_triangle_support.h.

References _holes.

Referenced by MeshTools::Generation::build_delaunay_square().

{_holes = holes;}
 

Real& TriangleInterface::desired_area () [inline]Sets and/or gets the desired triangle area.

Definition at line 187 of file mesh_triangle_support.h.

References _desired_area.

Referenced by MeshTools::Generation::build_delaunay_square().

{return _desired_area;}
 

ElemType& TriangleInterface::elem_type () [inline]Sets and/or gets the desired element type.

Definition at line 182 of file mesh_triangle_support.h.

References _elem_type.

Referenced by MeshTools::Generation::build_delaunay_square().

{return _elem_type;}
 

bool& TriangleInterface::insert_extra_points () [inline]Sets and/or gets the flag for inserting add'l points.

Definition at line 197 of file mesh_triangle_support.h.

References _insert_extra_points.

{return _insert_extra_points;}
 

bool& TriangleInterface::smooth_after_generating () [inline]Sets/gets flag which tells whether to do Delaunay mesh smoothing after generating the grid.

Definition at line 203 of file mesh_triangle_support.h.

References _smooth_after_generating.

{return _smooth_after_generating;}
 

void TriangleInterface::triangulate ()This is the main public interface for this function. Internally, it calls Triangle's triangulate routine.

Definition at line 287 of file mesh_triangle_support.C.

References _desired_area, _elem_type, _holes, _insert_extra_points, _mesh, _smooth_after_generating, _triangulation_type, MeshBase::add_point(), MeshBase::clear(), Triangle::copy_tri_to_mesh(), Triangle::destroy(), GENERATE_CONVEX_HULL, libMesh::init(), Triangle::INPUT, INVALID_TRIANGULATION_TYPE, MeshBase::n_nodes(), MeshBase::nodes_begin(), MeshBase::nodes_end(), Triangle::OUTPUT, MeshBase::point(), PSLG, segments, MeshBase::set_mesh_dimension(), LaplaceMeshSmoother::smooth(), libMeshEnums::TRI3, and libMeshEnums::TRI6.

Referenced by MeshTools::Generation::build_delaunay_square().

{
  // Will the triangulation have holes?
  const bool have_holes = ((_holes != NULL) && (!_holes->empty()));
  
  // If the initial PSLG is really simple, e.g. an L-shaped domain or
  // a square/rectangle, the resulting triangulation may be very
  // 'structured' looking.  Sometimes this is a problem if your
  // intention is to work with an 'unstructured' looking grid.  We can
  // attempt to work around this limitation by inserting midpoints
  // into the original PSLG.  Inserting additional points into a
  // set of points meant to be a convex hull usually makes less sense.

  // May or may not need to insert new points ...
  if ((_triangulation_type==PSLG) && (_insert_extra_points))
    {
      // Make a copy of the original points from the Mesh
      std::vector<Point> original_points (_mesh.n_nodes());
      
      MeshBase::node_iterator       node_it  = _mesh.nodes_begin();
      const MeshBase::node_iterator node_end = _mesh.nodes_end();
      
      for (unsigned int ctr=0; node_it != node_end; ++node_it)
        original_points[ctr++] = **node_it;
      
      // Clear out the mesh
      _mesh.clear();
      
      // Make sure the new Mesh will be 2D
      _mesh.set_mesh_dimension(2);
  
      // Insert a new point on each PSLG at some random location
      // np=index into new points vector
      // n =index into original points vector
      for (unsigned int np=0, n=0; np<2*original_points.size(); ++np)
        {
          // the even entries are the original points
          if (np%2==0)
            _mesh.add_point(original_points[n++]);

          else // the odd entries are the midpoints of the original PSLG segments
            _mesh.add_point (0.5*(original_points[n] + original_points[n-1]));
        }
    }
  
  // Regardless of whether we added additional points, the set of points to
  // triangulate is now sitting in the mesh.

  // If the holes vector is non-NULL (and non-empty) we need to determine
  // the number of additional points which the holes will add to the
  // triangulation.
  unsigned int n_hole_points = 0;

  if (have_holes)
    {
      for (unsigned int i=0; i<_holes->size(); ++i)
        n_hole_points += (*_holes)[i]->n_points();
    }
  
  // Triangle data structure for the mesh
  Triangle::triangulateio initial;
  Triangle::triangulateio final;

  // Pseudo-Constructor for the triangle io structs
  Triangle::init(initial);
  Triangle::init(final);
    
  initial.numberofpoints = _mesh.n_nodes() + n_hole_points;
  initial.pointlist      = static_cast<REAL*>(std::malloc(initial.numberofpoints * 2 * sizeof(REAL)));

  if (_triangulation_type==PSLG)
    {
      // Implicit segment ordering: One segment per point, including hole points
      if (this->segments.empty())
        initial.numberofsegments = initial.numberofpoints; 

      // User-defined segment ordering: One segment per entry in the segments vector
      else
        initial.numberofsegments = this->segments.size();
    }
  
  else if (_triangulation_type==GENERATE_CONVEX_HULL)
    initial.numberofsegments = n_hole_points; // One segment for each hole point

  // Debugging
  // std::cout << 'Number of segments set to: ' << initial.numberofsegments << std::endl;
  
  // Allocate space for the segments (2 int per segment)
  if (initial.numberofsegments > 0)
    {
      initial.segmentlist = static_cast<int*> (std::malloc(initial.numberofsegments * 2 * sizeof(int)));
    }  


  // Copy all the holes' points and segments into the triangle struct.

  // The hole_offset is a constant offset into the points vector which points
  // past the end of the last hole point added.
  unsigned int hole_offset=0;
  
  if (have_holes)
    for (unsigned int i=0; i<_holes->size(); ++i)
      {
        for (unsigned int ctr=0, h=0; h<(*_holes)[i]->n_points(); ctr+=2, ++h)
          {
            Point p = (*_holes)[i]->point(h);

            const unsigned int index0 = 2*hole_offset+ctr;
            const unsigned int index1 = 2*hole_offset+ctr+1;

            // Save the x,y locations in the triangle struct.
            initial.pointlist[index0] = p(0);
            initial.pointlist[index1] = p(1);

            // Set the points which define the segments
            initial.segmentlist[index0] = hole_offset+h;
            initial.segmentlist[index1] = (h==(*_holes)[i]->n_points()-1) ? hole_offset : hole_offset+h+1; // wrap around
          }
        
        // Update the hole_offset for the next hole
        hole_offset += (*_holes)[i]->n_points();
      }

  
  // Copy all the non-hole points and segments into the triangle struct.
  for (unsigned int ctr=0, n=0; n<_mesh.n_nodes(); ctr+=2, ++n)
    {
      const unsigned int index0 = 2*hole_offset+ctr;
      const unsigned int index1 = 2*hole_offset+ctr+1;
      
      initial.pointlist[index0] = _mesh.point(n)(0);
      initial.pointlist[index1] = _mesh.point(n)(1);

      // If the user requested a PSLG, the non-hole points are also segments
      if (_triangulation_type==PSLG)
        {
          // Use implicit ordering to define segments
          if (this->segments.empty())
            {
              initial.segmentlist[index0] = hole_offset+n;
              initial.segmentlist[index1] = (n==_mesh.n_nodes()-1) ? hole_offset : hole_offset+n+1; // wrap around
            }
        }
    }

  
  // If the user provided it, use his ordering to define the segments
  for (unsigned int ctr=0, s=0; s<this->segments.size(); ctr+=2, ++s)
    {
      const unsigned int index0 = 2*hole_offset+ctr;
      const unsigned int index1 = 2*hole_offset+ctr+1;
      
      initial.segmentlist[index0] = hole_offset + this->segments[s].first;
      initial.segmentlist[index1] = hole_offset + this->segments[s].second;
    }



  // Tell the input struct about the holes
  if (have_holes)
    {
      initial.numberofholes = _holes->size();
      initial.holelist      = static_cast<REAL*>(std::malloc(initial.numberofholes * 2 * sizeof(REAL)));
      for (unsigned int i=0, ctr=0; i<_holes->size(); ++i, ctr+=2)
        {
          Point inside_point = (*_holes)[i]->inside();
          initial.holelist[ctr]   = inside_point(0);
          initial.holelist[ctr+1] = inside_point(1);
        }
    }
  
  // Set the triangulation flags.
  // c ~ enclose convex hull with segments
  // z ~ use zero indexing
  // B ~ Suppresses boundary markers in the output
  // Q ~ run in 'quiet' mode
  // p ~ Triangulates a Planar Straight Line Graph
  //     If the `p' switch is used, `segmentlist' must point to a list of     
  //     segments, `numberofsegments' must be properly set, and               
  //     `segmentmarkerlist' must either be set to NULL (in which case all    
  //     markers default to zero), or must point to a list of markers.
  // D ~ Conforming Delaunay: use this switch if you want all triangles
  //     in the mesh to be Delaunay, and not just constrained Delaunay
  // q ~  Quality mesh generation with no angles smaller than 20 degrees.
  //      An alternate minimum angle may be specified after the q
  // a ~ Imposes a maximum triangle area constraint.
  // -P  Suppresses the output .poly file. Saves disk space, but you lose the ability to maintain
  //     constraining segments on later refinements of the mesh.
  // Create the flag strings, depends on element type
std::ostringstream flags;

// Default flags always used
flags << 'zBPQq';

  // Flags which are specific to the type of triangulation
  switch (_triangulation_type)
    {
    case GENERATE_CONVEX_HULL:
      {
        flags << 'c';
        break;
      }

    case PSLG:
      {
        flags << 'p';
        break;
      }
      
    case INVALID_TRIANGULATION_TYPE:
      {
        libmesh_error();
        break;
      }
      
    default:
      {
        libmesh_error();
      }
    }


  // Flags specific to the type of element
  switch (_elem_type)
    {
    case TRI3:
      {
        // do nothing.
        break;
      }

    case TRI6:
      {
        flags << 'o2';
        break;
      }
      
    default:
      {
        std::cerr << 'ERROR: Unrecognized triangular element type.' << std::endl;
        libmesh_error();
      }
    }


  // If we do have holes and the user asked to GENERATE_CONVEX_HULL,
  // need to add the p flag so the triangulation respects those segments.
  if ((_triangulation_type==GENERATE_CONVEX_HULL) && (have_holes))
    flags << 'p';
  
  // Finally, add the area constraint
  flags << 'a' << std::fixed << _desired_area;

  // Refine the initial output to conform to the area constraint
  Triangle::triangulate(const_cast<char*>(flags.str().c_str()),
                        &initial,
                        &final,
                        NULL); // voronoi ouput -- not used
  
  
  // Send the information computed by Triangle to the Mesh.
  Triangle::copy_tri_to_mesh(final,
                             _mesh,
                             _elem_type);
      
  // To the naked eye, a few smoothing iterations usually looks better,
  // so we do this by default unless the user says not to.
  if (this->_smooth_after_generating)
    LaplaceMeshSmoother(_mesh).smooth(2);

    
  // Clean up.    
  Triangle::destroy(initial,      Triangle::INPUT);
  Triangle::destroy(final,        Triangle::OUTPUT);
  
}
 

TriangulationType& TriangleInterface::triangulation_type () [inline]Sets and/or gets the desired triangulation type.

Definition at line 192 of file mesh_triangle_support.h.

References _triangulation_type.

Referenced by MeshTools::Generation::build_delaunay_square().

{return _triangulation_type;}
 

Member Data Documentation

 

Real TriangleInterface::_desired_area [private]The desired area for the elements in the resulting mesh.

Definition at line 245 of file mesh_triangle_support.h.

Referenced by desired_area(), and triangulate().  

ElemType TriangleInterface::_elem_type [private]The type of elements to generate. (Defaults to TRI3).

Definition at line 240 of file mesh_triangle_support.h.

Referenced by elem_type(), and triangulate().  

const std::vector<Hole*>* TriangleInterface::_holes [private]A pointer to a vector of Hole*s. If this is NULL, there are no holes!

Definition at line 234 of file mesh_triangle_support.h.

Referenced by attach_hole_list(), and triangulate().  

bool TriangleInterface::_insert_extra_points [private]Flag which tells whether or not to insert additional nodes before triangulation. This can sometimes be used to 'de-regularize' the resulting triangulation.

Definition at line 259 of file mesh_triangle_support.h.

Referenced by insert_extra_points(), and triangulate().  

UnstructuredMesh& TriangleInterface::_mesh [private]Reference to the mesh which is to be created by triangle.

Definition at line 228 of file mesh_triangle_support.h.

Referenced by triangulate().  

bool TriangleInterface::_smooth_after_generating [private]Flag which tells whether we should smooth the mesh after it is generated. True by default.

Definition at line 265 of file mesh_triangle_support.h.

Referenced by smooth_after_generating(), and triangulate().  

TriangulationType TriangleInterface::_triangulation_type [private]The type of triangulation to perform: choices are: convex hull PSLG

Definition at line 252 of file mesh_triangle_support.h.

Referenced by triangulate(), and triangulation_type().  

std::vector<std::pair<unsigned int, unsigned int> > TriangleInterface::segmentsWhen constructing a PSLG, it is often not possible to do so implicitly through the ordering of the points. You can use the segments vector to specify the segments explicitly, Ex: unit square numbered counter-clockwise starting from origin segments[0] = (0,1) segments[1] = (1,2) segments[2] = (2,3) segments[3] = (3,0) Note: for this case, you could use the implicit ordering!

Definition at line 222 of file mesh_triangle_support.h.

Referenced by triangulate().

 

Author

Generated automatically by Doxygen for libMesh from the source code.


 

Index

NAME
SYNOPSIS
Classes
Public Types
Public Member Functions
Public Attributes
Private Attributes
Detailed Description
Member Enumeration Documentation
enum TriangleInterface::TriangulationTypeThe TriangulationType is used with the general triangulate function defind below.
Constructor & Destructor Documentation
TriangleInterface::TriangleInterface (UnstructuredMesh &mesh) [inline]The constructor. A reference to the mesh containing the points which are to be triangulated must be provided. Unless otherwise specified, a convex hull will be computed for set of input points and the convex hull will be meshed.
TriangleInterface::~TriangleInterface () [inline]Empty destructor.
Member Function Documentation
void TriangleInterface::attach_hole_list (const std::vector< Hole * > *holes) [inline]Attaches a vector of Hole* pointers which will be meshed around.
Real& TriangleInterface::desired_area () [inline]Sets and/or gets the desired triangle area.
ElemType& TriangleInterface::elem_type () [inline]Sets and/or gets the desired element type.
bool& TriangleInterface::insert_extra_points () [inline]Sets and/or gets the flag for inserting add'l points.
bool& TriangleInterface::smooth_after_generating () [inline]Sets/gets flag which tells whether to do Delaunay mesh smoothing after generating the grid.
void TriangleInterface::triangulate ()This is the main public interface for this function. Internally, it calls Triangle's triangulate routine.
TriangulationType& TriangleInterface::triangulation_type () [inline]Sets and/or gets the desired triangulation type.
Member Data Documentation
Real TriangleInterface::_desired_area [private]The desired area for the elements in the resulting mesh.
ElemType TriangleInterface::_elem_type [private]The type of elements to generate. (Defaults to TRI3).
const std::vector<Hole*>* TriangleInterface::_holes [private]A pointer to a vector of Hole*s. If this is NULL, there are no holes!
bool TriangleInterface::_insert_extra_points [private]Flag which tells whether or not to insert additional nodes before triangulation. This can sometimes be used to 'de-regularize' the resulting triangulation.
UnstructuredMesh& TriangleInterface::_mesh [private]Reference to the mesh which is to be created by triangle.
bool TriangleInterface::_smooth_after_generating [private]Flag which tells whether we should smooth the mesh after it is generated. True by default.
TriangulationType TriangleInterface::_triangulation_type [private]The type of triangulation to perform: choices are: convex hull PSLG
std::vector<std::pair<unsigned int, unsigned int> > TriangleInterface::segmentsWhen constructing a PSLG, it is often not possible to do so implicitly through the ordering of the points. You can use the segments vector to specify the segments explicitly, Ex: unit square numbered counter-clockwise starting from origin segments[0] = (0,1) segments[1] = (1,2) segments[2] = (2,3) segments[3] = (3,0) Note: for this case, you could use the implicit ordering!
Author

This document was created by man2html, using the manual pages.
Time: 21:57:47 GMT, April 16, 2011