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:   if (transpose) {
 35:     PetscCall(MatTranspose(A, MAT_INITIAL_MATRIX, &B));
 36:     PetscCall(MatAXPY(A, 1.0, B, DIFFERENT_NONZERO_PATTERN));
 37:     PetscCall(MatDestroy(&B));
 38:   }
 39:   PetscCall(MatShift(A, 10.0));
 40:   PetscCall(MatCreateDense(PETSC_COMM_WORLD, m, PETSC_DECIDE, PETSC_DECIDE, m, NULL, &B));
 41:   PetscCall(MatCreateDense(PETSC_COMM_WORLD, m, PETSC_DECIDE, PETSC_DECIDE, m, NULL, &X));
 42:   PetscCall(MatSetRandom(B, NULL));
 43:   PetscCall(MatSetFromOptions(A));
 44:   PetscCall(PetscObjectTypeCompareAny((PetscObject)A, &flg, MATSEQAIJCUSPARSE, MATMPIAIJCUSPARSE, ""));
 45:   if (flg) {
 46:     PetscCall(MatConvert(B, MATDENSECUDA, MAT_INPLACE_MATRIX, &B));
 47:     PetscCall(MatConvert(X, MATDENSECUDA, MAT_INPLACE_MATRIX, &X));
 48:   } else {
 49:     PetscCall(PetscObjectTypeCompareAny((PetscObject)A, &flg, MATSEQAIJHIPSPARSE, MATMPIAIJHIPSPARSE, ""));
 50:     if (flg) {
 51:       PetscCall(MatConvert(B, MATDENSEHIP, MAT_INPLACE_MATRIX, &B));
 52:       PetscCall(MatConvert(X, MATDENSEHIP, MAT_INPLACE_MATRIX, &X));
 53:     }
 54:   }
 55:   PetscCall(KSPCreate(PETSC_COMM_WORLD, &ksp));
 56:   PetscCall(KSPSetOperators(ksp, A, A));
 57:   PetscCall(KSPSetFromOptions(ksp));
 58:   PetscCall(KSPGetPC(ksp, &pc));
 59:   PetscCall(PetscObjectTypeCompare((PetscObject)pc, PCHYPRE, &flg));
 60:   if (flg && PetscDefined(HAVE_HYPRE_DEVICE)) {
 61: #if defined(HYPRE_USING_HIP)
 62:     PetscCall(MatConvert(A, MATAIJHIPSPARSE, MAT_INPLACE_MATRIX, &A));
 63:     PetscCall(MatConvert(B, MATDENSEHIP, MAT_INPLACE_MATRIX, &B));
 64:     PetscCall(MatConvert(X, MATDENSEHIP, MAT_INPLACE_MATRIX, &X));
 65: #elif defined(HYPRE_USING_CUDA)
 66:     PetscCall(MatConvert(A, MATAIJCUSPARSE, MAT_INPLACE_MATRIX, &A));
 67:     PetscCall(MatConvert(B, MATDENSECUDA, MAT_INPLACE_MATRIX, &B));
 68:     PetscCall(MatConvert(X, MATDENSECUDA, MAT_INPLACE_MATRIX, &X));
 69: #endif
 70:   } else PetscCall(PCShellSetMatApply(pc, MatApply));
 71:   PetscCall(KSPMatSolve(ksp, B, X));
 72:   PetscCall(PCMatApply(pc, B, X));
 73:   if (transpose) {
 74:     PetscCall(KSPMatSolveTranspose(ksp, B, X));
 75:     PetscCall(PCMatApply(pc, B, X));
 76:     PetscCall(KSPMatSolve(ksp, B, X));
 77:   }
 78:   PetscCall(MatDestroy(&X));
 79:   PetscCall(MatDestroy(&B));
 80:   PetscCall(MatDestroy(&A));
 81:   PetscCall(KSPDestroy(&ksp));
 82:   PetscCall(PetscLogEventRegister("PCApply", PC_CLASSID, &event));
 83:   PetscCall(PetscLogEventGetPerfInfo(PETSC_DETERMINE, event, &info));
 84:   PetscCheck(!PetscDefined(USE_LOG) || m <= 1 || !info.count, PetscObjectComm((PetscObject)ksp), PETSC_ERR_PLIB, "PCApply() called %d times", info.count);
 85:   PetscCall(PetscLogEventRegister("PCMatApply", PC_CLASSID, &event));
 86:   PetscCall(PetscLogEventGetPerfInfo(PETSC_DETERMINE, event, &info));
 87:   PetscCheck(!PetscDefined(USE_LOG) || m <= 1 || info.count, PetscObjectComm((PetscObject)ksp), PETSC_ERR_PLIB, "PCMatApply() never called");
 88:   PetscCall(PetscFinalize());
 89:   return 0;
 90: }

 92: /*TEST
 93:    # KSPHPDDM does either pseudo-blocking or "true" blocking, all tests should succeed with other -ksp_hpddm_type
 94:    testset:
 95:       nsize: 1
 96:       args: -pc_type {{bjacobi lu ilu mat cholesky icc none shell asm gasm}shared output}
 97:       test:
 98:          suffix: 1
 99:          output_file: output/ex77_preonly.out
100:          args: -ksp_type preonly
101:       test:
102:          suffix: 1_hpddm
103:          output_file: output/ex77_preonly.out
104:          requires: hpddm
105:          args: -ksp_type hpddm -ksp_hpddm_type preonly

107:    testset:
108:       nsize: 1
109:       args: -pc_type ksp
110:       test:
111:          suffix: 2
112:          output_file: output/ex77_preonly.out
113:          args: -ksp_ksp_type preonly -ksp_type preonly
114:       test:
115:          suffix: 2_hpddm
116:          output_file: output/ex77_preonly.out
117:          requires: hpddm
118:          args: -ksp_ksp_type hpddm -ksp_type hpddm -ksp_hpddm_type preonly -ksp_ksp_hpddm_type preonly

120:    testset:
121:       nsize: 1
122:       requires: h2opus
123:       args: -pc_type h2opus -pc_h2opus_init_mat_h2opus_leafsize 10
124:       test:
125:          suffix: 3
126:          output_file: output/ex77_preonly.out
127:          args: -ksp_type preonly
128:       test:
129:          suffix: 3_hpddm
130:          output_file: output/ex77_preonly.out
131:          requires: hpddm
132:          args: -ksp_type hpddm -ksp_hpddm_type preonly

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

174:    testset:
175:       nsize: 1
176:       requires: suitesparse
177:       args: -pc_type qr
178:       test:
179:          suffix: 7
180:          output_file: output/ex77_preonly.out
181:          args: -ksp_type preonly
182:       test:
183:          suffix: 7_hpddm
184:          output_file: output/ex77_preonly.out
185:          requires: hpddm
186:          args: -ksp_type hpddm -ksp_hpddm_type preonly

188:    testset:
189:       nsize: 1
190:       requires: hpddm cuda
191:       args: -mat_type aijcusparse -ksp_type hpddm
192:       test:
193:          suffix: 8_hpddm
194:          output_file: output/ex77_preonly.out
195:       test:
196:          suffix: 8_hpddm_transpose
197:          output_file: output/ex77_preonly.out
198:          args: -pc_type icc -transpose

200:    testset:
201:       nsize: 1
202:       args: -pc_type {{cholesky icc none}shared output} -transpose
203:       test:
204:          suffix: 1_transpose
205:          output_file: output/ex77_preonly.out
206:          args: -ksp_type preonly
207:       test:
208:          suffix: 1_hpddm_transpose
209:          output_file: output/ex77_preonly.out
210:          requires: hpddm
211:          args: -ksp_type hpddm -ksp_hpddm_type preonly

213:    testset:
214:       requires: hypre !complex
215:       args: -pc_type hypre -pc_hypre_boomeramg_relax_type_all l1scaled-Jacobi -pc_hypre_boomeramg_no_CF
216:       test:
217:          suffix: 9
218:          output_file: output/ex77_preonly.out
219:          args: -ksp_type preonly
220:       test:
221:          suffix: 9_hpddm
222:          output_file: output/ex77_preonly.out
223:          requires: hpddm !hip
224:          args: -ksp_type hpddm -ksp_max_it 15 -ksp_error_if_not_converged

226: TEST*/