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       A, X, B; /* computed solutions and RHS */
 15:   KSP       ksp;     /* linear solver context */
 16:   PC        pc;      /* preconditioner context */
 17:   PetscInt  m = 10;
 18:   PetscBool flg, transpose = PETSC_FALSE;
 19: #if defined(PETSC_USE_LOG)
 20:   PetscLogEvent event;
 21: #endif
 22:   PetscEventPerfInfo info;

 24:   PetscFunctionBeginUser;
 25:   PetscCall(PetscInitialize(&argc, &args, NULL, help));
 26:   PetscCall(PetscLogIsActive(&flg));
 27:   if (!flg) 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(MatSetRandom(A, NULL));
 31:   PetscCall(MatSetOption(A, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_FALSE));
 32:   PetscCall(PetscOptionsGetBool(NULL, NULL, "-transpose", &transpose, NULL));
 33:   if (transpose) {
 34:     PetscCall(MatTranspose(A, MAT_INITIAL_MATRIX, &B));
 35:     PetscCall(MatAXPY(A, 1.0, B, DIFFERENT_NONZERO_PATTERN));
 36:     PetscCall(MatDestroy(&B));
 37:   }
 38:   PetscCall(MatShift(A, 10.0));
 39: #if defined(HYPRE_USING_HIP)
 40:   PetscCall(MatConvert(A, MATAIJHIPSPARSE, MAT_INPLACE_MATRIX, &A));
 41: #elif defined(HYPRE_USING_CUDA)
 42:   PetscCall(MatConvert(A, MATAIJCUSPARSE, MAT_INPLACE_MATRIX, &A));
 43: #endif
 44:   PetscCall(MatCreateDense(PETSC_COMM_WORLD, m, PETSC_DECIDE, PETSC_DECIDE, m, NULL, &B));
 45:   PetscCall(MatCreateDense(PETSC_COMM_WORLD, m, PETSC_DECIDE, PETSC_DECIDE, m, NULL, &X));
 46:   PetscCall(MatSetRandom(B, NULL));
 47:   PetscCall(MatSetFromOptions(A));
 48:   PetscCall(PetscObjectTypeCompareAny((PetscObject)A, &flg, MATSEQAIJCUSPARSE, MATMPIAIJCUSPARSE, ""));
 49:   if (flg) {
 50:     PetscCall(MatConvert(B, MATDENSECUDA, MAT_INPLACE_MATRIX, &B));
 51:     PetscCall(MatConvert(X, MATDENSECUDA, MAT_INPLACE_MATRIX, &X));
 52:   }
 53:   PetscCall(KSPCreate(PETSC_COMM_WORLD, &ksp));
 54:   PetscCall(KSPSetOperators(ksp, A, A));
 55:   PetscCall(KSPSetFromOptions(ksp));
 56:   PetscCall(KSPGetPC(ksp, &pc));
 57:   PetscCall(PCShellSetMatApply(pc, MatApply));
 58:   PetscCall(KSPMatSolve(ksp, B, X));
 59:   PetscCall(PCMatApply(pc, B, X));
 60:   if (transpose) {
 61:     PetscCall(KSPMatSolveTranspose(ksp, B, X));
 62:     PetscCall(PCMatApply(pc, B, X));
 63:     PetscCall(KSPMatSolve(ksp, B, X));
 64:   }
 65:   PetscCall(MatDestroy(&X));
 66:   PetscCall(MatDestroy(&B));
 67:   PetscCall(MatDestroy(&A));
 68:   PetscCall(KSPDestroy(&ksp));
 69:   PetscCall(PetscLogEventRegister("PCApply", PC_CLASSID, &event));
 70:   PetscCall(PetscLogEventGetPerfInfo(PETSC_DETERMINE, event, &info));
 71:   PetscCheck(!PetscDefined(USE_LOG) || m <= 1 || !info.count, PetscObjectComm((PetscObject)ksp), PETSC_ERR_PLIB, "PCApply() called %d times", info.count);
 72:   PetscCall(PetscLogEventRegister("PCMatApply", PC_CLASSID, &event));
 73:   PetscCall(PetscLogEventGetPerfInfo(PETSC_DETERMINE, event, &info));
 74:   PetscCheck(!PetscDefined(USE_LOG) || m <= 1 || info.count, PetscObjectComm((PetscObject)ksp), PETSC_ERR_PLIB, "PCMatApply() never called");
 75:   PetscCall(PetscFinalize());
 76:   return 0;
 77: }

 79: /*TEST
 80:    # KSPHPDDM does either pseudo-blocking or "true" blocking, all tests should succeed with other -ksp_hpddm_type
 81:    testset:
 82:       nsize: 1
 83:       args: -pc_type {{bjacobi lu ilu mat cholesky icc none shell asm gasm}shared output}
 84:       test:
 85:          suffix: 1
 86:          output_file: output/ex77_preonly.out
 87:          args: -ksp_type preonly
 88:       test:
 89:          suffix: 1_hpddm
 90:          output_file: output/ex77_preonly.out
 91:          requires: hpddm
 92:          args: -ksp_type hpddm -ksp_hpddm_type preonly

 94:    testset:
 95:       nsize: 1
 96:       args: -pc_type ksp
 97:       test:
 98:          suffix: 2
 99:          output_file: output/ex77_preonly.out
100:          args: -ksp_ksp_type preonly -ksp_type preonly
101:       test:
102:          suffix: 2_hpddm
103:          output_file: output/ex77_preonly.out
104:          requires: hpddm
105:          args: -ksp_ksp_type hpddm -ksp_type hpddm -ksp_hpddm_type preonly -ksp_ksp_hpddm_type preonly

107:    testset:
108:       nsize: 1
109:       requires: h2opus
110:       args: -pc_type h2opus -pc_h2opus_init_mat_h2opus_leafsize 10
111:       test:
112:          suffix: 3
113:          output_file: output/ex77_preonly.out
114:          args: -ksp_type preonly
115:       test:
116:          suffix: 3_hpddm
117:          output_file: output/ex77_preonly.out
118:          requires: hpddm
119:          args: -ksp_type hpddm -ksp_hpddm_type preonly

121:    testset:
122:       nsize: 1
123:       requires: spai
124:       args: -pc_type spai
125:       test:
126:          suffix: 4
127:          output_file: output/ex77_preonly.out
128:          args: -ksp_type preonly
129:       test:
130:          suffix: 4_hpddm
131:          output_file: output/ex77_preonly.out
132:          requires: hpddm
133:          args: -ksp_type hpddm -ksp_hpddm_type preonly
134:    # special code path in PCMatApply() for PCBJACOBI when a block is shared by multiple processes
135:    testset:
136:       nsize: 2
137:       args: -pc_type bjacobi -pc_bjacobi_blocks 1 -sub_pc_type none
138:       test:
139:          suffix: 5
140:          output_file: output/ex77_preonly.out
141:          args: -ksp_type preonly -sub_ksp_type preonly
142:       test:
143:          suffix: 5_hpddm
144:          output_file: output/ex77_preonly.out
145:          requires: hpddm
146:          args: -ksp_type hpddm -ksp_hpddm_type preonly -sub_ksp_type hpddm
147:    # special code path in PCMatApply() for PCGASM when a block is shared by multiple processes
148:    testset:
149:       nsize: 2
150:       args: -pc_type gasm -pc_gasm_total_subdomains 1 -sub_pc_type none
151:       test:
152:          suffix: 6
153:          output_file: output/ex77_preonly.out
154:          args: -ksp_type preonly -sub_ksp_type preonly
155:       test:
156:          suffix: 6_hpddm
157:          output_file: output/ex77_preonly.out
158:          requires: hpddm
159:          args: -ksp_type hpddm -ksp_hpddm_type preonly -sub_ksp_type hpddm

161:    testset:
162:       nsize: 1
163:       requires: suitesparse
164:       args: -pc_type qr
165:       test:
166:          suffix: 7
167:          output_file: output/ex77_preonly.out
168:          args: -ksp_type preonly
169:       test:
170:          suffix: 7_hpddm
171:          output_file: output/ex77_preonly.out
172:          requires: hpddm
173:          args: -ksp_type hpddm -ksp_hpddm_type preonly

175:    testset:
176:       nsize: 1
177:       requires: hpddm cuda
178:       args: -mat_type aijcusparse -ksp_type hpddm
179:       test:
180:          suffix: 8_hpddm
181:          output_file: output/ex77_preonly.out
182:       test:
183:          suffix: 8_hpddm_transpose
184:          output_file: output/ex77_preonly.out
185:          args: -pc_type icc -transpose

187:    testset:
188:       nsize: 1
189:       args: -pc_type {{cholesky icc none}shared output} -transpose
190:       test:
191:          suffix: 1_transpose
192:          output_file: output/ex77_preonly.out
193:          args: -ksp_type preonly
194:       test:
195:          suffix: 1_hpddm_transpose
196:          output_file: output/ex77_preonly.out
197:          requires: hpddm
198:          args: -ksp_type hpddm -ksp_hpddm_type preonly

200: TEST*/