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

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

 52: /*TEST

 54:    test:
 55:       args: -pc_type fieldsplit

 57: TEST*/