Program test

#include "petsc/finclude/petscdmda.h"
  use petscdmda
  use petscdmcomposite

  DM                  :: dmComp, dm1
  PetscErrorCode      :: ierr
  PetscInt            :: i1, i2, i20, ndm
  IS, pointer         :: isComp(:)

  PetscCallA(PetscInitialize(PETSC_NULL_CHARACTER,ierr))

  ! Initialize constants
  i1=1
  i2=2
  i20=20

  ! Create composite DM
  PetscCallA(DMCompositeCreate(PETSC_COMM_WORLD,dmComp,ierr))

  ! Create DM1
  PetscCallA(DMDACreate3d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,DMDA_STENCIL_BOX,i20,i20,i20,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE,i1,i2,PETSC_NULL_INTEGER_ARRAY,PETSC_NULL_INTEGER_ARRAY,PETSC_NULL_INTEGER_ARRAY,dm1,ierr))
  PetscCallA(DMSetUp(dm1,ierr))

  ! Add two DM1s to the composite DM
  PetscCallA(DMCompositeAddDM(dmComp,dm1,ierr))
  PetscCallA(DMCompositeAddDM(dmComp,dm1,ierr))

  ! Finish setting up composite DM
  PetscCallA(DMSetUp(dmComp,ierr))

  ! Get the global IS sets
  PetscCallA(DMCompositeGetNumberDM(dmComp, ndm, ierr))
  !allocate (isComp(ndm))
  PetscCallA(DMCompositeGetGlobalISs(dmComp,isComp,ierr))
      PetscCallA(DMCompositeRestoreGlobalISs(dmComp,isComp,ierr))
      PetscCallA(DMDestroy(dmComp,ierr))
            PetscCallA(DMDestroy(dm1,ierr))


  PetscCallA(PetscFinalize(ierr))


end program test