532 {
536
537 CHKERR KSPSetFromOptions(ksp);
538 PC pc;
539 CHKERR KSPGetPC(ksp, &pc);
540
541 PetscBool is_pcfs = PETSC_FALSE;
542 PetscObjectTypeCompare((PetscObject)pc, PCFIELDSPLIT, &is_pcfs);
543 if (is_pcfs) {
544
545 MOFEM_LOG(
"AT", Sev::inform) <<
"Setup Schur pc";
546
547 auto create_sub_dm = [&]() {
549
550 auto create_dm = [&](
551
552 std::string problem_name,
553 std::vector<std::string> fe_names,
554 std::vector<std::string> fields,
555
556 auto dm_type
557
558 ) {
560 auto create_dm_imp = [&]() {
564 for (auto fe : fe_names) {
566 }
568 for (auto field : fields) {
571 }
574 };
576 create_dm_imp(),
577 "Error in creating schurDM. It is possible that schurDM is "
578 "already created");
579 return dm;
580 };
581
582 auto schur_dm = create_dm(
583
584 "SCHUR",
585
586 {
simple->getDomainFEName(),
simple->getSkeletonFEName()},
587
588 {"HYBRID"},
589
590 "DMMOFEM_MG");
591
592 auto block_dm = create_dm(
593
594 "BLOCK",
595
596 {
simple->getDomainFEName(),
simple->getSkeletonFEName()},
597
598 {"BROKEN", "U"},
599
600 "DMMOFEM");
601
602 return std::make_tuple(schur_dm, block_dm);
603 };
604
605 auto get_nested_mat_data = [&](auto schur_dm, auto block_dm) {
608
609 {
610
611 {
612
613 simple->getDomainFEName(),
614
615 {
616
617 {"BROKEN", "BROKEN"},
618 {"U", "U"},
619 {"BROKEN", "U"},
620 {"U", "BROKEN"}
621
622 }
623
624 },
625
626 {
627
628 simple->getSkeletonFEName(),
629
630 {
631
632 {"BROKEN", "HYBRID"}, {"HYBRID", "BROKEN"}
633
634 }
635
636 }
637
638 }
639
640 );
641
643
644 {schur_dm, block_dm}, block_mat_data,
645
646 {"BROKEN", "U"}, {nullptr, nullptr}, true
647
648 );
649 };
650
651 auto set_ops = [&](auto schur_dm) {
655
656 boost::shared_ptr<BlockStructure> block_data;
658
660 pip_mng->getOpDomainLhsPipeline().push_front(
662 pip_mng->getOpDomainLhsPipeline().push_back(
663
666
667 );
668 }
669
670 auto pre_proc_schur_lhs_ptr = boost::make_shared<FEMethod>();
671 auto post_proc_schur_lhs_ptr = boost::make_shared<FEMethod>();
672
673 pre_proc_schur_lhs_ptr->preProcessHook = [this]() {
676 MOFEM_LOG(
"AT", Sev::verbose) <<
"Lhs Assemble Begin";
678 };
679
680 post_proc_schur_lhs_ptr->postProcessHook = [this, ao_up,
681 post_proc_schur_lhs_ptr]() {
683 MOFEM_LOG(
"AT", Sev::verbose) <<
"Lhs Assemble End";
684
686 CHKERR MatAssemblyBegin(
S, MAT_FINAL_ASSEMBLY);
687 CHKERR MatAssemblyEnd(
S, MAT_FINAL_ASSEMBLY);
688 if (1) {
689 auto S_from_block =
matDuplicate(
S, MAT_SHARE_NONZERO_PATTERN);
690
692 S_from_block, {"BROKEN", "U"},
693 {nullptr, nullptr}, ao_up);
694 CHKERR MatAssemblyBegin(S_from_block, MAT_FINAL_ASSEMBLY);
695 CHKERR MatAssemblyEnd(S_from_block, MAT_FINAL_ASSEMBLY);
696 CHKERR MatAYPX(S_from_block, -1,
S, DIFFERENT_NONZERO_PATTERN);
697 double norm;
698 CHKERR MatNorm(S_from_block, NORM_FROBENIUS, &norm);
699 MOFEM_LOG(
"AT", Sev::inform) <<
"Norm of difference: " << norm;
700 if (norm > 1e-6)
701 SETERRQ(
703 "Norm of difference between Schur and block matrix is larger "
704 "than accepted");
705 }
706 } else {
708 {"BROKEN", "U"}, {nullptr, nullptr}, ao_up);
709 CHKERR MatAssemblyBegin(
S, MAT_FINAL_ASSEMBLY);
710 CHKERR MatAssemblyEnd(
S, MAT_FINAL_ASSEMBLY);
711 }
712
713 MOFEM_LOG(
"AT", Sev::verbose) <<
"Lhs Assemble Finish";
715 };
716
718 ksp_ctx_ptr->getPreProcSetOperators().push_front(pre_proc_schur_lhs_ptr);
719 ksp_ctx_ptr->getPostProcSetOperators().push_back(post_proc_schur_lhs_ptr);
720
722 };
723
724 auto set_pc = [&](auto pc, auto block_dm) {
726 auto block_is =
getDMSubData(block_dm)->getSmartRowIs();
727 CHKERR PCFieldSplitSetIS(pc, NULL, block_is);
728 CHKERR PCFieldSplitSetSchurPre(pc, PC_FIELDSPLIT_SCHUR_PRE_USER,
S);
730 };
731
732
733
734 auto [schur_dm, block_dm] = create_sub_dm();
736 auto nested_mat_data = get_nested_mat_data(schur_dm, block_dm);
738 }
740 CHKERR MatSetDM(
S, PETSC_NULLPTR);
741
745
747 CHKERR set_pc(pc, block_dm);
748 DM solver_dm;
749 CHKERR KSPGetDM(ksp, &solver_dm);
751 CHKERR DMSetMatType(solver_dm, MATSHELL);
752
753 auto get_pc = [](auto ksp) {
754 PC pc_raw;
755 CHKERR KSPGetPC(ksp, &pc_raw);
756 return pc_raw;
757 };
758
759 auto set_diagonal_pc = [&](auto pc) {
761
765 CHKERR PCSetOperators(pc, A, P);
766 }
767
768 KSP *subksp;
769 CHKERR PCFieldSplitSchurGetSubKSP(pc, PETSC_NULLPTR, &subksp);
771
774 };
775
776 auto set_mg_for_schur_complement = [&](
auto pc,
auto schur_dm,
auto S) {
778
779 KSP *subksp;
780 CHKERR PCFieldSplitSchurGetSubKSP(pc, PETSC_NULLPTR, &subksp);
781 auto subpc = get_pc(subksp[1]);
782
783 auto set_pc_p_mg = [](
auto dm,
auto pc,
auto S) {
786 PetscBool same = PETSC_FALSE;
787 PetscObjectTypeCompare((PetscObject)pc, PCMG, &same);
788 if (same) {
789 MOFEM_LOG(
"TIMER", Sev::inform) <<
"Set up MG";
792 CHKERR PCSetFromOptions(pc);
793 }
795 };
796
797 PetscBool same = PETSC_FALSE;
798 PetscObjectTypeCompare((PetscObject)subpc, PCKSP, &same);
799 if (same) {
800 CHKERR PCSetFromOptions(subpc);
801 KSP inner_ksp;
802 CHKERR PCKSPGetKSP(subpc, &inner_ksp);
803 CHKERR KSPSetFromOptions(inner_ksp);
804 PC inner_pc;
805 CHKERR KSPGetPC(inner_ksp, &inner_pc);
806 CHKERR PCSetFromOptions(inner_pc);
807 CHKERR set_pc_p_mg(schur_dm, inner_pc,
S);
808 }
809
812 };
813
816 CHKERR set_diagonal_pc(pc);
817 CHKERR set_mg_for_schur_complement(pc, schur_dm,
S);
818 }
819
820 } else {
822 "PC is not set to PCFIELDSPLIT");
823 }
825}
@ MOFEM_ATOM_TEST_INVALID
#define NBFACETRI_L2(P)
Number of base functions on triangle for L2 space.
#define NBEDGE_L2(P)
Number of base functions on edge from L2 space.
MoFEMErrorCode DMMoFEMGetBlocMatData(DM dm, boost::shared_ptr< BlockStructure > &)
Get data for block mat.
auto getDMKspCtx(DM dm)
Get KSP context data structure used by DM.
auto createDMHybridisedL2Matrix(DM dm)
Get smart hybridised L2 matrix from DM.
SmartPetscObj< Mat > matDuplicate(Mat mat, MatDuplicateOption op)
auto createDMNestSchurMat(DM dm)
MoFEMErrorCode assembleBlockMatSchur(MoFEM::Interface &m_field, Mat B, Mat S, std::vector< std::string > fields_name, std::vector< boost::shared_ptr< Range > > field_ents, SmartPetscObj< AO > ao)
Assemble Schur matrix.
auto createDMBlockMat(DM dm)
static constexpr int approx_order