Actual source code: ex97f90.F90

  1: program ex97f90
  2: #include "petsc/finclude/petsc.h"
  3:     use petsc
  4:     implicit none

  6:     ! Get the fortran kind associated with PetscInt and PetscReal so that we can use literal constants.
  7:     PetscInt                           :: dummyPetscInt
  8:     PetscReal                          :: dummyPetscreal
  9:     integer,parameter                  :: kPI = kind(dummyPetscInt)
 10:     integer,parameter                  :: kPR = kind(dummyPetscReal)

 12:     type(tDM)                          :: dm
 13:     type(tDMLabel)                     :: label
 14:     character(len=PETSC_MAX_PATH_LEN)  :: ifilename,iobuffer
 15:     DMPolytopeType                     :: cellType
 16:     PetscInt                           :: pStart,pEnd,p
 17:     PetscErrorCode                     :: ierr
 18:     PetscBool                          :: flg

 20:     PetscCallA(PetscInitialize(ierr))

 22:     PetscCallA(PetscOptionsGetString(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,"-i",ifilename,flg,ierr))
 23:     if (.not. flg) then
 24:         SETERRA(PETSC_COMM_WORLD,PETSC_ERR_ARG_OUTOFRANGE,"missing input file name -i <input file name>")
 25:     end if

 27:     PetscCallA(DMPlexCreateFromFile(PETSC_COMM_WORLD,ifilename,PETSC_NULL_CHARACTER,PETSC_TRUE,dm,ierr))
 28:     PetscCallA(DMPlexDistributeSetDefault(dm,PETSC_FALSE,ierr))
 29:     PetscCallA(PetscObjectSetName(dm,"ex97f90",ierr))
 30:     PetscCallA(DMSetFromOptions(dm,ierr))
 31:     PetscCallA(DMViewFromOptions(dm,PETSC_NULL_OPTIONS,"-dm_view",ierr))

 33:     PetscCallA(DMGetLabel(dm,'celltype',label,ierr))
 34:     PetscCallA(DMLabelView(label,PETSC_VIEWER_STDOUT_WORLD,ierr))
 35:     PetscCallA(DMPlexGetHeightStratum(dm,0_kPI,pStart,pEnd,ierr))
 36:     Do p = pStart,pEnd-1
 37:         PetscCallA(DMPlexGetCellType(dm,p,cellType,ierr))
 38:         Write(IOBuffer,'("cell: ",i3," type: ",i3,"\n")' ) p,cellType
 39:         PetscCallA(PetscPrintf(PETSC_COMM_SELF,IOBuffer,ierr))
 40:     End Do
 41:     PetscCallA(DMDestroy(dm,ierr))

 43:     PetscCallA(PetscFinalize(ierr))
 44: end program ex97f90

 46: ! /*TEST
 47: !   build:
 48: !     requires: !complex
 49: !   testset:
 50: !     args: -i ${wPETSC_DIR}/share/petsc/datafiles/meshes/doublet-tet.msh -dm_view
 51: !     nsize: 1
 52: !     test:
 53: !       suffix: 0
 54: !       args:
 55: ! TEST*/