Actual source code: ex79.c

  1: #include <petsc.h>

  3: #if PetscDefined(HAVE_HYPRE_DEVICE)
  4: #include <petsc/private/petschypre.h>
  5: #endif

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

  9: PetscErrorCode MatApply(PC pc, Mat X, Mat Y)
 10: {
 11:   PetscFunctionBeginUser;
 12:   PetscCall(MatCopy(X, Y, SAME_NONZERO_PATTERN));
 13:   PetscFunctionReturn(PETSC_SUCCESS);
 14: }

 16: int main(int argc, char **args)
 17: {
 18:   Mat                A, X, B; /* computed solutions and RHS */
 19:   KSP                ksp;     /* linear solver context */
 20:   PC                 pc;      /* preconditioner context */
 21:   PetscInt           m = 10;
 22:   PetscBool          flg, transpose = PETSC_FALSE;
 23:   PetscLogEvent      event;
 24:   PetscEventPerfInfo info;

 26:   PetscFunctionBeginUser;
 27:   PetscCall(PetscInitialize(&argc, &args, NULL, help));
 28:   PetscCall(PetscLogDefaultBegin());
 29:   PetscCall(PetscOptionsGetInt(NULL, NULL, "-m", &m, NULL));
 30:   PetscCall(MatCreateAIJ(PETSC_COMM_WORLD, m, m, PETSC_DECIDE, PETSC_DECIDE, m, NULL, m, NULL, &A));
 31:   PetscCall(MatSetRandom(A, NULL));
 32:   PetscCall(MatSetOption(A, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_FALSE));
 33:   PetscCall(PetscOptionsGetBool(NULL, NULL, "-transpose", &transpose, NULL));
 34:   PetscCall(MatShift(A, 10.0));
 35:   PetscCall(MatCreateDense(PETSC_COMM_WORLD, m, PETSC_DECIDE, PETSC_DECIDE, m, NULL, &B));
 36:   PetscCall(MatCreateDense(PETSC_COMM_WORLD, m, PETSC_DECIDE, PETSC_DECIDE, m, NULL, &X));
 37:   PetscCall(MatSetRandom(B, NULL));
 38:   PetscCall(MatSetFromOptions(A));
 39:   PetscCall(PetscObjectTypeCompareAny((PetscObject)A, &flg, MATSEQAIJCUSPARSE, MATMPIAIJCUSPARSE, ""));
 40:   if (flg) {
 41:     PetscCall(MatConvert(B, MATDENSECUDA, MAT_INPLACE_MATRIX, &B));
 42:     PetscCall(MatConvert(X, MATDENSECUDA, MAT_INPLACE_MATRIX, &X));
 43:   } else {
 44:     PetscCall(PetscObjectTypeCompareAny((PetscObject)A, &flg, MATSEQAIJHIPSPARSE, MATMPIAIJHIPSPARSE, ""));
 45:     if (flg) {
 46:       PetscCall(MatConvert(B, MATDENSEHIP, MAT_INPLACE_MATRIX, &B));
 47:       PetscCall(MatConvert(X, MATDENSEHIP, MAT_INPLACE_MATRIX, &X));
 48:     }
 49:   }
 50:   PetscCall(KSPCreate(PETSC_COMM_WORLD, &ksp));
 51:   PetscCall(KSPSetOperators(ksp, A, A));
 52:   PetscCall(KSPSetFromOptions(ksp));
 53:   PetscCall(KSPGetPC(ksp, &pc));
 54:   PetscCall(PetscObjectTypeCompare((PetscObject)pc, PCHYPRE, &flg));
 55:   if (flg && PetscDefined(HAVE_HYPRE_DEVICE)) {
 56: #if defined(HYPRE_USING_HIP)
 57:     PetscCall(MatConvert(A, MATAIJHIPSPARSE, MAT_INPLACE_MATRIX, &A));
 58:     PetscCall(MatConvert(B, MATDENSEHIP, MAT_INPLACE_MATRIX, &B));
 59:     PetscCall(MatConvert(X, MATDENSEHIP, MAT_INPLACE_MATRIX, &X));
 60: #elif defined(HYPRE_USING_CUDA)
 61:     PetscCall(MatConvert(A, MATAIJCUSPARSE, MAT_INPLACE_MATRIX, &A));
 62:     PetscCall(MatConvert(B, MATDENSECUDA, MAT_INPLACE_MATRIX, &B));
 63:     PetscCall(MatConvert(X, MATDENSECUDA, MAT_INPLACE_MATRIX, &X));
 64: #endif
 65:   } else PetscCall(PCShellSetMatApply(pc, MatApply));
 66:   PetscCall(KSPMatSolve(ksp, B, X));
 67:   PetscCall(PCMatApply(pc, B, X));
 68:   if (transpose) {
 69:     PetscCall(KSPMatSolveTranspose(ksp, B, X));
 70:     PetscCall(PCMatApply(pc, B, X));
 71:     PetscCall(KSPMatSolve(ksp, B, X));
 72:   }
 73:   PetscCall(MatDestroy(&X));
 74:   PetscCall(MatDestroy(&B));
 75:   PetscCall(MatDestroy(&A));
 76:   PetscCall(KSPDestroy(&ksp));
 77:   PetscCall(PetscLogEventRegister("PCApply", PC_CLASSID, &event)); // there is no PCApplyTranspose event
 78:   PetscCall(PetscLogEventGetPerfInfo(PETSC_DETERMINE, event, &info));
 79:   PetscCheck(!PetscDefined(USE_LOG) || m <= 1 || !info.count, PetscObjectComm((PetscObject)ksp), PETSC_ERR_PLIB, "PCApply() called %d times", info.count);
 80:   PetscCall(PetscLogEventRegister("PCMatApply", PC_CLASSID, &event)); // there is no PCMatApplyTranspose event
 81:   PetscCall(PetscLogEventGetPerfInfo(PETSC_DETERMINE, event, &info));
 82:   PetscCheck(!PetscDefined(USE_LOG) || m <= 1 || info.count, PetscObjectComm((PetscObject)ksp), PETSC_ERR_PLIB, "PCMatApply() never called");
 83:   PetscCall(PetscFinalize());
 84:   return 0;
 85: }

 87: /*TEST
 88:    # KSPHPDDM does either pseudo-blocking or "true" blocking, all tests should succeed with other -ksp_hpddm_type
 89:    testset:
 90:       nsize: 1
 91:       args: -pc_type {{mat none shell gasm}shared output}
 92:       test:
 93:          suffix: 1
 94:          output_file: output/ex77_preonly.out
 95:          args: -ksp_type preonly
 96:       test:
 97:          suffix: 1_hpddm
 98:          output_file: output/ex77_preonly.out
 99:          requires: hpddm
100:          args: -ksp_type hpddm -ksp_hpddm_type preonly

102:    testset:
103:       nsize: 1
104:       args: -pc_type ksp -transpose
105:       test:
106:          suffix: 2
107:          output_file: output/ex77_preonly.out
108:          args: -ksp_ksp_type preonly -ksp_type preonly
109:       test:
110:          suffix: 2_hpddm
111:          output_file: output/ex77_preonly.out
112:          requires: hpddm
113:          args: -ksp_ksp_type hpddm -ksp_type hpddm -ksp_hpddm_type preonly -ksp_ksp_hpddm_type preonly

115:    testset:
116:       nsize: 1
117:       requires: h2opus
118:       args: -pc_type h2opus -pc_h2opus_init_mat_h2opus_leafsize 10
119:       test:
120:          suffix: 3
121:          output_file: output/ex77_preonly.out
122:          args: -ksp_type preonly
123:       test:
124:          suffix: 3_hpddm
125:          output_file: output/ex77_preonly.out
126:          requires: hpddm
127:          args: -ksp_type hpddm -ksp_hpddm_type preonly

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

169:    testset:
170:       nsize: 1
171:       requires: suitesparse
172:       args: -pc_type qr
173:       test:
174:          suffix: 7
175:          output_file: output/ex77_preonly.out
176:          args: -ksp_type preonly
177:       test:
178:          suffix: 7_hpddm
179:          output_file: output/ex77_preonly.out
180:          requires: hpddm
181:          args: -ksp_type hpddm -ksp_hpddm_type preonly

183:    testset:
184:       nsize: 1
185:       requires: hpddm cuda
186:       args: -mat_type aijcusparse -ksp_type hpddm
187:       test:
188:          suffix: 8_hpddm
189:          output_file: output/ex77_preonly.out
190:       test:
191:          suffix: 8_hpddm_transpose
192:          output_file: output/ex77_preonly.out
193:          args: -pc_type icc -transpose

195:    testset:
196:       nsize: 1
197:       args: -pc_type {{bjacobi lu ilu cholesky icc none}shared output} -transpose
198:       test:
199:          suffix: 1_transpose
200:          output_file: output/ex77_preonly.out
201:          args: -ksp_type preonly
202:       test:
203:          suffix: 1_hpddm_transpose
204:          output_file: output/ex77_preonly.out
205:          requires: hpddm
206:          args: -ksp_type hpddm -ksp_hpddm_type preonly

208:    testset:
209:       nsize: 2
210:       args: -pc_type asm -transpose
211:       test:
212:          suffix: 2_transpose
213:          output_file: output/ex77_preonly.out
214:          args: -ksp_type preonly
215:       test:
216:          suffix: 2_hpddm_transpose
217:          output_file: output/ex77_preonly.out
218:          requires: hpddm
219:          args: -ksp_type hpddm -ksp_hpddm_type preonly

221:    testset:
222:       requires: hypre !complex
223:       args: -pc_type hypre -pc_hypre_boomeramg_relax_type_all l1scaled-Jacobi -pc_hypre_boomeramg_no_CF
224:       test:
225:          suffix: 9
226:          output_file: output/ex77_preonly.out
227:          args: -ksp_type preonly
228:       test:
229:          suffix: 9_hpddm
230:          output_file: output/ex77_preonly.out
231:          requires: hpddm !hip
232:          args: -ksp_type hpddm -ksp_max_it 15 -ksp_error_if_not_converged

234: TEST*/