Actual source code: ex42.c


  2: static char help[] = "Solves a linear system in parallel with MINRES. Modified from ../tutorials/ex2.c \n\n";

  4: #include <petscksp.h>

  6: int main(int argc, char **args)
  7: {
  8:   Vec         x, b; /* approx solution, RHS */
  9:   Mat         A;    /* linear system matrix */
 10:   KSP         ksp;  /* linear solver context */
 11:   PetscInt    Ii, Istart, Iend, m = 11;
 12:   PetscScalar v;

 15:   PetscInitialize(&argc, &args, (char *)0, help);
 16:   PetscOptionsGetInt(NULL, NULL, "-m", &m, NULL);

 18:   /* Create parallel diagonal matrix */
 19:   MatCreate(PETSC_COMM_WORLD, &A);
 20:   MatSetSizes(A, PETSC_DECIDE, PETSC_DECIDE, m, m);
 21:   MatSetFromOptions(A);
 22:   MatMPIAIJSetPreallocation(A, 1, NULL, 1, NULL);
 23:   MatSeqAIJSetPreallocation(A, 1, NULL);
 24:   MatSetUp(A);
 25:   MatGetOwnershipRange(A, &Istart, &Iend);

 27:   for (Ii = Istart; Ii < Iend; Ii++) {
 28:     v = (PetscReal)Ii + 1;
 29:     MatSetValues(A, 1, &Ii, 1, &Ii, &v, INSERT_VALUES);
 30:   }
 31:   /* Make A sigular */
 32:   Ii = m - 1; /* last diagonal entry */
 33:   v  = 0.0;
 34:   MatSetValues(A, 1, &Ii, 1, &Ii, &v, INSERT_VALUES);
 35:   MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
 36:   MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);

 38:   /* A is symmetric. Set symmetric flag to enable KSP_type = minres */
 39:   MatSetOption(A, MAT_SYMMETRIC, PETSC_TRUE);

 41:   VecCreate(PETSC_COMM_WORLD, &b);
 42:   VecSetSizes(b, PETSC_DECIDE, m);
 43:   VecSetFromOptions(b);
 44:   VecDuplicate(b, &x);
 45:   VecSet(x, 1.0);
 46:   MatMult(A, x, b);
 47:   VecSet(x, 0.0);

 49:   /* Create linear solver context */
 50:   KSPCreate(PETSC_COMM_WORLD, &ksp);
 51:   KSPSetOperators(ksp, A, A);
 52:   KSPSetFromOptions(ksp);
 53:   KSPSolve(ksp, b, x);

 55:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 56:                       Check solution and clean up
 57:      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
 58:   VecView(x, PETSC_VIEWER_STDOUT_WORLD);

 60:   /* Free work space. */
 61:   KSPDestroy(&ksp);
 62:   VecDestroy(&x);
 63:   VecDestroy(&b);
 64:   MatDestroy(&A);

 66:   PetscFinalize();
 67:   return 0;
 68: }

 70: /*TEST

 72:    test:
 73:       args: -ksp_type minres -pc_type none -ksp_converged_reason

 75:    test:
 76:       suffix: 2
 77:       nsize: 3
 78:       args: -ksp_type minres -pc_type none -ksp_converged_reason

 80: TEST*/