Actual source code: ex42.c
2: static char help[] = "Solves a linear system in parallel with MINRES.\n\n";
4: #include <petscksp.h>
6: int main(int argc, char **args)
7: {
8: Vec x, b; /* approx solution, RHS */
9: Mat A; /* linear system matrix */
10: KSP ksp; /* linear solver context */
11: PC pc; /* preconditioner */
12: PetscScalar v = 0.0;
13: PetscInt Ii, Istart, Iend, m = 11;
14: PetscBool consistent = PETSC_TRUE;
16: PetscFunctionBeginUser;
17: PetscCall(PetscInitialize(&argc, &args, (char *)0, help));
18: PetscCall(PetscOptionsGetInt(NULL, NULL, "-m", &m, NULL));
19: PetscCall(PetscOptionsGetScalar(NULL, NULL, "-vv", &v, NULL));
20: PetscCall(PetscOptionsGetBool(NULL, NULL, "-consistent", &consistent, NULL));
22: /* Create parallel diagonal matrix */
23: PetscCall(MatCreate(PETSC_COMM_WORLD, &A));
24: PetscCall(MatSetSizes(A, PETSC_DECIDE, PETSC_DECIDE, m, m));
25: PetscCall(MatSetFromOptions(A));
26: PetscCall(MatMPIAIJSetPreallocation(A, 1, NULL, 1, NULL));
27: PetscCall(MatSeqAIJSetPreallocation(A, 1, NULL));
28: PetscCall(MatSetUp(A));
29: PetscCall(MatGetOwnershipRange(A, &Istart, &Iend));
31: for (Ii = Istart; Ii < Iend; Ii++) {
32: PetscScalar vv = (PetscReal)Ii + 1;
33: PetscCall(MatSetValues(A, 1, &Ii, 1, &Ii, &vv, INSERT_VALUES));
34: }
36: /* Make A singular or indefinite */
37: Ii = m - 1; /* last diagonal entry */
38: PetscCall(MatSetValues(A, 1, &Ii, 1, &Ii, &v, INSERT_VALUES));
39: PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
40: PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
42: PetscCall(MatCreateVecs(A, &x, &b));
43: if (consistent) {
44: PetscCall(VecSet(x, 1.0));
45: PetscCall(MatMult(A, x, b));
46: PetscCall(VecSet(x, 0.0));
47: } else {
48: PetscCall(VecSet(b, 1.0));
49: }
51: /* Create linear solver context */
52: PetscCall(KSPCreate(PETSC_COMM_WORLD, &ksp));
53: PetscCall(KSPSetOperators(ksp, A, A));
54: PetscCall(KSPSetType(ksp, KSPMINRES));
55: PetscCall(KSPGetPC(ksp, &pc));
56: PetscCall(PCSetType(pc, PCNONE));
57: PetscCall(KSPSetFromOptions(ksp));
58: PetscCall(KSPSolve(ksp, b, x));
60: /* Test reuse */
61: PetscCall(KSPSetInitialGuessNonzero(ksp, PETSC_TRUE));
62: PetscCall(VecSet(x, 0.0));
63: PetscCall(KSPSolve(ksp, b, x));
65: PetscCall(VecView(x, PETSC_VIEWER_STDOUT_WORLD));
67: /* Free work space. */
68: PetscCall(KSPDestroy(&ksp));
69: PetscCall(VecDestroy(&x));
70: PetscCall(VecDestroy(&b));
71: PetscCall(MatDestroy(&A));
73: PetscCall(PetscFinalize());
74: return 0;
75: }
77: /*TEST
79: test:
80: args: -ksp_converged_reason
82: test:
83: suffix: 2
84: nsize: 3
85: args: -ksp_converged_reason
87: test:
88: suffix: minres_qlp
89: args: -ksp_converged_reason -ksp_minres_qlp -ksp_minres_monitor
91: test:
92: suffix: minres_qlp_nonconsistent
93: args: -ksp_converged_reason -ksp_minres_qlp -ksp_minres_monitor -consistent 0
95: test:
96: suffix: minres_neg_curve
97: args: -ksp_converged_neg_curve -vv -1 -ksp_converged_reason -ksp_minres_qlp {{0 1}}
99: test:
100: suffix: cg_neg_curve
101: args: -ksp_converged_neg_curve -vv -1 -ksp_converged_reason -ksp_type {{cg stcg}}
102: requires: !complex
104: TEST*/