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;
17: PetscInitialize(&argc,&args,(char*)0,help);
18: MPI_Comm_size(PETSC_COMM_WORLD,&size);
20: PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL);
21: PetscOptionsGetBool(NULL,NULL,"-reverse",&reverse,NULL);
23: sx = PetscSinReal(n*PETSC_PI/2/(n+1));
24: alpha = 4.0*sx*sx; /* alpha is the largest eigenvalue of the matrix */
25: beta = 4.0;
27: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
28: Create the matrix
29: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
30: MatCreate(PETSC_COMM_WORLD,&A);
31: MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,n,n);
32: MatSetFromOptions(A);
33: MatSetUp(A);
35: value[0] = -1.0; value[1] = 2.0; value[2] = -1.0;
36: for (i=1; i<n-1; i++) {
37: col[0] = i-1; col[1] = i; col[2] = i+1;
38: MatSetValues(A,1,&i,3,col,value,INSERT_VALUES);
39: }
40: i = n - 1; col[0] = n - 2; col[1] = n - 1;
41: MatSetValues(A,1,&i,2,col,value,INSERT_VALUES);
42: i = 0; col[0] = 0; col[1] = 1; value[0] = 2.0; value[1] = -1.0;
43: MatSetValues(A,1,&i,2,col,value,INSERT_VALUES);
44: MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
45: MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
47: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
48: Create the linear solver and set various options
49: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
50: KSPCreate(PETSC_COMM_WORLD,&ksp);
51: KSPSetOperators(ksp,A,A);
52: MatShift(A,reverse?-alpha:-beta);
53: KSPGetPC(ksp,&pc);
54: PCSetType(pc,PCLU);
55: KSPSetFromOptions(ksp);
57: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
58: Factorize first matrix
59: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
60: PetscPrintf(PETSC_COMM_WORLD,"First matrix\n");
61: KSPSetUp(ksp);
62: KSPGetConvergedReason(ksp,&reason);
63: if (reason < 0) {
64: PetscViewerPushFormat(PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)ksp)),PETSC_VIEWER_DEFAULT);
65: KSPConvergedReasonView(ksp,PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)ksp)));
66: PetscViewerPopFormat(PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)ksp)));
67: } else {
68: PetscPrintf(PETSC_COMM_WORLD,"Success!\n");
69: }
71: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
72: Factorize second matrix
73: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
74: MatShift(A,reverse?alpha-beta:beta-alpha);
75: KSPSetOperators(ksp,A,A);
77: PetscPrintf(PETSC_COMM_WORLD,"Second matrix\n");
78: KSPSetUp(ksp);
79: KSPGetConvergedReason(ksp,&reason);
80: if (reason < 0) {
81: PetscViewerPushFormat(PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)ksp)),PETSC_VIEWER_DEFAULT);
82: KSPConvergedReasonView(ksp,PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)ksp)));
83: PetscViewerPopFormat(PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)ksp)));
84: } else {
85: PetscPrintf(PETSC_COMM_WORLD,"Success!\n");
86: PCGetFailedReason(pc,&pcreason);
87: PetscPrintf(PETSC_COMM_WORLD,"PC failed reason is %s\n",PCFailedReasons[pcreason]);
88: }
90: /*
91: Free work space.
92: */
93: MatDestroy(&A);
94: KSPDestroy(&ksp);
96: PetscFinalize();
97: return 0;
98: }
100: /*TEST
102: test:
103: args: -reverse
105: test:
106: suffix: 2
107: args: -reverse -pc_type cholesky
109: TEST*/