v0.9.1
Public Types | Public Member Functions | Public Attributes | List of all members
MoFEM::PrismsFromSurfaceInterface Struct Reference

merge node from two bit levels More...

#include <src/interfaces/PrismsFromSurfaceInterface.hpp>

Inheritance diagram for MoFEM::PrismsFromSurfaceInterface:
[legend]
Collaboration diagram for MoFEM::PrismsFromSurfaceInterface:
[legend]

Public Types

enum  SwapType { NO_SWAP = 0, SWAP_TRI_NODE_ORDER = 1, SWAP_TOP_AND_BOT_TRI = 2 }
 List of types of node swapping performed on a created prism Node swapping is required to satisfy the canonical ordering for the prism in case the surface normal is pointing inwards rather than outwards Currently supported options (assuming canonical ordering of nodes 0-5): NO_SWAP : node swapping is not performed SWAP_TRI_NODE_ORDER : swap the order of nodes on prism's triangle faces (1 <-> 2, 4 <-> 5) SWAP_TOP_AND_BOT_TRI : swap nodes between the top and the bottom triangles (0 <-> 3, 1 <-> 4, 2 <-> 5) More...
 

Public Member Functions

MoFEMErrorCode query_interface (const MOFEMuuid &uuid, UnknownInterface **iface) const
 
 PrismsFromSurfaceInterface (const MoFEM::Core &core)
 
MoFEMErrorCode createPrisms (const Range &ents, const SwapType swap_type, Range &prisms, int verb=-1)
 Make prisms from triangles. More...
 
DEPRECATED MoFEMErrorCode createPrisms (const Range &ents, Range &prisms, int verb=-1)
 
MoFEMErrorCode seedPrismsEntities (Range &prisms, const BitRefLevel &bit, int verb=-1)
 Seed prism entities by bit level. More...
 
MoFEMErrorCode createPrismsFromPrisms (const Range &prisms, bool from_down, Range &out_prisms, int verb=-1)
 Make prisms by extruding top or bottom prisms. More...
 
MoFEMErrorCode setThickness (const Range &prisms, const double director3[], const double director4[])
 
MoFEMErrorCode setNormalThickness (const Range &prisms, double thickness3, double thickness4)
 
MoFEMErrorCode updateMeshestByEdgeBlock (const Range &prisms)
 Add quads to bockset. More...
 
MoFEMErrorCode updateMeshestByTriBlock (const Range &prisms)
 Add prism to bockset. More...
 
- Public Member Functions inherited from MoFEM::UnknownInterface
template<class IFACE >
MoFEMErrorCode registerInterface (const MOFEMuuid &uuid, bool error_if_registration_failed=true)
 Register interface. More...
 
template<class IFACE , bool VERIFY = false>
MoFEMErrorCode getInterface (const MOFEMuuid &uuid, IFACE *&iface) const
 Get interface by uuid and return reference to pointer of interface. More...
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE *&iface) const
 Get interface refernce 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
 
virtual MoFEMErrorCode getLibVersion (Version &version) const
 Get library version. More...
 
virtual const MoFEMErrorCode getFileVersion (moab::Interface &moab, Version &version) const
 Get database major version. More...
 
virtual MoFEMErrorCode getInterfaceVersion (Version &version) const
 Get database major version. More...
 
template<>
MoFEMErrorCode getInterface (const MOFEMuuid &uuid, UnknownInterface *&iface) const
 

Public Attributes

MoFEM::CorecOre
 
std::map< EntityHandle, EntityHandlecreatedVertices
 

Additional Inherited Members

- Protected Member Functions inherited from MoFEM::UnknownInterface
boost::typeindex::type_index getClassIdx (const MOFEMuuid &uid) const
 Get type name for interface Id. More...
 
MOFEMuuid getUId (const boost::typeindex::type_index &class_idx) const
 Get interface Id for class name. More...
 

Detailed Description

merge node from two bit levels

Examples
prism_elements_from_surface.cpp, and prism_polynomial_approximation.cpp.

Definition at line 27 of file PrismsFromSurfaceInterface.hpp.

Member Enumeration Documentation

◆ SwapType

