Actual source code: zsda.c

  1: /*
  2:      Fortran interface for SDA routines.
  3: */
 4:  #include src/fortran/custom/zpetsc.h
 5:  #include petscda.h

  7: #ifdef PETSC_HAVE_FORTRAN_CAPS
  8: #define sdadestroy_           SDADESTROY
  9: #define sdalocaltolocalbegin_ SDALOCALTOLOCALBEGIN
 10: #define sdalocaltolocalend_   SDALOCALTOLOCALEND
 11: #define sdacreate1d_          SDACREATE1D
 12: #define sdacreate2d_          SDACREATE2D
 13: #define sdacreate3d_          SDACREATE3D
 14: #define sdagetghostcorners_   SDAGETGHOSTCORNERS
 15: #define sdagetcorners_        SDAGETCORNERS
 16: #define sdaarrayview_         SDAARRAYVIEW
 17: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
 18: #define sdadestroy_           sdadestroy
 19: #define sdalocaltolocalbegin_ sdalocaltolocalbegin
 20: #define sdalocaltolocalend_   sdalocaltolocalend
 21: #define sdacreate1d_          sdacreate1d
 22: #define sdacreate2d_          sdacreate2d
 23: #define sdacreate3d_          sdacreate3d
 24: #define sdagetghostcorners_   sdagetghostcorners
 25: #define sdagetcorners_        sdagetcorners
 26: #define sdaarrayview_         sdaarrayview
 27: #endif

 29: EXTERN PetscErrorCode SDAArrayView(SDA,PetscScalar*,PetscViewer);

 32: void PETSC_STDCALL sdaarrayview_(SDA *da,PetscScalar *values,PetscViewer *vin,PetscErrorCode *ierr)
 33: {
 34:   PetscViewer v;
 35:   PetscPatchDefaultViewers_Fortran(vin,v);
 36:   *SDAArrayView(*da,values,v);
 37: }

 39: void PETSC_STDCALL sdagetghostcorners_(SDA *da,PetscInt *x,PetscInt *y,PetscInt *z,PetscInt *m,PetscInt *n,PetscInt *p,PetscErrorCode *ierr)
 40: {
 41:   CHKFORTRANNULLINTEGER(x);
 42:   CHKFORTRANNULLINTEGER(y);
 43:   CHKFORTRANNULLINTEGER(z);
 44:   CHKFORTRANNULLINTEGER(m);
 45:   CHKFORTRANNULLINTEGER(n);
 46:   CHKFORTRANNULLINTEGER(p);
 47:   *SDAGetGhostCorners(*da,x,y,z,m,n,p);
 48: }

 50: void PETSC_STDCALL sdagetcorners_(SDA *da,PetscInt *x,PetscInt *y,PetscInt *z,PetscInt *m,PetscInt *n,PetscInt *p,PetscErrorCode *ierr)
 51: {
 52:   CHKFORTRANNULLINTEGER(x);
 53:   CHKFORTRANNULLINTEGER(y);
 54:   CHKFORTRANNULLINTEGER(z);
 55:   CHKFORTRANNULLINTEGER(m);
 56:   CHKFORTRANNULLINTEGER(n);
 57:   CHKFORTRANNULLINTEGER(p);
 58:   *SDAGetCorners(*da,x,y,z,m,n,p);
 59: }

 61: void PETSC_STDCALL sdadestroy_(SDA *sda,PetscErrorCode *ierr)
 62: {
 63:   *SDADestroy((SDA)PetscToPointer(sda));
 64:   PetscRmPointer(sda);
 65: }

 67: void PETSC_STDCALL sdalocaltolocalbegin_(SDA *sda,PetscScalar *g,InsertMode *mode,PetscScalar *l,
 68:                            PetscErrorCode *ierr)
 69: {
 70:   *SDALocalToLocalBegin((SDA)PetscToPointer(sda),g,*mode,l);
 71: }

 73: void PETSC_STDCALL sdalocaltolocalend_(SDA *sda,PetscScalar *g,InsertMode *mode,PetscScalar *l,
 74:                          PetscErrorCode *ierr){
 75:   *SDALocalToLocalEnd((SDA)PetscToPointer(sda),g,*mode,l);
 76: }

 78: void PETSC_STDCALL sdacreate2d_(MPI_Comm *comm,DAPeriodicType *wrap,DAStencilType
 79:                   *stencil_type,PetscInt *M,PetscInt *N,PetscInt *m,PetscInt *n,PetscInt *w,
 80:                   PetscInt *s,PetscInt *lx,PetscInt *ly,SDA *inra,PetscErrorCode *ierr)
 81: {
 82:   CHKFORTRANNULLINTEGER(lx);
 83:   CHKFORTRANNULLINTEGER(ly);
 84:   *SDACreate2d(
 85:             (MPI_Comm)PetscToPointerComm(*comm),*wrap,
 86:             *stencil_type,*M,*N,*m,*n,*w,*s,lx,ly,inra);
 87: }

 89: void PETSC_STDCALL sdacreate1d_(MPI_Comm *comm,DAPeriodicType *wrap,PetscInt *M,PetscInt *w,PetscInt *s,
 90:                  PetscInt *lc,SDA *inra,PetscErrorCode *ierr)
 91: {
 92:   CHKFORTRANNULLINTEGER(lc);
 93:   *SDACreate1d(
 94:            (MPI_Comm)PetscToPointerComm(*comm),*wrap,*M,*w,*s,lc,inra);
 95: }

 97: void PETSC_STDCALL sdacreate3d_(MPI_Comm *comm,DAPeriodicType *wrap,DAStencilType 
 98:                  *stencil_type,PetscInt *M,PetscInt *N,PetscInt *P,PetscInt *m,PetscInt *n,PetscInt *p,
 99:                  PetscInt *w,PetscInt *s,PetscInt *lx,PetscInt *ly,PetscInt *lz,SDA *inra,PetscErrorCode *ierr)
100: {
101:   CHKFORTRANNULLINTEGER(lx);
102:   CHKFORTRANNULLINTEGER(ly);
103:   CHKFORTRANNULLINTEGER(lz);
104:   *SDACreate3d(
105:            (MPI_Comm)PetscToPointerComm(*comm),*wrap,*stencil_type,
106:            *M,*N,*P,*m,*n,*p,*w,*s,lx,ly,lz,inra);
107: }