MoFEM::TetGenInterface Struct Reference

TetGen interface. More...

#include <src/interfaces/TetGenInterface.hpp>

Inheritance diagram for MoFEM::TetGenInterface:

Collaboration diagram for MoFEM::TetGenInterface:

Public Types
enum	tetGenNodesTypes { RIDGEVERTEX = 0, FREESEGVERTEX = 1, FREEFACETVERTEX = 2, FREEVOLVERTEX = 3 }
typedef std::map< EntityHandle, unsigned long >	moabTetGen_Map
typedef std::map< unsigned long, EntityHandle >	tetGenMoab_Map
typedef std::map< int, Range >	idxRange_Map

Public Member Functions
MoFEMErrorCode	query_interface (boost::typeindex::type_index type_index, UnknownInterface **iface) const
	TetGenInterface (const MoFEM::Core &core)
MoFEMErrorCode	inData (Range &ents, tetgenio &in, moabTetGen_Map &moab_tetgen_map, tetGenMoab_Map &tetgen_moab_map, Tag th=NULL)
	create TetGen data structure form range of moab entities More...
MoFEMErrorCode	setGeomData (tetgenio &in, moabTetGen_Map &moab_tetgen_map, tetGenMoab_Map &tetgen_moab_map, std::map< int, Range > &type_ents)
	set point tags and type More...
MoFEMErrorCode	outData (tetgenio &in, tetgenio &out, moabTetGen_Map &moab_tetgen_map, tetGenMoab_Map &tetgen_moab_map, Range *ents=NULL, bool id_in_tags=false, bool error_if_created=false, bool assume_first_nodes_the_same=false, Tag th=nullptr)
	get entities for TetGen data structure More...
MoFEMErrorCode	outData (tetgenio &in, tetgenio &out, moabTetGen_Map &moab_tetgen_map, tetGenMoab_Map &tetgen_moab_map, BitRefLevel bit, bool id_in_tags=false, bool error_if_created=false, bool assume_first_nodes_the_same=false, Tag th=nullptr)
	get entities for TetGen data structure More...
MoFEMErrorCode	setFaceData (std::vector< std::pair< Range, int > > &markers, tetgenio &in, moabTetGen_Map &moab_tetgen_map, tetGenMoab_Map &tetgen_moab_map)
	set markers to faces More...
MoFEMErrorCode	getTriangleMarkers (tetGenMoab_Map &tetgen_moab_map, tetgenio &out, Range *ents=NULL, idxRange_Map *ents_map=NULL, bool only_non_zero=true)
	get markers to faces More...
MoFEMErrorCode	setRegionData (std::vector< std::pair< EntityHandle, int > > &regions, tetgenio &in, Tag th=NULL)
	set region data to tetrahedral More...
MoFEMErrorCode	getRegionData (tetGenMoab_Map &tetgen_moab_map, tetgenio &out, Range *ents=NULL, idxRange_Map *ents_map=NULL)
	get region data to tetrahedral More...
MoFEMErrorCode	tetRahedralize (char switches[], tetgenio &in, tetgenio &out)
	run tetgen More...
MoFEMErrorCode	loadPoly (char file_name[], tetgenio &in)
	load poly file More...
MoFEMErrorCode	checkPlanar_Trinagle (double coords[], bool *result, const double eps=1e-9)
MoFEMErrorCode	groupPlanar_Triangle (Range &tris, std::vector< Range > &sorted, const double eps=1e-9, Tag th=NULL)
MoFEMErrorCode	groupRegion_Triangle (Range &tris, std::vector< std::vector< Range > > &sorted, const double eps=1e-9)
	Group surface triangles in planar regions. More...
MoFEMErrorCode	makePolygonFacet (Range &ents, Range &polygons, bool reduce_edges=false, Range *not_reducable_nodes=NULL, const double eps=1e-9, Tag th=NULL)
Public Member Functions inherited from MoFEM::UnknownInterface
template<class IFACE >
MoFEMErrorCode	registerInterface (bool error_if_registration_failed=true)
	Register interface. More...
template<class IFACE >
MoFEMErrorCode	getInterface (IFACE *&iface) const
	Get interface reference to pointer of interface. More...
template<class IFACE >
MoFEMErrorCode	getInterface (IFACE **const iface) const
	Get interface pointer to pointer of interface. More...
template<class IFACE , typename boost::enable_if< boost::is_pointer< IFACE >, int >::type = 0>
IFACE	getInterface () const
	Get interface pointer to pointer of interface. More...
template<class IFACE , typename boost::enable_if< boost::is_reference< IFACE >, int >::type = 0>
IFACE	getInterface () const
	Get reference to interface. More...
template<class IFACE >
IFACE *	getInterface () const
	Function returning pointer to interface. More...
virtual	~UnknownInterface ()=default

Public Attributes
MoFEM::Core &	cOre

Additional Inherited Members
Static Public Member Functions inherited from MoFEM::UnknownInterface
static MoFEMErrorCode	getLibVersion (Version &version)
	Get library version. More...
static MoFEMErrorCode	getFileVersion (moab::Interface &moab, Version &version)
	Get database major version. More...
static MoFEMErrorCode	setFileVersion (moab::Interface &moab, Version version=Version(MoFEM_VERSION_MAJOR, MoFEM_VERSION_MINOR, MoFEM_VERSION_BUILD))
	Get database major version. More...
static MoFEMErrorCode	getInterfaceVersion (Version &version)
	Get database major version. More...

Detailed Description

TetGen interface.

Definition at line 23 of file TetGenInterface.hpp.

Member Typedef Documentation

idxRange_Map

typedef std::map<int, Range> MoFEM::TetGenInterface::idxRange_Map

Definition at line 34 of file TetGenInterface.hpp.

moabTetGen_Map

typedef std::map<EntityHandle, unsigned long> MoFEM::TetGenInterface::moabTetGen_Map

Definition at line 32 of file TetGenInterface.hpp.

tetGenMoab_Map

typedef std::map<unsigned long, EntityHandle> MoFEM::TetGenInterface::tetGenMoab_Map

Definition at line 33 of file TetGenInterface.hpp.

Member Enumeration Documentation

tetGenNodesTypes

enum MoFEM::TetGenInterface::tetGenNodesTypes

Enumerator
RIDGEVERTEX
FREESEGVERTEX
FREEFACETVERTEX
FREEVOLVERTEX

