1 #include <petsc/private/fortranimpl.h> 2 #include <petsc/private/dmdaimpl.h> 3 4 #if defined(PETSC_HAVE_FORTRAN_CAPS) 5 #define dmdagetownershipranges_ DMDAGETOWNERSHIPRANGES 6 #define dmdagetneighbors_ DMDAGETNEIGHBORS 7 #define dmdagetrefinementfactor_ DMDAGETREFINEMENTFACTOR 8 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE) 9 #define dmdagetownershipranges_ dmdagetownershipranges 10 #define dmdagetneighbors_ dmdagetneighbors 11 #define dmdagetrefinementfactor_ dmdagetrefinementfactor 12 #endif 13 14 PETSC_EXTERN void dmdagetneighbors_(DM *da,PetscMPIInt *ranks,PetscErrorCode *ierr) 15 { 16 const PetscMPIInt *r; 17 PetscInt n,dim; 18 19 *ierr = DMDAGetNeighbors(*da,&r);if (*ierr) return; 20 *ierr = DMGetDimension(*da,&dim);if (*ierr) return; 21 if (dim == 2) n = 9; 22 else n = 27; 23 *ierr = PetscArraycpy(ranks,r,n); 24 } 25 26 PETSC_EXTERN void dmdagetownershipranges_(DM *da,PetscInt lx[],PetscInt ly[],PetscInt lz[],PetscErrorCode *ierr) 27 { 28 const PetscInt *gx,*gy,*gz; 29 PetscInt M,N,P,i; 30 31 CHKFORTRANNULLINTEGER(lx); 32 CHKFORTRANNULLINTEGER(ly); 33 CHKFORTRANNULLINTEGER(lz); 34 *ierr = DMDAGetInfo(*da,0,0,0,0,&M,&N,&P,0,0,0,0,0,0);if (*ierr) return; 35 *ierr = DMDAGetOwnershipRanges(*da,&gx,&gy,&gz);if (*ierr) return; 36 if (lx) { 37 for (i=0; i<M; i++) lx[i] = gx[i]; 38 } 39 if (ly) { 40 for (i=0; i<N; i++) ly[i] = gy[i]; 41 } 42 if (lz) { 43 for (i=0; i<P; i++) lz[i] = gz[i]; 44 } 45 } 46 47 PETSC_EXTERN void dmdagetrefinementfactor_(DM *da,PetscInt *lx,PetscInt *ly,PetscInt *lz,PetscErrorCode *ierr) 48 { 49 CHKFORTRANNULLINTEGER(lx); 50 CHKFORTRANNULLINTEGER(ly); 51 CHKFORTRANNULLINTEGER(lz); 52 *ierr = DMDAGetRefinementFactor(*da,lx,ly,lz); 53 } 54 55