13#ifndef EXECUTABLE_DIMENSION
14#define EXECUTABLE_DIMENSION 3
18#define SCHUR_ASSEMBLE 0
26#include <GenericElementInterface.hpp>
28#ifdef ENABLE_PYTHON_BINDING
29#include <boost/python.hpp>
30#include <boost/python/def.hpp>
31#include <boost/python/numpy.hpp>
32namespace bp = boost::python;
33namespace np = boost::python::numpy;
40 IntegrationType::GAUSS;
75 IT>::OpBaseTimesVector<1, SPACE_DIM, 1>;
110#ifdef WITH_MODULE_MFRONT_INTERFACE
111#include <MFrontMoFEMInterface.hpp>
149#ifdef ENABLE_PYTHON_BINDING
150 boost::shared_ptr<SDFPython> sdfPythonPtr;
194 <<
"Selected material model: Hooke (small strain)";
197 <<
"Selected deafult material model: Hencky (finite strain)";
205 enum bases { AINSWORTH, DEMKOWICZ, LASBASETOPT };
206 const char *list_bases[LASBASETOPT] = {
"ainsworth",
"demkowicz"};
207 PetscInt choice_base_value = AINSWORTH;
209 LASBASETOPT, &choice_base_value, PETSC_NULLPTR);
212 switch (choice_base_value) {
216 <<
"Set AINSWORTH_LEGENDRE_BASE for displacements";
221 <<
"Set DEMKOWICZ_JACOBI_BASE for displacements";
241 auto get_skin = [&]() {
246 CHKERR skin.find_skin(0, body_ents,
false, skin_ents);
250 auto filter_blocks = [&](
auto skin) {
251 bool is_contact_block =
false;
256 (boost::format(
"%s(.*)") %
"CONTACT").str()
265 <<
"Find contact block set: " <<
m->getName();
266 auto meshset =
m->getMeshset();
267 Range contact_meshset_range;
269 meshset,
SPACE_DIM - 1, contact_meshset_range,
true);
272 contact_meshset_range);
273 contact_range.merge(contact_meshset_range);
275 if (is_contact_block) {
277 <<
"Nb entities in contact surface: " << contact_range.size();
279 skin = intersect(skin, contact_range);
286 ParallelComm *pcomm =
288 CHKERR pcomm->filter_pstatus(skin, PSTATUS_SHARED | PSTATUS_MULTISHARED,
289 PSTATUS_NOT, -1, &boundary_ents);
290 return boundary_ents;
301 moab::Interface::UNION);
310 auto project_ho_geometry = [&]() {
315 PetscBool project_geometry = PETSC_TRUE;
317 &project_geometry, PETSC_NULLPTR);
318 if (project_geometry) {
319 CHKERR project_ho_geometry();
329 PetscBool use_mfront = PETSC_FALSE;
355 MOFEM_LOG(
"CONTACT", Sev::inform) <<
"Using MFront for material model";
366 PetscBool use_scale = PETSC_FALSE;
380#ifdef ENABLE_PYTHON_BINDING
381 auto file_exists = [](std::string myfile) {
382 std::ifstream file(myfile.c_str());
388 char sdf_file_name[255] =
"sdf.py";
390 sdf_file_name, 255, PETSC_NULLPTR);
392 if (file_exists(sdf_file_name)) {
393 MOFEM_LOG(
"CONTACT", Sev::inform) << sdf_file_name <<
" file found";
394 sdfPythonPtr = boost::make_shared<SDFPython>();
395 CHKERR sdfPythonPtr->sdfInit(sdf_file_name);
396 sdfPythonWeakPtr = sdfPythonPtr;
398 MOFEM_LOG(
"CONTACT", Sev::warning) << sdf_file_name <<
" file NOT found";
405 "Use executable contact_2d with axisymmetric model");
409 "Axisymmetric model is only available with MFront (set "
412 MOFEM_LOG(
"CONTACT", Sev::inform) <<
"Using axisymmetric model";
417 MOFEM_LOG(
"CONTACT", Sev::inform) <<
"Using plane strain model";
422#ifndef WITH_MODULE_MFRONT_INTERFACE
425 "MFrontInterface module was not found while use_mfront was set to 1");
429 boost::make_shared<MFrontMoFEMInterface<TRIDIMENSIONAL>>(
434 boost::make_shared<MFrontMoFEMInterface<AXISYMMETRICAL>>(
437 mfrontInterface = boost::make_shared<MFrontMoFEMInterface<PLANESTRAIN>>(
464 for (
auto f : {
"U",
"SIGMA"}) {
465 CHKERR bc_mng->removeBlockDOFsOnEntities(
simple->getProblemName(),
466 "REMOVE_X",
f, 0, 0);
467 CHKERR bc_mng->removeBlockDOFsOnEntities(
simple->getProblemName(),
468 "REMOVE_Y",
f, 1, 1);
469 CHKERR bc_mng->removeBlockDOFsOnEntities(
simple->getProblemName(),
470 "REMOVE_Z",
f, 2, 2);
471 CHKERR bc_mng->removeBlockDOFsOnEntities(
simple->getProblemName(),
472 "REMOVE_ALL",
f, 0, 3);
475 CHKERR bc_mng->removeBlockDOFsOnEntities(
simple->getProblemName(),
"FIX_X",
476 "SIGMA", 0, 0,
false,
true);
477 CHKERR bc_mng->removeBlockDOFsOnEntities(
simple->getProblemName(),
"FIX_Y",
478 "SIGMA", 1, 1,
false,
true);
479 CHKERR bc_mng->removeBlockDOFsOnEntities(
simple->getProblemName(),
"FIX_Z",
480 "SIGMA", 2, 2,
false,
true);
481 CHKERR bc_mng->removeBlockDOFsOnEntities(
simple->getProblemName(),
"FIX_ALL",
482 "SIGMA", 0, 3,
false,
true);
483 CHKERR bc_mng->removeBlockDOFsOnEntities(
484 simple->getProblemName(),
"NO_CONTACT",
"SIGMA", 0, 3,
false,
true);
489 simple->getProblemName(),
"U");
500 auto time_scale = boost::make_shared<ScaledTimeScale>();
501 auto body_force_time_scale =
502 boost::make_shared<ScaledTimeScale>(
"body_force_hist.txt");
504 auto integration_rule_vol = [](int, int,
int approx_order) {
507 auto integration_rule_boundary = [](int, int,
int approx_order) {
511 auto add_domain_base_ops = [&](
auto &pip) {
518 auto add_domain_ops_lhs = [&](
auto &pip) {
528 auto fe_domain_lhs = pip_mng->getDomainLhsFE();
530 auto get_inertia_and_mass_damping =
532 return (
rho *
scale) * fe_domain_lhs->ts_aa +
535 pip.push_back(
new OpMass(
"U",
"U", get_inertia_and_mass_damping));
539 auto fe_domain_lhs = pip_mng->getDomainLhsFE();
541 auto get_inertia_and_mass_damping =
545 pip.push_back(
new OpMass(
"U",
"U", get_inertia_and_mass_damping));
551 mField, pip,
"U",
"MAT_ELASTIC", Sev::verbose,
scale);
554 mField, pip,
"U",
"MAT_ELASTIC", Sev::verbose,
scale);
563 auto add_domain_ops_rhs = [&](
auto &pip) {
567 pip,
mField,
"U", {body_force_time_scale}, Sev::inform);
576 auto mat_acceleration = boost::make_shared<MatrixDouble>();
578 "U", mat_acceleration));
580 new OpInertiaForce(
"U", mat_acceleration, [](
double,
double,
double) {
587 auto mat_velocity = boost::make_shared<MatrixDouble>();
591 new OpInertiaForce(
"U", mat_velocity, [](
double,
double,
double) {
599 mField, pip,
"U",
"MAT_ELASTIC", Sev::inform, 1,
scale);
602 mField, pip,
"U",
"MAT_ELASTIC", Sev::inform,
scale);
614 auto add_boundary_base_ops = [&](
auto &pip) {
624 auto add_boundary_ops_lhs = [&](
auto &pip) {
634 pip,
mField,
"U", Sev::inform);
639 auto fe_boundary_lhs = pip_mng->getBoundaryLhsFE();
644 [
this, fe_boundary_lhs](
double,
double,
double) {
657 mField, pip,
simple->getDomainFEName(),
"SIGMA",
"U",
"GEOMETRY",
663 auto add_boundary_ops_rhs = [&](
auto &pip) {
673 pip,
mField,
"U", {time_scale}, Sev::inform);
676 auto u_disp = boost::make_shared<MatrixDouble>();
677 auto dot_u_disp = boost::make_shared<MatrixDouble>();
682 new OpSpringRhs(
"U", u_disp, [
this](
double,
double,
double) {
686 new OpSpringRhs(
"U", dot_u_disp, [
this](
double,
double,
double) {
698 CHKERR add_domain_base_ops(pip_mng->getOpDomainLhsPipeline());
699 CHKERR add_domain_base_ops(pip_mng->getOpDomainRhsPipeline());
700 CHKERR add_domain_ops_lhs(pip_mng->getOpDomainLhsPipeline());
701 CHKERR add_domain_ops_rhs(pip_mng->getOpDomainRhsPipeline());
703 CHKERR add_boundary_base_ops(pip_mng->getOpBoundaryLhsPipeline());
704 CHKERR add_boundary_base_ops(pip_mng->getOpBoundaryRhsPipeline());
705 CHKERR add_boundary_ops_lhs(pip_mng->getOpBoundaryLhsPipeline());
706 CHKERR add_boundary_ops_rhs(pip_mng->getOpBoundaryRhsPipeline());
712 CHKERR pip_mng->setDomainRhsIntegrationRule(integration_rule_vol);
713 CHKERR pip_mng->setDomainLhsIntegrationRule(integration_rule_vol);
714 CHKERR pip_mng->setBoundaryRhsIntegrationRule(integration_rule_boundary);
715 CHKERR pip_mng->setBoundaryLhsIntegrationRule(integration_rule_boundary);
723 static boost::shared_ptr<SetUpSchur>
738 auto set_section_monitor = [&](
auto solver) {
741 CHKERR TSGetSNES(solver, &snes);
742 PetscViewerAndFormat *vf;
743 CHKERR PetscViewerAndFormatCreate(PETSC_VIEWER_STDOUT_WORLD,
744 PETSC_VIEWER_DEFAULT, &vf);
747 (
MoFEMErrorCode (*)(SNES, PetscInt, PetscReal,
void *))SNESMonitorFields,
752 auto scatter_create = [&](
auto D,
auto coeff) {
754 CHKERR is_manager->isCreateProblemFieldAndRank(
simple->getProblemName(),
755 ROW,
"U", coeff, coeff, is);
757 CHKERR ISGetLocalSize(is, &loc_size);
761 CHKERR VecScatterCreate(
D, is,
v, PETSC_NULLPTR, &scatter);
766 auto set_time_monitor = [&](
auto dm,
auto solver) {
769 boost::shared_ptr<ForcesAndSourcesCore> null;
775 auto set_essential_bc = [&]() {
779 auto pre_proc_ptr = boost::make_shared<FEMethod>();
780 auto post_proc_rhs_ptr = boost::make_shared<FEMethod>();
781 auto post_proc_lhs_ptr = boost::make_shared<FEMethod>();
784 auto time_scale = boost::make_shared<TimeScale>();
786 auto get_bc_hook_rhs = [&]() {
788 {time_scale},
false);
791 pre_proc_ptr->preProcessHook = get_bc_hook_rhs();
793 auto get_post_proc_hook_rhs = [&]() {
795 mField, post_proc_rhs_ptr, 1.);
797 auto get_post_proc_hook_lhs = [&]() {
799 mField, post_proc_lhs_ptr, 1.);
801 post_proc_rhs_ptr->postProcessHook = get_post_proc_hook_rhs();
803 ts_ctx_ptr->getPreProcessIFunction().push_front(pre_proc_ptr);
804 ts_ctx_ptr->getPreProcessIJacobian().push_front(pre_proc_ptr);
805 ts_ctx_ptr->getPostProcessIFunction().push_back(post_proc_rhs_ptr);
806 post_proc_lhs_ptr->postProcessHook = get_post_proc_hook_lhs();
807 ts_ctx_ptr->getPostProcessIJacobian().push_back(post_proc_lhs_ptr);
812 auto set_schur_pc = [&](
auto solver) {
813 boost::shared_ptr<SetUpSchur> schur_ptr;
814 if (
AT == AssemblyType::BLOCK_SCHUR) {
822 auto dm =
simple->getDM();
834 CHKERR set_essential_bc();
837 auto solver = pip_mng->createTSIM();
838 CHKERR TSSetFromOptions(solver);
841 CHKERR TSSetIJacobian(solver,
B,
B, PETSC_NULLPTR, PETSC_NULLPTR);
842 auto schur_pc_ptr = set_schur_pc(solver);
845 CHKERR set_section_monitor(solver);
846 CHKERR set_time_monitor(dm, solver);
847 CHKERR TSSetSolution(solver,
D);
849 CHKERR TSSolve(solver, NULL);
851 auto solver = pip_mng->createTSIM2();
852 CHKERR TSSetFromOptions(solver);
855 CHKERR TSSetI2Jacobian(solver,
B,
B, PETSC_NULLPTR, PETSC_NULLPTR);
856 auto schur_pc_ptr = set_schur_pc(solver);
860 CHKERR set_section_monitor(solver);
861 CHKERR set_time_monitor(dm, solver);
862 CHKERR TS2SetSolution(solver,
D, DD);
864 CHKERR TSSolve(solver, NULL);
879 double analytical_active_area = 1.0;
881 double tol_force = 1e-3;
882 double tol_norm = 7.5;
883 double tol_area = 3e-2;
884 double fem_active_area = t_ptr[3];
889 fem_force = t_ptr[1];
894 fem_force = t_ptr[1];
901 fem_force = t_ptr[2];
902 analytical_active_area = M_PI / 4;
912 hertz_force = 15.873;
914 fem_force = t_ptr[1];
916 analytical_active_area = M_PI;
921 fem_force = t_ptr[1];
925 hertz_force = 0.5289;
926 fem_force = t_ptr[2];
931 "atom test %d does not exist",
atom_test);
933 if (fabs(fem_force - hertz_force) / hertz_force > tol_force) {
935 "atom test %d failed: Wrong FORCE output: %3.4e != %3.4e",
938 if (norm > tol_norm) {
940 "atom test %d failed: Wrong NORM output: %3.4e > %3.4e",
943 if (fabs(fem_active_area - analytical_active_area) > tol_area) {
945 "atom test %d failed: AREA computed %3.4e but should be %3.4e",
946 atom_test, fem_active_area, analytical_active_area);
959int main(
int argc,
char *argv[]) {
961#ifdef ENABLE_PYTHON_BINDING
967 const char param_file[] =
"param_file.petsc";
971 auto core_log = logging::core::get();
980 DMType dm_name =
"DMMOFEM";
982 DMType dm_name_mg =
"DMMOFEM_MG";
987 moab::Core mb_instance;
988 moab::Interface &moab = mb_instance;
1011#ifdef ENABLE_PYTHON_BINDING
1012 if (Py_FinalizeEx() < 0) {
1047 CHKERR TSGetSNES(solver, &snes);
1049 CHKERR SNESGetKSP(snes, &ksp);
1050 CHKERR KSPSetFromOptions(ksp);
1053 CHKERR KSPGetPC(ksp, &pc);
1055 PetscBool is_pcfs = PETSC_FALSE;
1056 PetscObjectTypeCompare((PetscObject)pc, PCFIELDSPLIT, &is_pcfs);
1059 MOFEM_LOG(
"CONTACT", Sev::inform) <<
"Setup Schur pc";
1064 "It is expected that Schur matrix is not allocated. This is "
1065 "possible only if PC is set up twice");
1077 CHKERR TSGetDM(solver, &solver_dm);
1078 CHKERR DMSetMatType(solver_dm, MATSHELL);
1085 auto swap_assemble = [](TS ts, PetscReal
t, Vec u, Vec u_t, PetscReal
a,
1086 Mat
A, Mat
B,
void *ctx) {
1089 CHKERR TSSetIJacobian(solver, A, P, swap_assemble, ts_ctx_ptr.get());
1091 auto swap_assemble = [](TS ts, PetscReal
t, Vec u, Vec u_t, Vec utt,
1092 PetscReal
a, PetscReal aa, Mat
A, Mat
B,
1096 CHKERR TSSetI2Jacobian(solver, A, P, swap_assemble, ts_ctx_ptr.get());
1098 CHKERR KSPSetOperators(ksp, A, P);
1107 MOFEM_LOG(
"CONTACT", Sev::inform) <<
"No Schur PC";
1120 auto create_dm = [&](
const char *name,
const char *
field_name,
auto dm_type) {
1122 auto create_dm_imp = [&]() {
1134 "Error in creating schurDM. It is possible that schurDM is "
1142 schurDM = create_dm(
"SCHUR",
"U",
"DMMOFEM_MG");
1143 blockDM = create_dm(
"BLOCK",
"SIGMA",
"DMMOFEM");
1145 if constexpr (
AT == AssemblyType::BLOCK_SCHUR) {
1147 auto get_nested_mat_data = [&](
auto schur_dm,
auto block_dm) {
1153 simple->getDomainFEName(),
1157 {
"U",
"U"}, {
"SIGMA",
"U"}, {
"U",
"SIGMA"}, {
"SIGMA",
"SIGMA"}
1165 {schur_dm, block_dm}, block_mat_data,
1167 {
"SIGMA"}, {
nullptr}, true
1177 "Only BLOCK_SCHUR is implemented");
1186 double eps_stab = 1e-4;
1192 using OpMassStab = B::OpMass<3, SPACE_DIM * SPACE_DIM>;
1202 pip->getOpBoundaryLhsPipeline().push_back(
1203 new OpMassStab(
"SIGMA",
"SIGMA",
1204 [eps_stab](
double,
double,
double) {
return eps_stab; }));
1205 pip->getOpBoundaryLhsPipeline().push_back(
1213 pip->getOpDomainLhsPipeline().push_back(
1219 auto pre_proc_schur_lhs_ptr = boost::make_shared<FEMethod>();
1220 auto post_proc_schur_lhs_ptr = boost::make_shared<FEMethod>();
1222 pre_proc_schur_lhs_ptr->preProcessHook = [
this]() {
1225 MOFEM_LOG(
"CONTACT", Sev::verbose) <<
"Lhs Assemble Begin";
1229 post_proc_schur_lhs_ptr->postProcessHook = [
this, ao_up,
1230 post_proc_schur_lhs_ptr]() {
1232 MOFEM_LOG(
"CONTACT", Sev::verbose) <<
"Lhs Assemble End";
1233 auto print_mat_norm = [
this](
auto a, std::string prefix) {
1236 CHKERR MatNorm(
a, NORM_FROBENIUS, &nrm);
1237 MOFEM_LOG(
"CONTACT", Sev::noisy) << prefix <<
" norm = " << nrm;
1240 CHKERR MatAssemblyBegin(
S, MAT_FINAL_ASSEMBLY);
1241 CHKERR MatAssemblyEnd(
S, MAT_FINAL_ASSEMBLY);
1243 mField, post_proc_schur_lhs_ptr, 1,
S, ao_up)();
1245 CHKERR print_mat_norm(
S,
"S");
1247 MOFEM_LOG(
"CONTACT", Sev::verbose) <<
"Lhs Assemble Finish";
1252 ts_ctx_ptr->getPreProcessIJacobian().push_front(pre_proc_schur_lhs_ptr);
1253 ts_ctx_ptr->getPostProcessIJacobian().push_back(post_proc_schur_lhs_ptr);
1261 CHKERR PCFieldSplitSetIS(pc, NULL, block_is);
1262 CHKERR PCFieldSplitSetSchurPre(pc, PC_FIELDSPLIT_SCHUR_PRE_USER,
S);
1269 CHKERR PCFieldSplitSchurGetSubKSP(pc, PETSC_NULLPTR, &subksp);
1270 auto get_pc = [](
auto ksp) {
1272 CHKERR KSPGetPC(ksp, &pc_raw);
1277 auto set_pc_p_mg = [](
auto dm,
auto pc,
auto S) {
1280 PetscBool same = PETSC_FALSE;
1281 PetscObjectTypeCompare((PetscObject)pc, PCMG, &same);
1285 CHKERR PCSetFromOptions(pc);
1290 auto set_pc_ksp = [&](
auto dm,
auto pc,
auto S) {
1292 PetscBool same = PETSC_FALSE;
1293 PetscObjectTypeCompare((PetscObject)pc, PCKSP, &same);
1295 CHKERR PCSetFromOptions(pc);
1297 CHKERR PCKSPGetKSP(pc, &inner_ksp);
1298 CHKERR KSPSetFromOptions(inner_ksp);
1300 CHKERR KSPGetPC(inner_ksp, &inner_pc);
1301 CHKERR PCSetFromOptions(inner_pc);
1302 CHKERR set_pc_p_mg(dm, inner_pc,
S);
1310 CHKERR PetscFree(subksp);
1314boost::shared_ptr<SetUpSchur>
1316 return boost::shared_ptr<SetUpSchur>(
new SetUpSchurImpl(m_field));
static auto filter_true_skin(MoFEM::Interface &m_field, Range &&skin)
Implementation of Hookes operator Hookes for linear elastic problems in MoFEM.
#define MOFEM_LOG_SYNCHRONISE(comm)
Synchronise "SYNC" channel.
void simple(double P1[], double P2[], double P3[], double c[], const int N)
ElementsAndOps< SPACE_DIM >::BoundaryEle BoundaryEle
ElementsAndOps< SPACE_DIM >::DomainEle DomainEle
[Define dimension]
#define CATCH_ERRORS
Catch errors.
FieldApproximationBase
approximation base
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base nme:nme847.
#define CHK_THROW_MESSAGE(err, msg)
Check and throw MoFEM exception.
FieldSpace
approximation spaces
@ HCURL
field with continuous tangents
@ HDIV
field with continuous normal traction
#define MYPCOMM_INDEX
default communicator number PCOMM
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define CHK_MOAB_THROW(err, msg)
Check error code of MoAB function and throw MoFEM exception.
@ MOFEM_ATOM_TEST_INVALID
@ MOFEM_DATA_INCONSISTENCY
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpBaseTimesVector< 1, SPACE_DIM, 1 > OpInertiaForce
PetscErrorCode DMMoFEMCreateSubDM(DM subdm, DM dm, const char problem_name[])
Must be called by user to set Sub DM MoFEM data structures.
PetscErrorCode DMMoFEMAddElement(DM dm, std::string fe_name)
add element to dm
PetscErrorCode DMMoFEMSetSquareProblem(DM dm, PetscBool square_problem)
set squared problem
PetscErrorCode DMMoFEMAddSubFieldRow(DM dm, const char field_name[])
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
MoFEMErrorCode DMRegister_MGViaApproxOrders(const char sname[])
Register DM for Multi-Grid via approximation orders.
auto createDMVector(DM dm)
Get smart vector from DM.
PetscErrorCode DMMoFEMAddSubFieldCol(DM dm, const char field_name[])
auto createDMMatrix(DM dm)
Get smart matrix from DM.
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.
virtual MoFEMErrorCode loop_dofs(const Problem *problem_ptr, const std::string &field_name, RowColData rc, DofMethod &method, int lower_rank, int upper_rank, int verb=DEFAULT_VERBOSITY)=0
Make a loop over dofs.
MoFEMErrorCode getCubitMeshsetPtr(const int ms_id, const CubitBCType cubit_bc_type, const CubitMeshSets **cubit_meshset_ptr) const
get cubit meshset
const double v
phase velocity of light in medium (cm/ns)
MoFEMErrorCode opFactoryDomainLhs(MoFEM::Interface &m_field, boost::ptr_deque< ForcesAndSourcesCore::UserDataOperator > &pip, std::string field_name, boost::shared_ptr< HenckyOps::CommonData > common_ptr, Sev sev)
MoFEMErrorCode opFactoryDomainRhs(MoFEM::Interface &m_field, boost::ptr_deque< ForcesAndSourcesCore::UserDataOperator > &pip, std::string field_name, boost::shared_ptr< HenckyOps::CommonData > common_ptr, Sev sev)
MoFEMErrorCode opFactoryDomainLhs(MoFEM::Interface &m_field, boost::ptr_deque< ForcesAndSourcesCore::UserDataOperator > &pip, std::string field_name, boost::shared_ptr< HookeOps::CommonData > common_ptr, Sev sev)
Assemble domain LHS K factory (LHS first overload) Initializes and pushes operators to assemble the L...
MoFEMErrorCode opFactoryDomainRhs(MoFEM::Interface &m_field, boost::ptr_deque< ForcesAndSourcesCore::UserDataOperator > &pip, std::string field_name, boost::shared_ptr< HookeOps::CommonData > common_ptr, Sev sev, bool is_non_linear=false)
Factory function to create and push internal force operators for the domain RHS. (RHS first overload)...
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
implementation of Data Operators for Forces and Sources
PetscErrorCode DMMoFEMTSSetMonitor(DM dm, TS ts, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
Set Monitor To TS solver.
PetscErrorCode TsSetIJacobian(TS ts, PetscReal t, Vec u, Vec u_t, PetscReal a, Mat A, Mat B, void *ctx)
Set function evaluating jacobian in TS solver.
auto getDMTsCtx(DM dm)
Get TS context data structure used by DM.
OpSchurAssembleBase * createOpSchurAssembleEnd(std::vector< std::string > fields_name, std::vector< boost::shared_ptr< Range > > field_ents, SmartPetscObj< AO > ao, SmartPetscObj< Mat > schur, bool sym_schur, bool symm_op)
Construct a new Op Schur Assemble End object.
PetscErrorCode PetscOptionsGetInt(PetscOptions *, const char pre[], const char name[], PetscInt *ivalue, PetscBool *set)
MoFEMErrorCode setSchurA00MatSolvePC(SmartPetscObj< PC > pc)
Set PC for A00 block.
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.
boost::shared_ptr< PCMGSetUpViaApproxOrdersCtx > createPCMGSetUpViaApproxOrdersCtx(DM dm, Mat A, bool use_shell_mat)
createPCMGSetUpViaApproxOrdersCtx
auto getDMSubData(DM dm)
Get sub problem data structure.
PetscErrorCode TsSetI2Jacobian(TS ts, PetscReal t, Vec u, Vec u_t, Vec u_tt, PetscReal a, PetscReal aa, Mat A, Mat B, void *ctx)
Calculation Jacobian for second order PDE in time.
MoFEMErrorCode PCMGSetUpViaApproxOrders(PC pc, boost::shared_ptr< PCMGSetUpViaApproxOrdersCtx > ctx, int verb)
Function build MG structure.
boost::shared_ptr< BlockStructure > createBlockMatStructure(DM dm, SchurFEOpsFEandFields schur_fe_op_vec)
Create a Mat Diag Blocks object.
boost::shared_ptr< NestSchurData > createSchurNestedMatrixStruture(std::pair< SmartPetscObj< DM >, SmartPetscObj< DM > > dms, boost::shared_ptr< BlockStructure > block_mat_data_ptr, std::vector< std::string > fields_names, std::vector< boost::shared_ptr< Range > > field_ents, bool add_preconditioner_block)
Get the Schur Nest Mat Array object.
auto createAOMappingIS(IS isapp, IS ispetsc)
Creates an application mapping using two index sets.
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)
MoFEMErrorCode DMMoFEMSetNestSchurData(DM dm, boost::shared_ptr< NestSchurData >)
auto createDMNestSchurMat(DM dm)
auto createDM(MPI_Comm comm, const std::string dm_type_name)
Creates smart DM object.
OpSchurAssembleBase * createOpSchurAssembleBegin()
auto createDMBlockMat(DM dm)
PetscBool is_quasi_static
ElementsAndOps< SPACE_DIM >::SideEle SideEle
constexpr FieldSpace CONTACT_SPACE
constexpr double t
plate stiffness
constexpr auto field_name
static constexpr int approx_order
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpMass< 1, SPACE_DIM > OpMass
[Only used with Hooke equation (linear material model)]
FTensor::Index< 'm', 3 > m
Add operators pushing bases from local to physical configuration.
Simple interface for fast problem set-up.
virtual moab::Interface & get_moab()=0
virtual MPI_Comm & get_comm() const =0
virtual int get_comm_rank() 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.
Definition of the displacement bc data structure.
Data on single entity (This is passed as argument to DataOperator::doWork)
Class (Function) to enforce essential constrains on the left hand side diagonal.
Class (Function) to enforce essential constrains on the right hand side diagonal.
Class (Function) to enforce essential constrains.
Section manager is used to create indexes and sections.
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.
Interface for managing meshsets containing materials and boundary conditions.
Approximate field values for given petsc vector.
Get values at integration pts for tensor field rank 1, i.e. vector field.
PipelineManager interface.
Projection of edge entities with one mid-node on hierarchical basis.
Simple interface for fast problem set-up.
MoFEMErrorCode getOptions()
get options
MoFEMErrorCode getDM(DM *dm)
Get DM.
intrusive_ptr for managing petsc objects
Force scale operator for reading two columns.
double getScale(const double time)
Get scaling at a given time.
TimeScale(std::string file_name="", bool error_if_file_not_given=false, ScalingFun def_scaling_fun=[](double time) { return time;})
TimeScale constructor.
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.
MoFEMErrorCode createSubDM()
MoFEMErrorCode setDiagonalPC(PC pc)
SetUpSchurImpl(MoFEM::Interface &m_field)
SmartPetscObj< DM > schurDM
MoFEMErrorCode setUp(SmartPetscObj< KSP >)
SmartPetscObj< DM > blockDM
virtual ~SetUpSchurImpl()
MoFEMErrorCode setPC(PC pc)
MoFEMErrorCode setOperator()
MoFEM::Interface & mField
[Push operators to pipeline]
static boost::shared_ptr< SetUpSchur > createSetUpSchur(MoFEM::Interface &m_field)
virtual MoFEMErrorCode setUp(SmartPetscObj< TS > solver)=0
static boost::shared_ptr< SetUpSchur > createSetUpSchur(MoFEM::Interface &m_field)
constexpr AssemblyType AT