Definition at line 49 of file TetGenInterface.hpp.

                        {
    RIDGEVERTEX = 0,
    FREESEGVERTEX = 1,
    FREEFACETVERTEX = 2,
    FREEVOLVERTEX = 3
  };

Constructor & Destructor Documentation

TetGenInterface()

MoFEM::TetGenInterface::TetGenInterface ( const MoFEM::Core &  core )

inline

Definition at line 29 of file TetGenInterface.hpp.

      : cOre(const_cast<MoFEM::Core &>(core)) {}

Member Function Documentation

checkPlanar_Trinagle()

MoFEMErrorCode MoFEM::TetGenInterface::checkPlanar_Trinagle	(	double	coords[],
		bool *	result,
		const double	eps = 1e-9
	)

Definition at line 662 of file TetGenInterface.cpp.

                                                                       {
  MoFEMFunctionBeginHot;
  double *pa = &coords[0];
  double *pb = &coords[3];
  double *pc = &coords[6];
  double *pd = &coords[9];
  double adx = pa[0] - pd[0];
  double bdx = pb[0] - pd[0];
  double cdx = pc[0] - pd[0];
  double ady = pa[1] - pd[1];
  double bdy = pb[1] - pd[1];
  double cdy = pc[1] - pd[1];
  double adz = pa[2] - pd[2];
  double bdz = pb[2] - pd[2];
  double cdz = pc[2] - pd[2];
  double v = adx * (bdy * cdz - bdz * cdy) + bdx * (cdy * adz - cdz * ady) +
             cdx * (ady * bdz - adz * bdy);
  double l = sqrt(pow(pa[0] - pb[0], 2) + pow(pa[1] - pb[1], 2) +
                  pow(pa[2] - pb[2], 2)) +
             sqrt(pow(pa[0] - pc[0], 2) + pow(pa[1] - pc[1], 2) +
                  pow(pa[2] - pc[2], 2)) +
             sqrt(pow(pa[0] - pd[0], 2) + pow(pa[1] - pd[1], 2) +
                  pow(pa[2] - pd[2], 2)) +
             sqrt(pow(pb[0] - pc[0], 2) + pow(pb[1] - pc[1], 2) +
                  pow(pb[2] - pc[2], 2)) +
             sqrt(pow(pb[0] - pd[0], 2) + pow(pb[1] - pd[1], 2) +
                  pow(pb[2] - pd[2], 2)) +
             sqrt(pow(pc[0] - pd[0], 2) + pow(pc[1] - pd[1], 2) +
                  pow(pc[2] - pd[2], 2));
  // std::cerr << fabs(v/pow(l,3)) << " ";
  *result = fabs(v / pow(l, 3)) < eps ? true : false;
  MoFEMFunctionReturnHot(0);
}

getRegionData()

MoFEMErrorCode MoFEM::TetGenInterface::getRegionData	(	tetGenMoab_Map &	tetgen_moab_map,
		tetgenio &	out,
		Range *	ents = NULL,
		idxRange_Map *	ents_map = NULL
	)

get region data to tetrahedral

Definition at line 606 of file TetGenInterface.cpp.

                                         {
  MoFEMFunctionBeginHot;
  Interface &m_field = cOre;
  int nbattributes = out.numberoftetrahedronattributes;
  if (nbattributes == 0) {
    SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
            "tetgen has no regions attributes");
  }
  Tag th_region;
  rval = m_field.get_moab().tag_get_handle("TETGEN_REGION", th_region);
  if (rval == MB_SUCCESS) {
    rval = m_field.get_moab().tag_delete(th_region);
    CHKERRQ_MOAB(rval);
  }
  double def_marker = 0;
  CHKERR m_field.get_moab().tag_get_handle(
      "TETGEN_REGION", nbattributes, MB_TYPE_DOUBLE, th_region,
      MB_TAG_CREAT | MB_TAG_SPARSE, &def_marker);
  int ii = 0;
  for (; ii < out.numberoftetrahedra; ii++) {
    int jj = 0;
    for (; jj < nbattributes; jj++) {
      double id = out.tetrahedronattributelist[ii * nbattributes + jj];
      int iii = MBTET | (ii << MB_TYPE_WIDTH);
      if (tetgen_moab_map.find(iii) == tetgen_moab_map.end()) {
        SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
                "data inconsistency between TetGen and MoAB");
      }
      EntityHandle ent = tetgen_moab_map[iii];
      CHKERR m_field.get_moab().tag_set_data(th_region, &ent, 1, &id);
      if (ents != NULL)
        ents->insert(ent);
      if (ents_map != NULL)
        (*ents_map)[id].insert(ent);
    }
  }
  MoFEMFunctionReturnHot(0);
}

getTriangleMarkers()

MoFEMErrorCode MoFEM::TetGenInterface::getTriangleMarkers	(	tetGenMoab_Map &	tetgen_moab_map,
		tetgenio &	out,
		Range *	ents = NULL,
		idxRange_Map *	ents_map = NULL,
		bool	only_non_zero = true
	)

get markers to faces

Parameters

markers	data structure with markers
in	tetgen data structure (look to TetGen user manual)
moab_tetgen_map	mapping MoAB to TetGen entities
tetgen_moab_map	mapping TetGen to MoAB entities

Definition at line 509 of file TetGenInterface.cpp.

                                                             {
  MoFEMFunctionBeginHot;
  ErrorCode rval;
  Interface &m_field = cOre;
  Tag th_marker;
  int def_marker = 0;
  rval = m_field.get_moab().tag_get_handle(
      "TETGEN_MARKER", 1, MB_TYPE_INTEGER, th_marker,
      MB_TAG_CREAT | MB_TAG_SPARSE, &def_marker);
  CHKERRQ_MOAB(rval);
  int ii = 0;
  for (; ii < out.numberoftrifaces; ii++) {
    if (only_non_zero) {
      if (out.trifacemarkerlist[ii] == 0) {
        continue;
      }
    }
    EntityHandle conn[3];
    for (int nn = 0; nn < 3; nn++) {
      int iii = MBVERTEX | (out.trifacelist[3 * ii + nn] << MB_TYPE_WIDTH);
      if (tetgen_moab_map.find(iii) == tetgen_moab_map.end()) {
        SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
                "data inconsistency between TetGen and MoAB");
      } else {
        conn[nn] = tetgen_moab_map[iii];
      }
    }
    Range face;
    rval = m_field.get_moab().get_adjacencies(conn, 3, 2, false, face);
    CHKERRQ_MOAB(rval);
    face = face.subset_by_type(MBTRI);
    if (face.size() != 1) {
      SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
               "data inconsistency between TetGen and MoAB, %u", face.size());
    }
    if (ents != NULL)
      ents->merge(face);
    if (ents_map != NULL)
      (*ents_map)[out.trifacemarkerlist[ii]].merge(face);
    rval = m_field.get_moab().tag_set_data(th_marker, &*face.begin(), 1,
                                           &out.trifacemarkerlist[ii]);
    CHKERRQ_MOAB(rval);
  }
  MoFEMFunctionReturnHot(0);
}