List of types of node swapping performed on a created prism Node swapping is required to satisfy the canonical ordering for the prism in case the surface normal is pointing inwards rather than outwards Currently supported options (assuming canonical ordering of nodes 0-5): NO_SWAP : node swapping is not performed SWAP_TRI_NODE_ORDER : swap the order of nodes on prism's triangle faces (1 <-> 2, 4 <-> 5) SWAP_TOP_AND_BOT_TRI : swap nodes between the top and the bottom triangles (0 <-> 3, 1 <-> 4, 2 <-> 5)

Enumerator
NO_SWAP 
SWAP_TRI_NODE_ORDER 
SWAP_TOP_AND_BOT_TRI 

Definition at line 49 of file PrismsFromSurfaceInterface.hpp.

Constructor & Destructor Documentation

◆ PrismsFromSurfaceInterface()

MoFEM::PrismsFromSurfaceInterface::PrismsFromSurfaceInterface ( const MoFEM::Core core)

Definition at line 33 of file PrismsFromSurfaceInterface.hpp.

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

Member Function Documentation

◆ createPrisms() [1/2]

MoFEMErrorCode MoFEM::PrismsFromSurfaceInterface::createPrisms ( const Range &  ents,
const SwapType  swap_type,
Range &  prisms,
int  verb = -1 
)

Make prisms from triangles.

Parameters
entsRange of triangles
swap_typeDefines how the nodes of the created prism are swapped (required for canonical ordering if the surface normal is pointing inwards)
prismsReturned range of prisms
verbVerbosity level
Returns
Error code
Examples
prism_elements_from_surface.cpp.

Definition at line 40 of file PrismsFromSurfaceInterface.cpp.

41  {
43  Interface &m_field = cOre;
44  Range tris = ents.subset_by_type(MBTRI);
45  for (Range::iterator tit = tris.begin(); tit != tris.end(); tit++) {
46  const EntityHandle *conn;
47  int number_nodes = 0;
48  CHKERR m_field.get_moab().get_connectivity(*tit, conn, number_nodes, false);
49  double coords[3 * number_nodes];
50  CHKERR m_field.get_moab().get_coords(conn, number_nodes, coords);
51  EntityHandle prism_nodes[6];
52  for (int nn = 0; nn < 3; nn++) {
53  prism_nodes[nn] = conn[nn];
54  if (createdVertices.find(conn[nn]) != createdVertices.end()) {
55  prism_nodes[3 + nn] = createdVertices[prism_nodes[nn]];
56  } else {
57  CHKERR m_field.get_moab().create_vertex(&coords[3 * nn],
58  prism_nodes[3 + nn]);
59  createdVertices[conn[nn]] = prism_nodes[3 + nn];
60  CHKERR m_field.get_moab().tag_set_data(cOre.get_th_RefParentHandle(),
61  &prism_nodes[3 + nn], 1,
62  &prism_nodes[nn]);
63  }
64  }
65  switch (swap_type) {
67  std::swap(prism_nodes[1], prism_nodes[2]);
68  std::swap(prism_nodes[4], prism_nodes[5]);
69  break;
71  std::swap(prism_nodes[0], prism_nodes[3]);
72  std::swap(prism_nodes[1], prism_nodes[4]);
73  std::swap(prism_nodes[2], prism_nodes[5]);
74  break;
75  case NO_SWAP:
76  default:
77  break;
78  }
79  EntityHandle prism;
80  CHKERR m_field.get_moab().create_element(MBPRISM, prism_nodes, 6, prism);
81  Range edges;
82  CHKERR m_field.get_moab().get_adjacencies(&prism, 1, 1, true, edges,
83  moab::Interface::UNION);
84  Range faces;
85  CHKERR m_field.get_moab().get_adjacencies(&prism, 1, 2, true, faces,
86  moab::Interface::UNION);
87  prisms.insert(prism);
88  for (int ee = 0; ee <= 2; ee++) {
89  EntityHandle e1;
90  CHKERR m_field.get_moab().side_element(prism, 1, ee, e1);
91  EntityHandle e2;
92  CHKERR m_field.get_moab().side_element(prism, 1, ee + 6, e2);
93  CHKERR m_field.get_moab().tag_set_data(cOre.get_th_RefParentHandle(), &e2,
94  1, &e1);
95  }
96  EntityHandle f3, f4;
97  {
98  CHKERR m_field.get_moab().side_element(prism, 2, 3, f3);
99  CHKERR m_field.get_moab().side_element(prism, 2, 4, f4);
100  CHKERR m_field.get_moab().tag_set_data(cOre.get_th_RefParentHandle(), &f4,
101  1, &f3);
102  }
103  if (number_nodes > 3) {
104  EntityHandle meshset;
105  CHKERR m_field.get_moab().create_meshset(MESHSET_SET, meshset);
106  CHKERR m_field.get_moab().add_entities(meshset, &f4, 1);
107  for (int ee = 0; ee <= 2; ee++) {
108  EntityHandle e2;
109  CHKERR m_field.get_moab().side_element(prism, 1, ee + 6, e2);
110  CHKERR m_field.get_moab().add_entities(meshset, &e2, 1);
111  }
112  CHKERR m_field.get_moab().convert_entities(meshset, true, false, false);
113  CHKERR m_field.get_moab().delete_entities(&meshset, 1);
114  const EntityHandle *conn_f4;
115  int number_nodes_f4 = 0;
116  CHKERR m_field.get_moab().get_connectivity(f4, conn_f4, number_nodes_f4,
117  false);
118  if (number_nodes_f4 != number_nodes) {
119  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
120  "data inconsistency");
121  }
122  CHKERR m_field.get_moab().set_coords(&conn_f4[3], 3, &coords[9]);
123  }
124  }
126 }
Tag get_th_RefParentHandle() const
Definition: Core.hpp:150
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
std::map< EntityHandle, EntityHandle > createdVertices
#define CHKERR
Inline error check.
Definition: definitions.h:602
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1879
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ createPrisms() [2/2]

