565 {
569
570 CHKERR KSPSetFromOptions(ksp);
571 PC pc;
572 CHKERR KSPGetPC(ksp, &pc);
573
574 PetscBool is_pcfs = PETSC_FALSE;
575 PetscObjectTypeCompare((PetscObject)pc, PCFIELDSPLIT, &is_pcfs);
576 if (is_pcfs) {
577
578 MOFEM_LOG(
"AT", Sev::inform) <<
"Setup Schur pc";
579
580 auto create_sub_dm = [&]() {
582
583 auto create_dm = [&](
584
585 std::string problem_name,
586 std::vector<std::string> fe_names,
587 std::vector<std::string> fields,
588
589 auto dm_type
590
591 ) {
593 auto create_dm_imp = [&]() {
597 for (auto fe : fe_names) {
599 }
601 for (auto field : fields) {
604 }
607 };
609 create_dm_imp(),
610 "Error in creating schurDM. It is possible that schurDM is "
611 "already created");
612 return dm;
613 };
614
615 auto schur_dm = create_dm(
616
617 "SCHUR",
618
619 {
simple->getDomainFEName(),
simple->getSkeletonFEName()},
620
621 {"HYBRID"},
622
623 "DMMOFEM_MG");
624
625 auto block_dm = create_dm(
626
627 "BLOCK",
628
629 {
simple->getDomainFEName(),
simple->getSkeletonFEName()},
630
631 {"BROKEN", "U"},
632
633 "DMMOFEM");
634
635 return std::make_tuple(schur_dm, block_dm);
636 };
637
638 auto get_nested_mat_data = [&](auto schur_dm, auto block_dm) {
641
642 {
643
644 {
645
646 simple->getDomainFEName(),
647
648 {
649
650 {"BROKEN", "BROKEN"},
651 {"U", "U"},
652 {"BROKEN", "U"},
653 {"U", "BROKEN"}
654
655 }
656
657 },
658
659 {
660
661 simple->getSkeletonFEName(),
662
663 {
664
665 {"BROKEN", "HYBRID"}, {"HYBRID", "BROKEN"}
666
667 }
668
669 }
670
671 }
672
673 );
674
676
677 {schur_dm, block_dm}, block_mat_data,
678
679 {"BROKEN", "U"}, {nullptr, nullptr}, true
680
681 );
682 };
683
684 auto set_ops = [&](auto schur_dm) {
688
689 boost::shared_ptr<BlockStructure> block_data;
691
693 pip_mng->getOpDomainLhsPipeline().push_front(
695 pip_mng->getOpDomainLhsPipeline().push_back(
696
699
700 );
701 }
702
703 auto pre_proc_schur_lhs_ptr = boost::make_shared<FEMethod>();
704 auto post_proc_schur_lhs_ptr = boost::make_shared<FEMethod>();
705
706 pre_proc_schur_lhs_ptr->preProcessHook = [this]() {
709 MOFEM_LOG(
"AT", Sev::verbose) <<
"Lhs Assemble Begin";
711 };
712
713 post_proc_schur_lhs_ptr->postProcessHook = [this, ao_up,
714 post_proc_schur_lhs_ptr]() {
716 MOFEM_LOG(
"AT", Sev::verbose) <<
"Lhs Assemble End";
717
719 CHKERR MatAssemblyBegin(
S, MAT_FINAL_ASSEMBLY);
720 CHKERR MatAssemblyEnd(
S, MAT_FINAL_ASSEMBLY);
721 if (1) {
722 auto S_from_block =
matDuplicate(
S, MAT_SHARE_NONZERO_PATTERN);
723
725 S_from_block, {"BROKEN", "U"},
726 {nullptr, nullptr}, ao_up);
727 CHKERR MatAssemblyBegin(S_from_block, MAT_FINAL_ASSEMBLY);
728 CHKERR MatAssemblyEnd(S_from_block, MAT_FINAL_ASSEMBLY);
729 CHKERR MatAYPX(S_from_block, -1,
S, DIFFERENT_NONZERO_PATTERN);
730 double norm;
731 CHKERR MatNorm(S_from_block, NORM_FROBENIUS, &norm);
732 MOFEM_LOG(
"AT", Sev::inform) <<
"Norm of difference: " << norm;
733 if (norm > 1e-6)
734 SETERRQ(
736 "Norm of difference between Schur and block matrix is larger "
737 "than accepted");
738 }
739 } else {
741 {"BROKEN", "U"}, {nullptr, nullptr}, ao_up);
742 CHKERR MatAssemblyBegin(
S, MAT_FINAL_ASSEMBLY);
743 CHKERR MatAssemblyEnd(
S, MAT_FINAL_ASSEMBLY);
744 }
745
746 MOFEM_LOG(
"AT", Sev::verbose) <<
"Lhs Assemble Finish";
748 };
749
751 ksp_ctx_ptr->getPreProcSetOperators().push_front(pre_proc_schur_lhs_ptr);
752 ksp_ctx_ptr->getPostProcSetOperators().push_back(post_proc_schur_lhs_ptr);
753
755 };
756
757 auto set_pc = [&](auto pc, auto block_dm) {
759 auto block_is =
getDMSubData(block_dm)->getSmartRowIs();
760 CHKERR PCFieldSplitSetIS(pc, NULL, block_is);
761 CHKERR PCFieldSplitSetSchurPre(pc, PC_FIELDSPLIT_SCHUR_PRE_USER,
S);
763 };
764
765
766
767 auto [schur_dm, block_dm] = create_sub_dm();
769 auto nested_mat_data = get_nested_mat_data(schur_dm, block_dm);
771 }
773 CHKERR MatSetDM(
S, PETSC_NULLPTR);
774
778
780 CHKERR set_pc(pc, block_dm);
781 DM solver_dm;
782 CHKERR KSPGetDM(ksp, &solver_dm);
784 CHKERR DMSetMatType(solver_dm, MATSHELL);
785
786 auto get_pc = [](auto ksp) {
787 PC pc_raw;
788 CHKERR KSPGetPC(ksp, &pc_raw);
789 return pc_raw;
790 };
791
792 auto set_diagonal_pc = [&](auto pc) {
794
798 CHKERR PCSetOperators(pc,
A, P);
799 }
800
801 KSP *subksp;
802 CHKERR PCFieldSplitSchurGetSubKSP(pc, PETSC_NULLPTR, &subksp);
804
807 };
808
809 auto set_mg_for_schur_complement = [&](
auto pc,
auto schur_dm,
auto S) {
811
812 KSP *subksp;
813 CHKERR PCFieldSplitSchurGetSubKSP(pc, PETSC_NULLPTR, &subksp);
814 auto subpc = get_pc(subksp[1]);
815
816 auto set_pc_p_mg = [](
auto dm,
auto pc,
auto S) {
819 PetscBool same = PETSC_FALSE;
820 PetscObjectTypeCompare((PetscObject)pc, PCMG, &same);
821 if (same) {
822 MOFEM_LOG(
"TIMER", Sev::inform) <<
"Set up MG";
825 CHKERR PCSetFromOptions(pc);
826 }
828 };
829
830 PetscBool same = PETSC_FALSE;
831 PetscObjectTypeCompare((PetscObject)subpc, PCKSP, &same);
832 if (same) {
833 CHKERR PCSetFromOptions(subpc);
834 KSP inner_ksp;
835 CHKERR PCKSPGetKSP(subpc, &inner_ksp);
836 CHKERR KSPSetFromOptions(inner_ksp);
837 PC inner_pc;
838 CHKERR KSPGetPC(inner_ksp, &inner_pc);
839 CHKERR PCSetFromOptions(inner_pc);
840 CHKERR set_pc_p_mg(schur_dm, inner_pc,
S);
841 }
842
845 };
846
849 CHKERR set_diagonal_pc(pc);
850 CHKERR set_mg_for_schur_complement(pc, schur_dm,
S);
851 }
852
853 } else {
855 "PC is not set to PCFIELDSPLIT");
856 }
858}
@ 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)
constexpr AssemblyType A
[Define dimension]