groupPlanar_Triangle()

MoFEMErrorCode MoFEM::TetGenInterface::groupPlanar_Triangle	(	Range &	tris,
		std::vector< Range > &	sorted,
		const double	eps = 1e-9,
		Tag	th = NULL
	)

Definition at line 697 of file TetGenInterface.cpp.

                                                                               {
  MoFEMFunctionBegin;
 
  Interface &m_field = cOre;
  Skinner skin(&m_field.get_moab());
 
  for (;;) {
 
    Range noplanar_to_anyother;
    std::vector<Range>::iterator vit = sorted.begin();
 
    do {
 
      bool repeat = false;
 
      // get edges on vit skin
      Range skin_edges;
      CHKERR skin.find_skin(0, *vit, false, skin_edges);
 
      // skin edge nodes
      Range skin_edges_nodes;
      CHKERR m_field.get_moab().get_connectivity(skin_edges, skin_edges_nodes,
                                                 true);
 
      // get tris adjacent to vit skin edges
      Range skin_edges_tris;
      CHKERR m_field.get_moab().get_adjacencies(
          skin_edges, 2, false, skin_edges_tris, moab::Interface::UNION);
      // get tris which are part of facet
      Range inner_tris = intersect(skin_edges_tris, *vit);
      Range outer_tris = intersect(skin_edges_tris, tris);
 
      // tris coplanar with vit tris
      Range coplanar_tris;
 
      Range::iterator tit = outer_tris.begin();
      for (; tit != outer_tris.end(); tit++) {
        Range tit_conn;
        CHKERR m_field.get_moab().get_connectivity(&*tit, 1, tit_conn, true);
        tit_conn = subtract(tit_conn, skin_edges_nodes);
        if (tit_conn.empty()) {
          coplanar_tris.insert(*tit);
          noplanar_to_anyother.erase(*tit);
          repeat = true;
        } else {
          Range tit_edges;
          CHKERR m_field.get_moab().get_adjacencies(
              &*tit, 1, 1, false, tit_edges, moab::Interface::UNION);
          tit_edges = intersect(tit_edges, skin_edges);
          if (tit_edges.size() != 1) {
            SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
                    "data inconsistency");
          }
          Range inner_tri;
          CHKERR m_field.get_moab().get_adjacencies(
              tit_edges, 2, false, inner_tri, moab::Interface::UNION);
          inner_tri = intersect(inner_tri, inner_tris);
          if (inner_tri.size() != 1) {
            SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
                    "data inconsistency");
          }
          // get connectivity
          int num_nodes;
          const EntityHandle *inner_tri_conn;
          CHKERR m_field.get_moab().get_connectivity(
              *inner_tri.begin(), inner_tri_conn, num_nodes, true);
          double coords[12];
          if (th) {
            CHKERR m_field.get_moab().tag_get_data(th, inner_tri_conn, 3,
                                                   coords);
            CHKERR m_field.get_moab().tag_get_data(th, &*tit_conn.begin(), 1,
                                                   &coords[9]);
          } else {
            CHKERR m_field.get_moab().get_coords(inner_tri_conn, 3, coords);
            CHKERR m_field.get_moab().get_coords(&*tit_conn.begin(), 1,
                                                 &coords[9]);
          }
          bool coplanar;
          CHKERR checkPlanar_Trinagle(coords, &coplanar, eps);
          if (coplanar) {
            coplanar_tris.insert(*tit);
            noplanar_to_anyother.erase(*tit);
            repeat = true;
          } else {
            noplanar_to_anyother.insert(*tit);
          }
        }
      }
 
      vit->merge(coplanar_tris);
      tris = subtract(tris, *vit);
 
      if (repeat) {
        vit = sorted.begin();
      } else {
        vit++;
      }
 
    } while (vit != sorted.end());
 
    if (noplanar_to_anyother.empty()) {
      MoFEMFunctionReturnHot(0);
    } else {
      Range seed;
      seed.insert(noplanar_to_anyother[0]);
      tris.erase(noplanar_to_anyother[0]);
      sorted.push_back(seed);
    }
  }
 
  MoFEMFunctionReturn(0);
}

groupRegion_Triangle()

MoFEMErrorCode MoFEM::TetGenInterface::groupRegion_Triangle	(	Range &	tris,
		std::vector< std::vector< Range > > &	sorted,
		const double	eps = 1e-9
	)

Group surface triangles in planar regions.

Parameters

tris	input triangles
sorted	output sorted planar faces
eps	tolerance

Definition at line 812 of file TetGenInterface.cpp.

                                                                        {
  MoFEMFunctionBeginHot;
 
  // PetscAttachDebugger();
 
  Range seed;
  seed.insert(tris[0]);
  tris.erase(tris[0]);
  std::vector<Range> vec_seed;
  vec_seed.push_back(seed);
  sorted.push_back(vec_seed);
 
  for (;;) {
    std::vector<Range> &vec = sorted.back();
    CHKERR groupPlanar_Triangle(tris, vec, eps);
    if (tris.empty()) {
      MoFEMFunctionReturnHot(0);
    } else {
      Range seed;
      seed.insert(tris[0]);
      tris.erase(tris[0]);
      std::vector<Range> vec_seed;
      vec_seed.push_back(seed);
      sorted.push_back(vec_seed);
    }
  }
 
  MoFEMFunctionReturnHot(0);
}

inData()

MoFEMErrorCode MoFEM::TetGenInterface::inData	(	Range &	ents,
		tetgenio &	in,
		moabTetGen_Map &	moab_tetgen_map,
		tetGenMoab_Map &	tetgen_moab_map,
		Tag	th = NULL
	)

create TetGen data structure form range of moab entities

Move mesh to TetGen data structures

Parameters

ents	range of entities (tetrahedrons or nodes)
in	tetgen data structure (look to TetGen user manual)
moab_tetgen_map	mapping moab to TetGen entities
tetgen_moab_map	mapping tetgen to moab entities

