Actual source code: ex26f.F90

  1: !
  2: !  Test VecGetSubVector()
  3: !  Contributed-by: Adrian Croucher <gitlab@mg.gitlab.com>

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

 10:       PetscMPIInt :: rank
 11:       PetscErrorCode :: ierr
 12:       PetscInt :: num_cells, subsize, i
 13:       PetscInt, parameter :: blocksize = 3, field = 0
 14:       Vec :: v, subv
 15:       IS :: index_set
 16:       PetscInt, allocatable :: subindices(:)

 18:       PetscCallA(PetscInitialize(ierr))
 19:       PetscCallMPIA(MPI_COMM_RANK(PETSC_COMM_WORLD, rank, ierr))

 21:       if (rank .eq. 0) then
 22:          num_cells = 1
 23:       else
 24:          num_cells = 0
 25:       end if

 27:       PetscCallA(VecCreate(PETSC_COMM_WORLD, v, ierr))
 28:       PetscCallA(VecSetSizes(v, num_cells * blocksize, PETSC_DECIDE, ierr))
 29:       PetscCallA(VecSetBlockSize(v, blocksize, ierr))
 30:       PetscCallA(VecSetFromOptions(v, ierr))

 32:       subsize = num_cells
 33:       allocate(subindices(0: subsize - 1))
 34:       subindices = [(i, i = 0, subsize - 1)] * blocksize + field
 35:       PetscCallA(ISCreateGeneral(PETSC_COMM_WORLD, subsize, subindices,PETSC_COPY_VALUES, index_set, ierr))
 36:       deallocate(subindices)

 38:       PetscCallA(VecGetSubVector(v, index_set, subv, ierr))
 39:       PetscCallA(VecRestoreSubVector(v, index_set, subv, ierr))
 40:       PetscCallA(ISDestroy(index_set, ierr))

 42:       PetscCallA(VecDestroy(v, ierr))
 43:       PetscCallA(PetscFinalize(ierr))
 44:       end

 46: !/*TEST
 47: !
 48: !   test:
 49: !      nsize: 2
 50: !      filter: sort -b
 51: !      filter_output: sort -b
 52: !
 53: !TEST*/