v0.15.0
Loading...
Searching...
No Matches
VtkInterface Struct Reference
Collaboration diagram for VtkInterface:
[legend]

Public Member Functions

 VtkInterface (MoFEM::Interface &m_field, moab::Interface &moab)
 
MoFEMErrorCode readVtk ()
 

Private Member Functions

MoFEMErrorCode readMesh ()
 
MoFEMErrorCode getBoundaryConditions ()
 
MoFEMErrorCode setBoundaryConditions ()
 
MoFEMErrorCode getMaterialProperties ()
 
MoFEMErrorCode getBlockSetBCFromMesh (std::string)
 
MoFEMErrorCode writeOutput ()
 
MoFEMErrorCode checkResults ()
 

Private Attributes

MoFEM::InterfacemField
 
moab::Interface & mOab
 
SimplesimpleInterface
 

Detailed Description

Definition at line 11 of file read_vtk.cpp.

Constructor & Destructor Documentation

◆ VtkInterface()

VtkInterface::VtkInterface ( MoFEM::Interface & m_field,
moab::Interface & moab )
inline

Definition at line 13 of file read_vtk.cpp.

14 : mField(m_field), mOab(moab) {}
moab::Interface & mOab
Definition read_vtk.cpp:20
MoFEM::Interface & mField
Definition read_vtk.cpp:19

Member Function Documentation

◆ checkResults()

MoFEMErrorCode VtkInterface::checkResults ( )
private

Definition at line 227 of file read_vtk.cpp.

227 {
229 PetscBool atom_test = PETSC_FALSE;
230 PetscOptionsGetBool(PETSC_NULLPTR, "", "-atom_test", &atom_test, PETSC_NULLPTR);
231 if (atom_test) {
232 // check mesh has been read correctly
233 Range ents;
234 CHKERR mOab.get_entities_by_dimension(0, 3, ents);
235 if (ents.size() != 125) {
236 SETERRQ(
237 PETSC_COMM_SELF, MOFEM_ATOM_TEST_INVALID,
238 "atom test %d failed: Wrong number of elements %ld, should be 125",
239 atom_test, ents.size(), 125);
240 }
241 // check boundary conditions have been set correctly
242 MeshsetsManager *meshsets_manager_ptr;
243 CHKERR mField.getInterface(meshsets_manager_ptr);
244 // iterate meshsets
246 Range ents;
247 CHKERR mField.get_moab().get_entities_by_dimension(it->meshset, 2, ents,
248 true);
249 if (it->getMeshsetId() == 1) {
250 if (it->getName() != "FIX_ALL") {
251 SETERRQ(PETSC_COMM_SELF, MOFEM_ATOM_TEST_INVALID,
252 "atom test %d failed: Wrong name of meshset %d, should be "
253 "FIX_ALL is %s",
254 atom_test, it->getMeshsetId(), it->getName().c_str());
255 }
256 if (ents.size() != 25) {
257 SETERRQ(PETSC_COMM_SELF, MOFEM_ATOM_TEST_INVALID,
258 "atom test %d failed: Wrong number of entities in meshset "
259 "%d with %ld, should be 25",
260 atom_test, it->getMeshsetId(), ents.size());
261 }
262 } else if (it->getMeshsetId() == 2) {
263 if (it->getName() != "FIX_X") {
264 SETERRQ(PETSC_COMM_SELF, MOFEM_ATOM_TEST_INVALID,
265 "atom test %d failed: Wrong name of meshset %d, should be "
266 "FIX_X is %s",
267 atom_test, it->getMeshsetId(), it->getName().c_str());
268 }
269 if (ents.size() != 25) {
270 SETERRQ(PETSC_COMM_SELF, MOFEM_ATOM_TEST_INVALID,
271 "atom test %d failed: Wrong number of entities in meshset "
272 "%d with %ld, should be 25",
273 atom_test, it->getMeshsetId(), ents.size());
274 }
275 } else if (it->getMeshsetId() == 100) {
276 if (it->getName() != "ADOLCMAT") {
277 SETERRQ(PETSC_COMM_SELF, MOFEM_ATOM_TEST_INVALID,
278 "atom test %d failed: Wrong name of meshset %d, should be "
279 "ADOLCMAT is %s",
280 atom_test, it->getMeshsetId(), it->getName().c_str());
281 }
282 }
283 }
284 }
286}
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
@ BLOCKSET
@ MOFEM_ATOM_TEST_INVALID
Definition definitions.h:40
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
#define _IT_CUBITMESHSETS_BY_SET_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet having a particular BC meshset in a moFEM field.
PetscErrorCode PetscOptionsGetBool(PetscOptions *, const char pre[], const char name[], PetscBool *bval, PetscBool *set)
int atom_test
Atom test.
Definition plastic.cpp:121
virtual moab::Interface & get_moab()=0
Interface for managing meshsets containing materials and boundary conditions.
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.