Definition at line 38 of file TetGenInterface.cpp.

                                {
  MoFEMFunctionBegin;
 
  Interface &m_field = cOre;
  Range::iterator it;
 
  Tag th_marker;
  int def_marker = 0;
  CHKERR m_field.get_moab().tag_get_handle(
      "TETGEN_MARKER", 1, MB_TYPE_INTEGER, th_marker,
      MB_TAG_CREAT | MB_TAG_SPARSE, &def_marker);
 
  // All indices start from 0
  in.firstnumber = 0;
 
  Range points = ents.subset_by_dimension(0);
  in.numberofpoints = points.size();
  if (points.size() > 0) {
    in.pointlist = new double[in.numberofpoints * 3];
    in.pointmarkerlist = new int[in.numberofpoints];
    if (th) {
      CHKERR m_field.get_moab().tag_get_data(th, points, in.pointlist);
    } else {
      CHKERR m_field.get_moab().get_coords(points, in.pointlist);
    }
    CHKERR m_field.get_moab().tag_get_data(th_marker, points,
                                           in.pointmarkerlist);
    it = points.begin();
    for (int ii = 0; it != points.end(); it++, ii++) {
      unsigned long iii = MBVERTEX | (ii << MB_TYPE_WIDTH);
      tetgen_moab_map[iii] = *it;
      moab_tetgen_map[*it] = iii;
    }
  }
 
  Range tets = ents.subset_by_type(MBTET);
  in.numberoftetrahedra = tets.size();
  if (in.numberoftetrahedra > 0) {
    in.tetrahedronlist = new int[4 * ents.subset_by_type(MBTET).size()];
    it = tets.begin();
    for (int ii = 0; it != tets.end(); it++, ii++) {
      int num_nodes;
      const EntityHandle *conn;
      CHKERR m_field.get_moab().get_connectivity(*it, conn, num_nodes, true);
#ifdef DEBUG_TETGEN
      {
        double coords[12];
        if (th) {
          CHKERR m_field.tag_get_data(th, conn, num_nodes, coords);
        } else {
          CHKERR m_field.get_moab().get_coords(conn, num_nodes, coords);
        }
        double diff_n[12];
        ShapeDiffMBTET(diff_n);
        FTensor::Tensor1<double *, 3> t_diff_n(&diff_n[0], &diff_n[1],
                                               &diff_n[2], 3);
        FTensor::Tensor1<double *, 3> t_coords(&coords[0], &coords[1],
                                               &coords[2], 3);
        FTensor::Tensor2<double, 3, 3> jac;
        FTensor::Index<'i', 3> i;
        FTensor::Index<'j', 3> j;
        jac(i, j) = 0;
        for (int nn = 0; nn != 4; nn++) {
          jac(i, j) += t_coords(i) * t_diff_n(j);
          ++t_coords;
          ++t_diff_n;
        }
        double det;
        determinantTensor3by3(jac, det);
        if (det <= 0) {
          SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
                   "Negative volume det = %6.4e", det);
        }
      }
#endif
      tetgen_moab_map[MBTET | (ii << MB_TYPE_WIDTH)] = *it;
      moab_tetgen_map[*it] = MBTET | (ii << MB_TYPE_WIDTH);
      for (int nn = 0; nn != 4; nn++) {
        if (moab_tetgen_map.find(conn[nn]) == moab_tetgen_map.end()) {
          SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
                  "data inconsistency between TetGen and MoAB");
        }
        in.tetrahedronlist[4 * ii + nn] =
            moab_tetgen_map[conn[nn]] >> MB_TYPE_WIDTH;
      }
    }
  }
 
  Range tris = ents.subset_by_type(MBTRI);
  in.numberoftrifaces = tris.size();
  if (in.numberoftrifaces) {
    in.trifacelist = new int[3 * in.numberoftrifaces];
    in.trifacemarkerlist = new int[in.numberoftrifaces];
    // std::fill(&in.trifacemarkerlist[0],&in.trifacemarkerlist[in.numberoftrifaces],1);
    CHKERR m_field.get_moab().tag_get_data(th_marker, tris,
                                           in.trifacemarkerlist);
    it = tris.begin();
    for (int ii = 0; it != tris.end(); it++, ii++) {
      int num_nodes;
      const EntityHandle *conn;
      CHKERR m_field.get_moab().get_connectivity(*it, conn, num_nodes, true);
      int order[] = {0, 1, 2};
      Range adj_tets;
      CHKERR m_field.get_moab().get_adjacencies(&*it, 1, 3, true, adj_tets);
      adj_tets = intersect(adj_tets, tets);
      if (adj_tets.size() == 1) {
        int side_number;
        int sense;
        int offset;
        CHKERR m_field.get_moab().side_number(adj_tets[0], *it, side_number,
                                              sense, offset);
        if (sense == -1) {
          order[0] = 1;
          order[1] = 0;
        }
      }
      tetgen_moab_map[MBTRI | (ii << MB_TYPE_WIDTH)] = *it;
      moab_tetgen_map[*it] = MBTRI | (ii << MB_TYPE_WIDTH);
      for (int nn = 0; nn < 3; nn++) {
        if (moab_tetgen_map.find(conn[order[nn]]) == moab_tetgen_map.end()) {
          SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
                  "data inconsistency between TetGen and MoAB");
        }
        in.trifacelist[3 * ii + nn] =
            moab_tetgen_map[conn[order[nn]]] >> MB_TYPE_WIDTH;
      }
    }
  }
 
  Range edges = ents.subset_by_type(MBEDGE);
  in.numberofedges = edges.size();
  if (in.numberofedges > 0) {
    in.edgelist = new int[2 * in.numberofedges];
    in.edgemarkerlist = new int[in.numberofedges];
    // std::fill(&in.edgemarkerlist[0],&in.edgemarkerlist[in.numberofedges],1);
    CHKERR m_field.get_moab().tag_get_data(th_marker, edges, in.edgemarkerlist);
    it = edges.begin();
    for (int ii = 0; it != edges.end(); it++, ii++) {
      int num_nodes;
      const EntityHandle *conn;
      CHKERR m_field.get_moab().get_connectivity(*it, conn, num_nodes, true);
      tetgen_moab_map[MBEDGE | (ii << MB_TYPE_WIDTH)] = *it;
      moab_tetgen_map[*it] = MBEDGE | (ii << MB_TYPE_WIDTH);
      for (int nn = 0; nn < 2; nn++) {
        if (moab_tetgen_map.find(conn[nn]) == moab_tetgen_map.end()) {
          SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
                  "data inconsistency between TetGen and MoAB");
        }
        in.edgelist[2 * ii + nn] = moab_tetgen_map[conn[nn]] >> MB_TYPE_WIDTH;
      }
    }
  }
 
  MoFEMFunctionReturn(0);
}

