799 auto dm =
simple->getDM();
800 auto solver = pipeline_mng->createTSIM();
803 auto set_section_monitor = [&](
auto solver) {
806 CHKERR TSGetSNES(solver, &snes);
807 CHKERR SNESMonitorSet(snes,
810 (
void *)(snes_ctx_ptr.get()),
nullptr);
814 auto create_post_process_element = [&]() {
815 auto post_proc_fe = boost::make_shared<PostProcEle>(
mField);
817 auto block_params = boost::make_shared<BlockedParameters>();
818 auto mDPtr = block_params->getDPtr();
820 "MAT_FLUID", block_params, Sev::verbose);
822 post_proc_fe->getOpPtrVector(), {H1, HDIV});
824 auto mat_grad_ptr = boost::make_shared<MatrixDouble>();
825 auto mat_strain_ptr = boost::make_shared<MatrixDouble>();
826 auto mat_stress_ptr = boost::make_shared<MatrixDouble>();
828 auto h_ptr = boost::make_shared<VectorDouble>();
829 auto mat_flux_ptr = boost::make_shared<MatrixDouble>();
831 post_proc_fe->getOpPtrVector().push_back(
833 post_proc_fe->getOpPtrVector().push_back(
836 auto u_ptr = boost::make_shared<MatrixDouble>();
837 post_proc_fe->getOpPtrVector().push_back(
839 post_proc_fe->getOpPtrVector().push_back(
842 post_proc_fe->getOpPtrVector().push_back(
844 post_proc_fe->getOpPtrVector().push_back(
846 "U", mat_strain_ptr, mat_stress_ptr, mDPtr));
850 post_proc_fe->getOpPtrVector().push_back(
854 post_proc_fe->getPostProcMesh(), post_proc_fe->getMapGaussPts(),
858 {{
"U", u_ptr}, {
"FLUX", mat_flux_ptr}},
862 {{
"STRAIN", mat_strain_ptr}, {
"STRESS", mat_stress_ptr}}
871 auto create_creaction_fe = [&]() {
872 auto fe_ptr = boost::make_shared<DomainEle>(
mField);
873 fe_ptr->getRuleHook = [](
int,
int,
int o) {
return 2 * o; };
875 auto &pip = fe_ptr->getOpPtrVector();
877 auto block_params = boost::make_shared<BlockedParameters>();
878 auto mDPtr = block_params->getDPtr();
883 auto u_grad_ptr = boost::make_shared<MatrixDouble>();
884 auto strain_ptr = boost::make_shared<MatrixDouble>();
885 auto stress_ptr = boost::make_shared<MatrixDouble>();
893 "U", strain_ptr, stress_ptr, mDPtr));
896 fe_ptr->postProcessHook =
902 auto monitor_ptr = boost::make_shared<FEMethod>();
903 auto post_proc_fe = create_post_process_element();
905 auto rections_fe = create_creaction_fe();
907 auto set_time_monitor = [&](
auto dm,
auto solver) {
909 monitor_ptr->preProcessHook = [&]() {
914 monitor_ptr->getCacheWeakPtr());
915 CHKERR post_proc_fe->writeFile(
916 "out_" + boost::lexical_cast<std::string>(monitor_ptr->ts_step) +
919 rections_fe->f = res;
922 monitor_ptr->getCacheWeakPtr());
925 CHKERR VecNorm(res, NORM_2, &nrm);
927 <<
"Residual norm " << nrm <<
" at time step "
928 << monitor_ptr->ts_step;
932 auto scalar_field_ptr = boost::make_shared<VectorDouble>();
933 auto vector_field_ptr = boost::make_shared<MatrixDouble>();
934 auto tensor_field_ptr = boost::make_shared<MatrixDouble>();
937 ->getData<DomainEle>();
941 field_eval_data,
simple->getDomainFEName());
944 field_eval_data,
simple->getDomainFEName());
948 auto no_rule = [](
int,
int,
int) {
return -1; };
950 auto field_eval_ptr = field_eval_data->feMethodPtr.lock();
951 field_eval_ptr->getRuleHook = no_rule;
953 auto u_grad_ptr = boost::make_shared<MatrixDouble>();
954 auto strain_ptr = boost::make_shared<MatrixDouble>();
955 auto stress_ptr = boost::make_shared<MatrixDouble>();
956 auto h_ptr = boost::make_shared<VectorDouble>();
958 auto block_params = boost::make_shared<BlockedParameters>();
959 auto mDPtr = block_params->getDPtr();
961 "MAT_FLUID", block_params, Sev::noisy);
963 field_eval_ptr->getOpPtrVector(), {H1, HDIV});
964 field_eval_ptr->getOpPtrVector().push_back(
967 field_eval_ptr->getOpPtrVector().push_back(
969 field_eval_ptr->getOpPtrVector().push_back(
971 field_eval_ptr->getOpPtrVector().push_back(
973 "U", strain_ptr, stress_ptr, mDPtr));
977 ->evalFEAtThePoint3D(
979 simple->getDomainFEName(), field_eval_data,
984 ->evalFEAtThePoint2D(
986 simple->getDomainFEName(), field_eval_data,
992 <<
"Eval point hydrostatic hight: " << *h_ptr;
994 <<
"Eval point skeleton stress pressure: " << *stress_ptr;
1000 auto null = boost::shared_ptr<FEMethod>();
1006 auto set_fieldsplit_preconditioner = [&](
auto solver) {
1010 CHKERR TSGetSNES(solver, &snes);
1012 CHKERR SNESGetKSP(snes, &ksp);
1014 CHKERR KSPGetPC(ksp, &pc);
1015 PetscBool is_pcfs = PETSC_FALSE;
1016 PetscObjectTypeCompare((PetscObject)pc, PCFIELDSPLIT, &is_pcfs);
1019 if (is_pcfs == PETSC_TRUE) {
1022 auto name_prb =
simple->getProblemName();
1025 CHKERR is_mng->isCreateProblemFieldAndRank(name_prb,
ROW,
"U", 0,
1028 CHKERR is_mng->isCreateProblemFieldAndRank(name_prb,
ROW,
"FLUX", 0, 0,
1031 CHKERR is_mng->isCreateProblemFieldAndRank(name_prb,
ROW,
"H", 0, 0,
1034 CHKERR ISExpand(is_H, is_flux, &is_tmp);
1037 auto is_all_bc = bc_mng->getBlockIS(name_prb,
"FLUIDFLUX",
"FLUX", 0, 0);
1039 CHKERR ISGetSize(is_all_bc, &is_all_bc_size);
1041 <<
"Field split block size " << is_all_bc_size;
1042 if (is_all_bc_size) {
1044 CHKERR ISDifference(is_Flux, is_all_bc, &is_tmp2);
1046 CHKERR PCFieldSplitSetIS(pc, PETSC_NULL,
1052 CHKERR PCFieldSplitSetIS(pc, PETSC_NULL, is_Flux);
1053 CHKERR PCFieldSplitSetIS(pc, PETSC_NULL, is_u);
1059 auto pre_proc_ptr = boost::make_shared<FEMethod>();
1060 auto post_proc_rhs_ptr = boost::make_shared<FEMethod>();
1061 auto post_proc_lhs_ptr = boost::make_shared<FEMethod>();
1062 auto time_scale = boost::make_shared<TimeScale>();
1064 auto get_bc_hook_rhs = [
this, pre_proc_ptr, time_scale]() {
1066 {time_scale},
false);
1070 auto get_post_proc_hook_rhs = [
this, post_proc_rhs_ptr]() {
1073 mField, post_proc_rhs_ptr,
nullptr, Sev::verbose)();
1075 mField, post_proc_rhs_ptr, 1.)();
1078 auto get_post_proc_hook_lhs = [
this, post_proc_lhs_ptr]() {
1080 post_proc_lhs_ptr, 1.);
1083 pre_proc_ptr->preProcessHook = get_bc_hook_rhs();
1084 post_proc_rhs_ptr->postProcessHook = get_post_proc_hook_rhs;
1085 post_proc_lhs_ptr->postProcessHook = get_post_proc_hook_lhs();
1088 ts_ctx_ptr->getPreProcessIFunction().push_front(pre_proc_ptr);
1089 ts_ctx_ptr->getPreProcessIJacobian().push_front(pre_proc_ptr);
1090 ts_ctx_ptr->getPostProcessIFunction().push_back(post_proc_rhs_ptr);
1091 ts_ctx_ptr->getPostProcessIJacobian().push_back(post_proc_lhs_ptr);
1094 CHKERR TSSetSolution(solver,
D);
1095 CHKERR TSSetFromOptions(solver);
1097 CHKERR set_section_monitor(solver);
1098 CHKERR set_fieldsplit_preconditioner(solver);
1099 CHKERR set_time_monitor(dm, solver);
1102 CHKERR TSSolve(solver, NULL);