Actual source code: ex87.c
1: #include <petscksp.h>
2: #include <petsc/private/petscimpl.h>
4: static char help[] = "Solves a saddle-point linear system using PCHPDDM.\n\n";
6: static PetscErrorCode MatAndISLoad(const char *prefix, const char *identifier, Mat A, IS is, Mat N, PetscMPIInt size);
8: int main(int argc, char **args)
9: {
10: Vec b, x; /* computed solution and RHS */
11: Mat A[4], aux[2], S; /* linear system matrix */
12: KSP ksp, *subksp; /* linear solver context */
13: PC pc;
14: IS is[2];
15: PetscMPIInt size;
16: PetscInt m, M, n, N, id = 0;
17: PetscViewer viewer;
18: const char *const system[] = {"elasticity", "stokes", "diffusion"};
19: /* "elasticity":
20: * 2D linear elasticity with rubber-like and steel-like material coefficients, i.e., Poisson's ratio \in {0.4999, 0.35} and Young's modulus \in {0.01 GPa, 200.0 GPa}
21: * discretized by order 2 (resp. 0) Lagrange finite elements in displacements (resp. pressure) on a triangle mesh
22: * "stokes":
23: * 2D lid-driven cavity with constant viscosity
24: * discretized by order 2 (resp. 1) Lagrange finite elements, i.e., lowest-order Taylor--Hood finite elements, in velocities (resp. pressure) on a triangle mesh
25: * if the option -empty_A11 is not set (or set to false), a pressure with a zero mean-value is computed
26: * "diffusion":
27: * 2D primal-dual nonsymmetric diffusion equation
28: * discretized by order 2 (resp. 1) Lagrange finite elements in primal (resp. dual) unknowns on a triangle mesh
29: */
30: char dir[PETSC_MAX_PATH_LEN], prefix[PETSC_MAX_PATH_LEN];
31: PetscBool flg[4] = {PETSC_FALSE, PETSC_FALSE, PETSC_FALSE, PETSC_FALSE};
33: PetscFunctionBeginUser;
34: PetscCall(PetscInitialize(&argc, &args, NULL, help));
35: PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
36: PetscCheck(size == 4, PETSC_COMM_WORLD, PETSC_ERR_USER, "This example requires 4 processes");
37: PetscCall(PetscOptionsGetEList(NULL, NULL, "-system", system, PETSC_STATIC_ARRAY_LENGTH(system), &id, NULL));
38: if (id == 1) PetscCall(PetscOptionsGetBool(NULL, NULL, "-empty_A11", flg, NULL));
39: for (PetscInt i = 0; i < 2; ++i) {
40: PetscCall(MatCreate(PETSC_COMM_WORLD, A + (i ? 3 : 0)));
41: if (id < 2) {
42: PetscCall(ISCreate(PETSC_COMM_SELF, is + i));
43: PetscCall(MatCreate(PETSC_COMM_SELF, aux + i));
44: } else {
45: is[i] = NULL;
46: aux[i] = NULL;
47: }
48: }
49: PetscCall(PetscStrncpy(dir, ".", sizeof(dir)));
50: PetscCall(PetscOptionsGetString(NULL, NULL, "-load_dir", dir, sizeof(dir), NULL));
51: /* loading matrices and auxiliary data for the diagonal blocks */
52: PetscCall(PetscSNPrintf(prefix, sizeof(prefix), "%s/%s", dir, id == 2 ? "C" : (id == 1 ? "B" : "A")));
53: PetscCall(MatAndISLoad(prefix, "00", A[0], is[0], aux[0], size));
54: PetscCall(MatAndISLoad(prefix, "11", A[3], is[1], aux[1], size));
55: /* loading the off-diagonal block with a coherent row/column layout */
56: PetscCall(MatCreate(PETSC_COMM_WORLD, A + 2));
57: PetscCall(MatGetLocalSize(A[0], &n, NULL));
58: PetscCall(MatGetSize(A[0], &N, NULL));
59: PetscCall(MatGetLocalSize(A[3], &m, NULL));
60: PetscCall(MatGetSize(A[3], &M, NULL));
61: PetscCall(MatSetSizes(A[2], m, n, M, N));
62: PetscCall(PetscSNPrintf(prefix, sizeof(prefix), "%s/%s10.dat", dir, id == 2 ? "C" : (id == 1 ? "B" : "A")));
63: PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, prefix, FILE_MODE_READ, &viewer));
64: PetscCall(MatLoad(A[2], viewer));
65: PetscCall(PetscViewerDestroy(&viewer));
66: if (id < 2) {
67: /* transposing the off-diagonal block */
68: PetscCall(PetscOptionsGetBool(NULL, NULL, "-transpose", flg + 1, NULL));
69: PetscCall(PetscOptionsGetBool(NULL, NULL, "-permute", flg + 2, NULL));
70: PetscCall(PetscOptionsGetBool(NULL, NULL, "-explicit", flg + 3, NULL));
71: if (flg[1]) {
72: if (flg[2]) {
73: PetscCall(MatTranspose(A[2], MAT_INITIAL_MATRIX, A + 1));
74: PetscCall(MatDestroy(A + 2));
75: }
76: if (!flg[3]) PetscCall(MatCreateTranspose(A[2 - flg[2]], A + 1 + flg[2]));
77: else PetscCall(MatTranspose(A[2 - flg[2]], MAT_INITIAL_MATRIX, A + 1 + flg[2]));
78: } else {
79: if (flg[2]) {
80: PetscCall(MatHermitianTranspose(A[2], MAT_INITIAL_MATRIX, A + 1));
81: PetscCall(MatDestroy(A + 2));
82: }
83: if (!flg[3]) PetscCall(MatCreateHermitianTranspose(A[2 - flg[2]], A + 1 + flg[2]));
84: else PetscCall(MatHermitianTranspose(A[2 - flg[2]], MAT_INITIAL_MATRIX, A + 1 + flg[2]));
85: }
86: } else {
87: PetscCall(MatCreate(PETSC_COMM_WORLD, A + 1));
88: PetscCall(MatSetSizes(A[1], n, m, N, M));
89: PetscCall(PetscSNPrintf(prefix, sizeof(prefix), "%s/C01.dat", dir));
90: PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, prefix, FILE_MODE_READ, &viewer));
91: PetscCall(MatLoad(A[1], viewer));
92: PetscCall(PetscViewerDestroy(&viewer));
93: }
94: if (flg[0]) PetscCall(MatDestroy(A + 3));
95: else {
96: PetscCall(PetscOptionsGetBool(NULL, NULL, "-diagonal_A11", flg, NULL));
97: if (flg[0]) {
98: PetscCall(MatDestroy(A + 3));
99: PetscCall(MatCreateConstantDiagonal(PETSC_COMM_WORLD, m, m, M, M, PETSC_SMALL, A + 3));
100: }
101: }
102: flg[1] = PETSC_FALSE;
103: PetscCall(PetscOptionsGetBool(NULL, NULL, "-all_transpose", flg + 1, NULL));
104: if (flg[1] && flg[2]) {
105: PetscCall(MatTranspose(A[1], MAT_INITIAL_MATRIX, &S));
106: PetscCall(MatDestroy(A + 1));
107: PetscCall(MatCreateHermitianTranspose(S, A + 1));
108: PetscCall(MatDestroy(&S));
109: }
110: /* global coefficient matrix */
111: PetscCall(MatCreateNest(PETSC_COMM_WORLD, 2, NULL, 2, NULL, A, &S));
112: PetscCall(KSPCreate(PETSC_COMM_WORLD, &ksp));
113: PetscCall(KSPSetOperators(ksp, S, S));
114: PetscCall(KSPGetPC(ksp, &pc));
115: /* outer preconditioner */
116: PetscCall(PCSetType(pc, PCFIELDSPLIT));
117: PetscCall(PCFieldSplitSetType(pc, PC_COMPOSITE_SCHUR));
118: PetscCall(PCFieldSplitSetSchurPre(pc, PC_FIELDSPLIT_SCHUR_PRE_SELF, NULL));
119: PetscCall(PCSetFromOptions(pc));
120: PetscCall(PCSetUp(pc));
121: PetscCall(PCFieldSplitGetSubKSP(pc, &n, &subksp));
122: PetscCall(KSPGetPC(subksp[0], &pc));
123: /* inner preconditioner associated to top-left block */
124: #if defined(PETSC_HAVE_HPDDM) && defined(PETSC_HAVE_DYNAMIC_LIBRARIES) && defined(PETSC_USE_SHARED_LIBRARIES)
125: PetscCall(PCSetType(pc, PCHPDDM));
126: PetscCall(PCHPDDMSetAuxiliaryMat(pc, is[0], aux[0], NULL, NULL));
127: #endif
128: PetscCall(PCSetFromOptions(pc));
129: PetscCall(KSPGetPC(subksp[1], &pc));
130: /* inner preconditioner associated to Schur complement, which will be set internally to a PCKSP */
131: #if defined(PETSC_HAVE_HPDDM) && defined(PETSC_HAVE_DYNAMIC_LIBRARIES) && defined(PETSC_USE_SHARED_LIBRARIES)
132: PetscCall(PCSetType(pc, PCHPDDM));
133: if (!flg[0]) PetscCall(PCHPDDMSetAuxiliaryMat(pc, is[1], aux[1], NULL, NULL));
134: #endif
135: PetscCall(PCSetFromOptions(pc));
136: PetscCall(PetscFree(subksp));
137: PetscCall(KSPSetFromOptions(ksp));
138: PetscCall(MatCreateVecs(S, &b, &x));
139: PetscCall(PetscSNPrintf(prefix, sizeof(prefix), "%s/rhs_%s.dat", dir, id == 2 ? "C" : (id == 1 ? "B" : "A")));
140: PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, prefix, FILE_MODE_READ, &viewer));
141: PetscCall(VecLoad(b, viewer));
142: PetscCall(PetscViewerDestroy(&viewer));
143: PetscCall(KSPSolve(ksp, b, x));
144: flg[1] = PETSC_FALSE;
145: PetscCall(PetscOptionsGetBool(NULL, NULL, "-viewer", flg + 1, NULL));
146: if (flg[1]) PetscCall(PCView(pc, PETSC_VIEWER_STDOUT_WORLD));
147: flg[1] = PETSC_FALSE;
148: PetscCall(PetscOptionsGetBool(NULL, NULL, "-successive_solves", flg + 1, NULL));
149: if (flg[1]) {
150: KSPConvergedReason reason[2];
151: PetscInt iterations[2];
152: PetscCall(KSPGetConvergedReason(ksp, reason));
153: PetscCall(KSPGetTotalIterations(ksp, iterations));
154: PetscCall(KSPMonitorCancel(ksp));
155: PetscCall(PetscOptionsClearValue(NULL, "-ksp_monitor"));
156: PetscCall(PetscObjectStateIncrease((PetscObject)S));
157: PetscCall(KSPSetUp(ksp));
158: PetscCall(KSPGetPC(ksp, &pc));
159: PetscCall(PCFieldSplitGetSubKSP(pc, &n, &subksp));
160: PetscCall(KSPGetPC(subksp[0], &pc));
161: #if defined(PETSC_HAVE_HPDDM) && defined(PETSC_HAVE_DYNAMIC_LIBRARIES) && defined(PETSC_USE_SHARED_LIBRARIES)
162: PetscCall(PCHPDDMSetAuxiliaryMat(pc, is[0], aux[0], NULL, NULL));
163: #endif
164: PetscCall(PCSetFromOptions(pc));
165: PetscCall(KSPGetPC(subksp[1], &pc));
166: #if defined(PETSC_HAVE_HPDDM) && defined(PETSC_HAVE_DYNAMIC_LIBRARIES) && defined(PETSC_USE_SHARED_LIBRARIES)
167: if (!flg[0]) PetscCall(PCHPDDMSetAuxiliaryMat(pc, is[1], aux[1], NULL, NULL));
168: #endif
169: PetscCall(PCSetFromOptions(pc));
170: PetscCall(PetscFree(subksp));
171: PetscCall(KSPSolve(ksp, b, x));
172: PetscCall(KSPGetConvergedReason(ksp, reason + 1));
173: PetscCall(KSPGetTotalIterations(ksp, iterations + 1));
174: iterations[1] -= iterations[0];
175: PetscCheck(reason[0] == reason[1] && PetscAbs(iterations[0] - iterations[1]) <= 3, PetscObjectComm((PetscObject)ksp), PETSC_ERR_PLIB, "Successive calls to KSPSolve() did not converge for the same reason (%s v. %s) or with the same number of iterations (+/- 3, %" PetscInt_FMT " v. %" PetscInt_FMT ")", KSPConvergedReasons[reason[0]], KSPConvergedReasons[reason[1]], iterations[0], iterations[1]);
176: }
177: PetscCall(VecDestroy(&x));
178: PetscCall(VecDestroy(&b));
179: PetscCall(KSPDestroy(&ksp));
180: PetscCall(MatDestroy(&S));
181: PetscCall(MatDestroy(A + 1));
182: PetscCall(MatDestroy(A + 2));
183: for (PetscInt i = 0; i < 2; ++i) {
184: PetscCall(MatDestroy(A + (i ? 3 : 0)));
185: PetscCall(MatDestroy(aux + i));
186: PetscCall(ISDestroy(is + i));
187: }
188: PetscCall(PetscFinalize());
189: return 0;
190: }
192: PetscErrorCode MatAndISLoad(const char *prefix, const char *identifier, Mat A, IS is, Mat aux, PetscMPIInt size)
193: {
194: Mat tmp[3];
195: IS sizes;
196: const PetscInt *idx;
197: PetscInt m;
198: PetscLayout map;
199: PetscViewer viewer;
200: char name[PETSC_MAX_PATH_LEN];
202: PetscFunctionBeginUser;
203: PetscCall(PetscSNPrintf(name, sizeof(name), "%s%s_sizes_%d.dat", prefix, identifier, size));
204: PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, name, FILE_MODE_READ, &viewer));
205: PetscCall(ISCreate(PETSC_COMM_WORLD, &sizes));
206: PetscCall(ISLoad(sizes, viewer));
207: PetscCall(ISGetIndices(sizes, &idx));
208: PetscCall(MatSetSizes(A, idx[0], idx[1], idx[2], idx[3]));
209: if (is && aux) {
210: PetscCall(MatCreate(PETSC_COMM_WORLD, tmp));
211: PetscCall(MatSetSizes(tmp[0], idx[4], idx[4], PETSC_DETERMINE, PETSC_DETERMINE));
212: PetscCall(MatSetUp(tmp[0]));
213: }
214: PetscCall(ISRestoreIndices(sizes, &idx));
215: PetscCall(ISDestroy(&sizes));
216: PetscCall(PetscViewerDestroy(&viewer));
217: PetscCall(PetscSNPrintf(name, sizeof(name), "%s%s.dat", prefix, identifier));
218: PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, name, FILE_MODE_READ, &viewer));
219: PetscCall(MatLoad(A, viewer));
220: PetscCall(PetscViewerDestroy(&viewer));
221: if (is && aux) {
222: PetscCall(ISCreate(PETSC_COMM_WORLD, &sizes));
223: PetscCall(MatGetLayouts(tmp[0], &map, NULL));
224: PetscCall(ISSetLayout(sizes, map));
225: PetscCall(PetscSNPrintf(name, sizeof(name), "%s%s_is_%d.dat", prefix, identifier, size));
226: PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, name, FILE_MODE_READ, &viewer));
227: PetscCall(ISLoad(sizes, viewer));
228: PetscCall(ISGetLocalSize(sizes, &m));
229: PetscCall(ISGetIndices(sizes, &idx));
230: PetscCall(ISSetType(is, ISGENERAL));
231: PetscCall(ISGeneralSetIndices(is, m, idx, PETSC_COPY_VALUES));
232: PetscCall(ISRestoreIndices(sizes, &idx));
233: PetscCall(ISDestroy(&sizes));
234: PetscCall(PetscViewerDestroy(&viewer));
235: PetscCall(PetscSNPrintf(name, sizeof(name), "%s%s_aux_%d.dat", prefix, identifier, size));
236: PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, name, FILE_MODE_READ, &viewer));
237: PetscCall(MatLoad(tmp[0], viewer));
238: PetscCall(PetscViewerDestroy(&viewer));
239: PetscCall(MatGetDiagonalBlock(tmp[0], tmp + 1));
240: PetscCall(MatDuplicate(tmp[1], MAT_COPY_VALUES, tmp + 2));
241: PetscCall(MatHeaderReplace(aux, tmp + 2));
242: PetscCall(MatDestroy(tmp));
243: }
244: PetscFunctionReturn(PETSC_SUCCESS);
245: }
247: /*TEST
249: testset:
250: requires: datafilespath hpddm slepc double !complex !defined(PETSC_USE_64BIT_INDICES) defined(PETSC_HAVE_DYNAMIC_LIBRARIES) defined(PETSC_USE_SHARED_LIBRARIES)
251: nsize: 4
252: args: -load_dir ${DATAFILESPATH}/matrices/hpddm/GENEO -ksp_monitor -ksp_rtol 1e-4 -fieldsplit_ksp_max_it 100 -fieldsplit_pc_hpddm_levels_1_eps_nev 10 -fieldsplit_pc_hpddm_levels_1_st_share_sub_ksp -fieldsplit_pc_hpddm_has_neumann -fieldsplit_pc_hpddm_define_subdomains -fieldsplit_1_pc_hpddm_schur_precondition geneo -fieldsplit_pc_hpddm_coarse_pc_type redundant -fieldsplit_pc_hpddm_coarse_redundant_pc_type cholesky -fieldsplit_pc_hpddm_levels_1_sub_pc_type lu -fieldsplit_ksp_type fgmres -ksp_type fgmres -ksp_max_it 10 -fieldsplit_1_pc_hpddm_coarse_correction balanced -fieldsplit_1_pc_hpddm_levels_1_eps_gen_non_hermitian -fieldsplit_1_pc_hpddm_coarse_p 2
253: test:
254: requires: mumps
255: suffix: 1
256: args: -viewer -system {{elasticity stokes}separate output} -fieldsplit_1_pc_hpddm_ksp_pc_side left -fieldsplit_1_pc_hpddm_levels_1_sub_mat_mumps_icntl_26 1
257: filter: grep -v -e "action of " -e " " -e "block size" -e "total: nonzeros=" -e "using I-node" -e "aij" -e "transpose" -e "diagonal" -e "total number of" -e " rows="
258: test:
259: requires: mumps
260: suffix: 2
261: output_file: output/ex87_1_system-stokes.out
262: args: -viewer -system stokes -empty_A11 -transpose {{false true}shared output} -permute {{false true}shared output} -fieldsplit_1_pc_hpddm_ksp_pc_side right -fieldsplit_1_pc_hpddm_coarse_mat_type baij -fieldsplit_1_pc_hpddm_levels_1_sub_mat_mumps_icntl_26 1 -explicit {{false true}shared output}
263: filter: grep -v -e "action of " -e " " -e "block size" -e "total: nonzeros=" -e "using I-node" -e "aij" -e "transpose" -e "diagonal" -e "total number of" -e " rows=" | sed -e "s/ right preconditioning/ left preconditioning/g" -e "s/ using UNPRECONDITIONED/ using PRECONDITIONED/g"
264: test:
265: suffix: 1_petsc
266: args: -system {{elasticity stokes}separate output} -fieldsplit_1_pc_hpddm_ksp_pc_side left -fieldsplit_1_pc_hpddm_levels_1_sub_pc_factor_mat_solver_type petsc -fieldsplit_1_pc_hpddm_levels_1_eps_threshold 0.3 -permute
267: test:
268: suffix: 2_petsc
269: output_file: output/ex87_1_petsc_system-stokes.out
270: args: -system stokes -empty_A11 -transpose -fieldsplit_1_pc_hpddm_ksp_pc_side right -fieldsplit_1_pc_hpddm_levels_1_sub_pc_factor_mat_solver_type petsc -fieldsplit_1_pc_hpddm_coarse_mat_type baij -fieldsplit_1_pc_hpddm_levels_1_eps_threshold 0.3 -fieldsplit_1_pc_hpddm_levels_1_sub_pc_factor_shift_type inblocks -successive_solves
271: filter: sed -e "s/type: transpose/type: hermitiantranspose/g"
272: test:
273: suffix: threshold
274: requires: !defined(PETSC_HAVE_MKL_SPARSE_SP2M_FEATURE)
275: output_file: output/ex87_1_petsc_system-elasticity.out
276: args: -fieldsplit_1_pc_hpddm_ksp_pc_side left -fieldsplit_1_pc_hpddm_levels_1_eps_threshold 0.2 -fieldsplit_1_pc_hpddm_coarse_mat_type {{baij sbaij}shared output} -successive_solves
277: testset:
278: requires: datafilespath hpddm slepc double !complex !defined(PETSC_USE_64BIT_INDICES) defined(PETSC_HAVE_DYNAMIC_LIBRARIES) defined(PETSC_USE_SHARED_LIBRARIES)
279: nsize: 4
280: args: -load_dir ${DATAFILESPATH}/matrices/hpddm/GENEO -ksp_monitor -ksp_rtol 1e-4 -fieldsplit_ksp_max_it 100 -fieldsplit_pc_hpddm_levels_1_st_share_sub_ksp -fieldsplit_pc_hpddm_define_subdomains -fieldsplit_1_pc_hpddm_schur_precondition geneo -fieldsplit_pc_hpddm_coarse_pc_type redundant -fieldsplit_pc_hpddm_coarse_redundant_pc_type cholesky -fieldsplit_pc_hpddm_levels_1_sub_pc_type lu -fieldsplit_ksp_type fgmres -ksp_type fgmres -ksp_max_it 10 -fieldsplit_1_pc_hpddm_coarse_correction balanced -fieldsplit_1_pc_hpddm_levels_1_eps_gen_non_hermitian -fieldsplit_1_pc_hpddm_coarse_p 2 -system stokes -fieldsplit_1_pc_hpddm_ksp_pc_side left -fieldsplit_1_pc_hpddm_levels_1_sub_pc_factor_mat_solver_type petsc -fieldsplit_1_pc_hpddm_levels_1_eps_threshold 0.3
281: test:
282: suffix: diagonal
283: output_file: output/ex87_1_petsc_system-stokes.out
284: args: -fieldsplit_pc_hpddm_levels_1_eps_nev 10 -fieldsplit_0_pc_hpddm_has_neumann -diagonal_A11 {{false true}shared output}
285: test:
286: suffix: harmonic_overlap_2
287: output_file: output/ex87_1_petsc_system-stokes.out
288: args: -fieldsplit_0_pc_hpddm_harmonic_overlap 2 -fieldsplit_0_pc_hpddm_levels_1_svd_nsv 20 -diagonal_A11 -permute {{false true}shared output} -all_transpose
290: test:
291: requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) !hpddm !memkind
292: nsize: 4
293: suffix: selfp
294: output_file: output/ex41_1.out
295: filter: grep -v "CONVERGED_RTOL iterations"
296: args: -load_dir ${DATAFILESPATH}/matrices/hpddm/GENEO -system stokes -ksp_rtol 1e-4 -ksp_converged_reason -ksp_max_it 30 -pc_type fieldsplit -pc_fieldsplit_type schur -fieldsplit_ksp_type preonly -pc_fieldsplit_schur_precondition selfp -fieldsplit_pc_type bjacobi -fieldsplit_sub_pc_type lu -transpose {{false true}shared output} -fieldsplit_1_mat_schur_complement_ainv_type lump
298: test:
299: requires: datafilespath hpddm slepc double !complex !defined(PETSC_USE_64BIT_INDICES) defined(PETSC_HAVE_DYNAMIC_LIBRARIES) defined(PETSC_USE_SHARED_LIBRARIES)
300: nsize: 4
301: suffix: nonsymmetric_least_squares
302: output_file: output/ex41_1.out
303: filter: grep -v "CONVERGED_RTOL iterations"
304: args: -load_dir ${DATAFILESPATH}/matrices/hpddm/GENEO -system diffusion -ksp_rtol 1e-4 -ksp_converged_reason -ksp_max_it 20 -pc_type fieldsplit -pc_fieldsplit_type schur -fieldsplit_ksp_type preonly -fieldsplit_0_pc_type jacobi -prefix_push fieldsplit_1_ -pc_hpddm_schur_precondition least_squares -pc_hpddm_define_subdomains -prefix_push pc_hpddm_levels_1_ -sub_pc_type lu -sub_pc_factor_shift_type nonzero -eps_nev 5 -st_share_sub_ksp -prefix_pop -prefix_pop
306: TEST*/