loadPoly()

MoFEMErrorCode MoFEM::TetGenInterface::loadPoly	(	char	file_name[],
		tetgenio &	in
	)

load poly file

Definition at line 654 of file TetGenInterface.cpp.

                                                                       {
  MoFEMFunctionBeginHot;
  if (!in.load_poly(file_name)) {
    SETERRQ(PETSC_COMM_SELF, MOFEM_OPERATION_UNSUCCESSFUL,
            "can not read TetGen poly file");
  }
  MoFEMFunctionReturnHot(0);
}

makePolygonFacet()

MoFEMErrorCode MoFEM::TetGenInterface::makePolygonFacet	(	Range &	ents,
		Range &	polygons,
		bool	reduce_edges = false,
		Range *	not_reducable_nodes = NULL,
		const double	eps = 1e-9,
		Tag	th = NULL
	)

make planar polygon facet

Parameters

ents	surface triangles
polygons	output list of polygons
reduce_edges	reduce edges if on the line
not_reducable_nodes	do not reduce node on edge if in this range
eps	tolerance

Bug:

assumes that are no holes

Definition at line 842 of file TetGenInterface.cpp.

                                                                           {
  MoFEMFunctionBegin;
  // FIXME: assumes that are no holes
 
  if (ents.empty()) {
    SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
            "no ents to build polygon");
  }
 
  Interface &m_field = cOre;
 
  Skinner skin(&m_field.get_moab());
 
  Range skin_edges;
  CHKERR skin.find_skin(0, ents, false, skin_edges);
 
  std::vector<EntityHandle> polygon_nodes;
  EntityHandle seed = skin_edges[0];
  Range seen_edges;
  seen_edges.insert(seed);
  skin_edges.erase(seed);
  int num_nodes;
  const EntityHandle *conn;
  CHKERR m_field.get_moab().get_connectivity(seed, conn, num_nodes, true);
  polygon_nodes.push_back(conn[0]);
  // std::cerr << std::endl;
  // std::cerr << conn[0] << " " << conn[1] << std::endl;
  do {
    EntityHandle last_node = polygon_nodes.back();
    Range adj_edges;
    CHKERR m_field.get_moab().get_adjacencies(&last_node, 1, 1, false,
                                              adj_edges);
    adj_edges = intersect(adj_edges, skin_edges);
    if (adj_edges.size() == 0) {
      break;
    }
    if (adj_edges.size() != 1) {
      SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
              "should be only one edge");
    }
    seen_edges.insert(adj_edges[0]);
    skin_edges.erase(adj_edges[0]);
    CHKERR m_field.get_moab().get_connectivity(adj_edges[0], conn, num_nodes,
                                               true);
    EntityHandle add_node = (last_node == conn[0]) ? conn[1] : conn[0];
    polygon_nodes.push_back(add_node);
    // std::cerr << "\t" << add_node << std::endl;
  } while (1);
 
  if (reduce_edges) {
    // std::cerr << "polygon " << polygon_nodes.size();
    std::vector<EntityHandle>::iterator pit = polygon_nodes.begin();
    // std::cerr << std::endl;
    for (; pit != polygon_nodes.end();) {
      if (not_reducable_nodes != NULL) {
        if (not_reducable_nodes->find(*pit) != not_reducable_nodes->end()) {
          pit++;
          continue;
        }
      }
      EntityHandle mm;
      if (pit == polygon_nodes.begin()) {
        mm = polygon_nodes.back();
      } else {
        mm = *(pit - 1);
      }
      EntityHandle mc = *pit;
      EntityHandle mp;
      if (polygon_nodes.back() == mc) {
        mp = polygon_nodes[0];
      } else {
        mp = *(pit + 1);
      }
      double coords[9];
      if (th) {
        CHKERR m_field.get_moab().tag_get_data(th, &mm, 1, &coords[3 * 0]);
        CHKERR m_field.get_moab().tag_get_data(th, &mc, 1, &coords[3 * 1]);
        CHKERR m_field.get_moab().tag_get_data(th, &mp, 1, &coords[3 * 2]);
      } else {
        CHKERR m_field.get_moab().get_coords(&mm, 1, &coords[3 * 0]);
        CHKERR m_field.get_moab().get_coords(&mc, 1, &coords[3 * 1]);
        CHKERR m_field.get_moab().get_coords(&mp, 1, &coords[3 * 2]);
      }
      cblas_daxpy(3, -1, &coords[3 * 1], 1, &coords[3 * 0], 1); // mm = mm -
                                                                // mc
      cblas_daxpy(3, -1, &coords[3 * 1], 1, &coords[3 * 2], 1); // mp = mp -
                                                                // mc
      double spin[9];
      CHKERR Spin(spin, &coords[3 * 0]);
      double l0 = cblas_dnrm2(3, &coords[3 * 0], 1);
      cblas_dgemv(CblasRowMajor, CblasNoTrans, 3, 3, 1. / l0, spin, 3,
                  &coords[3 * 2], 1, 0., &coords[3 * 1], 1);
      double dot = cblas_dnrm2(3, &coords[3 * 1], 1);
      // std::cerr << mm << " " << mc << " " << mp << " " << dot << std::endl;
      if (dot < eps) {
        polygon_nodes.erase(pit);
        pit = polygon_nodes.begin();
        // std::cerr << std::endl;
      } else {
        pit++;
      }
    }
  }
 
  Range existing_polygon;
  CHKERR m_field.get_moab().get_adjacencies(
      &polygon_nodes[0], polygon_nodes.size(), 2, true, existing_polygon);
  if (existing_polygon.empty()) {
    EntityHandle polygon;
    CHKERR m_field.get_moab().create_element(MBPOLYGON, &polygon_nodes[0],
                                             polygon_nodes.size(), polygon);
    polygons.insert(polygon);
  } else {
    polygons.merge(existing_polygon);
  }
 
  MoFEMFunctionReturn(0);
}

outData() [1/2]

MoFEMErrorCode MoFEM::TetGenInterface::outData	(	tetgenio &	in,
		tetgenio &	out,
		moabTetGen_Map &	moab_tetgen_map,
		tetGenMoab_Map &	tetgen_moab_map,
		BitRefLevel	bit,
		bool	id_in_tags = false,
		bool	error_if_created = false,
		bool	assume_first_nodes_the_same = false,
		Tag	th = nullptr
	)

