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*/