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