get entities for TetGen data structure

Parameters

ents	range of entities (tetrahedrons or nodes)
in	tetgen data structure (look to TetGen user manual)
moab_tetgen_map	mapping MoAB to TetGen entities
tetgen_moab_map	mapping TetGen to MoAB entities
ents	rerun entities which are in TetGen data structure
bit	set level to created entities
error_if_created	throw error if node need to be created

Definition at line 439 of file TetGenInterface.cpp.

                                                                   {
  MoFEMFunctionBegin;
  Interface &m_field = cOre;
  Range ents;
  CHKERR outData(in, out, moab_tetgen_map, tetgen_moab_map, &ents, id_in_tags,
                 error_if_created, assume_first_nodes_the_same,th);
  CHKERR m_field.getInterface<BitRefManager>()->setBitRefLevel(
      ents.subset_by_type(MBTET), bit);
  MoFEMFunctionReturn(0);
}

outData() [2/2]

MoFEMErrorCode MoFEM::TetGenInterface::outData	(	tetgenio &	in,
		tetgenio &	out,
		moabTetGen_Map &	moab_tetgen_map,
		tetGenMoab_Map &	tetgen_moab_map,
		Range *	ents = NULL,
		bool	id_in_tags = false,
		bool	error_if_created = false,
		bool	assume_first_nodes_the_same = false,
		Tag	th = nullptr
	)

get entities for TetGen data structure

Parameters

ents	range of entities (tetrahedrons or nodes)
in	tetgen data structure (look to TetGen user manual)
moab_tetgen_map	mapping MoAB to TetGen entities
tetgen_moab_map	mapping TetGen to moab entities
ents	rerun entities which are in TetGen Data structure
id_in_tags	use tags as entity handles, if that is a case use tag to find MoAB vertex id
error_if_created	throw error if node need to be created

Definition at line 246 of file TetGenInterface.cpp.

                                                                   {
  MoFEMFunctionBegin;
 
  Interface &m_field = cOre;
 
  Tag th_marker;
  int def_marker = 0;
  CHKERR m_field.get_moab().tag_get_handle(
      "TETGEN_MARKER", 1, MB_TYPE_INTEGER, th_marker,
      MB_TAG_CREAT | MB_TAG_SPARSE, &def_marker);
 
  int num_nodes = 0;
 
  int ii = 0;
  for (; ii < out.numberofpoints; ii++) {
    if (ii < in.numberofpoints) {
      bool mem_the_same = memcmp(&in.pointlist[3 * ii], &out.pointlist[3 * ii],
                                 3 * sizeof(double)) == 0;
      if (assume_first_nodes_the_same || mem_the_same) {
        unsigned long iii = MBVERTEX | (ii << MB_TYPE_WIDTH);
        if (tetgen_moab_map.find(iii) != tetgen_moab_map.end()) {
          if (!mem_the_same) {
            if (th)
              CHKERR m_field.get_moab().tag_set_data(th, &tetgen_moab_map[iii],
                                                     1, &out.pointlist[3 * ii]);
            else
              CHKERR m_field.get_moab().set_coords(&tetgen_moab_map[iii], 1,
                                                   &out.pointlist[3 * ii]);
          }
          CHKERR m_field.get_moab().tag_set_data(
              th_marker, &tetgen_moab_map[iii], 1, &out.pointmarkerlist[ii]);
          continue;
        } else {
          SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
                  "data inconsistency between TetGen and MoAB");
        }
      }
    }
    if (id_in_tags) {
      if (out.pointparamlist[ii].tag > 0) {
        EntityHandle node;
        CHKERR m_field.get_moab().handle_from_id(
            MBVERTEX, in.pointparamlist[ii].tag - 1, node);
        if (moab_tetgen_map.find(node) != moab_tetgen_map.end()) {
          continue;
        }
      }
    }
    if (error_if_created) {
      SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
              "node should not be created");
    }
    num_nodes++;
  }
 
  ReadUtilIface *iface;
  CHKERR m_field.get_moab().query_interface(iface);
 
  if (num_nodes) {
    vector<double *> arrays;
    EntityHandle startv;
    CHKERR iface->get_node_coords(3, num_nodes, 0, startv, arrays);
    Range verts(startv, startv + num_nodes - 1);
    int ii = in.numberofpoints;
    for (Range::iterator vit = verts.begin(); vit != verts.end(); vit++, ii++) {
      for (int nn = 0; nn != 3; nn++)
        arrays[nn][ii - in.numberofpoints] = out.pointlist[3 * ii + nn];
      moab_tetgen_map[*vit] = MBVERTEX | (ii << MB_TYPE_WIDTH);
      tetgen_moab_map[MBVERTEX | (ii << MB_TYPE_WIDTH)] = *vit;
      if (ents != NULL)
        ents->insert(*vit);
    }
    CHKERR m_field.get_moab().tag_set_data(
        th_marker, verts, &out.pointmarkerlist[in.numberofpoints]);
    if(th)
      CHKERR m_field.get_moab().tag_set_data(
          th, verts, &out.pointlist[3 * in.numberofpoints]);
  }
 
  std::vector<int> tetgen_ii;
 
  // Build tets
  std::vector<EntityHandle> conn_seq_tet;
  conn_seq_tet.reserve(4 * out.numberoftetrahedra);
  tetgen_ii.reserve(out.numberoftetrahedra);
  conn_seq_tet.clear();
  tetgen_ii.clear();
  ii = 0;
  for (; ii < out.numberoftetrahedra; ii++) {
    unsigned long iii = MBTET | (ii << MB_TYPE_WIDTH);
    if (ii < in.numberoftetrahedra) {
      if (memcmp(&in.tetrahedronlist[4 * ii], &out.tetrahedronlist[4 * ii],
                 4 * sizeof(int)) == 0) {
        if (tetgen_moab_map.find(iii) != tetgen_moab_map.end()) {
          if (ents != NULL)
            ents->insert(tetgen_moab_map[iii]);
          continue;
        } else {
          SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
                  "data inconsistency between TetGen and MoAB");
        }
      }
    }
    EntityHandle conn[4];
    for (int nn = 0; nn < 4; nn++) {
      int nnn = out.tetrahedronlist[4 * ii + nn];
      if (tetgen_moab_map.find(MBVERTEX | (nnn << MB_TYPE_WIDTH)) ==
          tetgen_moab_map.end()) {
        SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
                "data inconsistency between TetGen and MoAB");
      }
      conn[nn] = tetgen_moab_map.at(MBVERTEX | (nnn << MB_TYPE_WIDTH));
    }
    // auto get_coords = [&m_field](const EntityHandle *conn) {
    //   VectorDouble12 coords(12);
    //   CHKERR m_field.get_moab().get_coords(conn, 4, &coords[0]);
    //   return coords;
    // };
    // auto get_volume = [](VectorDouble12 &coords) {
    //   return Tools::tetVolume(&coords[0]);
    // };
    // if (get_volume(get_coords(conn)) < 0) {
    //   EntityHandle n0 = conn[0];
    //   conn[0] = conn[1];
    //   conn[1] = n0;
    // }
    Range tets;
    CHKERR m_field.get_moab().get_adjacencies(conn, 4, 3, false, tets);
    bool tet_found = false;
    for (auto tet : tets) {
      const EntityHandle *tet_conn;
      int num_nodes;
      CHKERR m_field.get_moab().get_connectivity(tet, tet_conn, num_nodes,
                                                 true);
      const EntityHandle *p = std::find(tet_conn, &tet_conn[4], conn[0]);
      if (p != &tet_conn[4]) {
        int s = std::distance(tet_conn, p);
        int n = 0;
        for (; n != 4; ++n) {
          const int idx[] = {0, 1, 2, 3, 0, 1, 2, 3};
          if (tet_conn[idx[s + n]] != conn[n])
            break;
        }
        if (n == 4 && !tet_found) {
          moab_tetgen_map[tet] = iii;
          tetgen_moab_map[iii] = tet;
          if (ents != NULL)
            ents->insert(tet);
          tet_found = true;
        } else if (n == 4) {
          SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY,
                  "More that one tet with the same connectivity");
        }
      } else {
        SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY, "Impossible case");
      }
    }
    if (!tet_found) {
      for (int nn = 0; nn != 4; nn++) {
        conn_seq_tet.push_back(conn[nn]);
      }
      tetgen_ii.push_back(ii);
    }
  }
 
  Range new_tets;
  if (!conn_seq_tet.empty()) {
    EntityHandle starte; // Connectivity
    EntityHandle *conn;
    int num_el = tetgen_ii.size();
    CHKERR iface->get_element_connect(num_el, 4, MBTET, 0, starte, conn);
    std::copy(conn_seq_tet.begin(), conn_seq_tet.end(), conn);
    CHKERR iface->update_adjacencies(starte, num_el, 4, conn);
    new_tets = Range(starte, starte + num_el - 1);
    std::vector<int>::iterator ii_it = tetgen_ii.begin();
    int ii = 0;
    for (Range::iterator tit = new_tets.begin(); tit != new_tets.end();
         tit++, ii_it++, ii++) {
      unsigned long iii = MBTET | ((*ii_it) << MB_TYPE_WIDTH);
      moab_tetgen_map[*tit] = iii;
      tetgen_moab_map[iii] = *tit;
    }
    if (ents != NULL)
      ents->merge(new_tets);
  }
 
  MoFEMFunctionReturn(0);
}

