Actual source code: zplexf90.c

  1: #include <petsc/private/fortranimpl.h>
  2: #include <petscdmplex.h>
  3: #include <petsc/private/f90impl.h>

  5: #if defined(PETSC_HAVE_FORTRAN_CAPS)
  6:   #define dmplexgetcone_                  DMPLEXGETCONE
  7:   #define dmplexrestorecone_              DMPLEXRESTORECONE
  8:   #define dmplexgetconeorientation_       DMPLEXGETCONEORIENTATION
  9:   #define dmplexrestoreconeorientation_   DMPLEXRESTORECONEORIENTATION
 10:   #define dmplexgetsupport_               DMPLEXGETSUPPORT
 11:   #define dmplexrestoresupport_           DMPLEXRESTORESUPPORT
 12:   #define dmplexgettransitiveclosure_     DMPLEXGETTRANSITIVECLOSURE
 13:   #define dmplexrestoretransitiveclosure_ DMPLEXRESTORETRANSITIVECLOSURE
 14:   #define dmplexvecgetclosure_            DMPLEXVECGETCLOSURE
 15:   #define dmplexvecrestoreclosure_        DMPLEXVECRESTORECLOSURE
 16:   #define dmplexvecsetclosure_            DMPLEXVECSETCLOSURE
 17:   #define dmplexmatsetclosure_            DMPLEXMATSETCLOSURE
 18:   #define dmplexgetjoin_                  DMPLEXGETJOIN
 19:   #define dmplexgetfulljoin_              DMPLEXGETFULLJOIN
 20:   #define dmplexrestorejoin_              DMPLEXRESTOREJOIN
 21:   #define dmplexgetmeet_                  DMPLEXGETMEET
 22:   #define dmplexgetfullmeet_              DMPLEXGETFULLMEET
 23:   #define dmplexrestoremeet_              DMPLEXRESTOREMEET
 24: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE) && !defined(FORTRANDOUBLEUNDERSCORE)
 25:   #define dmplexgetcone_                  dmplexgetcone
 26:   #define dmplexrestorecone_              dmplexrestorecone
 27:   #define dmplexgetconeorientation_       dmplexgetconeorientation
 28:   #define dmplexrestoreconeorientation_   dmplexrestoreconeorientation
 29:   #define dmplexgetsupport_               dmplexgetsupport
 30:   #define dmplexrestoresupport_           dmplexrestoresupport
 31:   #define dmplexgettransitiveclosure_     dmplexgettransitiveclosure
 32:   #define dmplexrestoretransitiveclosure_ dmplexrestoretransitiveclosure
 33:   #define dmplexvecgetclosure_            dmplexvecgetclosure
 34:   #define dmplexvecrestoreclosure_        dmplexvecrestoreclosure
 35:   #define dmplexvecsetclosure_            dmplexvecsetclosure
 36:   #define dmplexmatsetclosure_            dmplexmatsetclosure
 37:   #define dmplexgetjoin_                  dmplexgetjoin
 38:   #define dmplexgetfulljoin_              dmplexgetfulljoin
 39:   #define dmplexrestorejoin_              dmplexrestorejoin
 40:   #define dmplexgetmeet_                  dmplexgetmeet
 41:   #define dmplexgetfullmeet_              dmplexgetfullmeet
 42:   #define dmplexrestoremeet_              dmplexrestoremeet
 43: #endif

 45: /* Definitions of Fortran Wrapper routines */

 47: PETSC_EXTERN void dmplexgetcone_(DM *dm, PetscInt *p, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
 48: {
 49:   const PetscInt *v;
 50:   PetscInt        n;

 52:   *ierr = DMPlexGetConeSize(*dm, *p, &n);
 53:   if (*ierr) return;
 54:   *ierr = DMPlexGetCone(*dm, *p, &v);
 55:   if (*ierr) return;
 56:   *ierr = F90Array1dCreate((void *)v, MPIU_INT, 1, n, ptr PETSC_F90_2PTR_PARAM(ptrd));
 57: }

 59: PETSC_EXTERN void dmplexrestorecone_(DM *dm, PetscInt *p, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
 60: {
 61:   *ierr = F90Array1dDestroy(ptr, MPIU_INT PETSC_F90_2PTR_PARAM(ptrd));
 62:   if (*ierr) return;
 63: }

 65: PETSC_EXTERN void dmplexgetconeorientation_(DM *dm, PetscInt *p, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
 66: {
 67:   const PetscInt *v;
 68:   PetscInt        n;

 70:   *ierr = DMPlexGetConeSize(*dm, *p, &n);
 71:   if (*ierr) return;
 72:   *ierr = DMPlexGetConeOrientation(*dm, *p, &v);
 73:   if (*ierr) return;
 74:   *ierr = F90Array1dCreate((void *)v, MPIU_INT, 1, n, ptr PETSC_F90_2PTR_PARAM(ptrd));
 75: }

 77: PETSC_EXTERN void dmplexrestoreconeorientation_(DM *dm, PetscInt *p, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
 78: {
 79:   *ierr = F90Array1dDestroy(ptr, MPIU_INT PETSC_F90_2PTR_PARAM(ptrd));
 80:   if (*ierr) return;
 81: }

 83: PETSC_EXTERN void dmplexgetsupport_(DM *dm, PetscInt *p, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
 84: {
 85:   const PetscInt *v;
 86:   PetscInt        n;

 88:   *ierr = DMPlexGetSupportSize(*dm, *p, &n);
 89:   if (*ierr) return;
 90:   *ierr = DMPlexGetSupport(*dm, *p, &v);
 91:   if (*ierr) return;
 92:   *ierr = F90Array1dCreate((void *)v, MPIU_INT, 1, n, ptr PETSC_F90_2PTR_PARAM(ptrd));
 93: }

 95: PETSC_EXTERN void dmplexrestoresupport_(DM *dm, PetscInt *p, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
 96: {
 97:   *ierr = F90Array1dDestroy(ptr, MPIU_INT PETSC_F90_2PTR_PARAM(ptrd));
 98:   if (*ierr) return;
 99: }

101: PETSC_EXTERN void dmplexgettransitiveclosure_(DM *dm, PetscInt *p, PetscBool *useCone, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
102: {
103:   PetscInt *v = NULL;
104:   PetscInt  n;

106:   *ierr = DMPlexGetTransitiveClosure(*dm, *p, *useCone, &n, &v);
107:   if (*ierr) return;
108:   *ierr = F90Array1dCreate((void *)v, MPIU_INT, 1, n * 2, ptr PETSC_F90_2PTR_PARAM(ptrd));
109: }

111: PETSC_EXTERN void dmplexrestoretransitiveclosure_(DM *dm, PetscInt *p, PetscBool *useCone, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
112: {
113:   PetscInt *array;

115:   *ierr = F90Array1dAccess(ptr, MPIU_INT, (void **)&array PETSC_F90_2PTR_PARAM(ptrd));
116:   if (*ierr) return;
117:   *ierr = DMPlexRestoreTransitiveClosure(*dm, *p, *useCone, NULL, &array);
118:   if (*ierr) return;
119:   *ierr = F90Array1dDestroy(ptr, MPIU_INT PETSC_F90_2PTR_PARAM(ptrd));
120:   if (*ierr) return;
121: }

123: PETSC_EXTERN void dmplexvecgetclosure_(DM *dm, PetscSection *section, Vec *x, PetscInt *point, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
124: {
125:   PetscScalar *v = NULL;
126:   PetscInt     n;

128:   *ierr = DMPlexVecGetClosure(*dm, *section, *x, *point, &n, &v);
129:   if (*ierr) return;
130:   *ierr = F90Array1dCreate((void *)v, MPIU_SCALAR, 1, n, ptr PETSC_F90_2PTR_PARAM(ptrd));
131: }

133: PETSC_EXTERN void dmplexvecrestoreclosure_(DM *dm, PetscSection *section, Vec *v, PetscInt *point, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
134: {
135:   PetscScalar *array;

137:   *ierr = F90Array1dAccess(ptr, MPIU_SCALAR, (void **)&array PETSC_F90_2PTR_PARAM(ptrd));
138:   if (*ierr) return;
139:   *ierr = DMPlexVecRestoreClosure(*dm, *section, *v, *point, NULL, &array);
140:   if (*ierr) return;
141:   *ierr = F90Array1dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
142:   if (*ierr) return;
143: }

145: PETSC_EXTERN void dmplexvecsetclosure_(DM *dm, PetscSection *section, Vec *v, PetscInt *point, F90Array1d *ptr, InsertMode *mode, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
146: {
147:   PetscScalar *array;

149:   *ierr = F90Array1dAccess(ptr, MPIU_SCALAR, (void **)&array PETSC_F90_2PTR_PARAM(ptrd));
150:   if (*ierr) return;
151:   *ierr = DMPlexVecSetClosure(*dm, *section, *v, *point, array, *mode);
152: }

154: PETSC_EXTERN void dmplexmatsetclosure_(DM *dm, PetscSection *section, PetscSection *globalSection, Mat *A, PetscInt *point, F90Array1d *ptr, InsertMode *mode, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
155: {
156:   PetscScalar *array;

158:   *ierr = F90Array1dAccess(ptr, MPIU_SCALAR, (void **)&array PETSC_F90_2PTR_PARAM(ptrd));
159:   if (*ierr) return;
160:   *ierr = DMPlexMatSetClosure(*dm, *section, *globalSection, *A, *point, array, *mode);
161: }

163: PETSC_EXTERN void dmplexgetjoin_(DM *dm, PetscInt *numPoints, F90Array1d *pptr, F90Array1d *cptr, int *ierr PETSC_F90_2PTR_PROTO(pptrd) PETSC_F90_2PTR_PROTO(cptrd))
164: {
165:   PetscInt       *points;
166:   const PetscInt *coveredPoints;
167:   PetscInt        numCoveredPoints;

169:   *ierr = F90Array1dAccess(pptr, MPIU_INT, (void **)&points PETSC_F90_2PTR_PARAM(pptrd));
170:   if (*ierr) return;
171:   *ierr = DMPlexGetJoin(*dm, *numPoints, points, &numCoveredPoints, &coveredPoints);
172:   if (*ierr) return;
173:   *ierr = F90Array1dCreate((void *)coveredPoints, MPIU_INT, 1, numCoveredPoints, cptr PETSC_F90_2PTR_PARAM(cptrd));
174: }

176: PETSC_EXTERN void dmplexgetfulljoin_(DM *dm, PetscInt *numPoints, F90Array1d *pptr, F90Array1d *cptr, int *ierr PETSC_F90_2PTR_PROTO(pptrd) PETSC_F90_2PTR_PROTO(cptrd))
177: {
178:   PetscInt       *points;
179:   const PetscInt *coveredPoints;
180:   PetscInt        numCoveredPoints;

182:   *ierr = F90Array1dAccess(pptr, MPIU_INT, (void **)&points PETSC_F90_2PTR_PARAM(pptrd));
183:   if (*ierr) return;
184:   *ierr = DMPlexGetFullJoin(*dm, *numPoints, points, &numCoveredPoints, &coveredPoints);
185:   if (*ierr) return;
186:   *ierr = F90Array1dCreate((void *)coveredPoints, MPIU_INT, 1, numCoveredPoints, cptr PETSC_F90_2PTR_PARAM(cptrd));
187: }

189: PETSC_EXTERN void dmplexrestorejoin_(DM *dm, PetscInt *numPoints, F90Array1d *pptr, F90Array1d *cptr, int *ierr PETSC_F90_2PTR_PROTO(pptrd) PETSC_F90_2PTR_PROTO(cptrd))
190: {
191:   PetscInt *coveredPoints;

193:   *ierr = F90Array1dAccess(cptr, MPIU_INT, (void **)&coveredPoints PETSC_F90_2PTR_PARAM(cptrd));
194:   if (*ierr) return;
195:   *ierr = DMPlexRestoreJoin(*dm, 0, NULL, NULL, (const PetscInt **)&coveredPoints);
196:   if (*ierr) return;
197:   *ierr = F90Array1dDestroy(cptr, MPIU_INT PETSC_F90_2PTR_PARAM(cptrd));
198:   if (*ierr) return;
199: }

201: PETSC_EXTERN void dmplexgetmeet_(DM *dm, PetscInt *numPoints, F90Array1d *pptr, F90Array1d *cptr, int *ierr PETSC_F90_2PTR_PROTO(pptrd) PETSC_F90_2PTR_PROTO(cptrd))
202: {
203:   PetscInt       *points;
204:   const PetscInt *coveredPoints;
205:   PetscInt        numCoveredPoints;

207:   *ierr = F90Array1dAccess(pptr, MPIU_INT, (void **)&points PETSC_F90_2PTR_PARAM(pptrd));
208:   if (*ierr) return;
209:   *ierr = DMPlexGetMeet(*dm, *numPoints, points, &numCoveredPoints, &coveredPoints);
210:   if (*ierr) return;
211:   *ierr = F90Array1dCreate((void *)coveredPoints, MPIU_INT, 1, numCoveredPoints, cptr PETSC_F90_2PTR_PARAM(cptrd));
212: }

214: PETSC_EXTERN void dmplexgetfullmeet_(DM *dm, PetscInt *numPoints, F90Array1d *pptr, F90Array1d *cptr, int *ierr PETSC_F90_2PTR_PROTO(pptrd) PETSC_F90_2PTR_PROTO(cptrd))
215: {
216:   PetscInt       *points;
217:   const PetscInt *coveredPoints;
218:   PetscInt        numCoveredPoints;

220:   *ierr = F90Array1dAccess(pptr, MPIU_INT, (void **)&points PETSC_F90_2PTR_PARAM(pptrd));
221:   if (*ierr) return;
222:   *ierr = DMPlexGetFullMeet(*dm, *numPoints, points, &numCoveredPoints, &coveredPoints);
223:   if (*ierr) return;
224:   *ierr = F90Array1dCreate((void *)coveredPoints, MPIU_INT, 1, numCoveredPoints, cptr PETSC_F90_2PTR_PARAM(cptrd));
225: }

227: PETSC_EXTERN void dmplexrestoremeet_(DM *dm, PetscInt *numPoints, F90Array1d *pptr, F90Array1d *cptr, int *ierr PETSC_F90_2PTR_PROTO(pptrd) PETSC_F90_2PTR_PROTO(cptrd))
228: {
229:   PetscInt *coveredPoints;

231:   *ierr = F90Array1dAccess(cptr, MPIU_INT, (void **)&coveredPoints PETSC_F90_2PTR_PARAM(cptrd));
232:   if (*ierr) return;
233:   *ierr = DMPlexRestoreMeet(*dm, 0, NULL, NULL, (const PetscInt **)&coveredPoints);
234:   if (*ierr) return;
235:   *ierr = F90Array1dDestroy(cptr, MPIU_INT PETSC_F90_2PTR_PARAM(cptrd));
236:   if (*ierr) return;
237: }