Actual source code: ex1f.F90

  1: !
  2: !  Tests VecScatterCreateToAll Fortran stub
  3:       program main
  4: #include <petsc/finclude/petscvec.h>
  5:       use petscvec
  6:       implicit none

  8:       PetscErrorCode ierr
  9:       PetscInt  nlocal, row
 10:       PetscScalar num
 11:       PetscMPIInt rank
 12:       Vec v1, v2
 13:       VecScatter toall

 15:       PetscCallA(PetscInitialize(ierr))
 16:       PetscCallMPIA(MPI_Comm_rank(PETSC_COMM_WORLD,rank,ierr))

 18:       nlocal = 1
 19:       PetscCallA(VecCreateMPI(PETSC_COMM_WORLD,nlocal,PETSC_DECIDE,v1,ierr))

 21:       row = rank
 22:       num = rank
 23:       PetscCallA(VecSetValue(v1,row,num,INSERT_VALUES,ierr))
 24:       PetscCallA(VecAssemblyBegin(v1,ierr))
 25:       PetscCallA(VecAssemblyEnd(v1,ierr))

 27:       PetscCallA(VecScatterCreateToAll(v1,toall,v2,ierr))

 29:       PetscCallA(VecScatterBegin(toall,v1,v2,INSERT_VALUES,SCATTER_FORWARD,ierr))
 30:       PetscCallA(VecScatterEnd(toall,v1,v2,INSERT_VALUES,SCATTER_FORWARD,ierr))

 32:       PetscCallA(VecScatterDestroy(toall,ierr))
 33: ! Destroy v2 and then re-create it in VecScatterCreateToAll() to test if petsc can differentiate NULL projects with destroyed objects
 34:       PetscCallA(VecDestroy(v2,ierr))

 36:       PetscCallA(VecScatterCreateToAll(v1,toall,v2,ierr))
 37:       PetscCallA(VecScatterBegin(toall,v1,v2,INSERT_VALUES,SCATTER_FORWARD,ierr))
 38:       PetscCallA(VecScatterEnd(toall,v1,v2,INSERT_VALUES,SCATTER_FORWARD,ierr))

 40:       if (rank.eq.2) then
 41:          PetscCallA(PetscObjectSetName(v2, 'v2',ierr))
 42:          PetscCallA(VecView(v2,PETSC_VIEWER_STDOUT_SELF,ierr))
 43:       end if

 45:       PetscCallA(VecScatterDestroy(toall,ierr))
 46:       PetscCallA(VecDestroy(v1,ierr))
 47:       PetscCallA(VecDestroy(v2,ierr))
 48: ! It is OK to destroy again
 49:       PetscCallA(VecDestroy(v2,ierr))

 51:       PetscCallA(PetscFinalize(ierr))
 52:       end

 54: !/*TEST
 55: !
 56: !     test:
 57: !       nsize: 4
 58: !
 59: !TEST*/