};
};
private:
};
}
char meshFileName[255];
meshFileName, 255, PETSC_NULL);
}
enum bases { AINSWORTH, DEMKOWICZ, LASBASETOPT };
const char *list_bases[LASBASETOPT] = {"ainsworth", "demkowicz"};
PetscInt choice_base_value = AINSWORTH;
LASBASETOPT, &choice_base_value, PETSC_NULL);
switch (choice_base_value) {
case AINSWORTH:
<< "Set AINSWORTH_LEGENDRE_BASE for displacements";
break;
case DEMKOWICZ:
<< "Set DEMKOWICZ_JACOBI_BASE for displacements";
break;
default:
break;
}
}
auto time_scale = boost::make_shared<TimeScale>();
};
pipeline_mng->getOpBoundaryRhsPipeline(), mField, "U", {time_scale},
"FORCE", Sev::inform);
pipeline_mng->getOpDomainRhsPipeline(), mField, "U", {time_scale},
"BODY_FORCE", Sev::inform);
CHKERR bc_mng->removeBlockDOFsOnEntities(
simple->getProblemName(),
"REMOVE_X",
"U", 0, 0);
CHKERR bc_mng->removeBlockDOFsOnEntities(
simple->getProblemName(),
"REMOVE_Y",
"U", 1, 1);
CHKERR bc_mng->removeBlockDOFsOnEntities(
simple->getProblemName(),
"REMOVE_Z",
"U", 2, 2);
simple->getProblemName(),
"U");
CHKERR bc_mng->addBlockDOFsToMPCs(
simple->getProblemName(),
"U");
}
auto add_domain_ops_lhs = [&](auto &pip) {
CHKERR HenckyOps::opFactoryDomainLhs<SPACE_DIM, PETSC, GAUSS, DomainEleOp>(
mField, pip, "U", "MAT_ELASTIC", Sev::inform);
};
auto add_domain_ops_rhs = [&](auto &pip) {
CHKERR HenckyOps::opFactoryDomainRhs<SPACE_DIM, PETSC, GAUSS, DomainEleOp>(
mField, pip, "U", "MAT_ELASTIC", Sev::inform);
};
CHKERR add_domain_ops_lhs(pipeline_mng->getOpDomainLhsPipeline());
CHKERR add_domain_ops_rhs(pipeline_mng->getOpDomainRhsPipeline());
}
: dM(dm), postProc(post_proc){};
postProc->mField, postProc->getPostProcMesh(), {"U"});
"out_step_" + boost::lexical_cast<std::string>(ts_step) + ".h5m");
}
}
private:
boost::shared_ptr<PostProcEleBdy> postProc;
};
auto ts = pipeline_mng->createTSIM();
auto create_post_proc_fe = [dm,
this,
simple]() {
auto post_proc_ele_domain = [dm, this](auto &pip_domain) {
auto common_ptr = commonDataFactory<SPACE_DIM, GAUSS, DomainEleOp>(
mField, pip_domain, "U", "MAT_ELASTIC", Sev::inform);
return common_ptr;
};
auto post_proc_fe_bdy = boost::make_shared<PostProcEleBdy>(mField);
auto u_ptr = boost::make_shared<MatrixDouble>();
post_proc_fe_bdy->getOpPtrVector().push_back(
auto op_loop_side =
auto common_ptr = post_proc_ele_domain(op_loop_side->getOpPtrVector());
post_proc_fe_bdy->getOpPtrVector().push_back(op_loop_side);
post_proc_fe_bdy->getOpPtrVector().push_back(
post_proc_fe_bdy->getPostProcMesh(),
post_proc_fe_bdy->getMapGaussPts(),
{},
{{"U", u_ptr}},
{{"GRAD", common_ptr->matGradPtr},
{"FIRST_PIOLA", common_ptr->getMatFirstPiolaStress()}},
{}
)
);
return post_proc_fe_bdy;
};
auto add_extra_finite_elements_to_ksp_solver_pipelines = [&]() {
auto pre_proc_ptr = boost::make_shared<FEMethod>();
auto post_proc_rhs_ptr = boost::make_shared<FEMethod>();
auto post_proc_lhs_ptr = boost::make_shared<FEMethod>();
auto time_scale = boost::make_shared<TimeScale>();
auto get_bc_hook_rhs = [this, pre_proc_ptr, time_scale]() {
{time_scale}, false)();
};
pre_proc_ptr->preProcessHook = get_bc_hook_rhs;
auto get_post_proc_hook_rhs = [this, post_proc_rhs_ptr]() {
mField, post_proc_rhs_ptr, nullptr, Sev::verbose)();
mField, post_proc_rhs_ptr, 1.)();
};
auto get_post_proc_hook_lhs = [this, post_proc_lhs_ptr]() {
mField, post_proc_lhs_ptr, 1.)();
};
post_proc_rhs_ptr->postProcessHook = get_post_proc_hook_rhs;
post_proc_lhs_ptr->postProcessHook = get_post_proc_hook_lhs;
ts_ctx_ptr->getPreProcessIFunction().push_front(pre_proc_ptr);
ts_ctx_ptr->getPreProcessIJacobian().push_front(pre_proc_ptr);
ts_ctx_ptr->getPostProcessIFunction().push_back(post_proc_rhs_ptr);
ts_ctx_ptr->getPostProcessIJacobian().push_back(post_proc_lhs_ptr);
};
CHKERR add_extra_finite_elements_to_ksp_solver_pipelines();
auto create_monitor_fe = [dm](auto &&post_proc_fe) {
return boost::make_shared<Monitor>(dm, post_proc_fe);
};
boost::shared_ptr<FEMethod> null_fe;
auto monitor_ptr = create_monitor_fe(create_post_proc_fe());
null_fe, monitor_ptr);
double ftime = 1;
CHKERR TSSetDuration(ts, PETSC_DEFAULT, ftime);
CHKERR TSSetExactFinalTime(ts, TS_EXACTFINALTIME_MATCHSTEP);
PetscInt steps, snesfails, rejects, nonlinits, linits;
CHKERR TSGetStepNumber(ts, &steps);
CHKERR TSGetSNESFailures(ts, &snesfails);
CHKERR TSGetStepRejections(ts, &rejects);
CHKERR TSGetSNESIterations(ts, &nonlinits);
CHKERR TSGetKSPIterations(ts, &linits);
"steps %d (%d rejected, %d SNES fails), ftime %g, nonlinits "
"%d, linits %d",
steps, rejects, snesfails, ftime, nonlinits, linits);
}
SCATTER_FORWARD);
double nrm2;
CHKERR VecNorm(T, NORM_2, &nrm2);
MOFEM_LOG(
"EXAMPLE", Sev::inform) <<
"Solution norm " << nrm2;
auto post_proc_norm_fe = boost::make_shared<DomainEle>(mField);
auto post_proc_norm_rule_hook = [](
int,
int,
int p) ->
int {
return 2 * p; };
post_proc_norm_fe->getRuleHook = post_proc_norm_rule_hook;
post_proc_norm_fe->getOpPtrVector(), {H1});
enum NORMS { U_NORM_L2 = 0, PIOLA_NORM, LAST_NORM };
auto norms_vec =
(mField.get_comm_rank() == 0) ? LAST_NORM : 0, LAST_NORM);
CHKERR VecZeroEntries(norms_vec);
auto u_ptr = boost::make_shared<MatrixDouble>();
post_proc_norm_fe->getOpPtrVector().push_back(
post_proc_norm_fe->getOpPtrVector().push_back(
auto common_ptr = commonDataFactory<SPACE_DIM, GAUSS, DomainEleOp>(
mField, post_proc_norm_fe->getOpPtrVector(), "U", "MAT_ELASTIC",
Sev::inform);
post_proc_norm_fe->getOpPtrVector().push_back(
common_ptr->getMatFirstPiolaStress(), norms_vec, PIOLA_NORM));
post_proc_norm_fe);
CHKERR VecAssemblyBegin(norms_vec);
CHKERR VecAssemblyEnd(norms_vec);
if (mField.get_comm_rank() == 0) {
const double *norms;
CHKERR VecGetArrayRead(norms_vec, &norms);
<< "norm_u: " << std::scientific << std::sqrt(norms[U_NORM_L2]);
<< "norm_piola: " << std::scientific << std::sqrt(norms[PIOLA_NORM]);
CHKERR VecRestoreArrayRead(norms_vec, &norms);
}
}
PetscInt test_nb = 0;
PetscBool test_flg = PETSC_FALSE;
if (test_flg) {
SCATTER_FORWARD);
double nrm2;
CHKERR VecNorm(T, NORM_2, &nrm2);
MOFEM_LOG(
"EXAMPLE", Sev::verbose) <<
"Regression norm " << nrm2;
double regression_value = 0;
switch (test_nb) {
case 1:
regression_value = 1.02789;
break;
case 2:
regression_value = 1.8841e+00;
break;
case 3:
regression_value = 1.8841e+00;
break;
case 4:
regression_value = 4.9625e+00;
break;
case 5:
regression_value = 6.6394e+00;
break;
case 6:
regression_value = 4.98764e+00;
break;
case 7:
regression_value = 4.9473e+00;
break;
case 8:
regression_value = 2.5749e-01;
break;
default:
break;
}
if (fabs(nrm2 - regression_value) > 1e-2)
"Regression test field; wrong norm value. %6.4e != %6.4e", nrm2,
regression_value);
}
}
}
static char help[] =
"...\n\n";
int main(
int argc,
char *argv[]) {
const char param_file[] = "param_file.petsc";
auto core_log = logging::core::get();
core_log->add_sink(
try {
DMType dm_name = "DMMOFEM";
}
}