Test and example setting cubit meshsets.
static char help[] =
"...\n\n";
int main(
int argc,
char *argv[]) {
auto core_log = logging::core::get();
core_log->add_sink(
LogManager::createSink(LogManager::getStrmSelf(), "ATOM_TEST"));
LogManager::setLog("ATOM_TEST");
try {
moab::Core mb_instance;
moab::Interface &moab = mb_instance;
int rank;
MPI_Comm_rank(PETSC_COMM_WORLD, &rank);
LogManager::BitLineID | LogManager::BitScope);
MOFEM_LOG(
"ATOM_TEST", Sev::verbose) <<
"<<<< SIDESETs >>>>>";
bool add_block_is_there = false;
{
strncpy(mybc.
data.name,
"Pressure", 8);
}
if (it->getMeshsetId() != 1002)
continue;
add_block_is_there = true;
CHKERR it->getBcDataStructure(mydata);
MOFEM_LOG(
"ATOM_TEST", Sev::inform) << mydata;
}
if (!add_block_is_there)
"no added block set");
MOFEM_LOG(
"ATOM_TEST", Sev::inform) <<
"<<<< BLOCKSETs >>>>>";
add_block_is_there = false;
std::vector<double> attr(3);
attr[0] = 0;
attr[1] = 1;
attr[2] = 2;
std::string name = it->getName();
if (name.compare(0, 13, "ADD_BLOCK_SET") == 0) {
add_block_is_there = true;
std::vector<double> attributes;
it->getAttributes(attributes);
if (attributes.size() != 3) {
"should be 3 attributes but is %d", attributes.size());
}
if (attributes[0] != 0 || attributes[1] != 1 || attributes[2] != 2) {
"wrong values of attributes");
}
}
}
if (!add_block_is_there) {
"no added block set");
}
add_block_is_there = false;
{
mydata.
data.Poisson = 0.25;
mydata);
}
if (it->getMeshsetId() != 1001)
continue;
std::string name = it->getName();
if (name.compare(0, 13, "MAT_ELASTIC") == 0 &&
add_block_is_there = true;
CHKERR it->getAttributeDataStructure(mydata);
MOFEM_LOG(
"ATOM_TEST", Sev::inform) << mydata;
if (mydata.
data.Young != 1 || mydata.
data.Poisson != 0.25) {
"wrong values of attributes");
}
}
}
if (!add_block_is_there) {
"no added block set");
}
MOFEM_LOG(
"ATOM_TEST", Sev::inform) <<
"<<<< NODESET >>>>>";
std::memcpy(disp_bc.
data.name,
"Displacement", 12);
CHKERR it->getBcDataStructure(disp_data);
MOFEM_LOG(
"ATOM_TEST", Sev::inform) << disp_data;
}
<< "<<<< ADD BLOCKSETs FROM CONFIG FILE >>>>>";
);
MOFEM_LOG(
"ATOM_TEST", Sev::inform) <<
"Iterate blocksets";
bool mat_elastic_trans_is_found = true;
std::vector<double> attributes;
it->getAttributes(attributes);
std::ostringstream ss;
ss << "Attr: ";
for (unsigned int ii = 0; ii != attributes.size(); ii++) {
ss << attributes[ii] << " ";
}
MOFEM_LOG(
"ATOM_TEST", Sev::inform) << ss.str();
std::string block_name = it->getName();
if (block_name.compare(0, block_name.size(), "MAT_ELASTIC_TRANS_ISO") ==
0) {
MOFEM_LOG(
"ATOM_TEST", Sev::inform) <<
"Mat Trans Iso block ";
mat_elastic_trans_is_found = true;
CHKERR it->getAttributeDataStructure(mydata);
MOFEM_LOG(
"ATOM_TEST", Sev::inform) << mydata;
if (mydata.
data.Youngp != 1)
if (mydata.
data.Youngz != 2)
if (mydata.
data.Poissonp != 3)
if (mydata.
data.Poissonpz != 4)
if (mydata.
data.Shearzp != 5)
}
}
if(!mat_elastic_trans_is_found)
CHKERR it->getAttributeDataStructure(mydata);
MOFEM_LOG(
"ATOM_TEST", Sev::inform) <<
"Mat elastic found ";
MOFEM_LOG(
"ATOM_TEST", Sev::inform) << mydata;
}
CHKERR it->getAttributeDataStructure(mydata);
MOFEM_LOG(
"ATOM_TEST", Sev::inform) <<
"Mat thermal found ";
MOFEM_LOG(
"ATOM_TEST", Sev::inform) << mydata;
}
CHKERR it->getBcDataStructure(mydata);
MOFEM_LOG(
"ATOM_TEST", Sev::inform) << mydata;
}
CHKERR it->getBcDataStructure(mydata);
MOFEM_LOG(
"ATOM_TEST", Sev::inform) << mydata;
}
CHKERR it->getBcDataStructure(mydata);
MOFEM_LOG(
"ATOM_TEST", Sev::inform) << mydata;
}
CHKERR it->getBcDataStructure(mydata);
MOFEM_LOG(
"ATOM_TEST", Sev::inform) << mydata;
}
CHKERR it->getBcDataStructure(mydata);
MOFEM_LOG(
"ATOM_TEST", Sev::inform) << mydata;
}
}
}
}
#define CATCH_ERRORS
Catch errors.
@ MAT_ELASTICSET
block name is "MAT_ELASTIC"
@ MAT_THERMALSET
block name is "MAT_THERMAL"
@ MOFEM_OPERATION_UNSUCCESSFUL
@ MOFEM_ATOM_TEST_INVALID
#define CHKERR
Inline error check.
#define MOFEM_LOG(channel, severity)
Log.
#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
implementation of Data Operators for Forces and Sources
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Deprecated interface functions.
Definition of the displacement bc data structure.
Definition of the force bc data structure.
Definition of the heat flux bc data structure.
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 of the temperature bc data structure.
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.