v0.8.13
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 Member Functions

MoFEMErrorCode query_interface (const MOFEMuuid &uuid, UnknownInterface **iface) const
 
 PrismsFromSurfaceInterface (const MoFEM::Core &core)
 
MoFEMErrorCode createPrisms (const Range &ents, Range &prisms, int verb=-1)
 Make prisms from triangles. More...
 
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[])
 
- 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 ()
 
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

Definition at line 26 of file PrismsFromSurfaceInterface.hpp.

Constructor & Destructor Documentation

◆ PrismsFromSurfaceInterface()

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

Definition at line 31 of file PrismsFromSurfaceInterface.hpp.

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

Member Function Documentation

◆ createPrisms()

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

Make prisms from triangles.

Parameters
entsRange of triangles
prismsReturned range of prisms
verbVerbosity level
Returns
Error code

Definition at line 31 of file PrismsFromSurfaceInterface.cpp.

31  {
33 
34  Interface& m_field = cOre;
35  Range tris = ents.subset_by_type(MBTRI);
36  for(Range::iterator tit = tris.begin();tit!=tris.end();tit++) {
37  const EntityHandle* conn;
38  int number_nodes = 0;
39  rval = m_field.get_moab().get_connectivity(*tit,conn,number_nodes,false); CHKERRQ_MOAB(rval);
40  double coords[3*number_nodes];
41  rval = m_field.get_moab().get_coords(conn,number_nodes,coords); CHKERRQ_MOAB(rval);
42  EntityHandle prism_nodes[6];
43  for(int nn = 0;nn<3;nn++) {
44  prism_nodes[nn] = conn[nn];
45  if(createdVertices.find(conn[nn])!=createdVertices.end()) {
46  prism_nodes[3+nn] = createdVertices[prism_nodes[nn]];
47  } else {
48  rval = m_field.get_moab().create_vertex(&coords[3*nn],prism_nodes[3+nn]); CHKERRQ_MOAB(rval);
49  createdVertices[conn[nn]] = prism_nodes[3+nn];
50  rval = m_field.get_moab().tag_set_data(
51  cOre.get_th_RefParentHandle(),&prism_nodes[3+nn],1,&prism_nodes[nn]
52  ); CHKERRQ_MOAB(rval);
53  }
54  }
55  EntityHandle prism;
56  rval = m_field.get_moab().create_element(MBPRISM,prism_nodes,6,prism); CHKERRG(rval);
57  Range edges;
58  rval = m_field.get_moab().get_adjacencies(&prism,1,1,true,edges,moab::Interface::UNION); CHKERRQ_MOAB(rval);
59  Range faces;
60  rval = m_field.get_moab().get_adjacencies(&prism,1,2,true,faces,moab::Interface::UNION); CHKERRQ_MOAB(rval);
61  prisms.insert(prism);
62  for(int ee = 0;ee<=2;ee++) {
63  EntityHandle e1;
64  rval = m_field.get_moab().side_element(prism,1,ee,e1);
65  EntityHandle e2;
66  rval = m_field.get_moab().side_element(prism,1,ee+6,e2); CHKERRQ_MOAB(rval);
67  rval = m_field.get_moab().tag_set_data(cOre.get_th_RefParentHandle(),&e2,1,&e1); CHKERRQ_MOAB(rval);
68  }
69  EntityHandle f3,f4;
70  {
71  rval = m_field.get_moab().side_element(prism,2,3,f3);
72  rval = m_field.get_moab().side_element(prism,2,4,f4);
73  rval = m_field.get_moab().tag_set_data(cOre.get_th_RefParentHandle(),&f4,1,&f3); CHKERRQ_MOAB(rval);
74  }
75  if(number_nodes>3) {
76  EntityHandle meshset;
77  rval = m_field.get_moab().create_meshset(MESHSET_SET|MESHSET_TRACK_OWNER,meshset); CHKERRQ_MOAB(rval);
78  rval = m_field.get_moab().add_entities(meshset,&f4,1); CHKERRQ_MOAB(rval);
79  for(int ee = 0;ee<=2;ee++) {
80  EntityHandle e2;
81  rval = m_field.get_moab().side_element(prism,1,ee+6,e2); CHKERRQ_MOAB(rval);
82  rval = m_field.get_moab().add_entities(meshset,&e2,1); CHKERRQ_MOAB(rval);
83  }
84  rval = m_field.get_moab().convert_entities(meshset,true,false,false); CHKERRQ_MOAB(rval);
85  rval = m_field.get_moab().delete_entities(&meshset,1); CHKERRQ_MOAB(rval);
86  const EntityHandle* conn_f4;
87  int number_nodes_f4 = 0;
88  rval = m_field.get_moab().get_connectivity(f4,conn_f4,number_nodes_f4,false); CHKERRQ_MOAB(rval);
89  if(number_nodes_f4 != number_nodes) {
90  SETERRQ(PETSC_COMM_SELF,MOFEM_DATA_INCONSISTENCY,"data inconsistency");
91  }
92  rval = m_field.get_moab().set_coords(&conn_f4[3],3,&coords[9]); CHKERRQ_MOAB(rval);
93  }
94  }
96 }
#define CHKERRQ_MOAB(a)
check error code of MoAB function
Definition: definitions.h:533
Tag get_th_RefParentHandle() const
Definition: Core.hpp:150
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:519
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:562
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:526
static MoFEMErrorCodeGeneric< moab::ErrorCode > rval
Definition: Common.hpp:78
std::map< EntityHandle, EntityHandle > createdVertices

◆ 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

