Actual source code: ex55.c
1: static const char help[] = "Example demonstrating PCCOMPOSITE where one of the inner PCs uses a different operator\n\
2: \n";
4: #include <petscksp.h>
6: int main(int argc, char **argv)
7: {
8: PetscInt n = 10, i, col[3];
9: Vec x, b;
10: Mat A, B;
11: KSP ksp;
12: PC pc, subpc;
13: PetscScalar value[3];
16: PetscInitialize(&argc, &argv, (char *)0, help);
18: /* Create a diagonal matrix with a given distribution of diagonal elements */
19: MatCreate(PETSC_COMM_WORLD, &A);
20: MatSetSizes(A, PETSC_DECIDE, PETSC_DECIDE, n, n);
21: MatSetFromOptions(A);
22: MatSetUp(A);
23: /*
24: Assemble matrix
25: */
26: value[0] = -1.0;
27: value[1] = 2.0;
28: value[2] = -1.0;
29: for (i = 1; i < n - 1; i++) {
30: col[0] = i - 1;
31: col[1] = i;
32: col[2] = i + 1;
33: MatSetValues(A, 1, &i, 3, col, value, INSERT_VALUES);
34: }
35: i = n - 1;
36: col[0] = n - 2;
37: col[1] = n - 1;
38: MatSetValues(A, 1, &i, 2, col, value, INSERT_VALUES);
39: i = 0;
40: col[0] = 0;
41: col[1] = 1;
42: value[0] = 2.0;
43: value[1] = -1.0;
44: MatSetValues(A, 1, &i, 2, col, value, INSERT_VALUES);
45: MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
46: MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);
48: MatCreateVecs(A, &x, &b);
50: /* Create a KSP object */
51: KSPCreate(PETSC_COMM_WORLD, &ksp);
52: KSPSetOperators(ksp, A, A);
54: /* Set up a composite preconditioner */
55: KSPGetPC(ksp, &pc);
56: PCSetType(pc, PCCOMPOSITE); /* default composite with single Identity PC */
57: PCCompositeSetType(pc, PC_COMPOSITE_ADDITIVE);
58: PCCompositeAddPCType(pc, PCLU);
59: PCCompositeGetPC(pc, 0, &subpc);
60: /* B is set to the diagonal of A; this demonstrates that setting the operator for a subpc changes the preconditioning */
61: MatDuplicate(A, MAT_DO_NOT_COPY_VALUES, &B);
62: MatGetDiagonal(A, b);
63: MatDiagonalSet(B, b, ADD_VALUES);
64: PCSetOperators(subpc, B, B);
65: PCCompositeAddPCType(pc, PCNONE);
67: KSPSetFromOptions(ksp);
68: KSPSolve(ksp, b, x);
70: KSPDestroy(&ksp);
71: MatDestroy(&A);
72: MatDestroy(&B);
73: VecDestroy(&x);
74: VecDestroy(&b);
75: PetscFinalize();
76: return 0;
77: }
79: /*TEST
81: test:
82: args: -ksp_monitor -pc_composite_type multiplicative
84: TEST*/