752 {
754
758
760
761 auto set_section_monitor = [&](auto solver) {
763 SNES snes;
764 CHKERR TSGetSNES(solver, &snes);
765 CHKERR SNESMonitorSet(snes,
768 (void *)(snes_ctx_ptr.get()), nullptr);
770 };
771
772 auto create_post_process_element = [&]() {
773 auto pp_fe = boost::make_shared<PostProcEle>(
mField);
775
777 pp_fe->getOpPtrVector(), {H1}, "GEOMETRY");
778
779 auto x_ptr = boost::make_shared<MatrixDouble>();
780 pp_fe->getOpPtrVector().push_back(
782 auto u_ptr = boost::make_shared<MatrixDouble>();
783 pp_fe->getOpPtrVector().push_back(
785
786 pp_fe->getOpPtrVector().push_back(
791 pp_fe->getOpPtrVector().push_back(
794
796
799
800 pp_fe->getOpPtrVector().push_back(
802 pp_fe->getOpPtrVector().push_back(
804 pp_fe->getOpPtrVector().push_back(
806 pp_fe->getOpPtrVector().push_back(
809 pp_fe->getOpPtrVector().push_back(
811
812 pp_fe->getOpPtrVector().push_back(
813
815
816 pp_fe->getPostProcMesh(), pp_fe->getMapGaussPts(),
817
818 {},
819
820 {{"U", u_ptr}, {"GEOMETRY", x_ptr}},
821
824
825 {}
826
827 )
828
829 );
830
831 } else {
832 pp_fe->getOpPtrVector().push_back(
837
838 pp_fe->getOpPtrVector().push_back(
839
841
842 pp_fe->getPostProcMesh(), pp_fe->getMapGaussPts(),
843
844 {},
845
846 {{"U", u_ptr}, {"GEOMETRY", x_ptr}},
847
848 {},
849
852
853 )
854
855 );
856 }
857
858 pp_fe->getOpPtrVector().push_back(
860
861 pp_fe->getOpPtrVector().push_back(
862
864
865 pp_fe->getPostProcMesh(), pp_fe->getMapGaussPts(),
866
867 {{"PLASTIC_SURFACE", commonPlasticDataPtr->getPlasticSurfacePtr()},
868 {"PLASTIC_MULTIPLIER", commonPlasticDataPtr->getPlasticTauPtr()}},
869
870 {},
871
872 {},
873
876
877 )
878
879 );
880
881 return pp_fe;
882 };
883
884 auto scatter_create = [&](
auto D,
auto coeff) {
886 CHKERR is_manager->isCreateProblemFieldAndRank(
simple->getProblemName(),
887 ROW,
"U", coeff, coeff, is);
888 int loc_size;
889 CHKERR ISGetLocalSize(is, &loc_size);
892 VecScatter scatter;
893 CHKERR VecScatterCreate(
D, is,
v, PETSC_NULL, &scatter);
896 };
897
898 auto set_time_monitor = [&](auto dm, auto solver) {
900 boost::shared_ptr<Monitor> monitor_ptr(
903 boost::shared_ptr<ForcesAndSourcesCore> null;
905 monitor_ptr, null, null);
907 };
908
909 auto set_fieldsplit_preconditioner = [&](auto solver,
910 boost::shared_ptr<SetUpSchur>
911 &schur_ptr) {
913
914 SNES snes;
915 CHKERR TSGetSNES(solver, &snes);
916 KSP ksp;
917 CHKERR SNESGetKSP(snes, &ksp);
918 PC pc;
919 CHKERR KSPGetPC(ksp, &pc);
920 PetscBool is_pcfs = PETSC_FALSE;
921 PetscObjectTypeCompare((PetscObject)pc, PCFIELDSPLIT, &is_pcfs);
922
923
924 if (is_pcfs == PETSC_TRUE) {
925
927 auto name_prb =
simple->getProblemName();
928
934
940 for (
auto f : {
"U",
"EP",
"TAU"}) {
943 }
945
946 CHKERR bc_mng->removeBlockDOFsOnEntities(
"SUB_BC",
"FIX_X",
"U", 0, 0);
947 CHKERR bc_mng->removeBlockDOFsOnEntities(
"SUB_BC",
"FIX_Y",
"U", 1, 1);
948 CHKERR bc_mng->removeBlockDOFsOnEntities(
"SUB_BC",
"FIX_Z",
"U", 2, 2);
949 CHKERR bc_mng->removeBlockDOFsOnEntities(
"SUB_BC",
"FIX_ALL",
"U", 0,
950 2);
951
953 if (auto sub_data = prb_ptr->getSubData()) {
954 is_sub = sub_data->getSmartRowIs();
955 ao_sub = sub_data->getSmartRowMap();
956 int is_sub_size;
957 CHKERR ISGetSize(is_sub, &is_sub_size);
959 << "Field split second block size " << is_sub_size;
960
961 } else {
963 }
964
966 };
967
978 for (
auto f : {
"U"}) {
981 }
983
985 if (auto sub_data = prb_ptr->getSubData()) {
986 is_sub = sub_data->getSmartRowIs();
987 ao_sub = sub_data->getSmartRowMap();
988
989 int is_sub_size;
990 CHKERR ISGetSize(is_sub, &is_sub_size);
992 << "Field split second block size " << is_sub_size;
993
994 } else {
996 }
997
1000 };
1001
1004 is_all_bc = bc_mng->getBlockIS(name_prb, "FIX_X", "U", 0, 0);
1005 is_all_bc = bc_mng->getBlockIS(name_prb, "FIX_Y", "U", 1, 1, is_all_bc);
1006 is_all_bc = bc_mng->getBlockIS(name_prb, "FIX_Z", "U", 2, 2, is_all_bc);
1007 is_all_bc =
1008 bc_mng->getBlockIS(name_prb, "FIX_ALL", "U", 0, 2, is_all_bc);
1009 int is_all_bc_size;
1010 CHKERR ISGetSize(is_all_bc, &is_all_bc_size);
1012 << "Field split first block size " << is_all_bc_size;
1014 };
1015
1020
1021 CHKERR create_all_bc_is(is_all_bc);
1022 CHKERR create_sub_bc_dm(
simple->getDM(), dm_bc_sub, is_bc_sub, ao_bc_sub);
1023
1027
1028 CHKERR PCFieldSplitSetIS(pc, PETSC_NULL,
1029 is_all_bc);
1030 CHKERR PCFieldSplitSetIS(pc, PETSC_NULL, is_bc_sub);
1031
1032 if constexpr (
A == AssemblyType::SCHUR) {
1033
1040 IS is_union_raw;
1041 CHKERR ISExpand(is_epp, is_tau, &is_union_raw);
1043
1047 CHKERR create_sub_u_dm(dm_bc_sub, dm_u_sub, is_u_sub, ao_u_sub);
1048
1049
1050
1051 auto is_up = is_u_sub;
1052 CHKERR AOPetscToApplicationIS(ao_bc_sub, is_up);
1054 schur_ptr =
1057 KSP *ksps;
1058 CHKERR PCFieldSplitGetSubKSP(pc, &
n, &ksps);
1059 CHKERR schur_ptr->setUp(ksps[1]);
1060 }
1061 }
1062
1064 };
1065
1066 auto dm =
simple->getDM();
1072
1073 auto solver = pip->createTSIM();
1074
1075 auto active_pre_lhs = [&]() {
1080 };
1081
1082 auto active_post_lhs = [&]() {
1084 auto get_iter = [&]() {
1085 SNES snes;
1087 int iter;
1089 "Can not get iter");
1090 return iter;
1091 };
1092
1093 auto iter = get_iter();
1094 if (iter >= 0) {
1095
1096 std::array<int, 5> activity_data;
1097 std::fill(activity_data.begin(), activity_data.end(), 0);
1099 activity_data.data(), activity_data.size(), MPI_INT,
1101
1102 int &active_points = activity_data[0];
1103 int &avtive_full_elems = activity_data[1];
1104 int &avtive_elems = activity_data[2];
1105 int &nb_points = activity_data[3];
1106 int &nb_elements = activity_data[4];
1107
1108 if (nb_points) {
1109
1110 double proc_nb_points =
1111 100 * static_cast<double>(active_points) / nb_points;
1112 double proc_nb_active =
1113 100 * static_cast<double>(avtive_elems) / nb_elements;
1114 double proc_nb_full_active = 100;
1115 if (avtive_elems)
1116 proc_nb_full_active =
1117 100 * static_cast<double>(avtive_full_elems) / avtive_elems;
1118
1120 "EXAMPLE", Sev::inform,
1121 "Iter %d nb pts %d nb avtive pts %d (%3.3f\%) nb active elements %d "
1122 "(%3.3f\%) nb full active elems %d (%3.3f\%)",
1123 iter, nb_points, active_points, proc_nb_points, avtive_elems,
1124 proc_nb_active, avtive_full_elems, proc_nb_full_active, iter);
1125 }
1126 }
1127
1129 };
1130
1131 CHKERR TSSetSolution(solver,
D);
1132 CHKERR set_section_monitor(solver);
1133 CHKERR set_time_monitor(dm, solver);
1134 CHKERR TSSetSolution(solver,
D);
1135 CHKERR TSSetFromOptions(solver);
1136
1137 boost::shared_ptr<SetUpSchur> schur_ptr;
1138 CHKERR set_fieldsplit_preconditioner(solver, schur_ptr);
1139
1141 [&]() {
1143 if (schur_ptr)
1144 CHKERR schur_ptr->preProc();
1147 };
1149 [&]() {
1151
1152
1153 CHKERR active_post_lhs();
1155 };
1156
1158 [&]() {
1160 if (schur_ptr)
1161 CHKERR schur_ptr->postProc();
1163 };
1164
1168 BOOST_LOG_SCOPED_THREAD_ATTR("Timeline", attrs::timer());
1169 MOFEM_LOG(
"TIMER", Sev::verbose) <<
"TSSetUp";
1171 MOFEM_LOG(
"TIMER", Sev::verbose) <<
"TSSetUp <= done";
1172 MOFEM_LOG(
"TIMER", Sev::verbose) <<
"TSSolve";
1173 CHKERR TSSolve(solver, NULL);
1174 MOFEM_LOG(
"TIMER", Sev::verbose) <<
"TSSolve <= done";
1175
1176 CHKERR VecGhostUpdateBegin(
D, INSERT_VALUES, SCATTER_FORWARD);
1177 CHKERR VecGhostUpdateEnd(
D, INSERT_VALUES, SCATTER_FORWARD);
1179
1181}
#define MAX_DOFS_ON_ENTITY
Maximal number of DOFs on entity.
FTensor::Index< 'n', SPACE_DIM > n
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 DMMoFEMAddSubFieldCol(DM dm, const char field_name[])
#define MOFEM_LOG_TAG(channel, tag)
Tag channel.
const FTensor::Tensor2< T, Dim, Dim > Vec
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.
MoFEMErrorCode MoFEMSNESMonitorFields(SNES snes, PetscInt its, PetscReal fgnorm, SnesCtx *snes_ctx)
Sens monitor printing residual field by field.
auto createSmartDM(MPI_Comm comm, const std::string dm_type_name)
Creates smart DM object.
auto createAOMappingIS(IS isapp, IS ispetsc)
Creates an application mapping using two index sets.
auto getProblemPtr(DM dm)
get problem pointer from DM
auto smartGetDMSnesCtx(DM dm)
Get SNES context data structure used by DM.
std::tuple< SmartPetscObj< Vec >, SmartPetscObj< VecScatter > > uZScatter
std::tuple< SmartPetscObj< Vec >, SmartPetscObj< VecScatter > > uXScatter
std::tuple< SmartPetscObj< Vec >, SmartPetscObj< VecScatter > > uYScatter
Section manager is used to create indexes and sections.
intrusive_ptr for managing petsc objects
static boost::shared_ptr< SetUpSchur > createSetUpSchur(MoFEM::Interface &m_field, SmartPetscObj< DM >, SmartPetscObj< IS >, SmartPetscObj< AO >)
