xref: /petsc/src/ksp/ksp/tests/ex51.c (revision 19a666eebac44a514bcd9e8c725f1a1c73ecdf5f)
1 static char help[] = "Test PCFailedReason.\n\n";
2 
3 #include <petscksp.h>
4 
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   PetscErrorCode     ierr;
11   PetscInt           i,n = 10,col[3];
12   PetscMPIInt        size;
13   PetscScalar        value[3],alpha,beta,sx;
14   PetscBool          reverse=PETSC_FALSE;
15   KSPConvergedReason reason;
16   PCFailedReason     pcreason;
17 
18   ierr = PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr;
19   ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRQ(ierr);
20   if (size != 1) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_WRONG_MPI_SIZE,"This is a uniprocessor example only!");
21   ierr = PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL);CHKERRQ(ierr);
22   ierr = PetscOptionsGetBool(NULL,NULL,"-reverse",&reverse,NULL);CHKERRQ(ierr);
23 
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;
27 
28   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
29          Create the matrix
30      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
31   ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr);
32   ierr = MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,n,n);CHKERRQ(ierr);
33   ierr = MatSetFromOptions(A);CHKERRQ(ierr);
34   ierr = MatSetUp(A);CHKERRQ(ierr);
35 
36   value[0] = -1.0; value[1] = 2.0; value[2] = -1.0;
37   for (i=1; i<n-1; i++) {
38     col[0] = i-1; col[1] = i; col[2] = i+1;
39     ierr   = MatSetValues(A,1,&i,3,col,value,INSERT_VALUES);CHKERRQ(ierr);
40   }
41   i    = n - 1; col[0] = n - 2; col[1] = n - 1;
42   ierr = MatSetValues(A,1,&i,2,col,value,INSERT_VALUES);CHKERRQ(ierr);
43   i    = 0; col[0] = 0; col[1] = 1; value[0] = 2.0; value[1] = -1.0;
44   ierr = MatSetValues(A,1,&i,2,col,value,INSERT_VALUES);CHKERRQ(ierr);
45   ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
46   ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
47 
48   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
49                 Create the linear solver and set various options
50      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
51   ierr = KSPCreate(PETSC_COMM_WORLD,&ksp);CHKERRQ(ierr);
52   ierr = KSPSetOperators(ksp,A,A);CHKERRQ(ierr);
53   ierr = MatShift(A,reverse?-alpha:-beta);CHKERRQ(ierr);
54   ierr = KSPGetPC(ksp,&pc);CHKERRQ(ierr);
55   ierr = PCSetType(pc,PCLU);CHKERRQ(ierr);
56   ierr = KSPSetFromOptions(ksp);CHKERRQ(ierr);
57 
58   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
59                       Factorize first matrix
60      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
61   ierr = PetscPrintf(PETSC_COMM_WORLD,"First matrix\n");CHKERRQ(ierr);
62   ierr = KSPSetUp(ksp);CHKERRQ(ierr);
63   ierr = KSPGetConvergedReason(ksp,&reason);CHKERRQ(ierr);
64   if (reason < 0) {
65     ierr = PetscViewerPushFormat(PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)ksp)),PETSC_VIEWER_DEFAULT);CHKERRQ(ierr);
66     ierr = KSPConvergedReasonView(ksp,PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)ksp)));CHKERRQ(ierr);
67     ierr = PetscViewerPopFormat(PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)ksp)));CHKERRQ(ierr);
68   } else {
69     ierr = PetscPrintf(PETSC_COMM_WORLD,"Success!\n");CHKERRQ(ierr);
70   }
71 
72   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
73                       Factorize second matrix
74      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
75   ierr = MatShift(A,reverse?alpha-beta:beta-alpha);CHKERRQ(ierr);
76   ierr = KSPSetOperators(ksp,A,A);CHKERRQ(ierr);
77 
78   ierr = PetscPrintf(PETSC_COMM_WORLD,"Second matrix\n");CHKERRQ(ierr);
79   ierr = KSPSetUp(ksp);CHKERRQ(ierr);
80   ierr = KSPGetConvergedReason(ksp,&reason);CHKERRQ(ierr);
81   if (reason < 0) {
82     ierr = PetscViewerPushFormat(PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)ksp)),PETSC_VIEWER_DEFAULT);CHKERRQ(ierr);
83     ierr = KSPConvergedReasonView(ksp,PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)ksp)));CHKERRQ(ierr);
84     ierr = PetscViewerPopFormat(PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)ksp)));CHKERRQ(ierr);
85   } else {
86     ierr = PetscPrintf(PETSC_COMM_WORLD,"Success!\n");CHKERRQ(ierr);
87     ierr = PCGetFailedReason(pc,&pcreason);CHKERRQ(ierr);
88     ierr = PetscPrintf(PETSC_COMM_WORLD,"PC failed reason is %s\n",PCFailedReasons[pcreason]);CHKERRQ(ierr);
89   }
90 
91   /*
92      Free work space.
93   */
94   ierr = MatDestroy(&A);CHKERRQ(ierr);
95   ierr = KSPDestroy(&ksp);CHKERRQ(ierr);
96 
97   ierr = PetscFinalize();
98   return ierr;
99 }
100 
101 
102 /*TEST
103 
104    test:
105       args: -reverse
106 
107    test:
108       suffix: 2
109       args: -reverse -pc_type cholesky
110 
111 TEST*/
112