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: PETSC_EXTERN void dmplexgetcone_(DM *dm, PetscInt *p, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
46: {
47: const PetscInt *v;
48: PetscInt n;
50: *ierr = DMPlexGetConeSize(*dm, *p, &n);
51: if (*ierr) return;
52: *ierr = DMPlexGetCone(*dm, *p, &v);
53: if (*ierr) return;
54: *ierr = F90Array1dCreate((void *)v, MPIU_INT, 1, n, ptr PETSC_F90_2PTR_PARAM(ptrd));
55: }
57: PETSC_EXTERN void dmplexrestorecone_(DM *dm, PetscInt *p, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
58: {
59: *ierr = F90Array1dDestroy(ptr, MPIU_INT PETSC_F90_2PTR_PARAM(ptrd));
60: if (*ierr) return;
61: }
63: PETSC_EXTERN void dmplexgetconeorientation_(DM *dm, PetscInt *p, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
64: {
65: const PetscInt *v;
66: PetscInt n;
68: *ierr = DMPlexGetConeSize(*dm, *p, &n);
69: if (*ierr) return;
70: *ierr = DMPlexGetConeOrientation(*dm, *p, &v);
71: if (*ierr) return;
72: *ierr = F90Array1dCreate((void *)v, MPIU_INT, 1, n, ptr PETSC_F90_2PTR_PARAM(ptrd));
73: }
75: PETSC_EXTERN void dmplexrestoreconeorientation_(DM *dm, PetscInt *p, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
76: {
77: *ierr = F90Array1dDestroy(ptr, MPIU_INT PETSC_F90_2PTR_PARAM(ptrd));
78: if (*ierr) return;
79: }
81: PETSC_EXTERN void dmplexgetsupport_(DM *dm, PetscInt *p, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
82: {
83: const PetscInt *v;
84: PetscInt n;
86: *ierr = DMPlexGetSupportSize(*dm, *p, &n);
87: if (*ierr) return;
88: *ierr = DMPlexGetSupport(*dm, *p, &v);
89: if (*ierr) return;
90: *ierr = F90Array1dCreate((void *)v, MPIU_INT, 1, n, ptr PETSC_F90_2PTR_PARAM(ptrd));
91: }
93: PETSC_EXTERN void dmplexrestoresupport_(DM *dm, PetscInt *p, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
94: {
95: *ierr = F90Array1dDestroy(ptr, MPIU_INT PETSC_F90_2PTR_PARAM(ptrd));
96: if (*ierr) return;
97: }
99: PETSC_EXTERN void dmplexgettransitiveclosure_(DM *dm, PetscInt *p, PetscBool *useCone, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
100: {
101: PetscInt *v = NULL;
102: PetscInt n;
104: *ierr = DMPlexGetTransitiveClosure(*dm, *p, *useCone, &n, &v);
105: if (*ierr) return;
106: *ierr = F90Array1dCreate((void *)v, MPIU_INT, 1, n * 2, ptr PETSC_F90_2PTR_PARAM(ptrd));
107: }
109: PETSC_EXTERN void dmplexrestoretransitiveclosure_(DM *dm, PetscInt *p, PetscBool *useCone, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
110: {
111: PetscInt *array;
113: *ierr = F90Array1dAccess(ptr, MPIU_INT, (void **)&array PETSC_F90_2PTR_PARAM(ptrd));
114: if (*ierr) return;
115: *ierr = DMPlexRestoreTransitiveClosure(*dm, *p, *useCone, NULL, &array);
116: if (*ierr) return;
117: *ierr = F90Array1dDestroy(ptr, MPIU_INT PETSC_F90_2PTR_PARAM(ptrd));
118: if (*ierr) return;
119: }
121: PETSC_EXTERN void dmplexvecgetclosure_(DM *dm, PetscSection *section, Vec *x, PetscInt *point, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
122: {
123: PetscScalar *v = NULL;
124: PetscInt n;
126: *ierr = DMPlexVecGetClosure(*dm, *section, *x, *point, &n, &v);
127: if (*ierr) return;
128: *ierr = F90Array1dCreate((void *)v, MPIU_SCALAR, 1, n, ptr PETSC_F90_2PTR_PARAM(ptrd));
129: }
131: PETSC_EXTERN void dmplexvecrestoreclosure_(DM *dm, PetscSection *section, Vec *v, PetscInt *point, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
132: {
133: PetscScalar *array;
135: *ierr = F90Array1dAccess(ptr, MPIU_SCALAR, (void **)&array PETSC_F90_2PTR_PARAM(ptrd));
136: if (*ierr) return;
137: *ierr = DMPlexVecRestoreClosure(*dm, *section, *v, *point, NULL, &array);
138: if (*ierr) return;
139: *ierr = F90Array1dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
140: if (*ierr) return;
141: }
143: PETSC_EXTERN void dmplexvecsetclosure_(DM *dm, PetscSection *section, Vec *v, PetscInt *point, F90Array1d *ptr, InsertMode *mode, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
144: {
145: PetscScalar *array;
147: *ierr = F90Array1dAccess(ptr, MPIU_SCALAR, (void **)&array PETSC_F90_2PTR_PARAM(ptrd));
148: if (*ierr) return;
149: *ierr = DMPlexVecSetClosure(*dm, *section, *v, *point, array, *mode);
150: }
152: 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))
153: {
154: PetscScalar *array;
156: *ierr = F90Array1dAccess(ptr, MPIU_SCALAR, (void **)&array PETSC_F90_2PTR_PARAM(ptrd));
157: if (*ierr) return;
158: *ierr = DMPlexMatSetClosure(*dm, *section, *globalSection, *A, *point, array, *mode);
159: }
161: PETSC_EXTERN void dmplexgetjoin_(DM *dm, PetscInt *numPoints, F90Array1d *pptr, F90Array1d *cptr, int *ierr PETSC_F90_2PTR_PROTO(pptrd) PETSC_F90_2PTR_PROTO(cptrd))
162: {
163: PetscInt *points;
164: const PetscInt *coveredPoints;
165: PetscInt numCoveredPoints;
167: *ierr = F90Array1dAccess(pptr, MPIU_INT, (void **)&points PETSC_F90_2PTR_PARAM(pptrd));
168: if (*ierr) return;
169: *ierr = DMPlexGetJoin(*dm, *numPoints, points, &numCoveredPoints, &coveredPoints);
170: if (*ierr) return;
171: *ierr = F90Array1dCreate((void *)coveredPoints, MPIU_INT, 1, numCoveredPoints, cptr PETSC_F90_2PTR_PARAM(cptrd));
172: }
174: PETSC_EXTERN void dmplexgetfulljoin_(DM *dm, PetscInt *numPoints, F90Array1d *pptr, F90Array1d *cptr, int *ierr PETSC_F90_2PTR_PROTO(pptrd) PETSC_F90_2PTR_PROTO(cptrd))
175: {
176: PetscInt *points;
177: const PetscInt *coveredPoints;
178: PetscInt numCoveredPoints;
180: *ierr = F90Array1dAccess(pptr, MPIU_INT, (void **)&points PETSC_F90_2PTR_PARAM(pptrd));
181: if (*ierr) return;
182: *ierr = DMPlexGetFullJoin(*dm, *numPoints, points, &numCoveredPoints, &coveredPoints);
183: if (*ierr) return;
184: *ierr = F90Array1dCreate((void *)coveredPoints, MPIU_INT, 1, numCoveredPoints, cptr PETSC_F90_2PTR_PARAM(cptrd));
185: }
187: PETSC_EXTERN void dmplexrestorejoin_(DM *dm, PetscInt *numPoints, F90Array1d *pptr, F90Array1d *cptr, int *ierr PETSC_F90_2PTR_PROTO(pptrd) PETSC_F90_2PTR_PROTO(cptrd))
188: {
189: PetscInt *coveredPoints;
191: *ierr = F90Array1dAccess(cptr, MPIU_INT, (void **)&coveredPoints PETSC_F90_2PTR_PARAM(cptrd));
192: if (*ierr) return;
193: *ierr = DMPlexRestoreJoin(*dm, 0, NULL, NULL, (const PetscInt **)&coveredPoints);
194: if (*ierr) return;
195: *ierr = F90Array1dDestroy(cptr, MPIU_INT PETSC_F90_2PTR_PARAM(cptrd));
196: if (*ierr) return;
197: }
199: PETSC_EXTERN void dmplexgetmeet_(DM *dm, PetscInt *numPoints, F90Array1d *pptr, F90Array1d *cptr, int *ierr PETSC_F90_2PTR_PROTO(pptrd) PETSC_F90_2PTR_PROTO(cptrd))
200: {
201: PetscInt *points;
202: const PetscInt *coveredPoints;
203: PetscInt numCoveredPoints;
205: *ierr = F90Array1dAccess(pptr, MPIU_INT, (void **)&points PETSC_F90_2PTR_PARAM(pptrd));
206: if (*ierr) return;
207: *ierr = DMPlexGetMeet(*dm, *numPoints, points, &numCoveredPoints, &coveredPoints);
208: if (*ierr) return;
209: *ierr = F90Array1dCreate((void *)coveredPoints, MPIU_INT, 1, numCoveredPoints, cptr PETSC_F90_2PTR_PARAM(cptrd));
210: }
212: PETSC_EXTERN void dmplexgetfullmeet_(DM *dm, PetscInt *numPoints, F90Array1d *pptr, F90Array1d *cptr, int *ierr PETSC_F90_2PTR_PROTO(pptrd) PETSC_F90_2PTR_PROTO(cptrd))
213: {
214: PetscInt *points;
215: const PetscInt *coveredPoints;
216: PetscInt numCoveredPoints;
218: *ierr = F90Array1dAccess(pptr, MPIU_INT, (void **)&points PETSC_F90_2PTR_PARAM(pptrd));
219: if (*ierr) return;
220: *ierr = DMPlexGetFullMeet(*dm, *numPoints, points, &numCoveredPoints, &coveredPoints);
221: if (*ierr) return;
222: *ierr = F90Array1dCreate((void *)coveredPoints, MPIU_INT, 1, numCoveredPoints, cptr PETSC_F90_2PTR_PARAM(cptrd));
223: }
225: PETSC_EXTERN void dmplexrestoremeet_(DM *dm, PetscInt *numPoints, F90Array1d *pptr, F90Array1d *cptr, int *ierr PETSC_F90_2PTR_PROTO(pptrd) PETSC_F90_2PTR_PROTO(cptrd))
226: {
227: PetscInt *coveredPoints;
229: *ierr = F90Array1dAccess(cptr, MPIU_INT, (void **)&coveredPoints PETSC_F90_2PTR_PARAM(cptrd));
230: if (*ierr) return;
231: *ierr = DMPlexRestoreMeet(*dm, 0, NULL, NULL, (const PetscInt **)&coveredPoints);
232: if (*ierr) return;
233: *ierr = F90Array1dDestroy(cptr, MPIU_INT PETSC_F90_2PTR_PARAM(cptrd));
234: if (*ierr) return;
235: }