Actual source code: ex79.c

  1: #include <petsc.h>

  3: static char help[] = "Solves a linear system with a block of right-hand sides, apply a preconditioner to the same block.\n\n";

  5: PetscErrorCode MatApply(PC pc, Mat X, Mat Y)
  6: {
  7:   PetscFunctionBeginUser;
  8:   PetscCall(MatCopy(X, Y, SAME_NONZERO_PATTERN));
  9:   PetscFunctionReturn(PETSC_SUCCESS);
 10: }

 12: int main(int argc, char **args)
 13: {
 14:   Mat       X, B; /* computed solutions and RHS */
 15:   Mat       A;    /* linear system matrix */
 16:   KSP       ksp;  /* linear solver context */
 17:   PC        pc;   /* preconditioner context */
 18:   PetscInt  m = 10;
 19:   PetscBool flg;
 20: #if defined(PETSC_USE_LOG)
 21:   PetscLogEvent event;
 22: #endif
 23:   PetscEventPerfInfo info;

 25:   PetscFunctionBeginUser;
 26:   PetscCall(PetscInitialize(&argc, &args, NULL, help));
 27:   PetscCall(PetscLogDefaultBegin());
 28:   PetscCall(PetscOptionsGetInt(NULL, NULL, "-m", &m, NULL));
 29:   PetscCall(MatCreateAIJ(PETSC_COMM_WORLD, m, m, PETSC_DECIDE, PETSC_DECIDE, m, NULL, m, NULL, &A));
 30:   PetscCall(MatCreateDense(PETSC_COMM_WORLD, m, PETSC_DECIDE, PETSC_DECIDE, m, NULL, &B));
 31:   PetscCall(MatCreateDense(PETSC_COMM_WORLD, m, PETSC_DECIDE, PETSC_DECIDE, m, NULL, &X));
 32:   PetscCall(MatSetRandom(A, NULL));
 33:   PetscCall(MatSetOption(A, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_FALSE));
 34:   PetscCall(MatShift(A, 10.0));
 35:   PetscCall(MatSetRandom(B, NULL));
 36:   PetscCall(MatSetFromOptions(A));
 37:   PetscCall(PetscObjectTypeCompareAny((PetscObject)A, &flg, MATSEQAIJCUSPARSE, MATMPIAIJCUSPARSE, ""));
 38:   if (flg) {
 39:     PetscCall(MatConvert(B, MATDENSECUDA, MAT_INPLACE_MATRIX, &B));
 40:     PetscCall(MatConvert(X, MATDENSECUDA, MAT_INPLACE_MATRIX, &X));
 41:   }
 42:   PetscCall(KSPCreate(PETSC_COMM_WORLD, &ksp));
 43:   PetscCall(KSPSetOperators(ksp, A, A));
 44:   PetscCall(KSPSetFromOptions(ksp));
 45:   PetscCall(KSPGetPC(ksp, &pc));
 46:   PetscCall(PCShellSetMatApply(pc, MatApply));
 47:   PetscCall(KSPMatSolve(ksp, B, X));
 48:   PetscCall(PCMatApply(pc, B, X));
 49:   PetscCall(MatDestroy(&X));
 50:   PetscCall(MatDestroy(&B));
 51:   PetscCall(MatDestroy(&A));
 52:   PetscCall(KSPDestroy(&ksp));
 53:   PetscCall(PetscLogEventRegister("PCApply", PC_CLASSID, &event));
 54:   PetscCall(PetscLogEventGetPerfInfo(PETSC_DETERMINE, event, &info));
 55:   PetscCheck(!PetscDefined(USE_LOG) || m <= 1 || !info.count, PetscObjectComm((PetscObject)ksp), PETSC_ERR_PLIB, "PCApply() called %d times", info.count);
 56:   PetscCall(PetscLogEventRegister("PCMatApply", PC_CLASSID, &event));
 57:   PetscCall(PetscLogEventGetPerfInfo(PETSC_DETERMINE, event, &info));
 58:   PetscCheck(!PetscDefined(USE_LOG) || m <= 1 || info.count, PetscObjectComm((PetscObject)ksp), PETSC_ERR_PLIB, "PCMatApply() never called");
 59:   PetscCall(PetscFinalize());
 60:   return 0;
 61: }

 63: /*TEST
 64:    # KSPHPDDM does either pseudo-blocking or "true" blocking, all tests should succeed with other -ksp_hpddm_type
 65:    testset:
 66:       nsize: 1
 67:       args: -pc_type {{bjacobi lu ilu mat cholesky icc none shell asm gasm}shared output}
 68:       test:
 69:          suffix: 1
 70:          output_file: output/ex77_preonly.out
 71:          args: -ksp_type preonly
 72:       test:
 73:          suffix: 1_hpddm
 74:          output_file: output/ex77_preonly.out
 75:          requires: hpddm
 76:          args: -ksp_type hpddm -ksp_hpddm_type preonly

 78:    testset:
 79:       nsize: 1
 80:       args: -pc_type ksp
 81:       test:
 82:          suffix: 2
 83:          output_file: output/ex77_preonly.out
 84:          args: -ksp_ksp_type preonly -ksp_type preonly
 85:       test:
 86:          suffix: 2_hpddm
 87:          output_file: output/ex77_preonly.out
 88:          requires: hpddm
 89:          args: -ksp_ksp_type hpddm -ksp_type hpddm -ksp_hpddm_type preonly -ksp_ksp_hpddm_type preonly

 91:    testset:
 92:       nsize: 1
 93:       requires: h2opus
 94:       args: -pc_type h2opus -pc_h2opus_init_mat_h2opus_leafsize 10
 95:       test:
 96:          suffix: 3
 97:          output_file: output/ex77_preonly.out
 98:          args: -ksp_type preonly
 99:       test:
100:          suffix: 3_hpddm
101:          output_file: output/ex77_preonly.out
102:          requires: hpddm
103:          args: -ksp_type hpddm -ksp_hpddm_type preonly

105:    testset:
106:       nsize: 1
107:       requires: spai
108:       args: -pc_type spai
109:       test:
110:          suffix: 4
111:          output_file: output/ex77_preonly.out
112:          args: -ksp_type preonly
113:       test:
114:          suffix: 4_hpddm
115:          output_file: output/ex77_preonly.out
116:          requires: hpddm
117:          args: -ksp_type hpddm -ksp_hpddm_type preonly
118:    # special code path in PCMatApply() for PCBJACOBI when a block is shared by multiple processes
119:    testset:
120:       nsize: 2
121:       args: -pc_type bjacobi -pc_bjacobi_blocks 1 -sub_pc_type none
122:       test:
123:          suffix: 5
124:          output_file: output/ex77_preonly.out
125:          args: -ksp_type preonly -sub_ksp_type preonly
126:       test:
127:          suffix: 5_hpddm
128:          output_file: output/ex77_preonly.out
129:          requires: hpddm
130:          args: -ksp_type hpddm -ksp_hpddm_type preonly -sub_ksp_type hpddm
131:    # special code path in PCMatApply() for PCGASM when a block is shared by multiple processes
132:    testset:
133:       nsize: 2
134:       args: -pc_type gasm -pc_gasm_total_subdomains 1 -sub_pc_type none
135:       test:
136:          suffix: 6
137:          output_file: output/ex77_preonly.out
138:          args: -ksp_type preonly -sub_ksp_type preonly
139:       test:
140:          suffix: 6_hpddm
141:          output_file: output/ex77_preonly.out
142:          requires: hpddm
143:          args: -ksp_type hpddm -ksp_hpddm_type preonly -sub_ksp_type hpddm

145:    testset:
146:       nsize: 1
147:       requires: suitesparse
148:       args: -pc_type qr
149:       test:
150:          suffix: 7
151:          output_file: output/ex77_preonly.out
152:          args: -ksp_type preonly
153:       test:
154:          suffix: 7_hpddm
155:          output_file: output/ex77_preonly.out
156:          requires: hpddm
157:          args: -ksp_type hpddm -ksp_hpddm_type preonly

159:    testset:
160:       nsize: 1
161:       requires: hpddm cuda
162:       args: -mat_type aijcusparse -ksp_type hpddm
163:       test:
164:          suffix: 8_hpddm
165:          output_file: output/ex77_preonly.out

167: TEST*/