v0.14.0
Loading...
Searching...
No Matches
Functions | Variables
add_cubit_meshsets.cpp File Reference
#include <MoFEM.hpp>

Go to the source code of this file.

Functions

int main (int argc, char *argv[])
 

Variables

static char help [] = "...\n\n"
 

Function Documentation

◆ main()

int main ( int argc,
char * argv[] )

Definition at line 14 of file add_cubit_meshsets.cpp.

14 {
15
16 MoFEM::Core::Initialize(&argc, &argv, (char *)0, help);
17
18 auto core_log = logging::core::get();
19 core_log->add_sink(
21 LogManager::setLog("ATOM_TEST");
22
23 try {
24
25 moab::Core mb_instance;
26 moab::Interface &moab = mb_instance;
27 int rank;
28 MPI_Comm_rank(PETSC_COMM_WORLD, &rank);
29
30 // Create MoFEM (Joseph) database
31 MoFEM::Core core(moab);
32 MoFEM::Interface &m_field = core;
33
34 MeshsetsManager *meshsets_manager_ptr;
35 CHKERR m_field.getInterface(meshsets_manager_ptr);
36
37 MOFEM_LOG_CHANNEL("ATOM_TEST")
38 MOFEM_LOG_ATTRIBUTES("ATOM_TEST",
39 LogManager::BitLineID | LogManager::BitScope);
40 MOFEM_LOG_TAG("ATOM_TEST", "atom test");
41
42 MOFEM_LOG("ATOM_TEST", Sev::verbose) << "<<<< SIDESETs >>>>>";
43
44 bool add_block_is_there = false;
45 CHKERR meshsets_manager_ptr->addMeshset(SIDESET, 1002);
46 {
48 strncpy(mybc.data.name, "Pressure", 8);
49 mybc.data.flag1 = 0;
50 mybc.data.flag2 = 0;
51 mybc.data.value1 = 1;
52 CHKERR meshsets_manager_ptr->setBcData(SIDESET, 1002, mybc);
53 }
55 if (it->getMeshsetId() != 1002)
56 continue;
57 add_block_is_there = true;
59 CHKERR it->getBcDataStructure(mydata);
60 // Print data
61 MOFEM_LOG("ATOM_TEST", Sev::inform) << mydata;
62 }
63 if (!add_block_is_there)
64 SETERRQ(PETSC_COMM_WORLD, MOFEM_OPERATION_UNSUCCESSFUL,
65 "no added block set");
66
67 MOFEM_LOG("ATOM_TEST", Sev::inform) << "<<<< BLOCKSETs >>>>>";
68
69 add_block_is_there = false;
70 CHKERR meshsets_manager_ptr->addMeshset(BLOCKSET, 1000, "ADD_BLOCK_SET");
71 std::vector<double> attr(3);
72 attr[0] = 0;
73 attr[1] = 1;
74 attr[2] = 2;
75 CHKERR meshsets_manager_ptr->setAtributes(BLOCKSET, 1000, attr);
77 // Get block name
78 std::string name = it->getName();
79 if (name.compare(0, 13, "ADD_BLOCK_SET") == 0) {
80 add_block_is_there = true;
81 std::vector<double> attributes;
82 it->getAttributes(attributes);
83 if (attributes.size() != 3) {
84 SETERRQ1(PETSC_COMM_WORLD, MOFEM_ATOM_TEST_INVALID,
85 "should be 3 attributes but is %d", attributes.size());
86 }
87 if (attributes[0] != 0 || attributes[1] != 1 || attributes[2] != 2) {
88 SETERRQ(PETSC_COMM_WORLD, MOFEM_ATOM_TEST_INVALID,
89 "wrong values of attributes");
90 }
91 }
92 }
93 if (!add_block_is_there) {
94 SETERRQ(PETSC_COMM_WORLD, MOFEM_OPERATION_UNSUCCESSFUL,
95 "no added block set");
96 }
97 add_block_is_there = false;
98 CHKERR meshsets_manager_ptr->addMeshset(BLOCKSET, 1001, "MAT_ELASTIC");
99 {
100 Mat_Elastic mydata;
101 mydata.data.Young = 1;
102 mydata.data.Poisson = 0.25;
103 CHKERR meshsets_manager_ptr->setAtributesByDataStructure(BLOCKSET, 1001,
104 mydata);
105 }
107 if (it->getMeshsetId() != 1001)
108 continue;
109 // Get block name
110 std::string name = it->getName();
111 if (name.compare(0, 13, "MAT_ELASTIC") == 0 &&
112 (it->getBcType() & CubitBCType(MAT_ELASTICSET)).any()) {
113 add_block_is_there = true;
114 Mat_Elastic mydata;
115 CHKERR it->getAttributeDataStructure(mydata);
116 // Print data
117 MOFEM_LOG("ATOM_TEST", Sev::inform) << mydata;
118 if (mydata.data.Young != 1 || mydata.data.Poisson != 0.25) {
119 SETERRQ(PETSC_COMM_WORLD, MOFEM_ATOM_TEST_INVALID,
120 "wrong values of attributes");
121 }
122 }
123 }
124 if (!add_block_is_there) {
125 SETERRQ(PETSC_COMM_WORLD, MOFEM_OPERATION_UNSUCCESSFUL,
126 "no added block set");
127 }
128
129 MOFEM_LOG("ATOM_TEST", Sev::inform) << "<<<< NODESET >>>>>";
130
131 CHKERR meshsets_manager_ptr->addMeshset(NODESET, 1010);
133 std::memcpy(disp_bc.data.name, "Displacement", 12);
134 disp_bc.data.flag1 = 1;
135 disp_bc.data.flag2 = 1;
136 disp_bc.data.flag3 = 1;
137 disp_bc.data.flag4 = 0;
138 disp_bc.data.flag5 = 0;
139 disp_bc.data.flag6 = 0;
140 disp_bc.data.value1 = 0;
141 disp_bc.data.value2 = 0;
142 disp_bc.data.value3 = 0;
143 disp_bc.data.value4 = 0;
144 disp_bc.data.value5 = 0;
145 disp_bc.data.value6 = 0;
146
147 CHKERR meshsets_manager_ptr->setBcData(NODESET, 1010, disp_bc);
148
150 m_field, NODESET | DISPLACEMENTSET, it)) {
151 DisplacementCubitBcData disp_data;
152 CHKERR it->getBcDataStructure(disp_data);
153 MOFEM_LOG("ATOM_TEST", Sev::inform) << disp_data;
154 }
155
156 MOFEM_LOG("ATOM_TEST", Sev::inform)
157 << "<<<< ADD BLOCKSETs FROM CONFIG FILE >>>>>";
158
159 CHKERR meshsets_manager_ptr->addMeshset(BLOCKSET, 1002, "ADD_BLOCK_SET");
160 CHKERR meshsets_manager_ptr->addMeshset(BLOCKSET, 1003, "ADD_BLOCK_SET");
161 CHKERR meshsets_manager_ptr->addMeshset(BLOCKSET, 1004, "ADD_BLOCK_SET");
162 CHKERR meshsets_manager_ptr->addMeshset(BLOCKSET, 1005, "ADD_BLOCK_SET");
163 CHKERR meshsets_manager_ptr->addMeshset(BLOCKSET, 1006, "ADD_BLOCK_SET");
164 CHKERR meshsets_manager_ptr->addMeshset(BLOCKSET, 1007, "ADD_BLOCK_SET");
165 CHKERR meshsets_manager_ptr->addMeshset(BLOCKSET, 1008, "ADD_BLOCK_SET");
166 CHKERR meshsets_manager_ptr->addMeshset(BLOCKSET, 1009, "ADD_BLOCK_SET");
167 CHKERR meshsets_manager_ptr->addMeshset(BLOCKSET, 1010, "foo");
168
169 CHKERR meshsets_manager_ptr->setMeshsetFromFile(
170 /*"add_cubit_meshsets.in"*/);
171
172 // List all meshsets
173 MOFEM_LOG("ATOM_TEST", Sev::inform) << "Iterate blocksets";
174
175 bool mat_elastic_trans_is_found = true;
176 for (_IT_CUBITMESHSETS_FOR_LOOP_(m_field, it)) {
177 MOFEM_LOG("ATOM_TEST", Sev::inform) << *it;
178 if ((it->getBcType() & CubitBCType(BLOCKSET)).any()) {
179 std::vector<double> attributes;
180 it->getAttributes(attributes);
181 std::ostringstream ss;
182 ss << "Attr: ";
183 for (unsigned int ii = 0; ii != attributes.size(); ii++) {
184 ss << attributes[ii] << " ";
185 }
186 MOFEM_LOG("ATOM_TEST", Sev::inform) << ss.str();
187
188 std::string block_name = it->getName();
189 if (block_name.compare(0, block_name.size(), "MAT_ELASTIC_TRANS_ISO") ==
190 0) {
191 MOFEM_LOG("ATOM_TEST", Sev::inform) << "Mat Trans Iso block ";
192 mat_elastic_trans_is_found = true;
194 CHKERR it->getAttributeDataStructure(mydata);
195 // Print data
196 MOFEM_LOG("ATOM_TEST", Sev::inform) << mydata;
197 if (mydata.data.Youngp != 1)
198 SETERRQ(PETSC_COMM_SELF, MOFEM_ATOM_TEST_INVALID, "Wrong value");
199 if (mydata.data.Youngz != 2)
200 SETERRQ(PETSC_COMM_SELF, MOFEM_ATOM_TEST_INVALID, "Wrong value");
201 if (mydata.data.Poissonp != 3)
202 SETERRQ(PETSC_COMM_SELF, MOFEM_ATOM_TEST_INVALID, "Wrong value");
203 if (mydata.data.Poissonpz != 4)
204 SETERRQ(PETSC_COMM_SELF, MOFEM_ATOM_TEST_INVALID, "Wrong value");
205 if (mydata.data.Shearzp != 5)
206 SETERRQ(PETSC_COMM_SELF, MOFEM_ATOM_TEST_INVALID, "Wrong value");
207 }
208 }
209 if(!mat_elastic_trans_is_found)
210 SETERRQ(PETSC_COMM_SELF, MOFEM_ATOM_TEST_INVALID, "Block not found");
211
212 if ((it->getBcType() & CubitBCType(MAT_ELASTICSET)).any()) {
213 Mat_Elastic mydata;
214 CHKERR it->getAttributeDataStructure(mydata);
215 MOFEM_LOG("ATOM_TEST", Sev::inform) << "Mat elastic found ";
216 MOFEM_LOG("ATOM_TEST", Sev::inform) << mydata;
217 }
218 if ((it->getBcType() & CubitBCType(MAT_THERMALSET)).any()) {
219 Mat_Thermal mydata;
220 CHKERR it->getAttributeDataStructure(mydata);
221 MOFEM_LOG("ATOM_TEST", Sev::inform) << "Mat thermal found ";
222 MOFEM_LOG("ATOM_TEST", Sev::inform) << mydata;
223 }
224 if ((it->getBcType() & CubitBCType(DISPLACEMENTSET)).any()) {
226 CHKERR it->getBcDataStructure(mydata);
227 MOFEM_LOG("ATOM_TEST", Sev::inform) << mydata;
228 }
229 if ((it->getBcType() & CubitBCType(FORCESET)).any()) {
230 ForceCubitBcData mydata;
231 CHKERR it->getBcDataStructure(mydata);
232 MOFEM_LOG("ATOM_TEST", Sev::inform) << mydata;
233 }
234 if ((it->getBcType() & CubitBCType(PRESSURESET)).any()) {
235 PressureCubitBcData mydata;
236 CHKERR it->getBcDataStructure(mydata);
237 MOFEM_LOG("ATOM_TEST", Sev::inform) << mydata;
238 }
239 if ((it->getBcType() & CubitBCType(TEMPERATURESET)).any()) {
241 CHKERR it->getBcDataStructure(mydata);
242 MOFEM_LOG("ATOM_TEST", Sev::inform) << mydata;
243 }
244 if ((it->getBcType() & CubitBCType(HEATFLUXSET)).any()) {
245 HeatFluxCubitBcData mydata;
246 CHKERR it->getBcDataStructure(mydata);
247 MOFEM_LOG("ATOM_TEST", Sev::inform) << mydata;
248 }
249 }
250 }
252
254}
static char help[]
#define CATCH_ERRORS
Catch errors.
@ TEMPERATURESET
@ PRESSURESET
@ MAT_ELASTICSET
block name is "MAT_ELASTIC"
@ FORCESET
@ HEATFLUXSET
@ NODESET
@ SIDESET
@ DISPLACEMENTSET
@ MAT_THERMALSET
block name is "MAT_THERMAL"
@ BLOCKSET
@ MOFEM_OPERATION_UNSUCCESSFUL
Definition definitions.h:34
@ MOFEM_ATOM_TEST_INVALID
Definition definitions.h:40
#define CHKERR
Inline error check.
static LoggerType & setLog(const std::string channel)
Set ans resset chanel logger.
#define MOFEM_LOG(channel, severity)
Log.
SeverityLevel
Severity levels.
#define MOFEM_LOG_TAG(channel, tag)
Tag channel.
#define MOFEM_LOG_CHANNEL(channel)
Set and reset channel.
#define MOFEM_LOG_ATTRIBUTES(channel, bit)
Add attributes to channel.
#define _IT_CUBITMESHSETS_BY_BCDATA_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet in a moFEM field.
MoFEMErrorCode addMeshset(const CubitBCType cubit_bc_type, const int ms_id, const std::string name="")
add cubit meshset
#define _IT_CUBITMESHSETS_FOR_LOOP_(MESHSET_MANAGER, IT)
Iterator that loops over all the Cubit MeshSets in a moFEM field.
#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.
std::bitset< 32 > CubitBCType
Definition Types.hpp:52
Core (interface) class.
Definition Core.hpp:82
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
Definition Core.cpp:72
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Definition Core.cpp:112
Deprecated interface functions.
Definition of the displacement bc data structure.
Definition BCData.hpp:76
Definition of the force bc data structure.
Definition BCData.hpp:139
Definition of the heat flux bc data structure.
Definition BCData.hpp:427
Log manager is used to build and partition problems.
static boost::shared_ptr< SinkType > createSink(boost::shared_ptr< std::ostream > stream_ptr, std::string comm_filter)
Create a sink object.
static boost::shared_ptr< std::ostream > getStrmSelf()
Get the strm self object.
Transverse Isotropic material data structure.
Elastic material data structure.
Thermal material data structure.
Interface for managing meshsets containing materials and boundary conditions.
MoFEMErrorCode setBcData(const CubitBCType cubit_bc_type, const int ms_id, const GenericCubitBcData &data)
set boundary data structure to meshset
MoFEMErrorCode setAtributes(const CubitBCType cubit_bc_type, const int ms_id, const std::vector< double > &attributes, const std::string name="")
set attributes to cubit meshset
MoFEMErrorCode setMeshsetFromFile(const string file_name, const bool clean_file_options=true)
add blocksets reading config file
MoFEMErrorCode setAtributesByDataStructure(const CubitBCType cubit_bc_type, const int ms_id, const GenericAttributeData &data, const std::string name="")
set (material) data structure to cubit meshset
Definition of the pressure bc data structure.
Definition BCData.hpp:379
Definition of the temperature bc data structure.
Definition BCData.hpp:306
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.

Variable Documentation

◆ help

char help[] = "...\n\n"
static

Definition at line 12 of file add_cubit_meshsets.cpp.