static char help[] =
"...\n\n";
int main(
int argc,
char *argv[]) {
auto core_log = logging::core::get();
core_log->add_sink(
try {
int rank;
MPI_Comm_rank(PETSC_COMM_WORLD, &rank);
MeshsetsManager *meshsets_manager_ptr;
MOFEM_LOG(
"ATOM_TEST", Sev::verbose) <<
"<<<< SIDESETs >>>>>";
bool add_block_is_there = false;
{
PressureCubitBcData mybc;
strncpy(mybc.data.name, "Pressure", 8);
mybc.data.flag1 = 0;
mybc.data.flag2 = 0;
mybc.data.value1 = 1;
}
if (it->getMeshsetId() != 1002)
continue;
add_block_is_there = true;
PressureCubitBcData mydata;
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;
CHKERR meshsets_manager_ptr->addMeshset(
BLOCKSET, 1000,
"ADD_BLOCK_SET");
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;
CHKERR meshsets_manager_ptr->addMeshset(
BLOCKSET, 1001,
"MAT_ELASTIC");
{
Mat_Elastic mydata;
mydata.data.Young = 1;
mydata.data.Poisson = 0.25;
CHKERR meshsets_manager_ptr->setAtributesByDataStructure(
BLOCKSET, 1001,
mydata);
}
if (it->getMeshsetId() != 1001)
continue;
std::string name = it->getName();
if (name.compare(0, 13, "MAT_ELASTIC") == 0 &&
add_block_is_there = true;
Mat_Elastic mydata;
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 >>>>>";
DisplacementCubitBcData disp_bc;
std::memcpy(disp_bc.data.name, "Displacement", 12);
disp_bc.data.flag1 = 1;
disp_bc.data.flag2 = 1;
disp_bc.data.flag3 = 1;
disp_bc.data.flag4 = 0;
disp_bc.data.flag5 = 0;
disp_bc.data.flag6 = 0;
disp_bc.data.value1 = 0;
disp_bc.data.value2 = 0;
disp_bc.data.value3 = 0;
disp_bc.data.value4 = 0;
disp_bc.data.value5 = 0;
disp_bc.data.value6 = 0;
DisplacementCubitBcData disp_data;
CHKERR it->getBcDataStructure(disp_data);
MOFEM_LOG(
"ATOM_TEST", Sev::inform) << disp_data;
}
<< "<<<< ADD BLOCKSETs FROM CONFIG FILE >>>>>";
CHKERR meshsets_manager_ptr->addMeshset(
BLOCKSET, 1002,
"ADD_BLOCK_SET");
CHKERR meshsets_manager_ptr->addMeshset(
BLOCKSET, 1003,
"ADD_BLOCK_SET");
CHKERR meshsets_manager_ptr->addMeshset(
BLOCKSET, 1004,
"ADD_BLOCK_SET");
CHKERR meshsets_manager_ptr->addMeshset(
BLOCKSET, 1005,
"ADD_BLOCK_SET");
CHKERR meshsets_manager_ptr->addMeshset(
BLOCKSET, 1006,
"ADD_BLOCK_SET");
CHKERR meshsets_manager_ptr->addMeshset(
BLOCKSET, 1007,
"ADD_BLOCK_SET");
CHKERR meshsets_manager_ptr->addMeshset(
BLOCKSET, 1008,
"ADD_BLOCK_SET");
CHKERR meshsets_manager_ptr->addMeshset(
BLOCKSET, 1009,
"ADD_BLOCK_SET");
CHKERR meshsets_manager_ptr->setMeshsetFromFile(
);
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;
Mat_Elastic_TransIso mydata;
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)
Mat_Elastic mydata;
CHKERR it->getAttributeDataStructure(mydata);
MOFEM_LOG(
"ATOM_TEST", Sev::inform) <<
"Mat elastic found ";
MOFEM_LOG(
"ATOM_TEST", Sev::inform) << mydata;
}
Mat_Thermal mydata;
CHKERR it->getAttributeDataStructure(mydata);
MOFEM_LOG(
"ATOM_TEST", Sev::inform) <<
"Mat thermal found ";
MOFEM_LOG(
"ATOM_TEST", Sev::inform) << mydata;
}
DisplacementCubitBcData mydata;
CHKERR it->getBcDataStructure(mydata);
MOFEM_LOG(
"ATOM_TEST", Sev::inform) << mydata;
}
ForceCubitBcData mydata;
CHKERR it->getBcDataStructure(mydata);
MOFEM_LOG(
"ATOM_TEST", Sev::inform) << mydata;
}
PressureCubitBcData mydata;
CHKERR it->getBcDataStructure(mydata);
MOFEM_LOG(
"ATOM_TEST", Sev::inform) << mydata;
}
TemperatureCubitBcData mydata;
CHKERR it->getBcDataStructure(mydata);
MOFEM_LOG(
"ATOM_TEST", Sev::inform) << mydata;
}
HeatFluxCubitBcData mydata;
CHKERR it->getBcDataStructure(mydata);
MOFEM_LOG(
"ATOM_TEST", Sev::inform) << mydata;
}
}
}
}
int main(int argc, char *argv[])
#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.
static LoggerType & setLog(const std::string channel)
Set ans resset chanel logger.
#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.
#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
DeprecatedCoreInterface Interface
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.
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.
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.