[petsc-users] about Vecgetarray()

Satish Balay balay at mcs.anl.gov
Sat Apr 17 10:00:11 CDT 2010

> On Fri, 16 Apr 2010, 陈乐平(Leping Chen) wrote:

> >        PetscOffset xx_i
> >        Vec         x
> >        double precision u(6) 
> >        
> >         call VecGetArray(x,u,xx_i,ierr)
> >         do 30 i=1,n
> >           u(i) = 1000.0*i

>           u(i+xx_i) = 1000.0*i

> >   30   continue
> >        call VecRestoreArray(x,u,xx_i,ierr)

Hopefully the following clears up the confusion.

You are expecting the following to work:
integer n=6
double precision u(6)
call VecGetArray(x,u,ierr)
  do 30 i=1,n
  u(i) = 1000.0*i
30   continue
call VecRestoreArray(x,u,ierr)

Here you are expecting VecGetArray to *copy* the values of 'x' into
'u'.  And VecRestoreArray() to *copy* values from 'u' to 'x'. However
these copies are inefficient - so we do not do that.

Here VecGetArray() tries to return the *pointer* to the array stored
in the Vec. This is possible in C - but not fortran77. Hence the concept
of using offset in the f77 interface.

>           u(i+xx_i) = 1000.0*i

Note: the above is out-of-bounds access of 'u' wrt F77 - and it is a
workarround against the limitation of F77 language. [so its not pure
F77 code].

So if you need pure language compliant code - check VecGetArrayF90()
as indicated before. It returns a F90 pointer from a Vec.


More information about the petsc-users mailing list