Actual source code: ex1f.F90

petsc-3.13.4 2020-08-01
Report Typos and Errors
  1: !
  2: ! Test the workaround for a bug in OpenMPI-2.1.1 on Ubuntu 18.04.2
  3: ! See https://lists.mcs.anl.gov/pipermail/petsc-dev/2019-July/024803.html
  4: !
  5: ! Contributed-by:         Fabian Jakub  <Fabian.Jakub@physik.uni-muenchen.de>
  6: program main
  7: #include "petsc/finclude/petsc.h"

  9:   use petsc
 10:   implicit none

 12:   PetscInt, parameter :: Ndof=1, stencil_size=1
 13:   PetscInt, parameter :: Nx=3, Ny=3
 14:   PetscErrorCode :: myid, commsize, ierr
 15:   PetscScalar, pointer :: xv1d(:)

 17:   type(tDM) :: da
 18:   type(tVec) :: gVec!, naturalVec


 21:   call PetscInitialize(PETSC_NULL_CHARACTER, ierr)
 22:   call mpi_comm_rank(PETSC_COMM_WORLD, myid, ierr)
 23:   call mpi_comm_size(PETSC_COMM_WORLD, commsize, ierr)

 25:   call DMDACreate2d(PETSC_COMM_WORLD, &
 26:     DM_BOUNDARY_PERIODIC, DM_BOUNDARY_PERIODIC, &
 27:     DMDA_STENCIL_STAR, &
 28:     Nx, Ny, PETSC_DECIDE, PETSC_DECIDE, Ndof, stencil_size, &
 29:     PETSC_NULL_INTEGER, PETSC_NULL_INTEGER, da, ierr)
 30:   call DMSetup(da, ierr)
 31:   call DMSetFromOptions(da, ierr)

 33:   call DMCreateGlobalVector(da, gVec, ierr)
 34:   call VecGetArrayF90(gVec, xv1d, ierr)
 35:   xv1d(:) = real(myid, kind(xv1d))
 36:   !print *,myid, 'xv1d', xv1d, ':', xv1d
 37:   call VecRestoreArrayF90(gVec, xv1d, ierr)

 39:   call PetscObjectViewFromOptions(gVec, PETSC_NULL_VEC, "-show_gVec", ierr)

 41:   call VecDestroy(gVec, ierr)
 42:   call DMDestroy(da, ierr)
 43:   call PetscFinalize(ierr)
 44: end program

 46: !/*TEST
 47: !
 48: !   test:
 49: !      nsize: 9
 50: !      args: -show_gVec
 51: !TEST*/