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;

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

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

 53: /*TEST

 55:    test:
 56:       args: -pc_type fieldsplit

 58: TEST*/