1*6dd63270SBarry Smith #include <petsc/private/ftnimpl.h> 2*6dd63270SBarry Smith #include <petscdmplex.h> 3*6dd63270SBarry Smith #include <petscds.h> 4*6dd63270SBarry Smith 5*6dd63270SBarry Smith #if defined(PETSC_HAVE_FORTRAN_CAPS) 6*6dd63270SBarry Smith #define dmplexgetcellfields_ DMPLEXGETCELLFIELDS 7*6dd63270SBarry Smith #define dmplexrestorecellfields_ DMPLEXRESTORECELLFIELDS 8*6dd63270SBarry Smith #define dmplexgetfacefields_ DMPLEXGETFACEFIELDS 9*6dd63270SBarry Smith #define dmplexrestorefacefields_ DMPLEXRESTOREFACEFIELDS 10*6dd63270SBarry Smith #define dmplexgetfacegeometry_ DMPLEXGETFACEGEOMETRY 11*6dd63270SBarry Smith #define dmplexrestorefacegeometry_ DMPLEXRESTOREFACEGEOMETRY 12*6dd63270SBarry Smith #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE) 13*6dd63270SBarry Smith #define dmplexgetcellfields_ dmplexgetcellfields 14*6dd63270SBarry Smith #define dmplexrestorecellfields_ dmplexrestorecellfields 15*6dd63270SBarry Smith #define dmplexgetfacefields_ dmplexgetfacefields 16*6dd63270SBarry Smith #define dmplexrestorefacefields_ dmplexrestorefacefields 17*6dd63270SBarry Smith #define dmplexgetfacegeometry_ dmplexgetfacegeometry 18*6dd63270SBarry Smith #define dmplexrestorefacegeometry_ dmplexrestorefacegeometry 19*6dd63270SBarry Smith #endif 20*6dd63270SBarry Smith 21*6dd63270SBarry Smith PETSC_EXTERN void dmplexgetcellfields_(DM *dm, IS *cellIS, Vec *locX, Vec *locX_t, Vec *locA, F90Array1d *uPtr, F90Array1d *utPtr, F90Array1d *aPtr, int *ierr PETSC_F90_2PTR_PROTO(uPtrd) PETSC_F90_2PTR_PROTO(utPtrd) PETSC_F90_2PTR_PROTO(aPtrd)) 22*6dd63270SBarry Smith { 23*6dd63270SBarry Smith PetscDS prob; 24*6dd63270SBarry Smith PetscScalar *u, *u_t, *a; 25*6dd63270SBarry Smith PetscInt numCells, totDim, totDimAux = 0; 26*6dd63270SBarry Smith 27*6dd63270SBarry Smith *ierr = ISGetLocalSize(*cellIS, &numCells); 28*6dd63270SBarry Smith if (*ierr) return; 29*6dd63270SBarry Smith *ierr = DMPlexGetCellFields(*dm, *cellIS, *locX, *locX_t, *locA, &u, &u_t, &a); 30*6dd63270SBarry Smith if (*ierr) return; 31*6dd63270SBarry Smith *ierr = DMGetDS(*dm, &prob); 32*6dd63270SBarry Smith if (*ierr) return; 33*6dd63270SBarry Smith *ierr = PetscDSGetTotalDimension(prob, &totDim); 34*6dd63270SBarry Smith if (*ierr) return; 35*6dd63270SBarry Smith if (locA) { 36*6dd63270SBarry Smith DM dmAux; 37*6dd63270SBarry Smith PetscDS probAux; 38*6dd63270SBarry Smith 39*6dd63270SBarry Smith *ierr = VecGetDM(*locA, &dmAux); 40*6dd63270SBarry Smith if (*ierr) return; 41*6dd63270SBarry Smith *ierr = DMGetDS(dmAux, &probAux); 42*6dd63270SBarry Smith if (*ierr) return; 43*6dd63270SBarry Smith *ierr = PetscDSGetTotalDimension(probAux, &totDimAux); 44*6dd63270SBarry Smith if (*ierr) return; 45*6dd63270SBarry Smith } 46*6dd63270SBarry Smith *ierr = F90Array1dCreate((void *)u, MPIU_SCALAR, 1, numCells * totDim, uPtr PETSC_F90_2PTR_PARAM(uPtrd)); 47*6dd63270SBarry Smith if (*ierr) return; 48*6dd63270SBarry Smith *ierr = F90Array1dCreate((void *)u_t, MPIU_SCALAR, 1, locX_t ? numCells * totDim : 0, utPtr PETSC_F90_2PTR_PARAM(utPtrd)); 49*6dd63270SBarry Smith if (*ierr) return; 50*6dd63270SBarry Smith *ierr = F90Array1dCreate((void *)a, MPIU_SCALAR, 1, locA ? numCells * totDimAux : 0, aPtr PETSC_F90_2PTR_PARAM(aPtrd)); 51*6dd63270SBarry Smith } 52*6dd63270SBarry Smith 53*6dd63270SBarry Smith PETSC_EXTERN void dmplexrestorecellfields_(DM *dm, IS *cellIS, Vec *locX, Vec *locX_t, Vec *locA, F90Array1d *uPtr, F90Array1d *utPtr, F90Array1d *aPtr, int *ierr PETSC_F90_2PTR_PROTO(uPtrd) PETSC_F90_2PTR_PROTO(utPtrd) PETSC_F90_2PTR_PROTO(aPtrd)) 54*6dd63270SBarry Smith { 55*6dd63270SBarry Smith PetscScalar *u, *u_t, *a; 56*6dd63270SBarry Smith 57*6dd63270SBarry Smith *ierr = F90Array1dAccess(uPtr, MPIU_SCALAR, (void **)&u PETSC_F90_2PTR_PARAM(uPtrd)); 58*6dd63270SBarry Smith if (*ierr) return; 59*6dd63270SBarry Smith *ierr = F90Array1dAccess(utPtr, MPIU_SCALAR, (void **)&u_t PETSC_F90_2PTR_PARAM(utPtrd)); 60*6dd63270SBarry Smith if (*ierr) return; 61*6dd63270SBarry Smith *ierr = F90Array1dAccess(aPtr, MPIU_SCALAR, (void **)&a PETSC_F90_2PTR_PARAM(aPtrd)); 62*6dd63270SBarry Smith if (*ierr) return; 63*6dd63270SBarry Smith *ierr = DMPlexRestoreCellFields(*dm, *cellIS, *locX, NULL, NULL, &u, u_t ? &u_t : NULL, a ? &a : NULL); 64*6dd63270SBarry Smith if (*ierr) return; 65*6dd63270SBarry Smith *ierr = F90Array1dDestroy(uPtr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(uPtrd)); 66*6dd63270SBarry Smith if (*ierr) return; 67*6dd63270SBarry Smith *ierr = F90Array1dDestroy(utPtr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(utPtrd)); 68*6dd63270SBarry Smith if (*ierr) return; 69*6dd63270SBarry Smith *ierr = F90Array1dDestroy(aPtr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(aPtrd)); 70*6dd63270SBarry Smith if (*ierr) return; 71*6dd63270SBarry Smith } 72*6dd63270SBarry Smith 73*6dd63270SBarry Smith PETSC_EXTERN void dmplexgetfacefields_(DM *dm, PetscInt *fStart, PetscInt *fEnd, Vec *locX, Vec *locX_t, Vec *faceGeometry, Vec *cellGeometry, Vec *locGrad, PetscInt *Nface, F90Array1d *uLPtr, F90Array1d *uRPtr, int *ierr PETSC_F90_2PTR_PROTO(uLPtrd) PETSC_F90_2PTR_PROTO(uRPtrd)) 74*6dd63270SBarry Smith { 75*6dd63270SBarry Smith PetscDS prob; 76*6dd63270SBarry Smith PetscScalar *uL, *uR; 77*6dd63270SBarry Smith PetscInt numFaces = *fEnd - *fStart, totDim; 78*6dd63270SBarry Smith 79*6dd63270SBarry Smith *ierr = DMPlexGetFaceFields(*dm, *fStart, *fEnd, *locX, *locX_t, *faceGeometry, *cellGeometry, *locGrad, Nface, &uL, &uR); 80*6dd63270SBarry Smith if (*ierr) return; 81*6dd63270SBarry Smith *ierr = DMGetDS(*dm, &prob); 82*6dd63270SBarry Smith if (*ierr) return; 83*6dd63270SBarry Smith *ierr = PetscDSGetTotalDimension(prob, &totDim); 84*6dd63270SBarry Smith if (*ierr) return; 85*6dd63270SBarry Smith *ierr = F90Array1dCreate((void *)uL, MPIU_SCALAR, 1, numFaces * totDim, uLPtr PETSC_F90_2PTR_PARAM(uLPtrd)); 86*6dd63270SBarry Smith if (*ierr) return; 87*6dd63270SBarry Smith *ierr = F90Array1dCreate((void *)uR, MPIU_SCALAR, 1, numFaces * totDim, uRPtr PETSC_F90_2PTR_PARAM(uRPtrd)); 88*6dd63270SBarry Smith if (*ierr) return; 89*6dd63270SBarry Smith } 90*6dd63270SBarry Smith 91*6dd63270SBarry Smith PETSC_EXTERN void dmplexrestorefacefields_(DM *dm, PetscInt *fStart, PetscInt *fEnd, Vec *locX, Vec *locX_t, Vec *faceGeometry, Vec *cellGeometry, Vec *locGrad, PetscInt *Nface, F90Array1d *uLPtr, F90Array1d *uRPtr, int *ierr PETSC_F90_2PTR_PROTO(uLPtrd) PETSC_F90_2PTR_PROTO(uRPtrd)) 92*6dd63270SBarry Smith { 93*6dd63270SBarry Smith PetscScalar *uL, *uR; 94*6dd63270SBarry Smith 95*6dd63270SBarry Smith *ierr = F90Array1dAccess(uLPtr, MPIU_SCALAR, (void **)&uL PETSC_F90_2PTR_PARAM(uLPtrd)); 96*6dd63270SBarry Smith if (*ierr) return; 97*6dd63270SBarry Smith *ierr = F90Array1dAccess(uRPtr, MPIU_SCALAR, (void **)&uR PETSC_F90_2PTR_PARAM(uRPtrd)); 98*6dd63270SBarry Smith if (*ierr) return; 99*6dd63270SBarry Smith *ierr = DMPlexRestoreFaceFields(*dm, *fStart, *fEnd, *locX, NULL, *faceGeometry, *cellGeometry, NULL, Nface, &uL, &uR); 100*6dd63270SBarry Smith if (*ierr) return; 101*6dd63270SBarry Smith *ierr = F90Array1dDestroy(uLPtr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(uLPtrd)); 102*6dd63270SBarry Smith if (*ierr) return; 103*6dd63270SBarry Smith *ierr = F90Array1dDestroy(uRPtr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(uRPtrd)); 104*6dd63270SBarry Smith if (*ierr) return; 105*6dd63270SBarry Smith } 106*6dd63270SBarry Smith 107*6dd63270SBarry Smith PETSC_EXTERN void dmplexgetfacegeometry_(DM *dm, PetscInt *fStart, PetscInt *fEnd, Vec *faceGeometry, Vec *cellGeometry, PetscInt *Nface, F90Array1d *gPtr, F90Array1d *vPtr, int *ierr PETSC_F90_2PTR_PROTO(gPtrd) PETSC_F90_2PTR_PROTO(vPtrd)) 108*6dd63270SBarry Smith { 109*6dd63270SBarry Smith PetscFVFaceGeom *g; 110*6dd63270SBarry Smith PetscReal *v; 111*6dd63270SBarry Smith PetscInt numFaces = *fEnd - *fStart, structSize = sizeof(PetscFVFaceGeom) / sizeof(PetscScalar); 112*6dd63270SBarry Smith 113*6dd63270SBarry Smith *ierr = DMPlexGetFaceGeometry(*dm, *fStart, *fEnd, *faceGeometry, *cellGeometry, Nface, &g, &v); 114*6dd63270SBarry Smith if (*ierr) return; 115*6dd63270SBarry Smith *ierr = F90Array1dCreate((void *)g, MPIU_SCALAR, 1, numFaces * structSize, gPtr PETSC_F90_2PTR_PARAM(gPtrd)); 116*6dd63270SBarry Smith if (*ierr) return; 117*6dd63270SBarry Smith *ierr = F90Array1dCreate((void *)v, MPIU_REAL, 1, numFaces * 2, vPtr PETSC_F90_2PTR_PARAM(vPtrd)); 118*6dd63270SBarry Smith if (*ierr) return; 119*6dd63270SBarry Smith } 120*6dd63270SBarry Smith 121*6dd63270SBarry Smith PETSC_EXTERN void dmplexrestorefacegeometry_(DM *dm, PetscInt *fStart, PetscInt *fEnd, Vec *faceGeometry, Vec *cellGeometry, PetscInt *Nface, F90Array1d *gPtr, F90Array1d *vPtr, int *ierr PETSC_F90_2PTR_PROTO(gPtrd) PETSC_F90_2PTR_PROTO(vPtrd)) 122*6dd63270SBarry Smith { 123*6dd63270SBarry Smith PetscFVFaceGeom *g; 124*6dd63270SBarry Smith PetscReal *v; 125*6dd63270SBarry Smith 126*6dd63270SBarry Smith *ierr = F90Array1dAccess(gPtr, MPIU_SCALAR, (void **)&g PETSC_F90_2PTR_PARAM(gPtrd)); 127*6dd63270SBarry Smith if (*ierr) return; 128*6dd63270SBarry Smith *ierr = F90Array1dAccess(vPtr, MPIU_REAL, (void **)&v PETSC_F90_2PTR_PARAM(vPtrd)); 129*6dd63270SBarry Smith if (*ierr) return; 130*6dd63270SBarry Smith *ierr = DMPlexRestoreFaceGeometry(*dm, *fStart, *fEnd, *faceGeometry, *cellGeometry, Nface, &g, &v); 131*6dd63270SBarry Smith if (*ierr) return; 132*6dd63270SBarry Smith *ierr = F90Array1dDestroy(gPtr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(vPtrd)); 133*6dd63270SBarry Smith if (*ierr) return; 134*6dd63270SBarry Smith *ierr = F90Array1dDestroy(vPtr, MPIU_REAL PETSC_F90_2PTR_PARAM(gPtrd)); 135*6dd63270SBarry Smith if (*ierr) return; 136*6dd63270SBarry Smith } 137