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: PetscCall(MatCreateDense(PETSC_COMM_WORLD, m, PETSC_DECIDE, PETSC_DECIDE, m, NULL, &B));
40: PetscCall(MatCreateDense(PETSC_COMM_WORLD, m, PETSC_DECIDE, PETSC_DECIDE, m, NULL, &X));
41: PetscCall(MatSetRandom(B, NULL));
42: PetscCall(MatSetFromOptions(A));
43: PetscCall(PetscObjectTypeCompareAny((PetscObject)A, &flg, MATSEQAIJCUSPARSE, MATMPIAIJCUSPARSE, ""));
44: if (flg) {
45: PetscCall(MatConvert(B, MATDENSECUDA, MAT_INPLACE_MATRIX, &B));
46: PetscCall(MatConvert(X, MATDENSECUDA, MAT_INPLACE_MATRIX, &X));
47: }
48: PetscCall(KSPCreate(PETSC_COMM_WORLD, &ksp));
49: PetscCall(KSPSetOperators(ksp, A, A));
50: PetscCall(KSPSetFromOptions(ksp));
51: PetscCall(KSPGetPC(ksp, &pc));
52: PetscCall(PCShellSetMatApply(pc, MatApply));
53: PetscCall(KSPMatSolve(ksp, B, X));
54: PetscCall(PCMatApply(pc, B, X));
55: if (transpose) {
56: PetscCall(KSPMatSolveTranspose(ksp, B, X));
57: PetscCall(PCMatApply(pc, B, X));
58: PetscCall(KSPMatSolve(ksp, B, X));
59: }
60: PetscCall(MatDestroy(&X));
61: PetscCall(MatDestroy(&B));
62: PetscCall(MatDestroy(&A));
63: PetscCall(KSPDestroy(&ksp));
64: PetscCall(PetscLogEventRegister("PCApply", PC_CLASSID, &event));
65: PetscCall(PetscLogEventGetPerfInfo(PETSC_DETERMINE, event, &info));
66: PetscCheck(!PetscDefined(USE_LOG) || m <= 1 || !info.count, PetscObjectComm((PetscObject)ksp), PETSC_ERR_PLIB, "PCApply() called %d times", info.count);
67: PetscCall(PetscLogEventRegister("PCMatApply", PC_CLASSID, &event));
68: PetscCall(PetscLogEventGetPerfInfo(PETSC_DETERMINE, event, &info));
69: PetscCheck(!PetscDefined(USE_LOG) || m <= 1 || info.count, PetscObjectComm((PetscObject)ksp), PETSC_ERR_PLIB, "PCMatApply() never called");
70: PetscCall(PetscFinalize());
71: return 0;
72: }
74: /*TEST
75: # KSPHPDDM does either pseudo-blocking or "true" blocking, all tests should succeed with other -ksp_hpddm_type
76: testset:
77: nsize: 1
78: args: -pc_type {{bjacobi lu ilu mat cholesky icc none shell asm gasm}shared output}
79: test:
80: suffix: 1
81: output_file: output/ex77_preonly.out
82: args: -ksp_type preonly
83: test:
84: suffix: 1_hpddm
85: output_file: output/ex77_preonly.out
86: requires: hpddm
87: args: -ksp_type hpddm -ksp_hpddm_type preonly
89: testset:
90: nsize: 1
91: args: -pc_type ksp
92: test:
93: suffix: 2
94: output_file: output/ex77_preonly.out
95: args: -ksp_ksp_type preonly -ksp_type preonly
96: test:
97: suffix: 2_hpddm
98: output_file: output/ex77_preonly.out
99: requires: hpddm
100: args: -ksp_ksp_type hpddm -ksp_type hpddm -ksp_hpddm_type preonly -ksp_ksp_hpddm_type preonly
102: testset:
103: nsize: 1
104: requires: h2opus
105: args: -pc_type h2opus -pc_h2opus_init_mat_h2opus_leafsize 10
106: test:
107: suffix: 3
108: output_file: output/ex77_preonly.out
109: args: -ksp_type preonly
110: test:
111: suffix: 3_hpddm
112: output_file: output/ex77_preonly.out
113: requires: hpddm
114: args: -ksp_type hpddm -ksp_hpddm_type preonly
116: testset:
117: nsize: 1
118: requires: spai
119: args: -pc_type spai
120: test:
121: suffix: 4
122: output_file: output/ex77_preonly.out
123: args: -ksp_type preonly
124: test:
125: suffix: 4_hpddm
126: output_file: output/ex77_preonly.out
127: requires: hpddm
128: args: -ksp_type hpddm -ksp_hpddm_type preonly
129: # special code path in PCMatApply() for PCBJACOBI when a block is shared by multiple processes
130: testset:
131: nsize: 2
132: args: -pc_type bjacobi -pc_bjacobi_blocks 1 -sub_pc_type none
133: test:
134: suffix: 5
135: output_file: output/ex77_preonly.out
136: args: -ksp_type preonly -sub_ksp_type preonly
137: test:
138: suffix: 5_hpddm
139: output_file: output/ex77_preonly.out
140: requires: hpddm
141: args: -ksp_type hpddm -ksp_hpddm_type preonly -sub_ksp_type hpddm
142: # special code path in PCMatApply() for PCGASM when a block is shared by multiple processes
143: testset:
144: nsize: 2
145: args: -pc_type gasm -pc_gasm_total_subdomains 1 -sub_pc_type none
146: test:
147: suffix: 6
148: output_file: output/ex77_preonly.out
149: args: -ksp_type preonly -sub_ksp_type preonly
150: test:
151: suffix: 6_hpddm
152: output_file: output/ex77_preonly.out
153: requires: hpddm
154: args: -ksp_type hpddm -ksp_hpddm_type preonly -sub_ksp_type hpddm
156: testset:
157: nsize: 1
158: requires: suitesparse
159: args: -pc_type qr
160: test:
161: suffix: 7
162: output_file: output/ex77_preonly.out
163: args: -ksp_type preonly
164: test:
165: suffix: 7_hpddm
166: output_file: output/ex77_preonly.out
167: requires: hpddm
168: args: -ksp_type hpddm -ksp_hpddm_type preonly
170: testset:
171: nsize: 1
172: requires: hpddm cuda
173: args: -mat_type aijcusparse -ksp_type hpddm
174: test:
175: suffix: 8_hpddm
176: output_file: output/ex77_preonly.out
177: test:
178: suffix: 8_hpddm_transpose
179: output_file: output/ex77_preonly.out
180: args: -pc_type icc -transpose
182: testset:
183: nsize: 1
184: args: -pc_type {{cholesky icc none}shared output} -transpose
185: test:
186: suffix: 1_transpose
187: output_file: output/ex77_preonly.out
188: args: -ksp_type preonly
189: test:
190: suffix: 1_hpddm_transpose
191: output_file: output/ex77_preonly.out
192: requires: hpddm
193: args: -ksp_type hpddm -ksp_hpddm_type preonly
195: TEST*/