Actual source code: zdaf.c

  1: #include <petsc/private/ftnimpl.h>
  2: #include <petsc/private/dmdaimpl.h>

  4: #if defined(PETSC_HAVE_FORTRAN_CAPS)
  5:   #define dmdagetownershipranges_     DMDAGETOWNERSHIPRANGES
  6:   #define dmdarestoreownershipranges_ DMDARESTOREOWNERSHIPRANGES
  7:   #define dmdagetneighbors_           DMDAGETNEIGHBORS
  8:   #define dmdarestoreneighbors_       DMDARESTORENEIGHBORS
  9: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
 10:   #define dmdagetownershipranges_     dmdagetownershipranges
 11:   #define dmdarestoreownershipranges_ dmdarestoreownershipranges
 12:   #define dmdagetneighbors_           dmdagetneighbors
 13:   #define dmdarestoreneighbors_       dmdarestoreneighbors
 14: #endif

 16: PETSC_EXTERN void dmdagetneighbors_(DM *da, F90Array1d *ptr, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
 17: {
 18:   const PetscMPIInt *r;
 19:   PetscInt           n, dim;

 21:   *ierr = DMDAGetNeighbors(*da, &r);
 22:   if (*ierr) return;
 23:   *ierr = DMGetDimension(*da, &dim);
 24:   if (*ierr) return;
 25:   if (dim == 2) n = 9;
 26:   else n = 27;
 27:   *ierr = F90Array1dCreate((PetscInt *)r, MPI_INT, 1, n, ptr PETSC_F90_2PTR_PARAM(ptrd));
 28: }

 30: PETSC_EXTERN void dmdarestoreneighbors_(DM *da, F90Array1d *ptr, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
 31: {
 32:   *ierr = F90Array1dDestroy(ptr, MPI_INT PETSC_F90_2PTR_PARAM(ptrd));
 33: }

 35: PETSC_EXTERN void dmdagetownershipranges_(DM *da, F90Array1d *lx, F90Array1d *ly, F90Array1d *lz, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(lxd) PETSC_F90_2PTR_PROTO(lyd) PETSC_F90_2PTR_PROTO(lzd))
 36: {
 37:   const PetscInt *gx, *gy, *gz;
 38:   PetscInt        M, N, P;

 40:   *ierr = DMDAGetInfo(*da, NULL, NULL, NULL, NULL, &M, &N, &P, NULL, NULL, NULL, NULL, NULL, NULL);
 41:   if (*ierr) return;
 42:   *ierr = DMDAGetOwnershipRanges(*da, &gx, &gy, &gz);
 43:   if (*ierr) return;
 44:   if ((void *)lx != PETSC_NULL_INTEGER_POINTER_Fortran) {
 45:     *ierr = F90Array1dCreate((PetscInt *)gx, MPIU_INT, 1, M, lx PETSC_F90_2PTR_PARAM(lxd));
 46:     if (*ierr) return;
 47:   }
 48:   if ((void *)ly != PETSC_NULL_INTEGER_POINTER_Fortran) {
 49:     *ierr = F90Array1dCreate((PetscInt *)gy, MPIU_INT, 1, N, ly PETSC_F90_2PTR_PARAM(lyd));
 50:     if (*ierr) return;
 51:   }
 52:   if ((void *)lz != PETSC_NULL_INTEGER_POINTER_Fortran) *ierr = F90Array1dCreate((PetscInt *)gz, MPIU_INT, 1, P, lz PETSC_F90_2PTR_PARAM(lzd));
 53: }

 55: PETSC_EXTERN void dmdarestoreownershipranges_(DM *da, F90Array1d *lx, F90Array1d *ly, F90Array1d *lz, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(lxd) PETSC_F90_2PTR_PROTO(lyd) PETSC_F90_2PTR_PROTO(lzd))
 56: {
 57:   if ((void *)lx != PETSC_NULL_INTEGER_POINTER_Fortran) {
 58:     *ierr = F90Array1dDestroy(lx, MPIU_INT PETSC_F90_2PTR_PARAM(lxd));
 59:     if (*ierr) return;
 60:   }
 61:   if ((void *)ly != PETSC_NULL_INTEGER_POINTER_Fortran) {
 62:     *ierr = F90Array1dDestroy(ly, MPIU_INT PETSC_F90_2PTR_PARAM(lyd));
 63:     if (*ierr) return;
 64:   }
 65:   if ((void *)lz != PETSC_NULL_INTEGER_POINTER_Fortran) *ierr = F90Array1dDestroy(lz, MPIU_INT PETSC_F90_2PTR_PARAM(lzd));
 66: }