Actual source code: ex47.c

  1: /*
  2:     Tests attaching null space to IS for fieldsplit preconditioner
  3: */
  4: #include <petscksp.h>

  6: int main(int argc, char **argv)
  7: {
  8:   Mat          A;
  9:   KSP          ksp;
 10:   PC           pc;
 11:   IS           zero, one;
 12:   MatNullSpace nullsp;
 13:   Vec          x, b;
 14:   MPI_Comm     comm;

 17:   PetscInitialize(&argc, &argv, NULL, NULL);
 18:   comm = PETSC_COMM_WORLD;
 19:   MatCreate(comm, &A);
 20:   MatSetSizes(A, 4, 4, PETSC_DECIDE, PETSC_DECIDE);
 21:   MatSetUp(A);
 22:   MatSetFromOptions(A);
 23:   MatCreateVecs(A, &x, &b);
 24:   VecSet(x, 2.0);
 25:   VecSet(b, 12.0);
 26:   MatDiagonalSet(A, x, INSERT_VALUES);
 27:   MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
 28:   MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);
 29:   ISCreateStride(comm, 2, 0, 1, &zero);
 30:   ISCreateStride(comm, 2, 2, 1, &one);
 31:   MatNullSpaceCreate(comm, PETSC_TRUE, 0, NULL, &nullsp);
 32:   PetscObjectCompose((PetscObject)zero, "nullspace", (PetscObject)nullsp);
 33:   KSPCreate(comm, &ksp);
 34:   KSPSetOperators(ksp, A, A);
 35:   KSPSetUp(ksp);
 36:   KSPGetPC(ksp, &pc);
 37:   KSPSetFromOptions(ksp);
 38:   PCFieldSplitSetIS(pc, "0", zero);
 39:   PCFieldSplitSetIS(pc, "1", one);
 40:   KSPSolve(ksp, b, x);
 41:   KSPDestroy(&ksp);
 42:   MatNullSpaceDestroy(&nullsp);
 43:   ISDestroy(&zero);
 44:   ISDestroy(&one);
 45:   MatDestroy(&A);
 46:   VecDestroy(&x);
 47:   VecDestroy(&b);

 49:   PetscFinalize();
 50:   return 0;
 51: }

 53: /*TEST

 55:    test:
 56:       args: -pc_type fieldsplit

 58: TEST*/