diff --git a/include/petsc/finclude/ftn-custom/petscdmda.h90 b/include/petsc/finclude/ftn-custom/petscdmda.h90 index b16b866..45eed12 100644 --- a/include/petsc/finclude/ftn-custom/petscdmda.h90 +++ b/include/petsc/finclude/ftn-custom/petscdmda.h90 @@ -77,3 +77,65 @@ PetscErrorCode ierr End Subroutine End Interface DMDAVecRestoreArrayF90 + + Interface DMDAVecGetArrayReadF90 + Subroutine DMDAVecGetArrayReadF901(da1, v,d1,ierr) + USE_DM_HIDE + DM_HIDE da1 + VEC_HIDE v + PetscScalar,pointer :: d1(:) + PetscErrorCode ierr + End Subroutine + Subroutine DMDAVecGetArrayReadF902(da1, v,d1,ierr) + USE_DM_HIDE + DM_HIDE da1 + VEC_HIDE v + PetscScalar,pointer :: d1(:,:) + PetscErrorCode ierr + End Subroutine + Subroutine DMDAVecGetArrayReadF903(da1, v,d1,ierr) + USE_DM_HIDE + DM_HIDE da1 + VEC_HIDE v + PetscScalar,pointer :: d1(:,:,:) + PetscErrorCode ierr + End Subroutine + Subroutine DMDAVecGetArrayReadF904(da1, v,d1,ierr) + USE_DM_HIDE + DM_HIDE da1 + VEC_HIDE v + PetscScalar,pointer :: d1(:,:,:,:) + PetscErrorCode ierr + End Subroutine + End Interface DMDAVecGetArrayReadF90 + + Interface DMDAVecRestoreArrayReadF90 + Subroutine DMDAVecRestoreArrayReadF901(da1, v,d1,ierr) + USE_DM_HIDE + DM_HIDE da1 + VEC_HIDE v + PetscScalar,pointer :: d1(:) + PetscErrorCode ierr + End Subroutine + Subroutine DMDAVecRestoreArrayReadF902(da1, v,d1,ierr) + USE_DM_HIDE + DM_HIDE da1 + VEC_HIDE v + PetscScalar,pointer :: d1(:,:) + PetscErrorCode ierr + End Subroutine + Subroutine DMDAVecRestoreArrayReadF903(da1, v,d1,ierr) + USE_DM_HIDE + DM_HIDE da1 + VEC_HIDE v + PetscScalar,pointer :: d1(:,:,:) + PetscErrorCode ierr + End Subroutine + Subroutine DMDAVecRestoreArrayReadF904(da1, v,d1,ierr) + USE_DM_HIDE + DM_HIDE da1 + VEC_HIDE v + PetscScalar,pointer :: d1(:,:,:,:) + PetscErrorCode ierr + End Subroutine + End Interface DMDAVecRestoreArrayReadF90 diff --git a/src/dm/impls/da/dagetarray.c b/src/dm/impls/da/dagetarray.c index 392d364..1c47886 100644 --- a/src/dm/impls/da/dagetarray.c +++ b/src/dm/impls/da/dagetarray.c @@ -274,13 +274,13 @@ PetscErrorCode DMDAVecRestoreArrayDOF(DM da,Vec vec,void *array) appropriate DMGlobalToLocalBegin() and DMGlobalToLocalEnd() to have correct values in the ghost locations. - Fortran Notes: From Fortran use DMDAVecGetArrayF90() and pass for the array type PetscScalar,pointer :: array(:,...,:) of the appropriate + Fortran Notes: From Fortran use DMDAVecGetArrayReadF90() and pass for the array type PetscScalar,pointer :: array(:,...,:) of the appropriate dimension. For a DMDA created with a dof of 1 use the dimension of the DMDA, for a DMDA created with a dof greater than 1 use one more than the dimension of the DMDA. The order of the indices is array(xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1) (when dof is 1) otherwise array(0:dof-1,xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1) where the values are obtained from DMDAGetCorners() for a global array or DMDAGetGhostCorners() for a local array. Include petsc/finclude/petscdmda.h90 to access this routine. - Due to bugs in the compiler DMDAVecGetArrayF90() does not work with gfortran versions before 4.5 + Due to bugs in the compiler DMDAVecGetArrayReadF90() does not work with gfortran versions before 4.5 Level: intermediate @@ -342,7 +342,7 @@ PetscErrorCode DMDAVecGetArrayRead(DM da,Vec vec,void *array) Level: intermediate - Fortran Notes: From Fortran use DMDAVecRestoreArrayF90() + Fortran Notes: From Fortran use DMDAVecRestoreArrayReadF90() .keywords: distributed array, get, corners, nodes, local indices, coordinates @@ -407,7 +407,7 @@ PetscErrorCode DMDAVecRestoreArrayRead(DM da,Vec vec,void *array) In C, the indexing is "backwards" from what expects: array[k][j][i][DOF] NOT array[i][j][k][DOF]! - In Fortran 90 you do not need a version of DMDAVecRestoreArrayDOF() just use DMDAVecRestoreArrayF90() and declare your array with one higher dimension, + In Fortran 90 you do not need a version of DMDAVecRestoreArrayDOF() just use DMDAVecRestoreArrayReadF90() and declare your array with one higher dimension, see src/dm/examples/tutorials/ex11f90.F Level: intermediate diff --git a/src/dm/impls/da/f90-custom/zda1f90.c b/src/dm/impls/da/f90-custom/zda1f90.c index 6fb5b8e..fd9d3d5 100644 --- a/src/dm/impls/da/f90-custom/zda1f90.c +++ b/src/dm/impls/da/f90-custom/zda1f90.c @@ -11,6 +11,14 @@ #define dmdavecrestorearrayf903_ DMDAVECRESTOREARRAYF903 #define dmdavecgetarrayf904_ DMDAVECGETARRAYF904 #define dmdavecrestorearrayf904_ DMDAVECRESTOREARRAYF904 +#define dmdavecgetarrayreadf901_ DMDAVECGETARRAYREADF901 +#define dmdavecrestorearrayreadf901_ DMDAVECRESTOREARRAYREADF901 +#define dmdavecgetarrayreadf902_ DMDAVECGETARRAYREADF902 +#define dmdavecrestorearrayreadf902_ DMDAVECRESTOREARRAYREADF902 +#define dmdavecgetarrayreadf903_ DMDAVECGETARRAYREADF903 +#define dmdavecrestorearrayreadf903_ DMDAVECRESTOREARRAYREADF903 +#define dmdavecgetarrayreadf904_ DMDAVECGETARRAYREADF904 +#define dmdavecrestorearrayreadf904_ DMDAVECRESTOREARRAYREADF904 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE) #define dmdagetlocalinfof90_ dmdagetlocalinfof90 #define dmdavecgetarrayf901_ dmdavecgetarrayf901 @@ -21,6 +29,14 @@ #define dmdavecrestorearrayf903_ dmdavecrestorearrayf903 #define dmdavecgetarrayf904_ dmdavecgetarrayf904 #define dmdavecrestorearrayf904_ dmdavecrestorearrayf904 +#define dmdavecgetarrayreadf901_ dmdavecgetarrayreadf901 +#define dmdavecrestorearrayreadf901_ dmdavecrestorearrayreadf901 +#define dmdavecgetarrayreadf902_ dmdavecgetarrayreadf902 +#define dmdavecrestorearrayreadf902_ dmdavecrestorearrayreadf902 +#define dmdavecgetarrayreadf903_ dmdavecgetarrayreadf903 +#define dmdavecrestorearrayreadf903_ dmdavecrestorearrayreadf903 +#define dmdavecgetarrayreadf904_ dmdavecgetarrayreadf904 +#define dmdavecrestorearrayreadf904_ dmdavecrestorearrayreadf904 #endif PETSC_EXTERN void PETSC_STDCALL dmdagetlocalinfof90_(DM *da,DMDALocalInfo *info,PetscErrorCode *ierr) @@ -177,3 +193,151 @@ PETSC_EXTERN void PETSC_STDCALL dmdavecrestorearrayf904_(DM *da,Vec *v,F90Array4 *ierr = F90Array4dDestroy(&a,PETSC_SCALAR PETSC_F90_2PTR_PARAM(ptrd)); } +PETSC_EXTERN void PETSC_STDCALL dmdavecgetarrayreadf901_(DM *da,Vec *v,F90Array1d *a,PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) +{ + PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof; + const PetscScalar *aa; + + *ierr = DMDAGetCorners(*da,&xs,&ys,&zs,&xm,&ym,&zm);if (*ierr) return; + *ierr = DMDAGetGhostCorners(*da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);if (*ierr) return; + *ierr = DMDAGetInfo(*da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);if (*ierr) return; + + /* Handle case where user passes in global vector as opposed to local */ + *ierr = VecGetLocalSize(*v,&N);if (*ierr) return; + if (N == xm*ym*zm*dof) { + gxm = xm; + gym = ym; + gzm = zm; + gxs = xs; + gys = ys; + gzs = zs; + } else if (N != gxm*gym*gzm*dof) { + *ierr = PETSC_ERR_ARG_INCOMP; + } + *ierr = VecGetArrayRead(*v,&aa);if (*ierr) return; + *ierr = F90Array1dCreate((void*)aa,PETSC_SCALAR,gxs,gxm,a PETSC_F90_2PTR_PARAM(ptrd));if (*ierr) return; +} + +PETSC_EXTERN void PETSC_STDCALL dmdavecrestorearrayreadf901_(DM *da,Vec *v,F90Array1d *a,PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) +{ + const PetscScalar *fa; + *ierr = F90Array1dAccess(a,PETSC_SCALAR,(void**)&fa PETSC_F90_2PTR_PARAM(ptrd)); + *ierr = VecRestoreArrayRead(*v,&fa);if (*ierr) return; + *ierr = F90Array1dDestroy(&a,PETSC_SCALAR PETSC_F90_2PTR_PARAM(ptrd)); +} + +PETSC_EXTERN void PETSC_STDCALL dmdavecgetarrayreadf902_(DM *da,Vec *v,F90Array2d *a,PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) +{ + PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof; + const PetscScalar *aa; + + *ierr = DMDAGetCorners(*da,&xs,&ys,&zs,&xm,&ym,&zm);if (*ierr) return; + *ierr = DMDAGetGhostCorners(*da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);if (*ierr) return; + *ierr = DMDAGetInfo(*da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);if (*ierr) return; + + /* Handle case where user passes in global vector as opposed to local */ + *ierr = VecGetLocalSize(*v,&N);if (*ierr) return; + if (N == xm*ym*zm*dof) { + gxm = xm; + gym = ym; + gzm = zm; + gxs = xs; + gys = ys; + gzs = zs; + } else if (N != gxm*gym*gzm*dof) { + *ierr = PETSC_ERR_ARG_INCOMP; + } + if (dim == 1) { + gys = gxs; + gym = gxm; + gxs = 0; + gxm = dof; + } + *ierr = VecGetArrayRead(*v,&aa);if (*ierr) return; + *ierr = F90Array2dCreate((void*)aa,PETSC_SCALAR,gxs,gxm,gys,gym,a PETSC_F90_2PTR_PARAM(ptrd));if (*ierr) return; +} + +PETSC_EXTERN void PETSC_STDCALL dmdavecrestorearrayreadf902_(DM *da,Vec *v,F90Array2d *a,PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) +{ + const PetscScalar *fa; + *ierr = F90Array2dAccess(a,PETSC_SCALAR,(void**)&fa PETSC_F90_2PTR_PARAM(ptrd)); + *ierr = VecRestoreArrayRead(*v,&fa);if (*ierr) return; + *ierr = F90Array2dDestroy(&a,PETSC_SCALAR PETSC_F90_2PTR_PARAM(ptrd)); +} + +PETSC_EXTERN void PETSC_STDCALL dmdavecgetarrayreadf903_(DM *da,Vec *v,F90Array3d *a,PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) +{ + PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof; + const PetscScalar *aa; + + *ierr = DMDAGetCorners(*da,&xs,&ys,&zs,&xm,&ym,&zm);if (*ierr) return; + *ierr = DMDAGetGhostCorners(*da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);if (*ierr) return; + *ierr = DMDAGetInfo(*da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);if (*ierr) return; + + /* Handle case where user passes in global vector as opposed to local */ + *ierr = VecGetLocalSize(*v,&N);if (*ierr) return; + if (N == xm*ym*zm*dof) { + gxm = xm; + gym = ym; + gzm = zm; + gxs = xs; + gys = ys; + gzs = zs; + } else if (N != gxm*gym*gzm*dof) { + *ierr = PETSC_ERR_ARG_INCOMP; + } + if (dim == 2) { + gzs = gys; + gzm = gym; + gys = gxs; + gym = gxm; + gxs = 0; + gxm = dof; + } + *ierr = VecGetArrayRead(*v,&aa);if (*ierr) return; + *ierr = F90Array3dCreate((void*)aa,PETSC_SCALAR,gxs,gxm,gys,gym,gzs,gzm,a PETSC_F90_2PTR_PARAM(ptrd));if (*ierr) return; +} + +PETSC_EXTERN void PETSC_STDCALL dmdavecrestorearrayreadf903_(DM *da,Vec *v,F90Array3d *a,PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) +{ + const PetscScalar *fa; + *ierr = F90Array3dAccess(a,PETSC_SCALAR,(void**)&fa PETSC_F90_2PTR_PARAM(ptrd)); + *ierr = VecRestoreArrayRead(*v,&fa);if (*ierr) return; + *ierr = F90Array3dDestroy(&a,PETSC_SCALAR PETSC_F90_2PTR_PARAM(ptrd)); +} + +PETSC_EXTERN void PETSC_STDCALL dmdavecgetarrayreadf904_(DM *da,Vec *v,F90Array4d *a,PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) +{ + PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof,zero = 0; + const PetscScalar *aa; + + *ierr = DMDAGetCorners(*da,&xs,&ys,&zs,&xm,&ym,&zm);if (*ierr) return; + *ierr = DMDAGetGhostCorners(*da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);if (*ierr) return; + *ierr = DMDAGetInfo(*da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);if (*ierr) return; + + /* Handle case where user passes in global vector as opposed to local */ + *ierr = VecGetLocalSize(*v,&N);if (*ierr) return; + if (N == xm*ym*zm*dof) { + gxm = xm; + gym = ym; + gzm = zm; + gxs = xs; + gys = ys; + gzs = zs; + } else if (N != gxm*gym*gzm*dof) { + *ierr = PETSC_ERR_ARG_INCOMP; + } + *ierr = VecGetArrayRead(*v,&aa);if (*ierr) return; + *ierr = F90Array4dCreate((void*)aa,PETSC_SCALAR,zero,dof,gxs,gxm,gys,gym,gzs,gzm,a PETSC_F90_2PTR_PARAM(ptrd));if (*ierr) return; +} + +PETSC_EXTERN void PETSC_STDCALL dmdavecrestorearrayreadf904_(DM *da,Vec *v,F90Array4d *a,PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) +{ + const PetscScalar *fa; + /* + F90Array4dAccess is not implemented, so the following call would fail + */ + *ierr = F90Array4dAccess(a,PETSC_SCALAR,(void**)&fa PETSC_F90_2PTR_PARAM(ptrd)); + *ierr = VecRestoreArrayRead(*v,&fa);if (*ierr) return; + *ierr = F90Array4dDestroy(&a,PETSC_SCALAR PETSC_F90_2PTR_PARAM(ptrd)); +} diff --git a/src/docs/tex/manual/part2.tex b/src/docs/tex/manual/part2.tex index 3e1d4d3..838c3bf 100644 --- a/src/docs/tex/manual/part2.tex +++ b/src/docs/tex/manual/part2.tex @@ -5931,7 +5931,7 @@ Current routines include: \begin{tabbing} VecGetArrayF90(), VecRestoreArrayF90()\\ VecDuplicateVecsF90(), VecDestroyVecsF90()\\ - DMDAVecGetArrayF90(), ISLocalToGlobalMappingGetIndicesF90()\\ + DMDAVecGetArrayF90(), DMDAVecGetArrayReadF90(), ISLocalToGlobalMappingGetIndicesF90()\\ MatDenseGetArrayF90(), MatDenseRestoreArrayF90()\\ ISGetIndicesF90(), ISRestoreIndicesF90() \end{tabbing}