DEPRECATED MoFEMErrorCode MoFEM::PrismsFromSurfaceInterface::createPrisms ( const Range &  ents,
Range &  prisms,
int  verb = -1 
)
Deprecated:
Use the function with the same name and a parameter swap_type, permitting to swap the order of each triangle's nodes or alternatively swap nodes between top and bottom triangles, which is required for the canonical ordering if surface normal is pointing inwards

Definition at line 33 of file PrismsFromSurfaceInterface.cpp.

34  {
36  CHKERR createPrisms(ents, NO_SWAP, prisms, verb);
38 }
MoFEMErrorCode createPrisms(const Range &ents, const SwapType swap_type, Range &prisms, int verb=-1)
Make prisms from triangles.
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define CHKERR
Inline error check.
Definition: definitions.h:602
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ createPrismsFromPrisms()

MoFEMErrorCode MoFEM::PrismsFromSurfaceInterface::createPrismsFromPrisms ( const Range &  prisms,
bool  from_down,
Range &  out_prisms,
int  verb = -1 
)

Make prisms by extruding top or bottom prisms.

Parameters
prismsInput prisms
from_downUse top or down face, if true from f3
out_prismsReturned prisms entities
verbVerbosity level
Returns
Error code
Examples
prism_elements_from_surface.cpp.

Definition at line 161 of file PrismsFromSurfaceInterface.cpp.

162  {
164  Interface &m_field = cOre;
165  Range tris;
166  for (Range::iterator pit = prisms.begin(); pit != prisms.end(); pit++) {
167  EntityHandle face;
168  if (from_down) {
169  CHKERR m_field.get_moab().side_element(*pit, 2, 3, face);
170  } else {
171  CHKERR m_field.get_moab().side_element(*pit, 2, 4, face);
172  }
173  tris.insert(face);
174  }
175  CHKERR createPrisms(tris, out_prisms, verb);
177 }
MoFEMErrorCode createPrisms(const Range &ents, const SwapType swap_type, Range &prisms, int verb=-1)
Make prisms from triangles.
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define CHKERR
Inline error check.
Definition: definitions.h:602
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1879
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ query_interface()

MoFEMErrorCode MoFEM::PrismsFromSurfaceInterface::query_interface ( const MOFEMuuid uuid,
UnknownInterface **  iface 
) const
virtual

Implements MoFEM::UnknownInterface.

Definition at line 21 of file PrismsFromSurfaceInterface.cpp.

22  {
24  *iface = NULL;
25  if (uuid == IDD_MOFEMPrismsFromSurface) {
26  *iface = const_cast<PrismsFromSurfaceInterface *>(this);
28  }
29  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "unknown interface");
31 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
static const MOFEMuuid IDD_MOFEMPrismsFromSurface
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ seedPrismsEntities()