query_interface()

MoFEMErrorCode MoFEM::TetGenInterface::query_interface ( boost::typeindex::type_index  type_index, UnknownInterface **  iface  ) const

virtual

Implements MoFEM::UnknownInterface.

Definition at line 31 of file TetGenInterface.cpp.

                                                                 {
  *iface = const_cast<TetGenInterface *>(this);
  return 0;
}

setFaceData()

MoFEMErrorCode MoFEM::TetGenInterface::setFaceData	(	std::vector< std::pair< Range, int > > &	markers,
		tetgenio &	in,
		moabTetGen_Map &	moab_tetgen_map,
		tetGenMoab_Map &	tetgen_moab_map
	)

set markers to faces

Parameters

markers	data structure with markers
in	tetgen data structure (look to TetGen user manual)
moab_tetgen_map	mapping MoAB to TetGen entities
tetgen_moab_map	mapping TetGen to MoAB entities

Definition at line 454 of file TetGenInterface.cpp.

                                                          {
  MoFEMFunctionBeginHot;
  ErrorCode rval;
  Interface &m_field = cOre;
  in.numberoffacets = markers.size();
  in.facetlist = new tetgenio::facet[in.numberoffacets];
  in.facetmarkerlist = new int[in.numberoffacets];
  std::vector<std::pair<Range, int>>::iterator mit = markers.begin();
  for (int ii = 0; mit != markers.end(); mit++, ii++) {
    in.facetmarkerlist[ii] = mit->second;
    Range &faces = mit->first;
    tetgenio::facet *f = &(in.facetlist[ii]);
    f->numberofpolygons = faces.size();
    f->polygonlist = new tetgenio::polygon[f->numberofpolygons];
    int jj = 0;
    for (int dd = 3; dd >= 0; dd--) {
      Range dd_faces = faces.subset_by_dimension(dd);
      if (dd_faces.empty())
        continue;
      Range::iterator it = dd_faces.begin();
      for (; it != dd_faces.end(); it++, jj++) {
        int num_nodes;
        const EntityHandle *conn;
        tetgenio::polygon *p = &(f->polygonlist[jj]);
        switch (type_from_handle(*it)) {
        case MBVERTEX: {
          p->numberofvertices = 1;
          conn = &*it;
        } break;
        default: {
          rval =
              m_field.get_moab().get_connectivity(*it, conn, num_nodes, true);
          CHKERRQ_MOAB(rval);
          p->numberofvertices = num_nodes;
        }
        }
        p->vertexlist = new int[p->numberofvertices];
        for (int nn = 0; nn < p->numberofvertices; nn++) {
          if (moab_tetgen_map.find(conn[nn]) == moab_tetgen_map.end()) {
            SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
                    "data inconsistency between TetGen and MoAB");
          }
          p->vertexlist[nn] = moab_tetgen_map[conn[nn]] >> MB_TYPE_WIDTH;
        }
      }
    }
    // holes
    f->numberofholes = 0;
    f->holelist = NULL;
  }
  MoFEMFunctionReturnHot(0);
}

setGeomData()

MoFEMErrorCode MoFEM::TetGenInterface::setGeomData	(	tetgenio &	in,
		moabTetGen_Map &	moab_tetgen_map,
		tetGenMoab_Map &	tetgen_moab_map,
		std::map< int, Range > &	type_ents
	)

