Attaching the modified version that works with complex [using complex build of
petsc]
BTW: I'm using the following gfortran [but I don't think the version makes any
difference] with petsc-3.8 (maint branch):
$ gfortran --version
GNU Fortran (GCC) 8.0.1 20180222 (Red Hat 8.0.1-0.16)
Satish
--------
balay@asterix /home/balay/download-pine
$ make ex2f
mpif90 -c -Wall -ffree-line-length-0 -Wno-unused-dummy-argument -g
-I/home/balay/petsc/include -I/home/balay/petsc/arch-cmplx/include -o ex2f.o
ex2f.F90
mpif90 -Wall -ffree-line-length-0 -Wno-unused-dummy-argument -g -o ex2f
ex2f.o -Wl,-rpath,/home/balay/petsc/arch-cmplx/lib
-L/home/balay/petsc/arch-cmplx/lib -Wl,-rpath,/home/balay/soft/mpich-3.3b1/lib
-L/home/balay/soft/mpich-3.3b1/lib
-Wl,-rpath,/usr/lib/gcc/x86_64-redhat-linux/8
-L/usr/lib/gcc/x86_64-redhat-linux/8 -lpetsc -llapack -lblas -lX11 -lpthread
-lm -lmpifort -lgfortran -lm -lgfortran -lm -lquadmath -lmpicxx -lstdc++ -lm
-Wl,-rpath,/home/balay/soft/mpich-3.3b1/lib -L/home/balay/soft/mpich-3.3b1/lib
-Wl,-rpath,/usr/lib/gcc/x86_64-redhat-linux/8
-L/usr/lib/gcc/x86_64-redhat-linux/8 -ldl
-Wl,-rpath,/home/balay/soft/mpich-3.3b1/lib -lmpi -lgcc_s -ldl
/usr/bin/rm -f ex2f.o
balay@asterix /home/balay/download-pine
$ ./ex2f
setting values
assembling vector
copying vector to array
N = 10
m = 10
getting vecX size
getting ownership range
ownership range: 0 10
copy values to global array
vecX_pt:
SIZE = 10
(1.0000000000000000,0.0000000000000000)
(2.0000000000000000,0.0000000000000000)
(3.0000000000000000,0.0000000000000000)
(4.0000000000000000,0.0000000000000000)
(5.0000000000000000,0.0000000000000000)
(6.0000000000000000,0.0000000000000000)
(7.0000000000000000,0.0000000000000000)
(8.0000000000000000,0.0000000000000000)
(9.0000000000000000,0.0000000000000000)
(10.000000000000000,0.0000000000000000)
copy_buffer =
(1.0000000000000000,0.0000000000000000)
(2.0000000000000000,0.0000000000000000)
(3.0000000000000000,0.0000000000000000)
(4.0000000000000000,0.0000000000000000)
(5.0000000000000000,0.0000000000000000)
(6.0000000000000000,0.0000000000000000)
(7.0000000000000000,0.0000000000000000)
(8.0000000000000000,0.0000000000000000)
(9.0000000000000000,0.0000000000000000)
(10.000000000000000,0.0000000000000000)
calling allreduce
buffer =
(1.0000000000000000,0.0000000000000000)
(2.0000000000000000,0.0000000000000000)
(3.0000000000000000,0.0000000000000000)
(4.0000000000000000,0.0000000000000000)
(5.0000000000000000,0.0000000000000000)
(6.0000000000000000,0.0000000000000000)
(7.0000000000000000,0.0000000000000000)
(8.0000000000000000,0.0000000000000000)
(9.0000000000000000,0.0000000000000000)
(10.000000000000000,0.0000000000000000)
balay@asterix /home/balay/download-pine
$
On Fri, 23 Mar 2018, Satish Balay wrote:
> I'm attaching the modified code - that appears to work fine for me. My
> changes are:
>
> 1. use .F90 suffix - so the compiler does the f90 freeform/preprocessing
> automatically
> 2. remove slepc references and use petsc only
> 3. compile with a petsc formatted makefile
> 4. replace complex with real*8 - thats compatible with PetscScalar
>
> So perhaps you can modify your version of test code - to be closer to
> this one [or vice versa] to see what is triggering the difference in
> behavior.
>
> Satish
>
> On Fri, 23 Mar 2018, Samuel Lanthaler wrote:
>
> > Hi all,
> >
> > I am having trouble using the function VecGetArrayReadF90 under gfortran. I
> > have created a minimal example and put it in the attachment. Basically, it
> > appears that the input pointer is empty upon return from
> > VecGetArrayReadF90. The code runs fine, when compiled with ifort, so I
> > don't know what's going on.
> >
> > I am attaching the code and makefile that I used, as well as the command
> > line output that I can see. Though the errors are just to do with the fact
> > that the output pointer is empty, so the command line output is not all
> > that helpful...
> >
> > Thanks in advance for your help!
> > Sam
> >
>
PROGRAM allreduce
#include <petsc/finclude/petsc.h>
USE petsc
IMPLICIT NONE
! --- pure PETSc
PetscErrorCode :: ierr
PetscScalar, ALLOCATABLE :: vals(:)
PetscInt, ALLOCATABLE :: idxn(:)
PetscInt :: m,i
Vec :: vecX ! eigenvector (real/imaginary parts)
!
COMPLEX*16, ALLOCATABLE :: vecXcopy(:)
INTEGER :: N
! initialize PETSc
CALL PetscInitialize(PETSC_NULL_CHARACTER,ierr)
! Set up a new matrix
m = 10
ALLOCATE(vals(10),idxn(10))
!
CALL VecCreate(PETSC_COMM_WORLD,vecX,ierr);
CALL VecSetType(vecX,VECMPI,ierr);
CALL VecSetSizes(vecX,PETSC_DECIDE,m,ierr);
! set values of vector
DO i=1,m
vals(i) = i
idxn(i) = i-1
END DO
!
PRINT*,'setting values'
CALL VecSetValues(vecX,m,idxn,vals,INSERT_VALUES,ierr);
! assemble vector
PRINT*,'assembling vector'
CALL VecAssemblyBegin(vecX,ierr);
CALL VecAssemblyEnd(vecX,ierr);
!
PRINT*,'copying vector to array'
N = m
PRINT*,'N = ',N
PRINT*,'m = ',m
ALLOCATE(vecXcopy(N))
CALL PetscVec2Array(vecX,vecXcopy,N)
!
CALL VecDestroy(vecX,ierr);
CALL PetscFinalize(ierr);
END PROGRAM allreduce
!> Copy Petsc vector to fortran array
SUBROUTINE PetscVec2Array(vecX,buffer,N)
! ------
#include <petsc/finclude/petsc.h>
USE petsc
! ------
IMPLICIT NONE
Vec, INTENT(in) :: vecX
COMPLEX*16, INTENT(inout) :: buffer(0:N-1)
INTEGER, INTENT(in) :: N
!
INTEGER :: istart, iend, Np, i,j
INTEGER :: me2,nprocs
COMPLEX*16, ALLOCATABLE :: copy_buffer(:)
!
PetscScalar, POINTER :: vecX_pt(:) => NULL()
PetscInt :: petsc_N
PetscErrorCode :: ierr
!
CALL MPI_COMM_RANK(MPI_COMM_WORLD,me2,ierr);
CALL MPI_COMM_SIZE(MPI_COMM_WORLD,Nprocs,ierr);
! initailize buffer to 0.
buffer = 0.
ALLOCATE(copy_buffer(0:N-1))
! check length of petsc vector is N
PRINT*,'getting vecX size'
CALL VecGetSize(vecX,petsc_N,ierr); CHKERRQ(ierr)
Np = petsc_N
IF(N.NE.Np) THEN
PRINT*,'ERROR: In PetscVec2Array:'
PRINT*,'petsc_N = ',petsc_N, Np
PRINT*,'N = ',N
CALL ABORT
END IF
! vecX_pt will contain the element in the local range [istart,iend-1]
PRINT*,'getting ownership range'
CALL VecGetOwnershipRange(vecX,istart,iend,ierr); CHKERRQ(ierr)
PRINT*,'ownership range: ',istart,iend
! initialize the pointer
CALL VecGetArrayReadF90(vecX,vecX_pt,ierr); CHKERRQ(ierr)
! copy values to global array
PRINT*,'copy values to global array'
DO i=1,Nprocs
!
CALL MPI_BARRIER(MPI_COMM_WORLD,ierr)
!
IF(me2.EQ.i-1) THEN
PRINT*,'vecX_pt:'
PRINT*,'SIZE = ',SIZE(VecX_pt)
DO j=1,SIZE(vecX_pt)
PRINT*,' ',vecX_pt(j)
END DO
END IF
END DO
copy_buffer(istart:iend-1) = vecX_pt(:)
! free pointer
CALL VecRestoreArrayReadF90(vecX,vecX_pt,ierr); CHKERRQ(ierr)
CALL MPI_BARRIER(MPI_COMM_WORLD,ierr)
DO i=1,Nprocs
!
CALL MPI_BARRIER(MPI_COMM_WORLD,ierr)
!
IF(me2.EQ.i-1) THEN
PRINT*,'copy_buffer = '
DO j=istart,iend-1
PRINT*,' ',copy_buffer(j)
END DO
END IF
END DO
CALL MPI_BARRIER(MPI_COMM_WORLD,ierr)
! combine all components from all processes
PRINT*,'calling allreduce'
CALL
MPI_ALLREDUCE(copy_buffer,buffer,N,MPI_DOUBLE_COMPLEX,MPI_SUM,MPI_COMM_WORLD,ierr)
CALL MPI_BARRIER(MPI_COMM_WORLD,ierr)
DO i=1,Nprocs
!
CALL MPI_BARRIER(MPI_COMM_WORLD,ierr)
!
IF(me2.EQ.i-1) THEN
PRINT*,'buffer = '
DO j=0,N-1
PRINT*,' ',buffer(j)
END DO
END IF
END DO
! free local memory
DEALLOCATE(copy_buffer)
END SUBROUTINE PetscVec2Array
CFLAGS =
FFLAGS =
CPPFLAGS =
FPPFLAGS =
include ${PETSC_DIR}/lib/petsc/conf/variables
include ${PETSC_DIR}/lib/petsc/conf/rules
ex1f: ex1f.o chkopts
-${FLINKER} -o ex1f ex1f.o ${PETSC_KSP_LIB}
${RM} ex1f.o
ex2f: ex2f.o chkopts
-${FLINKER} -o ex2f ex2f.o ${PETSC_KSP_LIB}
${RM} ex2f.o
