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:   PetscFunctionBeginUser;
 18:   PetscCall(PetscInitialize(&argc, &args, (char *)0, help));
 19:   PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
 20:   PetscCheck(size == 1, PETSC_COMM_WORLD, PETSC_ERR_WRONG_MPI_SIZE, "This is a uniprocessor example only!");
 21:   PetscCall(PetscOptionsGetInt(NULL, NULL, "-n", &n, NULL));
 22:   PetscCall(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:   PetscCall(MatCreate(PETSC_COMM_WORLD, &A));
 32:   PetscCall(MatSetSizes(A, PETSC_DECIDE, PETSC_DECIDE, n, n));
 33:   PetscCall(MatSetFromOptions(A));
 34:   PetscCall(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:     PetscCall(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:   PetscCall(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:   PetscCall(MatSetValues(A, 1, &i, 2, col, value, INSERT_VALUES));
 55:   PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
 56:   PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));

 58:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 59:                 Create the linear solver and set various options
 60:      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
 61:   PetscCall(KSPCreate(PETSC_COMM_WORLD, &ksp));
 62:   PetscCall(KSPSetOperators(ksp, A, A));
 63:   PetscCall(MatShift(A, reverse ? -alpha : -beta));
 64:   PetscCall(KSPGetPC(ksp, &pc));
 65:   PetscCall(PCSetType(pc, PCLU));
 66:   PetscCall(KSPSetFromOptions(ksp));

 68:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 69:                       Factorize first matrix
 70:      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
 71:   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "First matrix\n"));
 72:   PetscCall(KSPSetUp(ksp));
 73:   PetscCall(KSPGetConvergedReason(ksp, &reason));
 74:   if (reason < 0) {
 75:     PetscCall(PetscViewerPushFormat(PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)ksp)), PETSC_VIEWER_DEFAULT));
 76:     PetscCall(KSPConvergedReasonView(ksp, PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)ksp))));
 77:     PetscCall(PetscViewerPopFormat(PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)ksp))));
 78:   } else {
 79:     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Success!\n"));
 80:   }

 82:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 83:                       Factorize second matrix
 84:      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
 85:   PetscCall(MatShift(A, reverse ? alpha - beta : beta - alpha));
 86:   PetscCall(KSPSetOperators(ksp, A, A));

 88:   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Second matrix\n"));
 89:   PetscCall(KSPSetUp(ksp));
 90:   PetscCall(KSPGetConvergedReason(ksp, &reason));
 91:   if (reason < 0) {
 92:     PetscCall(PetscViewerPushFormat(PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)ksp)), PETSC_VIEWER_DEFAULT));
 93:     PetscCall(KSPConvergedReasonView(ksp, PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)ksp))));
 94:     PetscCall(PetscViewerPopFormat(PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)ksp))));
 95:   } else {
 96:     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Success!\n"));
 97:     PetscCall(PCGetFailedReason(pc, &pcreason));
 98:     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "PC failed reason is %s\n", PCFailedReasons[pcreason]));
 99:   }

101:   /*
102:      Free work space.
103:   */
104:   PetscCall(MatDestroy(&A));
105:   PetscCall(KSPDestroy(&ksp));

107:   PetscCall(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*/