Actual source code: ex51.c
1: static char help[] = "Test PCFailedReason.\n\n";
3: #include <petscksp.h>
5: int main(int argc, char **args)
6: {
7: Mat A; /* linear system matrix */
8: KSP ksp; /* linear solver context */
9: PC pc; /* preconditioner context */
10: PetscInt i, n = 10, col[3];
11: PetscMPIInt size;
12: PetscScalar value[3], alpha, beta, sx;
13: PetscBool reverse = PETSC_FALSE;
14: KSPConvergedReason reason;
15: PCFailedReason pcreason;
18: PetscInitialize(&argc, &args, (char *)0, help);
19: MPI_Comm_size(PETSC_COMM_WORLD, &size);
21: PetscOptionsGetInt(NULL, NULL, "-n", &n, NULL);
22: PetscOptionsGetBool(NULL, NULL, "-reverse", &reverse, NULL);
24: sx = PetscSinReal(n * PETSC_PI / 2 / (n + 1));
25: alpha = 4.0 * sx * sx; /* alpha is the largest eigenvalue of the matrix */
26: beta = 4.0;
28: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
29: Create the matrix
30: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
31: MatCreate(PETSC_COMM_WORLD, &A);
32: MatSetSizes(A, PETSC_DECIDE, PETSC_DECIDE, n, n);
33: MatSetFromOptions(A);
34: MatSetUp(A);
36: value[0] = -1.0;
37: value[1] = 2.0;
38: value[2] = -1.0;
39: for (i = 1; i < n - 1; i++) {
40: col[0] = i - 1;
41: col[1] = i;
42: col[2] = i + 1;
43: MatSetValues(A, 1, &i, 3, col, value, INSERT_VALUES);
44: }
45: i = n - 1;
46: col[0] = n - 2;
47: col[1] = n - 1;
48: MatSetValues(A, 1, &i, 2, col, value, INSERT_VALUES);
49: i = 0;
50: col[0] = 0;
51: col[1] = 1;
52: value[0] = 2.0;
53: value[1] = -1.0;
54: MatSetValues(A, 1, &i, 2, col, value, INSERT_VALUES);
55: MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
56: MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);
58: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
59: Create the linear solver and set various options
60: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
61: KSPCreate(PETSC_COMM_WORLD, &ksp);
62: KSPSetOperators(ksp, A, A);
63: MatShift(A, reverse ? -alpha : -beta);
64: KSPGetPC(ksp, &pc);
65: PCSetType(pc, PCLU);
66: KSPSetFromOptions(ksp);
68: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
69: Factorize first matrix
70: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
71: PetscPrintf(PETSC_COMM_WORLD, "First matrix\n");
72: KSPSetUp(ksp);
73: KSPGetConvergedReason(ksp, &reason);
74: if (reason < 0) {
75: PetscViewerPushFormat(PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)ksp)), PETSC_VIEWER_DEFAULT);
76: KSPConvergedReasonView(ksp, PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)ksp)));
77: PetscViewerPopFormat(PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)ksp)));
78: } else {
79: PetscPrintf(PETSC_COMM_WORLD, "Success!\n");
80: }
82: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
83: Factorize second matrix
84: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
85: MatShift(A, reverse ? alpha - beta : beta - alpha);
86: KSPSetOperators(ksp, A, A);
88: PetscPrintf(PETSC_COMM_WORLD, "Second matrix\n");
89: KSPSetUp(ksp);
90: KSPGetConvergedReason(ksp, &reason);
91: if (reason < 0) {
92: PetscViewerPushFormat(PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)ksp)), PETSC_VIEWER_DEFAULT);
93: KSPConvergedReasonView(ksp, PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)ksp)));
94: PetscViewerPopFormat(PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)ksp)));
95: } else {
96: PetscPrintf(PETSC_COMM_WORLD, "Success!\n");
97: PCGetFailedReason(pc, &pcreason);
98: PetscPrintf(PETSC_COMM_WORLD, "PC failed reason is %s\n", PCFailedReasons[pcreason]);
99: }
101: /*
102: Free work space.
103: */
104: MatDestroy(&A);
105: KSPDestroy(&ksp);
107: PetscFinalize();
108: return 0;
109: }
111: /*TEST
113: test:
114: args: -reverse
116: test:
117: suffix: 2
118: args: -reverse -pc_type cholesky
120: TEST*/