◆ getBlockSetBCFromMesh()

MoFEMErrorCode VtkInterface::getBlockSetBCFromMesh ( std::string block_name)
private

Definition at line 53 of file read_vtk.cpp.

53 {
55
56 MeshsetsManager *meshsets_mng;
57 CHKERR mField.getInterface(meshsets_mng);
58
59 Tag mat_tag;
60 rval = mOab.tag_get_handle(block_name.c_str(), 1, MB_TYPE_INTEGER, mat_tag,
61 MB_TAG_EXCL | MB_TAG_SPARSE);
62 if (rval != MB_ALREADY_ALLOCATED) {
63 // std::cerr << "Error creating tag " << block_name << std::endl;
65 }
66 Range block_ents;
67 for (auto mbtype : {MBTET, MBTRI, MBHEX, MBQUAD, MBEDGE, MBVERTEX}) {
68 CHKERR mOab.get_entities_by_type_and_tag(0, mbtype, &mat_tag, NULL, 1,
69 block_ents, moab::Interface::UNION);
70 }
71
72 if (!block_ents.empty()) {
73 std::map<int, Range> mat_map;
74 for (auto &ent : block_ents) {
75 int mat_id;
76 CHKERR mOab.tag_get_data(mat_tag, &ent, 1, &mat_id);
77 mat_map[mat_id].insert(ent);
78 }
79
80 for (auto &[mat_id, mat_ents] : mat_map) {
81 std::string meshset_name = block_name + "_" + std::to_string(mat_id);
82 CHKERR meshsets_mng->addMeshset(BLOCKSET, mat_id, meshset_name);
83 CHKERR meshsets_mng->addEntitiesToMeshset(BLOCKSET, mat_id, mat_ents);
84 MOFEM_LOG("WORLD", Sev::inform)
85 << "Blockset " << meshset_name << " set added.";
86 }
87 }
88
90}
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define MOFEM_LOG(channel, severity)
Log.
MoFEMErrorCode addMeshset(const CubitBCType cubit_bc_type, const int ms_id, const std::string name="")
add cubit meshset
static MoFEMErrorCodeGeneric< moab::ErrorCode > rval
MoFEMErrorCode addEntitiesToMeshset(const CubitBCType cubit_bc_type, const int ms_id, const Range &ents)
add entities to cubit meshset

◆ getBoundaryConditions()

MoFEMErrorCode VtkInterface::getBoundaryConditions ( )
private

Definition at line 92 of file read_vtk.cpp.

