1: ! Test MatCreateNest() with NULL index sets
  2: #include <petsc/finclude/petscmat.h>
  3: program main
  4:   use petscmat
  5:   implicit none

  7:   Mat                   :: A, D, Id, Acopy
  8:   Mat, dimension(4)      :: mats
  9:   Vec                   :: v, w
 10:   PetscInt              :: i, rstart, rend
 11:   PetscInt, parameter    :: n = 6, nb = 2
 12:   PetscScalar, parameter :: one = 1.0, two = 2.0
 13:   PetscErrorCode        :: ierr

 15:   PetscCallA(PetscInitialize(ierr))

 17:   PetscCallA(MatCreateConstantDiagonal(PETSC_COMM_WORLD, PETSC_DETERMINE, PETSC_DETERMINE, n, n, one, Id, ierr))
 18:   PetscCallA(MatCreateVecs(Id, v, w, ierr))
 19:   PetscCallA(VecGetOwnershipRange(v, rstart, rend, ierr))
 20:   do i = rstart, rend - 1
 21:     PetscCallA(VecSetValue(v, i, two/(i + 1), INSERT_VALUES, ierr))
 22:   end do
 23:   PetscCallA(VecAssemblyBegin(v, ierr))
 24:   PetscCallA(VecAssemblyEnd(v, ierr))
 25:   PetscCallA(MatCreateDiagonal(v, D, ierr))

 27:   mats = [PETSC_NULL_MAT, D, Id, PETSC_NULL_MAT]
 28:   PetscCallA(MatCreateNest(PETSC_COMM_WORLD, nb, PETSC_NULL_IS_ARRAY, nb, PETSC_NULL_IS_ARRAY, mats, A, ierr))
 29:   PetscCallA(MatView(A, PETSC_VIEWER_STDOUT_WORLD, ierr))
 30:   PetscCallA(MatNestSetSubMats(A, nb, PETSC_NULL_IS_ARRAY, nb, PETSC_NULL_IS_ARRAY, mats, ierr))
 31:   PetscCallA(MatView(A, PETSC_VIEWER_STDOUT_WORLD, ierr))

 33:   ! test MatCopy()
 34:   PetscCallA(MatDuplicate(A, MAT_DO_NOT_COPY_VALUES, Acopy, ierr))
 35:   PetscCallA(MatCopy(A, Acopy, DIFFERENT_NONZERO_PATTERN, ierr))

 37:   PetscCallA(MatDestroy(Acopy, ierr))
 38:   PetscCallA(MatDestroy(Id, ierr))
 39:   PetscCallA(VecDestroy(v, ierr))
 40:   PetscCallA(VecDestroy(w, ierr))
 41:   PetscCallA(MatDestroy(D, ierr))
 42:   PetscCallA(MatDestroy(A, ierr))
 43:   PetscCallA(PetscFinalize(ierr))

 45: end program

 47: !/*TEST
 48: !
 49: !   test:
 50: !      nsize: 2
 51: !
 52: !TEST*/