Definition at line 132 of file PrismsFromSurfaceInterface.cpp.

132  {
134  Interface& m_field = cOre;
135  Range tris;
136  for(Range::iterator pit = prisms.begin();pit!=prisms.end();pit++) {
137  EntityHandle face;
138  if(from_down) {
139  rval = m_field.get_moab().side_element(*pit,2,3,face); CHKERRQ_MOAB(rval);
140  } else {
141  rval = m_field.get_moab().side_element(*pit,2,4,face); CHKERRQ_MOAB(rval);
142  }
143  tris.insert(face);
144  }
145  ierr = createPrisms(tris,out_prisms,verb); CHKERRG(ierr);
147 }
#define CHKERRQ_MOAB(a)
check error code of MoAB function
Definition: definitions.h:533
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:519
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:562
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:526
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Common.hpp:80
static MoFEMErrorCodeGeneric< moab::ErrorCode > rval
Definition: Common.hpp:78
MoFEMErrorCode createPrisms(const Range &ents, Range &prisms, int verb=-1)
Make prisms from triangles.

◆ query_interface()

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

Implements MoFEM::UnknownInterface.

Definition at line 20 of file PrismsFromSurfaceInterface.cpp.

20  {
22  *iface = NULL;
23  if(uuid == IDD_MOFEMPrismsFromSurface) {
24  *iface = const_cast<PrismsFromSurfaceInterface*>(this);
26  }
27  SETERRQ(PETSC_COMM_SELF,MOFEM_DATA_INCONSISTENCY,"unknown interface");
29 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:519
static const MOFEMuuid IDD_MOFEMPrismsFromSurface
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:526
PrismsFromSurfaceInterface(const MoFEM::Core &core)

◆ 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

Definition at line 98 of file PrismsFromSurfaceInterface.cpp.

99  {
101  Interface &m_field = cOre;
102  const RefEntity_multiIndex *const_refined_entities_ptr;
103  CHKERR m_field.get_ref_ents(&const_refined_entities_ptr);
104  MPI_Comm comm = m_field.get_comm();
105  RefEntity_multiIndex *refined_entities_ptr;
106  refined_entities_ptr =
107  const_cast<RefEntity_multiIndex *>(const_refined_entities_ptr);
108  if (!prisms.empty()) {
109  int dim = m_field.get_moab().dimension_from_handle(prisms[0]);
110  for (int dd = 0; dd <= dim; dd++) {
111  Range ents;
112  CHKERR m_field.get_moab().get_adjacencies(prisms, dd, true, ents,
113  moab::Interface::UNION);
114  Range::iterator eit = ents.begin();
115  for (; eit != ents.end(); eit++) {
116  std::pair<RefEntity_multiIndex::iterator, bool> p_ent =
117  refined_entities_ptr->insert(boost::shared_ptr<RefEntity>(
118  new RefEntity(m_field.get_basic_entity_data_ptr(), *eit)));
119  *(const_cast<RefEntity *>(p_ent.first->get())->getBitRefLevelPtr()) |=
120  bit;
121  if (verb >= VERY_VERBOSE) {
122  std::ostringstream ss;
123  ss << *(p_ent.first);
124  PetscSynchronizedPrintf(comm, "%s\n", ss.str().c_str());
125  }
126  }
127  }
128  }
130 }
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
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:495
#define CHKERR
Inline error check.
Definition: definitions.h:614
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
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:439

◆ 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

Definition at line 149 of file PrismsFromSurfaceInterface.cpp.

149  {
150 
152  Interface& m_field = cOre;
153  Range nodes_f3,nodes_f4;
154  for(Range::iterator pit = prisms.begin();pit!=prisms.end();pit++) {
155  for(int ff = 3;ff<=4;ff++) {
156  EntityHandle face;
157  rval = m_field.get_moab().side_element(*pit,2,ff,face); CHKERRQ_MOAB(rval);
158  const EntityHandle* conn;
159  int number_nodes = 0;
160  rval = m_field.get_moab().get_connectivity(face,conn,number_nodes,false); CHKERRQ_MOAB(rval);
161  if(ff == 3) {
162  nodes_f3.insert(&conn[0],&conn[number_nodes]);
163  } else {
164  nodes_f4.insert(&conn[0],&conn[number_nodes]);
165  }
166  }
167  }
168  double coords[3];
169  for(Range::iterator nit = nodes_f3.begin();nit!=nodes_f3.end();nit++) {
170  rval = m_field.get_moab().get_coords(&*nit,1,coords); CHKERRQ_MOAB(rval);
171  cblas_daxpy(3,1,director3,1,coords,1);
172  rval = m_field.get_moab().set_coords(&*nit,1,coords); CHKERRQ_MOAB(rval);
173  }
174  for(Range::iterator nit = nodes_f4.begin();nit!=nodes_f4.end();nit++) {
175  rval = m_field.get_moab().get_coords(&*nit,1,coords); CHKERRQ_MOAB(rval);
176  cblas_daxpy(3,1,director4,1,coords,1);
177  rval = m_field.get_moab().set_coords(&*nit,1,coords); CHKERRQ_MOAB(rval);
178  }
180 }
#define CHKERRQ_MOAB(a)
check error code of MoAB function
Definition: definitions.h:533
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 MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:519
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:526
static MoFEMErrorCodeGeneric< moab::ErrorCode > rval
Definition: Common.hpp:78

Member Data Documentation

◆ cOre

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

Definition at line 30 of file PrismsFromSurfaceInterface.hpp.

◆ createdVertices

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

Definition at line 34 of file PrismsFromSurfaceInterface.hpp.


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