92 {
94 MeshsetsManager *meshsets_manager_ptr;
95 CHKERR mField.getInterface(meshsets_manager_ptr);
96
97 Tag bc_tag;
98 CHKERR mOab.tag_get_handle("BOUNDARY_CONDITIONS", 1, MB_TYPE_INTEGER, bc_tag,
99 MB_TAG_CREAT | MB_TAG_SPARSE);
100
101 std::vector<std::tuple<int, std::string, std::string>> conditions = {
102 {1, "FIX_ALL", "FIX_ALL"},
103 {2, "FIX_X", "FIX_X"},
104 {3, "FIX_Y", "FIX_Y"},
105 {4, "FIX_Z", "FIX_Z"},
106 {5, "FIX_X_2", "DISP_X"},
107 {6, "FIX_Y_2", "DISP_Y"},
108 {7, "FIX_Z_2", "DISP_Z"},
109 {8, "FORCE_X", "FORCE_X"},
110 {9, "FORCE_Y", "FORCE_Y"},
111 {10, "FORCE_Z", "FORCE_Z"},
112 {11, "VEL_X", "VEL_X"},
113 {12, "VEL_Y", "VEL_Y"},
114 {13, "VEL_Z", "VEL_Z"},
115 {14, "ACCL_X", "ACCL_X"},
116 {15, "ACCL_Y", "ACCL_Y"},
117 {16, "ACCL_Z", "ACCL_Z"},
118 {17, "TEMP", "TEMP"},
119 {18, "PRESSURE", "PRESSURE"},
120 {19, "HEAT_FLUX", "HEAT_FLUX"},
121 {20, "CONTACT", "CONTACT"},
122 {21, "DISPLACEMENT", "DISPLACEMENT"},
123 {22, "ROTATE_ALL", "ROTATE_ALL"},
124 {23, "ROTATE_X", "ROTATE_X"},
125 {24, "ROTATE_Y", "ROTATE_Y"},
126 {25, "ROTATE_Z", "ROTATE_Z"},
127 {26, "TEMPERATURE", "TEMPERATURE"},
128 {50, "TIE_MATRIX", "TIE_MATRIX"}};
129
130 for (auto &condition : conditions) {
131 int desired_val = std::get<0>(condition);
132 const void *desired_val_ptr = &desired_val;
133 Range bc;
134 CHKERR mOab.get_entities_by_type_and_tag(0, MBVERTEX, &bc_tag,
135 &desired_val_ptr, 1, bc);
136
137 if (!bc.empty()) {
138 Range faces;
139 CHKERR mOab.get_adjacencies(bc, 2, true, faces, moab::Interface::UNION);
140 Range edges;
141 CHKERR mOab.get_adjacencies(bc, 1, true, edges, moab::Interface::UNION);
142
143 Range adj_nodes;
144 CHKERR mOab.get_adjacencies(faces, 0, false, adj_nodes,
145 moab::Interface::UNION);
146
147 Range non_BC_nodes;
148 non_BC_nodes = subtract(adj_nodes, bc);
149 Range non_BC_faces;
150 CHKERR mOab.get_adjacencies(non_BC_nodes, 2, false, non_BC_faces,
151 moab::Interface::UNION);
152 Range non_BC_edges;
153 CHKERR mOab.get_adjacencies(non_BC_nodes, 1, false, non_BC_edges,
154 moab::Interface::UNION);
155
156 faces = subtract(faces, non_BC_faces);
157 edges = subtract(edges, non_BC_edges);
158
159 MOFEM_LOG_C("WORLD", Sev::inform, "Number of nodes assigned to %s = %d",
160 std::get<2>(condition).c_str(), bc.size());
161 MOFEM_LOG_C("WORLD", Sev::inform, "Number of edges assigned to %s = %d",
162 std::get<2>(condition).c_str(), edges.size());
163 MOFEM_LOG_C("WORLD", Sev::inform, "Number of faces assigned to %s = %d",
164 std::get<2>(condition).c_str(), faces.size());
165
166 CHKERR meshsets_manager_ptr->addMeshset(BLOCKSET, desired_val,
167 std::get<1>(condition));
168 CHKERR meshsets_manager_ptr->addEntitiesToMeshset(BLOCKSET, desired_val,
169 bc);
170 CHKERR meshsets_manager_ptr->addEntitiesToMeshset(BLOCKSET, desired_val,
171 edges);
172 CHKERR meshsets_manager_ptr->addEntitiesToMeshset(BLOCKSET, desired_val,
173 faces);
174 } else {
175 CHKERR getBlockSetBCFromMesh(std::get<1>(condition));
176 }
177 }
178
180}
#define MOFEM_LOG_C(channel, severity, format,...)
MoFEMErrorCode getBlockSetBCFromMesh(std::string)
Definition read_vtk.cpp:53

