Actual source code: ex70.c


  2: static char help[] = "Solves an ill-conditioned tridiagonal linear system with KSP for testing GMRES breakdown tolerance.\n\n";

  4: #include <petscksp.h>

  6: int main(int argc, char **args)
  7: {
  8:   Vec         x, b, u; /* approx solution, RHS, exact solution */
  9:   Mat         A;       /* linear system matrix */
 10:   KSP         ksp;     /* linear solver context */
 11:   PetscInt    i, n = 10, col[3];
 12:   PetscMPIInt size;
 13:   PetscScalar value[3];

 16:   PetscInitialize(&argc, &args, (char *)0, help);
 17:   MPI_Comm_size(PETSC_COMM_WORLD, &size);

 20:   /*
 21:      Create vectors.  Note that we form 1 vector from scratch and
 22:      then duplicate as needed.
 23:   */
 24:   VecCreate(PETSC_COMM_WORLD, &x);
 25:   PetscObjectSetName((PetscObject)x, "Solution");
 26:   VecSetSizes(x, PETSC_DECIDE, n);
 27:   VecSetFromOptions(x);
 28:   VecDuplicate(x, &b);
 29:   VecDuplicate(x, &u);

 31:   MatCreate(PETSC_COMM_WORLD, &A);
 32:   MatSetSizes(A, PETSC_DECIDE, PETSC_DECIDE, n, n);
 33:   MatSetFromOptions(A);
 34:   MatSetUp(A);

 36:   /*
 37:      Set big off-diag values to make the system ill-conditioned
 38:   */
 39:   value[0] = 10.0;
 40:   value[1] = 2.0;
 41:   value[2] = 1.0;
 42:   for (i = 1; i < n - 1; i++) {
 43:     col[0] = i - 1;
 44:     col[1] = i;
 45:     col[2] = i + 1;
 46:     MatSetValues(A, 1, &i, 3, col, value, INSERT_VALUES);
 47:   }
 48:   i      = n - 1;
 49:   col[0] = n - 2;
 50:   col[1] = n - 1;
 51:   MatSetValues(A, 1, &i, 2, col, value, INSERT_VALUES);
 52:   i        = 0;
 53:   col[0]   = 0;
 54:   col[1]   = 1;
 55:   value[0] = 2.0;
 56:   value[1] = -1.0;
 57:   MatSetValues(A, 1, &i, 2, col, value, INSERT_VALUES);
 58:   MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
 59:   MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);

 61:   VecSet(u, 1.0);
 62:   MatMult(A, u, b);

 64:   KSPCreate(PETSC_COMM_WORLD, &ksp);
 65:   KSPSetOperators(ksp, A, A);
 66:   KSPSetFromOptions(ksp);
 67:   KSPSolve(ksp, b, x);

 69:   KSPSetInitialGuessNonzero(ksp, PETSC_TRUE);
 70:   PetscOptionsInsertString(NULL, "-ksp_type preonly -ksp_initial_guess_nonzero false");
 71:   PetscOptionsClearValue(NULL, "-ksp_converged_reason");
 72:   KSPSetFromOptions(ksp);
 73:   KSPSolve(ksp, b, x);

 75:   VecDestroy(&x);
 76:   VecDestroy(&u);
 77:   VecDestroy(&b);
 78:   MatDestroy(&A);
 79:   KSPDestroy(&ksp);

 81:   PetscFinalize();
 82:   return 0;
 83: }

 85: /*TEST

 87:    test:
 88:       requires: double !complex
 89:       args: -ksp_rtol  1e-18 -pc_type sor -ksp_converged_reason -ksp_gmres_breakdown_tolerance 1.e-9
 90:       output_file: output/ex70.out

 92: TEST*/