MoFEMErrorCode MoFEM::PrismsFromSurfaceInterface::seedPrismsEntities ( Range &  prisms,
const BitRefLevel bit,
int  verb = -1 
)

Seed prism entities by bit level.

Parameters
prismsRange of entities
bitBitRefLevel
verbVerbosity level
Returns
Error code
Examples
prism_elements_from_surface.cpp, and prism_polynomial_approximation.cpp.

Definition at line 128 of file PrismsFromSurfaceInterface.cpp.

129  {
130  Interface &m_field = cOre;
131  auto ref_ents_ptr = m_field.get_ref_ents();
133  MPI_Comm comm = m_field.get_comm();
134  RefEntity_multiIndex *refined_entities_ptr;
135  refined_entities_ptr =
136  const_cast<RefEntity_multiIndex *>(ref_ents_ptr);
137  if (!prisms.empty()) {
138  int dim = m_field.get_moab().dimension_from_handle(prisms[0]);
139  for (int dd = 0; dd <= dim; dd++) {
140  Range ents;
141  CHKERR m_field.get_moab().get_adjacencies(prisms, dd, true, ents,
142  moab::Interface::UNION);
143  Range::iterator eit = ents.begin();
144  for (; eit != ents.end(); eit++) {
145  std::pair<RefEntity_multiIndex::iterator, bool> p_ent =
146  refined_entities_ptr->insert(boost::shared_ptr<RefEntity>(
147  new RefEntity(m_field.get_basic_entity_data_ptr(), *eit)));
148  *(const_cast<RefEntity *>(p_ent.first->get())->getBitRefLevelPtr()) |=
149  bit;
150  if (verb >= VERY_VERBOSE) {
151  std::ostringstream ss;
152  ss << *(p_ent.first);
153  PetscSynchronizedPrintf(comm, "%s\n", ss.str().c_str());
154  }
155  }
156  }
157  }
159 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
const int dim
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T * > &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: ddTensor0.hpp:33
MoFEMErrorCode get_ref_ents(const RefEntity_multiIndex **refined_entities_ptr) const
Get ref entities multi-index from database.
Definition: Core.cpp:746
#define CHKERR
Inline error check.
Definition: definitions.h:602
multi_index_container< boost::shared_ptr< RefEntity >, indexed_by< ordered_unique< tag< Ent_mi_tag >, member< RefEntity::BasicEntity, EntityHandle, &RefEntity::ent > >, ordered_non_unique< tag< Ent_Ent_mi_tag >, const_mem_fun< RefEntity, EntityHandle, &RefEntity::getParentEnt > >, ordered_non_unique< tag< EntType_mi_tag >, const_mem_fun< RefEntity::BasicEntity, EntityType, &RefEntity::getEntType > >, ordered_non_unique< tag< ParentEntType_mi_tag >, const_mem_fun< RefEntity, EntityType, &RefEntity::getParentEntType > >, ordered_non_unique< tag< Composite_EntType_and_ParentEntType_mi_tag >, composite_key< RefEntity, const_mem_fun< RefEntity::BasicEntity, EntityType, &RefEntity::getEntType >, const_mem_fun< RefEntity, EntityType, &RefEntity::getParentEntType > > >, ordered_non_unique< tag< Composite_ParentEnt_And_EntType_mi_tag >, composite_key< RefEntity, const_mem_fun< RefEntity::BasicEntity, EntityType, &RefEntity::getEntType >, const_mem_fun< RefEntity, EntityHandle, &RefEntity::getParentEnt > > > > > RefEntity_multiIndex
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1879
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ setNormalThickness()

MoFEMErrorCode MoFEM::PrismsFromSurfaceInterface::setNormalThickness ( const Range &  prisms,
double  thickness3,
double  thickness4 
)

Set normal thickness

Parameters
prismsRange of prisms
thicknessnormal thickness
Returns

Definition at line 213 of file PrismsFromSurfaceInterface.cpp.

214  {
215  Interface &m_field = cOre;
217 
218  auto add_normal = [&](std::map<EntityHandle, std::array<double, 3>> &nodes,
219  EntityHandle face) {
221  const EntityHandle *conn;
222  int number_nodes;
223  CHKERR m_field.get_moab().get_connectivity(face, conn, number_nodes, false);
224  std::array<double, 9> coords;
225  CHKERR m_field.get_moab().get_coords(conn, number_nodes, coords.data());
226  std::array<double, 3> normal;
227  CHKERR Tools::getTriNormal(coords.data(), normal.data());
228  double a = sqrt(normal[0] * normal[0] + normal[1] * normal[1] +
229  normal[2] * normal[2]);
230  for (auto d : {0, 1, 2})
231  normal[d] /= a;
232  for (auto n : {0, 1, 2}) {
233  try {
234  for (auto d : {0, 1, 2})
235  nodes.at(conn[n])[d] += normal[d];
236  } catch (...) {
237  nodes.insert(
238  std::pair<EntityHandle, std::array<double, 3>>(conn[n], normal));
239  }
240  }
242  };
243 
244  auto apply_map = [&](auto &nodes, double t) {
246  for (auto &m : nodes) {
247  std::array<double, 3> coords;
248  CHKERR m_field.get_moab().get_coords(&m.first, 1, coords.data());
249  auto &normal = m.second;
250  double a = sqrt(normal[0] * normal[0] + normal[1] * normal[1] +
251  normal[2] * normal[2]);
252  for (auto d : {0, 1, 2})
253  coords[d] += (normal[d] / a) * t;
254  CHKERR m_field.get_moab().set_coords(&m.first, 1, coords.data());
255  }
257  };
258 
259  map<EntityHandle, std::array<double, 3>> nodes_f3, nodes_f4;
260  for (Range::iterator pit = prisms.begin(); pit != prisms.end(); pit++) {
261  for (int ff = 3; ff <= 4; ff++) {
262  EntityHandle face;
263  CHKERR m_field.get_moab().side_element(*pit, 2, ff, face);
264  if (ff == 3)
265  CHKERR add_normal(nodes_f3, face);
266  else
267  CHKERR add_normal(nodes_f4, face);
268  }
269  }
270 
271  CHKERR apply_map(nodes_f3, thickness3);
272  CHKERR apply_map(nodes_f4, thickness4);
273 
275 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
FTensor::Index< 'n', 2 > n
Definition: PlasticOps.hpp:68
FTensor::Tensor1< double, 2 > normal(FTensor::Tensor1< T1, 3 > &t_coords, FTensor::Tensor1< T2, 2 > &t_disp)
Definition: ContactOps.hpp:128
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514
const Tensor1_Expr< const dTensor0< T, Dim, i >, typename promote< T, double >::V, Dim, i > d(const Tensor0< T * > &a, const Index< i, Dim > index, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: dTensor0.hpp:27
FTensor::Index< 'm', 2 > m
Definition: PlasticOps.hpp:67
static MoFEMErrorCode getTriNormal(const double *coords, double *normal)
Get the Tri Normal objectGet triangle normal.
Definition: Tools.cpp:258
#define CHKERR
Inline error check.
Definition: definitions.h:602
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1879
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ setThickness()

MoFEMErrorCode MoFEM::PrismsFromSurfaceInterface::setThickness ( const Range &  prisms,
const double  director3[],
const double  director4[] 
)

Set uniform thickness

Parameters
prismsRange of prisms
director3Displacement of face 3
director4Displacement of face 4
Returns
Examples
prism_elements_from_surface.cpp.

Definition at line 179 of file PrismsFromSurfaceInterface.cpp.

180  {
182  Interface &m_field = cOre;
183  Range nodes_f3, nodes_f4;
184  for (Range::iterator pit = prisms.begin(); pit != prisms.end(); pit++) {
185  for (int ff = 3; ff <= 4; ff++) {
186  EntityHandle face;
187  CHKERR m_field.get_moab().side_element(*pit, 2, ff, face);
188  const EntityHandle *conn;
189  int number_nodes = 0;
190  CHKERR m_field.get_moab().get_connectivity(face, conn, number_nodes,
191  false);
192  if (ff == 3) {
193  nodes_f3.insert(&conn[0], &conn[number_nodes]);
194  } else {
195  nodes_f4.insert(&conn[0], &conn[number_nodes]);
196  }
197  }
198  }
199  double coords[3];
200  for (Range::iterator nit = nodes_f3.begin(); nit != nodes_f3.end(); nit++) {
201  CHKERR m_field.get_moab().get_coords(&*nit, 1, coords);
202  cblas_daxpy(3, 1, director3, 1, coords, 1);
203  CHKERR m_field.get_moab().set_coords(&*nit, 1, coords);
204  }
205  for (Range::iterator nit = nodes_f4.begin(); nit != nodes_f4.end(); nit++) {
206  CHKERR m_field.get_moab().get_coords(&*nit, 1, coords);
207  cblas_daxpy(3, 1, director4, 1, coords, 1);
208  CHKERR m_field.get_moab().set_coords(&*nit, 1, coords);
209  }
211 }
void cblas_daxpy(const int N, const double alpha, const double *X, const int incX, double *Y, const int incY)
Definition: cblas_daxpy.c:11
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define CHKERR
Inline error check.
Definition: definitions.h:602
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1879
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ updateMeshestByEdgeBlock()

MoFEMErrorCode MoFEM::PrismsFromSurfaceInterface::updateMeshestByEdgeBlock ( const Range &  prisms)

Add quads to bockset.

If quad is adjacent to extruded edge, is added to given blockset

Parameters
prisms
Returns
MoFEMErrorCode

Definition at line 278 of file PrismsFromSurfaceInterface.cpp.

278  {
279  Interface &m_field = cOre;
281  Range prisms_edges;
282  CHKERR m_field.get_moab().get_adjacencies(prisms, 1, true, prisms_edges,
283  moab::Interface::UNION);
284  Range prisms_faces;
285  CHKERR m_field.get_moab().get_adjacencies(prisms, 2, true, prisms_faces,
286  moab::Interface::UNION);
287  for (_IT_CUBITMESHSETS_FOR_LOOP_(m_field, it)) {
288  Range edges;
289  CHKERR m_field.get_moab().get_entities_by_type(it->meshset, MBEDGE, edges,
290  true);
291  edges = intersect(edges, prisms_edges);
292  if (!edges.empty()) {
293  Range edges_faces;
294  CHKERR m_field.get_moab().get_adjacencies(edges, 2, false, edges_faces,
295  moab::Interface::UNION);
296  edges_faces = intersect(edges_faces, prisms_faces.subset_by_type(MBQUAD));
297  EntityHandle meshset = it->getMeshset();
298  CHKERR m_field.get_moab().add_entities(meshset, edges_faces);
299  }
300  }
302 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define CHKERR
Inline error check.
Definition: definitions.h:602
#define _IT_CUBITMESHSETS_FOR_LOOP_(MESHSET_MANAGER, IT)
Iterator that loops over all the Cubit MeshSets in a moFEM field.
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1879
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ updateMeshestByTriBlock()

MoFEMErrorCode MoFEM::PrismsFromSurfaceInterface::updateMeshestByTriBlock ( const Range &  prisms)

Add prism to bockset.

If quad is adjacent to extruded triangle, is added to given blockset

Parameters
prisms
Returns
MoFEMErrorCode

Definition at line 305 of file PrismsFromSurfaceInterface.cpp.

305  {
306  Interface &m_field = cOre;
308  Range prisms_tris;
309  CHKERR m_field.get_moab().get_adjacencies(prisms, 2, true, prisms_tris,
310  moab::Interface::UNION);
311  prisms_tris = prisms_tris.subset_by_type(MBTRI);
312  for (_IT_CUBITMESHSETS_FOR_LOOP_(m_field, it)) {
313  Range tris;
314  CHKERR m_field.get_moab().get_entities_by_type(it->meshset, MBTRI, tris,
315  true);
316  tris = intersect(tris, prisms_tris);
317  if (!tris.empty()) {
318  Range tris_ents;
319  CHKERR m_field.get_moab().get_adjacencies(tris, 3, false, tris_ents,
320  moab::Interface::UNION);
321  tris_ents = intersect(tris_ents, prisms);
322  EntityHandle meshset = it->getMeshset();
323  CHKERR m_field.get_moab().add_entities(meshset, tris_ents);
324  }
325  }
327 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define CHKERR
Inline error check.
Definition: definitions.h:602
#define _IT_CUBITMESHSETS_FOR_LOOP_(MESHSET_MANAGER, IT)
Iterator that loops over all the Cubit MeshSets in a moFEM field.
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1879
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

Member Data Documentation

◆ cOre

MoFEM::Core& MoFEM::PrismsFromSurfaceInterface::cOre

Definition at line 32 of file PrismsFromSurfaceInterface.hpp.

◆ createdVertices

std::map<EntityHandle, EntityHandle> MoFEM::PrismsFromSurfaceInterface::createdVertices

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