Actual source code: ex88.c

  1: static char help[] = "Tests solving linear system on 0 by 0 matrix, and KSPLSQR and user convergence test handling.\n\n";

  3: #include <petscksp.h>

  5: typedef struct {
  6:   PetscBool converged;
  7: } ConvergedCtx;

  9: static PetscErrorCode TestConvergence(KSP ksp, PetscInt it, PetscReal rnorm, KSPConvergedReason *reason, void *ctx)
 10: {
 11:   ConvergedCtx *user = (ConvergedCtx *)ctx;

 13:   PetscFunctionBeginUser;
 14:   if (user->converged) *reason = KSP_CONVERGED_USER;
 15:   else *reason = KSP_DIVERGED_USER;
 16:   PetscFunctionReturn(PETSC_SUCCESS);
 17: }

 19: int main(int argc, char **args)
 20: {
 21:   Mat         C;
 22:   PetscInt    n = 1, m0, m1, i;
 23:   Vec         b, x;
 24:   KSP         ksp;
 25:   PetscScalar one = 1;
 26:   PetscRandom rctx;

 28:   ConvergedCtx user;
 29:   user.converged = PETSC_TRUE;

 31:   PetscFunctionBeginUser;
 32:   PetscCall(PetscInitialize(&argc, &args, NULL, help));

 34:   /* create stiffness matrix */
 35:   PetscCall(MatCreate(PETSC_COMM_WORLD, &C));
 36:   PetscCall(MatSetSizes(C, n, n, PETSC_DETERMINE, PETSC_DETERMINE));
 37:   PetscCall(MatSetFromOptions(C));
 38:   PetscCall(MatSetUp(C));
 39:   PetscCall(MatGetOwnershipRange(C, &m0, &m1));
 40:   for (i = m0; i < m1; i++) PetscCall(MatSetValues(C, 1, &i, 1, &i, &one, INSERT_VALUES));
 41:   PetscCall(MatAssemblyBegin(C, MAT_FINAL_ASSEMBLY));
 42:   PetscCall(MatAssemblyEnd(C, MAT_FINAL_ASSEMBLY));

 44:   PetscCall(PetscOptionsGetBool(NULL, NULL, "-mark_converged", &user.converged, NULL));
 45:   /* create right-hand side and solution */
 46:   PetscCall(VecCreate(PETSC_COMM_WORLD, &x));
 47:   PetscCall(VecSetSizes(x, n, PETSC_DETERMINE));
 48:   PetscCall(VecSetFromOptions(x));
 49:   PetscCall(VecDuplicate(x, &b));
 50:   PetscCall(VecSet(x, 0.0));
 51:   PetscCall(PetscRandomCreate(PETSC_COMM_WORLD, &rctx));
 52:   PetscCall(VecSetRandom(b, rctx));
 53:   PetscCall(PetscRandomDestroy(&rctx));

 55:   /* solve linear system */
 56:   PetscCall(KSPCreate(PETSC_COMM_WORLD, &ksp));
 57:   PetscCall(KSPSetOperators(ksp, C, C));
 58:   PetscCall(KSPSetConvergenceTest(ksp, TestConvergence, &user, NULL));
 59:   PetscCall(KSPSetFromOptions(ksp));
 60:   PetscCall(KSPSolve(ksp, b, x));
 61:   PetscCall(KSPConvergedReasonView(ksp, NULL));

 63:   PetscCall(KSPDestroy(&ksp));
 64:   PetscCall(VecDestroy(&x));
 65:   PetscCall(VecDestroy(&b));
 66:   PetscCall(MatDestroy(&C));
 67:   PetscCall(PetscFinalize());
 68:   return 0;
 69: }

 71: /*TEST

 73:     test:
 74:       args:

 76:     test:
 77:       suffix: 2
 78:       args: -mark_converged 0

 80: TEST*/