[petsc-users] DMCreateGlobalVector in fortran.
Barry Smith
bsmith at mcs.anl.gov
Sat Feb 4 19:46:00 CST 2017
For DMDAVecGetArrayF90 you need to declare the "array" arguments as Fortran pointers you don't declare them like
> u0(-1:IMax+2,-1:JMax+1,-1:KMax+1)
> On Feb 4, 2017, at 7:34 PM, Matthew Knepley <knepley at gmail.com> wrote:
>
> On Sat, Feb 4, 2017 at 7:30 PM, Manuel Valera <mvalera at mail.sdsu.edu> wrote:
> Thanks, that helped a lot,
>
> Now im a little confused with my next error, in my serial code i have several arrays with dimensions declared as:
>
> u0(-1:IMax+2,-1:JMax+1,-1:KMax+1)
>
> and i create my DMDA for this array as:
>
> call CDA3(daub,IMax+4,JMax+3,KMax+3)
> call DAs(daub,u0p,u0pl,u0)
>
> With the subroutines as:
>
> SUBROUTINE CDA3(da3d,dim3dx,dim3dy,dim3dz)
>
>
>
> DM, intent(inout) :: da3d
> PetscInt, intent(in) :: dim3dx,dim3dy,dim3dz
> DMBoundaryType :: bx,by,bz
>
>
> bx = DM_BOUNDARY_NONE
> by = DM_BOUNDARY_NONE
> bz = DM_BOUNDARY_NONE
>
> call DMDACreate3d(PETSC_COMM_WORLD,bx,by,bz,DMDA_STENCIL_BOX,dim3dx,dim3dy,dim3dz,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE,1,1,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,da3d,ierr)
>
> call DMSetFromOptions(da3d,ierr)
> call DMSetUp(da3d,ierr)
>
>
> END SUBROUTINE CDA3
>
> SUBROUTINE DAs(da,globalv,localv,array)
>
> ! use PetscObjects, only :: ierr
>
>
> ! Umbrella program to update and communicate the arrays in a
> ! distributed fashion using the DMDA objects from PETSc.
> ! Manuel Valera 1/20/17
>
> ! Arguments:
> ! da = DMDA array either 1d,2d or 3d, already created and setup
> ! globalv = global vector to be operated into
> ! localv = local chunk each processor works in.
> ! array = the array to be petscified. ONLY 3D ARRAYS rn.
> ! just in case.
>
> Vec,intent(inout) :: globalv,localv
>
>
> PetscScalar,dimension(:,:,:) :: array
> DM,intent(inout) :: da
>
>
>
>
>
>
> call DMCreateGlobalVector(da,globalv,ierr)
> call DMCreateLocalVector(da,localv,ierr)
>
> call DMGlobalToLocalBegin(da,globalv,INSERT_VALUES,localv,ierr)
> call DMGlobalToLocalEnd(da,globalv,INSERT_VALUES,localv,ierr)
>
>
>
> END SUBROUTINE DAs
>
> Which i later open into an array to operate into it, but when i try to access u0(:,j-2,:) it gives me the following error:
>
> Fortran runtime error: Index '-1' of dimension 2 of array 'u0' below lower bound of 0
>
> So i'm guessing the DMDA is creating the arrays with dimensions starting at 0 instead of -1 as it should, right? what can be done in this case?
>
> 1) The DMDA only takes sizes, not index bounds, in the constructor.
>
> 2) Its not the DMDA you must be talking about, but the array you get from DMDAVecGetArrayF90() or VecGetArrayF90()
>
> 3) The VecGetArrayF90() always start from 0 since it is indexed in local indices
>
> 4) The DMDVecGetArrayF90() arrays are indexed by global indices for the patch your process owns (including ghosts if its a local vector)
>
> Thanks,
>
> Matt
>
> Thanks,
> Manuel
>
> On Sat, Feb 4, 2017 at 4:18 PM, Barry Smith <bsmith at mcs.anl.gov> wrote:
>
> See: http://www.mcs.anl.gov/petsc/documentation/changes/dev.html search for DMDACreate and update your code.
>
> Barry
>
>
> > On Feb 4, 2017, at 6:02 PM, Manuel Valera <mvalera at mail.sdsu.edu> wrote:
> >
> > Hello devs,
> >
> > Ive been creating the framework to distribute the arrays from my model using the global/local vectors and the DMDAs, for this i created a small subroutine that is responsible for making the necessary calls:
> >
> > SUBROUTINE DAs(da,globalv,localv)
> >
> > ! use PetscObjects, only :: ierr
> >
> >
> > ! Umbrella program to update and communicate the arrays in a
> > ! distributed fashion using the DMDA objects from PETSc.
> > ! Manuel Valera 1/20/17
> >
> > ! Arguments:
> > ! da = DMDA array either 1d,2d or 3d, already created and setup
> > ! globalv = global vector to be operated into
> > ! localv = local chunk each processor works in.
> >
> >
> >
> > Vec,intent(inout) :: globalv,localv
> >
> >
> >
> > DM,intent(inout) :: da
> >
> > call DMCreateGlobalVector(da,globalv,ierr)
> > call DMCreateLocalVector(da,localv,ierr)
> >
> > call DMGlobalToLocalBegin(da,globalv,INSERT_VALUES,localv,ierr)
> > call DMGlobalToLocalEnd(da,globalv,INSERT_VALUES,localv,ierr)
> >
> >
> >
> > END SUBROUTINE DAs
> >
> >
> > I call this from another program where i have declared a vector as:
> >
> > #include "petsc/finclude/petscdm.h"
> > #include "petsc/finclude/petscdmda.h"
> >
> > use petsc
> > use petscdm
> > use petscdmda
> >
> > IMPLICIT NONE
> >
> > Vec :: xpdmda,yp,zp,xcp,ycp,zcp
> > DM :: da,dac,dau,dav,daw
> >
> > CONTAINS
> > SUBROUTINE InitGridP
> >
> > IMPLICIT NONE
> > !PetscInt :: ierr
> >
> > ! CREATE DMDA:
> >
> > call CDA3(da,IMax-1,JMax-1,KMax-1)
> >
> > !xpdmda,yp,zp have the same dimensions than x,y,z
> > !xpl,ypl,zpl are the local versions of xp,yp,zp
> >
> > !We create the global and local vectors and communicate their values:
> > call DAs(da,xpdmda,xpl,x)
> >
> > And at this call it crashes with the following message:
> >
> > [0]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------
> > [0]PETSC ERROR: Null argument, when expecting valid pointer
> > [0]PETSC ERROR: Null Object: Parameter # 2
> > [0]PETSC ERROR: See
> > http://www.mcs.anl.gov/petsc/documentation/faq.html
> > for trouble shooting.
> > [0]PETSC ERROR: Petsc Development GIT revision: v3.7.5-2949-gc3c607c GIT Date: 2017-01-20 17:34:16 -0600
> > [0]PETSC ERROR: ./ucmsLEP on a arch-linux2-c-debug named ocean by valera Sat Feb 4 15:53:53 2017
> > [0]PETSC ERROR: Configure options --with-cc=gcc --with-cxx=g++ --with-fc=gfortran --download-fblaslapack --download-mpich --download-hdf5 --download-netcdf --download-hypre --download-metis --download-parmetis --download-trillinos --with-debugging=1
> > [0]PETSC ERROR: #3 VecSetLocalToGlobalMapping() line 79 in /usr/dataC/home/valera/petsc/src/vec/vec/interface/vector.c
> > [0]PETSC ERROR: #4 DMCreateGlobalVector_DA() line 41 in /usr/dataC/home/valera/petsc/src/dm/impls/da/dadist.c
> > [0]PETSC ERROR: #5 DMCreateGlobalVector() line 840 in /usr/dataC/home/valera/petsc/src/dm/interface/dm.c
> >
> >
> >
> > So my question is, do i have to setup the vectors in any additional way? i followed the fortran DMDA examples closely and i didnt notice anything extra,
> >
> > Thanks,
> >
> > Manuel
>
>
>
>
>
> --
> What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.
> -- Norbert Wiener
More information about the petsc-users
mailing list