◆ getMaterialProperties()

MoFEMErrorCode VtkInterface::getMaterialProperties ( )
private

Definition at line 192 of file read_vtk.cpp.

192 {
194 MeshsetsManager *meshsets_manager_ptr;
195 CHKERR mField.getInterface(meshsets_manager_ptr);
196
197 Range ents;
198 rval = mOab.get_entities_by_dimension(0, 3, ents, true);
199 if (rval != MB_SUCCESS) {
200 MOFEM_LOG("WORLD", Sev::warning)
201 << "No hexes/tets in the mesh, no material block set. Not Implemented";
202 } else {
203 CHKERR getBlockSetBCFromMesh("MAT_ELASTIC");
204 // for backward compatibility
205 CHKERR meshsets_manager_ptr->addMeshset(BLOCKSET, 100, "ADOLCMAT");
206 CHKERR meshsets_manager_ptr->addEntitiesToMeshset(BLOCKSET, 100, ents);
207 MOFEM_LOG("WORLD", Sev::inform) << "Material block ADOLCMAT set added.";
208 }
210}

◆ readMesh()

MoFEMErrorCode VtkInterface::readMesh ( )
private

Definition at line 43 of file read_vtk.cpp.

43 {
45
49
51}
MoFEMErrorCode loadFile(const std::string options, const std::string mesh_file_name, LoadFileFunc loadFunc=defaultLoadFileFunc)
Load mesh file.
Definition Simple.cpp:191
MoFEMErrorCode getOptions()
get options
Definition Simple.cpp:180
Simple * simpleInterface
Definition read_vtk.cpp:21

◆ readVtk()

MoFEMErrorCode VtkInterface::readVtk ( )

Definition at line 32 of file read_vtk.cpp.

32 {
41}
MoFEMErrorCode setBoundaryConditions()
Definition read_vtk.cpp:182
MoFEMErrorCode getMaterialProperties()
Definition read_vtk.cpp:192
MoFEMErrorCode readMesh()
Definition read_vtk.cpp:43
MoFEMErrorCode checkResults()
Definition read_vtk.cpp:227
MoFEMErrorCode writeOutput()
Definition read_vtk.cpp:212
MoFEMErrorCode getBoundaryConditions()
Definition read_vtk.cpp:92

◆ setBoundaryConditions()

MoFEMErrorCode VtkInterface::setBoundaryConditions ( )
private

Definition at line 182 of file read_vtk.cpp.

182 {
184 // Add meshsets if config file provided
185 MeshsetsManager *meshsets_interface_ptr;
186 CHKERR mField.getInterface(meshsets_interface_ptr);
187 CHKERR meshsets_interface_ptr->setMeshsetFromFile();
188
190}
MoFEMErrorCode setMeshsetFromFile(const string file_name, const bool clean_file_options=true)
add blocksets reading config file

◆ writeOutput()

MoFEMErrorCode VtkInterface::writeOutput ( )
private

Definition at line 212 of file read_vtk.cpp.

212 {
214
215 char mesh_out_file[255] = "out_vtk.h5m";
216
217 PetscOptionsBegin(mField.get_comm(), "", "Read vtk tool", "none");
218 CHKERR PetscOptionsString("-output_file", "output mesh file name", "",
219 "out.h5m", mesh_out_file, 255, PETSC_NULLPTR);
220 PetscOptionsEnd();
221
222 CHKERR mOab.write_file(mesh_out_file);
223
225}
virtual MPI_Comm & get_comm() const =0

Member Data Documentation

◆ mField

MoFEM::Interface& VtkInterface::mField
private

Definition at line 19 of file read_vtk.cpp.

◆ mOab

moab::Interface& VtkInterface::mOab
private

Definition at line 20 of file read_vtk.cpp.

◆ simpleInterface

Simple* VtkInterface::simpleInterface
private

Definition at line 21 of file read_vtk.cpp.


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