31using namespace boost::multi_index;
 
   32using namespace boost::multiprecision;
 
   33using namespace boost::numeric;
 
   37#ifdef ENABLE_PYTHON_BINDING  
   38  #include <boost/python.hpp> 
   39  #include <boost/python/def.hpp> 
   40  #include <boost/python/numpy.hpp> 
   41namespace bp = boost::python;
 
   42namespace np = boost::python::numpy;
 
   48static char help[] = 
"...\n\n";
 
   50int main(
int argc, 
char *argv[]) {
 
   53  const char param_file[] = 
"param_file.petsc";
 
   57  auto core_log = logging::core::get();
 
   66  #ifdef ENABLE_PYTHON_BINDING 
   69  MOFEM_LOG(
"EP", Sev::inform) << 
"Python initialised";
 
   71  MOFEM_LOG(
"EP", Sev::inform) << 
"Python NOT initialised";
 
   82    PetscBool flg = PETSC_TRUE;
 
   88    char restart_file[255];
 
   89    PetscBool restart_flg = PETSC_TRUE;
 
   96    DMType dm_name = 
"DMMOFEM";
 
   98    DMType dm_name_mg = 
"DMMOFEM_MG";
 
  102    moab::Core moab_core;
 
  103    moab::Interface &moab = moab_core;
 
  105    ParallelComm *pcomm = ParallelComm::get_pcomm(&moab, 
MYPCOMM_INDEX);
 
  107      pcomm = 
new ParallelComm(&moab, PETSC_COMM_WORLD);
 
  109    PetscBool fully_distributed = PETSC_FALSE;
 
  111                               &fully_distributed, PETSC_NULLPTR);
 
  112    if (fully_distributed) {
 
  114      if (pcomm->proc_config().proc_size() == 1)
 
  117        option = 
"PARALLEL=READ_PART;" 
  118                 "PARALLEL_RESOLVE_SHARED_ENTS;" 
  119                 "PARTITION=PARALLEL_PARTITION";
 
  127    MOFEM_LOG(
"EP", Sev::inform) << 
"Initialise MoFEM database";
 
  130    MOFEM_LOG(
"EP", Sev::inform) << 
"Initialise MoFEM database <- done";
 
  147    auto get_adj = [&](
Range ents, 
int dim) {
 
  149      CHKERR moab.get_adjacencies(ents, dim, 
false, adj,
 
  150                                  moab::Interface::UNION);
 
  194    const char *list_materials[LastMaterial] = {
 
  195        "stvenant_kirchhoff", 
"mooney_rivlin", 
"hencky", 
"neo_hookean"};
 
  196    PetscInt choice_material = MooneyRivlin;
 
  198                                LastMaterial, &choice_material, PETSC_NULLPTR);
 
  200    switch (choice_material) {
 
  201    case StVenantKirchhoff:
 
  202      MOFEM_LOG(
"EP", Sev::inform) << 
"StVenantKirchhoff material model";
 
  207      MOFEM_LOG(
"EP", Sev::inform) << 
"MooneyRivlin material model";
 
  208      MOFEM_LOG(
"EP", Sev::inform) << 
"MU(1, 0)/2 = " << 
MU(1, 0) / 2;
 
  214      MOFEM_LOG(
"EP", Sev::inform) << 
"Hencky material model";
 
  218      MOFEM_LOG(
"EP", Sev::inform) << 
"Neo-Hookean material model";
 
  232    CHKERR VecSetOption(f_elastic, VEC_IGNORE_NEGATIVE_INDICES, PETSC_TRUE);
 
  234    PetscInt start_step = 0;
 
  237      CHKERR PetscViewerBinaryOpen(PETSC_COMM_WORLD, restart_file,
 
  238                                   FILE_MODE_READ, &viewer);
 
  239      CHKERR VecLoad(x_elastic, viewer);
 
  240      CHKERR PetscViewerDestroy(&viewer);
 
  244      std::string restart_file_str(restart_file);
 
  245      std::regex restart_pattern(R
"(restart_([1-9]\d*)\.dat)"); 
  247      if (std::regex_match(restart_file_str, match, restart_pattern)) {
 
  248        start_step = std::stoi(match[1]);
 
  251                "Restart file name must be in the format restart_##.dat");
 
  255    auto ts_elastic = 
createTS(PETSC_COMM_WORLD);
 
  256    CHKERR TSSetType(ts_elastic, TSBEULER);
 
  259    CHKERR TSGetAdapt(ts_elastic, &adapt);
 
  260    CHKERR TSAdaptSetType(adapt, TSADAPTNONE);
 
  261    CHKERR TSSetTime(ts_elastic, time);
 
  262    CHKERR TSSetStepNumber(ts_elastic, start_step);
 
  276#ifdef ENABLE_PYTHON_BINDING 
  277  if (Py_FinalizeEx() < 0) {
 
 
Eshelbian plasticity interface.
#define CATCH_ERRORS
Catch errors.
#define MYPCOMM_INDEX
default communicator number PCOMM
@ MOFEM_DATA_INCONSISTENCY
#define CHKERR
Inline error check.
PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode, ScatterMode scatter_mode)
set local (or ghosted) vector values on mesh for partition only
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
MoFEMErrorCode DMRegister_MGViaApproxOrders(const char sname[])
Register DM for Multi-Grid via approximation orders.
PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g)
DMShellSetCreateGlobalVector.
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.
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
implementation of Data Operators for Forces and Sources
PetscErrorCode PetscOptionsGetBool(PetscOptions *, const char pre[], const char name[], PetscBool *bval, PetscBool *set)
PetscErrorCode PetscOptionsGetScalar(PetscOptions *, const char pre[], const char name[], PetscScalar *dval, PetscBool *set)
SmartPetscObj< Vec > vectorDuplicate(Vec vec)
Create duplicate vector of smart vector.
auto createTS(MPI_Comm comm)
PetscErrorCode PetscOptionsGetEList(PetscOptions *, const char pre[], const char name[], const char *const *list, PetscInt next, PetscInt *value, PetscBool *set)
PetscErrorCode PetscOptionsGetString(PetscOptions *, const char pre[], const char name[], char str[], size_t size, PetscBool *set)
auto get_temp_meshset_ptr(moab::Interface &moab)
Create smart pointer to temporary meshset.
PetscErrorCode TSAdaptCreateMoFEM(TSAdapt adapt)
Craete MOFEM adapt.
MoFEMErrorCode setElasticElementOps(const int tag)
MoFEMErrorCode addVolumeFiniteElement(const EntityHandle meshset=0)
MoFEMErrorCode addFields(const EntityHandle meshset=0)
MoFEMErrorCode createCrackSurfaceMeshset()
MoFEMErrorCode addBoundaryFiniteElement(const EntityHandle meshset=0)
MoFEMErrorCode getSpatialRotationBc()
MoFEMErrorCode addMaterial_Hencky(double E, double nu)
MoFEMErrorCode projectGeometry(const EntityHandle meshset=0)
MoFEMErrorCode setBlockTagsOnSkin()
MoFEMErrorCode addMaterial_HMHHStVenantKirchhoff(const int tape, const double lambda, const double mu, const double sigma_y)
MoFEMErrorCode addMaterial_HMHMooneyRivlin(const int tape, const double alpha, const double beta, const double lambda, const double sigma_y)
static PetscBool dynamicRelaxation
MoFEMErrorCode solveElastic(TS ts, Vec x)
MoFEMErrorCode setElasticElementToTs(DM dm)
MoFEMErrorCode solveDynamicRelaxation(TS ts, Vec x)
MoFEMErrorCode getSpatialTractionFreeBc(const EntityHandle meshset=0)
MoFEMErrorCode getExternalStrain()
MoFEMErrorCode getSpatialTractionBc()
MoFEMErrorCode addDMs(const BitRefLevel bit=BitRefLevel().set(0), const EntityHandle meshset=0)
MoFEMErrorCode getSpatialDispBc()
[Getting norms]
MoFEMErrorCode createExchangeVectors(Sev sev)
MoFEMErrorCode addMaterial_HMHNeohookean(const int tape, const double c10, const double K)
SmartPetscObj< DM > dmElastic
Elastic problem.
static Range getPartEntities(moab::Interface &moab, int part)
static MoFEMErrorCode loadFileRootProcAllRestDistributed(moab::Interface &moab, const char *file_name, int dim, LoadFileFun proc_skin_fun=defaultProcSkinFun, const char *options="PARALLEL=BCAST;PARTITION=")
Root proc has whole mesh, other procs only part of it.
virtual MPI_Comm & get_comm() const =0
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 > getStrmWorld()
Get the strm world object.
static boost::shared_ptr< std::ostream > getStrmSync()
Get the strm sync object.
static boost::shared_ptr< std::ostream > getStrmSelf()
Get the strm self object.
Interface for managing meshsets containing materials and boundary conditions.
intrusive_ptr for managing petsc objects
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.