set point tags and type

Set type of entity, look in TetGen manual for details

std::map<int,Range> types_ents;
//RIDGEVERTEX
types_ents[TetGenInterface::RIDGEVERTEX].merge(region_tets_skin_without_boundary_nodes);
//FREESEGVERTEX
types_ents[TetGenInterface::FREESEGVERTEX].merge(crack_surface_tris_skin_nodes);
//FREEFACETVERTEX
types_ents[TetGenInterface::FREEFACETVERTEX].merge(region_tets_skin_nodes);
types_ents[TetGenInterface::FREEFACETVERTEX] =
subtract(types_ents[TetGenInterface::FREEFACETVERTEX],types_ents[TetGenInterface::RIDGEVERTEX]);
types_ents[TetGenInterface::FREEFACETVERTEX] =
subtract(types_ents[TetGenInterface::FREEFACETVERTEX],types_ents[TetGenInterface::FREESEGVERTEX]);
//FREEVOLVERTEX
types_ents[TetGenInterface::FREEVOLVERTEX].merge(region_nodes);
types_ents[TetGenInterface::FREEVOLVERTEX] =
subtract(types_ents[TetGenInterface::FREEVOLVERTEX],types_ents[TetGenInterface::RIDGEVERTEX]);
types_ents[TetGenInterface::FREEVOLVERTEX] =
subtract(types_ents[TetGenInterface::FREEVOLVERTEX],types_ents[TetGenInterface::FREESEGVERTEX]);
types_ents[TetGenInterface::FREEVOLVERTEX] =
subtract(types_ents[TetGenInterface::FREEVOLVERTEX],types_ents[TetGenInterface::FREEFACETVERTEX]);

Definition at line 197 of file TetGenInterface.cpp.

                                                                           {
  MoFEMFunctionBegin;
 
  Interface &m_field = cOre;
  //
  // ErrorCode rval;
  in.pointparamlist = new tetgenio::pointparam[in.numberofpoints];
  // std::vector<bool> points_is_set(in.numberofpoints,false);
  std::map<int, Range>::iterator mit = type_ents.begin();
  for (; mit != type_ents.end(); mit++) {
    if (mit->first < 0 && mit->first > 3) {
      SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
              "Wrong TetGen point type");
    }
    Range::iterator it = mit->second.begin();
    for (; it != mit->second.end(); it++) {
      moabTetGen_Map::iterator miit = moab_tetgen_map.find(*it);
      if (miit == moab_tetgen_map.end()) {
        SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
                "Data inconsistency between TetGen and MoAB");
        continue;
      }
      int id = miit->second >> MB_TYPE_WIDTH;
      in.pointparamlist[id].uv[0] = 0;
      in.pointparamlist[id].uv[1] = 0;
      in.pointparamlist[id].type = mit->first;
      in.pointparamlist[id].tag = m_field.get_moab().id_from_handle(*it) + 1;
      // points_is_set[id] = true;
    }
  }
 
  // // Check only if tag and type is set to all points
  // for(
  //   std::vector<bool>::iterator bit = points_is_set.begin();
  //   bit!=points_is_set.end();bit++
  // ) {
  //   if(!*bit) {
  //     SETERRQ(PETSC_COMM_SELF,MOFEM_DATA_INCONSISTENCY,"Point type for
  //     TetGen is not set");
  //   }
  // }
 
  MoFEMFunctionReturn(0);
}

setRegionData()

MoFEMErrorCode MoFEM::TetGenInterface::setRegionData	(	std::vector< std::pair< EntityHandle, int > > &	regions,
		tetgenio &	in,
		Tag	th = NULL
	)

set region data to tetrahedral

Definition at line 556 of file TetGenInterface.cpp.

                                                                          {
  MoFEMFunctionBeginHot;
  ErrorCode rval;
  Interface &m_field = cOre;
  in.numberofregions = regions.size();
  in.regionlist = new double[5 * in.numberofregions];
  int kk = 0;
  std::vector<std::pair<EntityHandle, int>>::iterator it = regions.begin();
  for (int ii = 0; it != regions.end(); it++, ii++) {
    double coords[3];
    switch (type_from_handle(it->first)) {
    case MBVERTEX: {
      if (th) {
        rval = m_field.get_moab().tag_get_data(th, &it->first, 1, coords);
        CHKERRQ_MOAB(rval);
      } else {
        rval = m_field.get_moab().get_coords(&it->first, 1, coords);
        CHKERRQ_MOAB(rval);
      }
    } break;
    case MBTET: {
      int num_nodes;
      const EntityHandle *conn;
      rval =
          m_field.get_moab().get_connectivity(it->first, conn, num_nodes, true);
      CHKERRQ_MOAB(rval);
      double _coords[12];
      if (th) {
        rval = m_field.get_moab().tag_get_data(th, conn, num_nodes, _coords);
        CHKERRQ_MOAB(rval);
      } else {
        rval = m_field.get_moab().get_coords(conn, num_nodes, _coords);
        CHKERRQ_MOAB(rval);
      }
      coords[0] = (_coords[0] + _coords[3] + _coords[6] + _coords[9]) / 4.;
      coords[1] = (_coords[1] + _coords[4] + _coords[7] + _coords[10]) / 4.;
      coords[2] = (_coords[2] + _coords[5] + _coords[8] + _coords[11]) / 4.;
    } break;
    default:
      SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
    }
    for (int nn = 0; nn < 3; nn++) {
      in.regionlist[kk++] = coords[nn];
    }
    in.regionlist[kk++] = it->second;
    in.regionlist[kk++] = it->second;
  }
  MoFEMFunctionReturnHot(0);
}

tetRahedralize()

MoFEMErrorCode MoFEM::TetGenInterface::tetRahedralize	(	char	switches[],
		tetgenio &	in,
		tetgenio &	out
	)

run tetgen

Definition at line 646 of file TetGenInterface.cpp.

                                                              {
  MoFEMFunctionBeginHot;
  tetgenbehavior a;
  a.parse_commandline(switches);
  tetrahedralize(&a, &in, &out);
  MoFEMFunctionReturnHot(0);
}

Member Data Documentation

cOre

MoFEM::Core& MoFEM::TetGenInterface::cOre

Definition at line 28 of file TetGenInterface.hpp.

---

The documentation for this struct was generated from the following files:

• src/interfaces/TetGenInterface.hpp
• src/interfaces/impl/TetGenInterface.cpp