1 #include <petsc/private/ftnimpl.h>
2 #include <petscdmplex.h>
3
4 #if defined(PETSC_HAVE_FORTRAN_CAPS)
5 #define dmplexgetcone_ DMPLEXGETCONE
6 #define dmplexrestorecone_ DMPLEXRESTORECONE
7 #define dmplexgetconeorientation_ DMPLEXGETCONEORIENTATION
8 #define dmplexrestoreconeorientation_ DMPLEXRESTORECONEORIENTATION
9 #define dmplexgetsupport_ DMPLEXGETSUPPORT
10 #define dmplexrestoresupport_ DMPLEXRESTORESUPPORT
11 #define dmplexgettransitiveclosure_ DMPLEXGETTRANSITIVECLOSURE
12 #define dmplexrestoretransitiveclosure_ DMPLEXRESTORETRANSITIVECLOSURE
13 #define dmplexvecgetclosure_ DMPLEXVECGETCLOSURE
14 #define dmplexvecrestoreclosure_ DMPLEXVECRESTORECLOSURE
15 #define dmplexvecsetclosure_ DMPLEXVECSETCLOSURE
16 #define dmplexmatsetclosure_ DMPLEXMATSETCLOSURE
17 #define dmplexgetclosureindices_ DMPLEXGETCLOSUREINDICES
18 #define dmplexrestoreclosureindices_ DMPLEXRESTORECLOSUREINDICES
19 #define dmplexgetjoin_ DMPLEXGETJOIN
20 #define dmplexgetfulljoin_ DMPLEXGETFULLJOIN
21 #define dmplexrestorejoin_ DMPLEXRESTOREJOIN
22 #define dmplexgetmeet_ DMPLEXGETMEET
23 #define dmplexgetfullmeet_ DMPLEXGETFULLMEET
24 #define dmplexrestoremeet_ DMPLEXRESTOREMEET
25 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
26 #define dmplexgetcone_ dmplexgetcone
27 #define dmplexrestorecone_ dmplexrestorecone
28 #define dmplexgetconeorientation_ dmplexgetconeorientation
29 #define dmplexrestoreconeorientation_ dmplexrestoreconeorientation
30 #define dmplexgetsupport_ dmplexgetsupport
31 #define dmplexrestoresupport_ dmplexrestoresupport
32 #define dmplexgettransitiveclosure_ dmplexgettransitiveclosure
33 #define dmplexrestoretransitiveclosure_ dmplexrestoretransitiveclosure
34 #define dmplexvecgetclosure_ dmplexvecgetclosure
35 #define dmplexvecrestoreclosure_ dmplexvecrestoreclosure
36 #define dmplexvecsetclosure_ dmplexvecsetclosure
37 #define dmplexmatsetclosure_ dmplexmatsetclosure
38 #define dmplexgetclosureindices_ dmplexgetclosureindices
39 #define dmplexrestoreclosureindices_ dmplexrestoreclosureindices
40 #define dmplexgetjoin_ dmplexgetjoin
41 #define dmplexgetfulljoin_ dmplexgetfulljoin
42 #define dmplexrestorejoin_ dmplexrestorejoin
43 #define dmplexgetmeet_ dmplexgetmeet
44 #define dmplexgetfullmeet_ dmplexgetfullmeet
45 #define dmplexrestoremeet_ dmplexrestoremeet
46 #endif
47
dmplexgetcone_(DM * dm,PetscInt * p,F90Array1d * ptr,int * ierr PETSC_F90_2PTR_PROTO (ptrd))48 PETSC_EXTERN void dmplexgetcone_(DM *dm, PetscInt *p, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
49 {
50 const PetscInt *v;
51 PetscInt n;
52
53 *ierr = DMPlexGetConeSize(*dm, *p, &n);
54 if (*ierr) return;
55 *ierr = DMPlexGetCone(*dm, *p, &v);
56 if (*ierr) return;
57 *ierr = F90Array1dCreate((void *)v, MPIU_INT, 1, n, ptr PETSC_F90_2PTR_PARAM(ptrd));
58 }
59
dmplexrestorecone_(DM * dm,PetscInt * p,F90Array1d * ptr,int * ierr PETSC_F90_2PTR_PROTO (ptrd))60 PETSC_EXTERN void dmplexrestorecone_(DM *dm, PetscInt *p, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
61 {
62 *ierr = F90Array1dDestroy(ptr, MPIU_INT PETSC_F90_2PTR_PARAM(ptrd));
63 if (*ierr) return;
64 }
65
dmplexgetconeorientation_(DM * dm,PetscInt * p,F90Array1d * ptr,int * ierr PETSC_F90_2PTR_PROTO (ptrd))66 PETSC_EXTERN void dmplexgetconeorientation_(DM *dm, PetscInt *p, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
67 {
68 const PetscInt *v;
69 PetscInt n;
70
71 *ierr = DMPlexGetConeSize(*dm, *p, &n);
72 if (*ierr) return;
73 *ierr = DMPlexGetConeOrientation(*dm, *p, &v);
74 if (*ierr) return;
75 *ierr = F90Array1dCreate((void *)v, MPIU_INT, 1, n, ptr PETSC_F90_2PTR_PARAM(ptrd));
76 }
77
dmplexrestoreconeorientation_(DM * dm,PetscInt * p,F90Array1d * ptr,int * ierr PETSC_F90_2PTR_PROTO (ptrd))78 PETSC_EXTERN void dmplexrestoreconeorientation_(DM *dm, PetscInt *p, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
79 {
80 *ierr = F90Array1dDestroy(ptr, MPIU_INT PETSC_F90_2PTR_PARAM(ptrd));
81 if (*ierr) return;
82 }
83
dmplexgetsupport_(DM * dm,PetscInt * p,F90Array1d * ptr,int * ierr PETSC_F90_2PTR_PROTO (ptrd))84 PETSC_EXTERN void dmplexgetsupport_(DM *dm, PetscInt *p, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
85 {
86 const PetscInt *v;
87 PetscInt n;
88
89 *ierr = DMPlexGetSupportSize(*dm, *p, &n);
90 if (*ierr) return;
91 *ierr = DMPlexGetSupport(*dm, *p, &v);
92 if (*ierr) return;
93 *ierr = F90Array1dCreate((void *)v, MPIU_INT, 1, n, ptr PETSC_F90_2PTR_PARAM(ptrd));
94 }
95
dmplexrestoresupport_(DM * dm,PetscInt * p,F90Array1d * ptr,int * ierr PETSC_F90_2PTR_PROTO (ptrd))96 PETSC_EXTERN void dmplexrestoresupport_(DM *dm, PetscInt *p, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
97 {
98 *ierr = F90Array1dDestroy(ptr, MPIU_INT PETSC_F90_2PTR_PARAM(ptrd));
99 if (*ierr) return;
100 }
101
dmplexgettransitiveclosure_(DM * dm,PetscInt * p,PetscBool * useCone,PetscInt * N,F90Array1d * ptr,int * ierr PETSC_F90_2PTR_PROTO (ptrd))102 PETSC_EXTERN void dmplexgettransitiveclosure_(DM *dm, PetscInt *p, PetscBool *useCone, PetscInt *N, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
103 {
104 PetscInt *v = NULL;
105 PetscInt n;
106
107 CHKFORTRANNULL(N);
108 *ierr = DMPlexGetTransitiveClosure(*dm, *p, *useCone, &n, &v);
109 if (*ierr) return;
110 *ierr = F90Array1dCreate((void *)v, MPIU_INT, 1, n * 2, ptr PETSC_F90_2PTR_PARAM(ptrd));
111 if (N) *N = n;
112 }
113
dmplexrestoretransitiveclosure_(DM * dm,PetscInt * p,PetscBool * useCone,PetscInt * N,F90Array1d * ptr,int * ierr PETSC_F90_2PTR_PROTO (ptrd))114 PETSC_EXTERN void dmplexrestoretransitiveclosure_(DM *dm, PetscInt *p, PetscBool *useCone, PetscInt *N, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
115 {
116 PetscInt *array;
117
118 *ierr = F90Array1dAccess(ptr, MPIU_INT, (void **)&array PETSC_F90_2PTR_PARAM(ptrd));
119 if (*ierr) return;
120 *ierr = DMPlexRestoreTransitiveClosure(*dm, *p, *useCone, NULL, &array);
121 if (*ierr) return;
122 *ierr = F90Array1dDestroy(ptr, MPIU_INT PETSC_F90_2PTR_PARAM(ptrd));
123 if (*ierr) return;
124 }
125
dmplexvecgetclosure_(DM * dm,PetscSection * section,Vec * x,PetscInt * point,PetscInt * N,F90Array1d * ptr,int * ierr PETSC_F90_2PTR_PROTO (ptrd))126 PETSC_EXTERN void dmplexvecgetclosure_(DM *dm, PetscSection *section, Vec *x, PetscInt *point, PetscInt *N, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
127 {
128 PetscScalar *v = NULL;
129 PetscInt n;
130
131 CHKFORTRANNULL(N);
132 *ierr = DMPlexVecGetClosure(*dm, *section, *x, *point, &n, &v);
133 if (*ierr) return;
134 *ierr = F90Array1dCreate((void *)v, MPIU_SCALAR, 1, n, ptr PETSC_F90_2PTR_PARAM(ptrd));
135 if (N) *N = n;
136 }
137
dmplexvecrestoreclosure_(DM * dm,PetscSection * section,Vec * v,PetscInt * point,PetscInt * N,F90Array1d * ptr,int * ierr PETSC_F90_2PTR_PROTO (ptrd))138 PETSC_EXTERN void dmplexvecrestoreclosure_(DM *dm, PetscSection *section, Vec *v, PetscInt *point, PetscInt *N, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
139 {
140 PetscScalar *array;
141
142 *ierr = F90Array1dAccess(ptr, MPIU_SCALAR, (void **)&array PETSC_F90_2PTR_PARAM(ptrd));
143 if (*ierr) return;
144 *ierr = DMPlexVecRestoreClosure(*dm, *section, *v, *point, NULL, &array);
145 if (*ierr) return;
146 *ierr = F90Array1dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
147 if (*ierr) return;
148 }
149
dmplexgetclosureindices_(DM * dm,PetscSection * section,PetscSection * idxSection,PetscInt * point,PetscBool * useClPerm,PetscInt * numIndices,F90Array1d * idxPtr,PetscInt * outOffsets,F90Array1d * valPtr,int * ierr PETSC_F90_2PTR_PROTO (idxPtrd)PETSC_F90_2PTR_PROTO (valPtrd))150 PETSC_EXTERN void dmplexgetclosureindices_(DM *dm, PetscSection *section, PetscSection *idxSection, PetscInt *point, PetscBool *useClPerm, PetscInt *numIndices, F90Array1d *idxPtr, PetscInt *outOffsets, F90Array1d *valPtr, int *ierr PETSC_F90_2PTR_PROTO(idxPtrd) PETSC_F90_2PTR_PROTO(valPtrd))
151 {
152 PetscInt *indices;
153 PetscScalar *values;
154
155 CHKFORTRANNULL(outOffsets);
156 if (FORTRANNULLSCALARPOINTER(valPtr)) *ierr = DMPlexGetClosureIndices(*dm, *section, *idxSection, *point, *useClPerm, numIndices, &indices, outOffsets, NULL);
157 else *ierr = DMPlexGetClosureIndices(*dm, *section, *idxSection, *point, *useClPerm, numIndices, &indices, outOffsets, &values);
158 if (*ierr) return;
159 *ierr = F90Array1dCreate((void *)indices, MPIU_INT, 1, *numIndices, idxPtr PETSC_F90_2PTR_PARAM(idxPtrd));
160 if (*ierr) return;
161 if (FORTRANNULLSCALARPOINTER(valPtr)) *ierr = F90Array1dCreate((void *)values, MPIU_SCALAR, 1, *numIndices, valPtr PETSC_F90_2PTR_PARAM(valPtrd));
162 }
163
dmplexrestoreclosureindices_(DM * dm,PetscSection * section,PetscSection * idxSection,PetscInt * point,PetscBool * useClPerm,PetscInt * numIndices,F90Array1d * idxPtr,PetscInt * outOffsets,F90Array1d * valPtr,int * ierr PETSC_F90_2PTR_PROTO (idxPtrd)PETSC_F90_2PTR_PROTO (valPtrd))164 PETSC_EXTERN void dmplexrestoreclosureindices_(DM *dm, PetscSection *section, PetscSection *idxSection, PetscInt *point, PetscBool *useClPerm, PetscInt *numIndices, F90Array1d *idxPtr, PetscInt *outOffsets, F90Array1d *valPtr, int *ierr PETSC_F90_2PTR_PROTO(idxPtrd) PETSC_F90_2PTR_PROTO(valPtrd))
165 {
166 PetscInt *indices;
167 PetscScalar *values = NULL;
168
169 CHKFORTRANNULL(outOffsets);
170 *ierr = F90Array1dAccess(idxPtr, MPIU_INT, (void **)&indices PETSC_F90_2PTR_PARAM(idxPtrd));
171 if (*ierr) return;
172 if (!FORTRANNULLSCALARPOINTER(valPtr)) {
173 *ierr = F90Array1dAccess(valPtr, MPIU_SCALAR, (void **)&values PETSC_F90_2PTR_PARAM(valPtrd));
174 if (*ierr) return;
175 *ierr = DMPlexRestoreClosureIndices(*dm, *section, *idxSection, *point, *useClPerm, numIndices, &indices, outOffsets, &values);
176 } else *ierr = DMPlexRestoreClosureIndices(*dm, *section, *idxSection, *point, *useClPerm, numIndices, &indices, outOffsets, NULL);
177 if (*ierr) return;
178 *ierr = F90Array1dDestroy(idxPtr, MPIU_INT PETSC_F90_2PTR_PARAM(idxPtrd));
179 if (*ierr) return;
180 if (!FORTRANNULLSCALARPOINTER(valPtr)) *ierr = F90Array1dDestroy(valPtr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(valPtrd));
181 }
182
dmplexgetjoin_(DM * dm,PetscInt * numPoints,PetscInt * points,PetscInt * N,F90Array1d * cptr,int * ierr PETSC_F90_2PTR_PROTO (cptrd))183 PETSC_EXTERN void dmplexgetjoin_(DM *dm, PetscInt *numPoints, PetscInt *points, PetscInt *N, F90Array1d *cptr, int *ierr PETSC_F90_2PTR_PROTO(cptrd))
184 {
185 const PetscInt *coveredPoints;
186 PetscInt n;
187
188 CHKFORTRANNULL(N);
189 *ierr = DMPlexGetJoin(*dm, *numPoints, points, &n, &coveredPoints);
190 if (*ierr) return;
191 *ierr = F90Array1dCreate((void *)coveredPoints, MPIU_INT, 1, n, cptr PETSC_F90_2PTR_PARAM(cptrd));
192 if (N) *N = n;
193 }
194
dmplexgetfulljoin_(DM * dm,PetscInt * numPoints,PetscInt * points,PetscInt * N,F90Array1d * cptr,int * ierr PETSC_F90_2PTR_PROTO (cptrd))195 PETSC_EXTERN void dmplexgetfulljoin_(DM *dm, PetscInt *numPoints, PetscInt *points, PetscInt *N, F90Array1d *cptr, int *ierr PETSC_F90_2PTR_PROTO(cptrd))
196 {
197 const PetscInt *coveredPoints;
198 PetscInt n;
199
200 CHKFORTRANNULL(N);
201 *ierr = DMPlexGetFullJoin(*dm, *numPoints, points, &n, &coveredPoints);
202 if (*ierr) return;
203 *ierr = F90Array1dCreate((void *)coveredPoints, MPIU_INT, 1, n, cptr PETSC_F90_2PTR_PARAM(cptrd));
204 if (N) *N = n;
205 }
206
dmplexrestorejoin_(DM * dm,PetscInt * numPoints,PetscInt * points,PetscInt * N,F90Array1d * cptr,int * ierr PETSC_F90_2PTR_PROTO (cptrd))207 PETSC_EXTERN void dmplexrestorejoin_(DM *dm, PetscInt *numPoints, PetscInt *points, PetscInt *N, F90Array1d *cptr, int *ierr PETSC_F90_2PTR_PROTO(cptrd))
208 {
209 PetscInt *coveredPoints;
210
211 *ierr = F90Array1dAccess(cptr, MPIU_INT, (void **)&coveredPoints PETSC_F90_2PTR_PARAM(cptrd));
212 if (*ierr) return;
213 *ierr = DMPlexRestoreJoin(*dm, 0, NULL, NULL, (const PetscInt **)&coveredPoints);
214 if (*ierr) return;
215 *ierr = F90Array1dDestroy(cptr, MPIU_INT PETSC_F90_2PTR_PARAM(cptrd));
216 if (*ierr) return;
217 }
218
dmplexgetmeet_(DM * dm,PetscInt * numPoints,PetscInt * points,PetscInt * N,F90Array1d * cptr,int * ierr PETSC_F90_2PTR_PROTO (cptrd))219 PETSC_EXTERN void dmplexgetmeet_(DM *dm, PetscInt *numPoints, PetscInt *points, PetscInt *N, F90Array1d *cptr, int *ierr PETSC_F90_2PTR_PROTO(cptrd))
220 {
221 const PetscInt *coveredPoints;
222 PetscInt n;
223
224 CHKFORTRANNULL(N);
225 *ierr = DMPlexGetMeet(*dm, *numPoints, points, &n, &coveredPoints);
226 if (*ierr) return;
227 *ierr = F90Array1dCreate((void *)coveredPoints, MPIU_INT, 1, n, cptr PETSC_F90_2PTR_PARAM(cptrd));
228 if (N) *N = n;
229 }
230
dmplexgetfullmeet_(DM * dm,PetscInt * numPoints,PetscInt * points,PetscInt * N,F90Array1d * cptr,int * ierr PETSC_F90_2PTR_PROTO (cptrd))231 PETSC_EXTERN void dmplexgetfullmeet_(DM *dm, PetscInt *numPoints, PetscInt *points, PetscInt *N, F90Array1d *cptr, int *ierr PETSC_F90_2PTR_PROTO(cptrd))
232 {
233 const PetscInt *coveredPoints;
234 PetscInt n;
235
236 CHKFORTRANNULL(N);
237 if (*ierr) return;
238 *ierr = DMPlexGetFullMeet(*dm, *numPoints, points, &n, &coveredPoints);
239 if (*ierr) return;
240 *ierr = F90Array1dCreate((void *)coveredPoints, MPIU_INT, 1, n, cptr PETSC_F90_2PTR_PARAM(cptrd));
241 if (N) *N = n;
242 }
243
dmplexrestoremeet_(DM * dm,PetscInt * numPoints,PetscInt * points,PetscInt * N,F90Array1d * cptr,int * ierr PETSC_F90_2PTR_PROTO (cptrd))244 PETSC_EXTERN void dmplexrestoremeet_(DM *dm, PetscInt *numPoints, PetscInt *points, PetscInt *N, F90Array1d *cptr, int *ierr PETSC_F90_2PTR_PROTO(cptrd))
245 {
246 PetscInt *coveredPoints;
247
248 *ierr = F90Array1dAccess(cptr, MPIU_INT, (void **)&coveredPoints PETSC_F90_2PTR_PARAM(cptrd));
249 if (*ierr) return;
250 *ierr = DMPlexRestoreMeet(*dm, 0, NULL, NULL, (const PetscInt **)&coveredPoints);
251 if (*ierr) return;
252 *ierr = F90Array1dDestroy(cptr, MPIU_INT PETSC_F90_2PTR_PARAM(cptrd));
253 if (*ierr) return;
254 }
255