1: !
  2: !  Tests MatGetNullSpace() and PetscObjectIsNull()
  3: !
  4: program main

  6: #include <petsc/finclude/petscmat.h>
  7:   use petscmat
  8:   implicit none

 10:   PetscErrorCode ierr
 11:   Mat A, B
 12:   Mat C, SC
 13:   MatNullSpace sp, sp1
 14:   PetscInt one, zero, rend
 15:   PetscScalar sone
 16:   Vec x, y

 18:   zero = 0
 19:   one = 1
 20:   sone = 1
 21:   PetscCallA(PetscInitialize(ierr))

 23:   PetscCallA(MatCreate(PETSC_COMM_WORLD, A, ierr))
 24:   PetscCallA(MatCreate(PETSC_COMM_WORLD, B, ierr))

 26:   PetscCallA(MatGetNullSpace(A, sp, ierr))
 27:   PetscCheckA(PetscObjectIsNull(sp), PETSC_COMM_SELF, PETSC_ERR_PLIB, 'Matrix null space should not exist')

 29:   PetscCallA(MatSetNullSpace(A, PETSC_NULL_MATNULLSPACE, ierr))
 30:   PetscCallA(MatGetNullSpace(A, sp, ierr))
 31:   PetscCheckA(PetscObjectIsNull(sp), PETSC_COMM_SELF, PETSC_ERR_PLIB, 'Matrix null space should not exist')

 33:   PetscCallA(MatNullSpaceCreate(PETSC_COMM_WORLD, PETSC_TRUE, zero, PETSC_NULL_VEC_ARRAY, sp, ierr))
 34:   PetscCallA(MatSetNullSpace(A, sp, ierr))
 35:   PetscCallA(MatGetNullSpace(A, sp1, ierr))
 36:   PetscCheckA(.not. PetscObjectIsNull(sp), PETSC_COMM_SELF, PETSC_ERR_PLIB, 'Matrix null space should not exist')
 37:   PetscCallA(MatNullSpaceDestroy(sp, ierr))

 39:   PetscCallA(MatCreateSeqDense(PETSC_COMM_WORLD, one, one, PETSC_NULL_SCALAR_ARRAY, C, ierr))
 40:   PetscCallA(MatSetValues(C, one, [zero], one, [zero], [sone], INSERT_VALUES, ierr))
 41:   PetscCallA(MatAssemblyBegin(C, MAT_FINAL_ASSEMBLY, ierr))
 42:   PetscCallA(MatAssemblyEnd(C, MAT_FINAL_ASSEMBLY, ierr))
 43:   PetscCallA(MatCreateSchurComplement(C, C, C, C, PETSC_NULL_MAT, SC, ierr))
 44:   PetscCallA(MatGetOwnershipRange(SC, PETSC_NULL_INTEGER, rend, ierr))
 45:   PetscCallA(VecCreateSeq(PETSC_COMM_SELF, one, x, ierr))
 46:   PetscCallA(VecDuplicate(x, y, ierr))
 47:   PetscCallA(VecSetValues(x, one, [zero], [sone], INSERT_VALUES, ierr))
 48:   PetscCallA(VecAssemblyBegin(x, ierr))
 49:   PetscCallA(VecAssemblyEnd(x, ierr))
 50:   PetscCallA(MatMult(SC, x, y, ierr))
 51:   PetscCallA(VecView(y, PETSC_VIEWER_STDOUT_SELF, ierr))
 52:   PetscCallA(VecSetRandom(x, PETSC_NULL_RANDOM, ierr))
 53:   PetscCallA(VecView(x, PETSC_VIEWER_STDOUT_SELF, ierr))

 55:   PetscCallA(MatDestroy(SC, ierr))
 56:   PetscCallA(MatDestroy(C, ierr))
 57:   PetscCallA(VecDestroy(x, ierr))
 58:   PetscCallA(VecDestroy(y, ierr))
 59:   PetscCallA(MatDestroy(A, ierr))
 60:   PetscCallA(MatDestroy(B, ierr))
 61:   PetscCallA(PetscFinalize(ierr))
 62: end

 64: !/*TEST
 65: !
 66: !   test:
 67: !      requires: !complex
 68: !
 69: !TEST*/