1 #include <petsc/private/f90impl.h>
2 #include <petscdmda.h>
3
4 #if defined(PETSC_HAVE_FORTRAN_CAPS)
5 #define dmdagetlocalinfof90_ DMDAGETLOCALINFOF90
6 #define dmdavecgetarrayf901_ DMDAVECGETARRAYF901
7 #define dmdavecrestorearrayf901_ DMDAVECRESTOREARRAYF901
8 #define dmdavecgetarrayf902_ DMDAVECGETARRAYF902
9 #define dmdavecrestorearrayf902_ DMDAVECRESTOREARRAYF902
10 #define dmdavecgetarrayf903_ DMDAVECGETARRAYF903
11 #define dmdavecrestorearrayf903_ DMDAVECRESTOREARRAYF903
12 #define dmdavecgetarrayf904_ DMDAVECGETARRAYF904
13 #define dmdavecrestorearrayf904_ DMDAVECRESTOREARRAYF904
14 #define dmdavecgetarrayreadf901_ DMDAVECGETARRAYREADF901
15 #define dmdavecrestorearrayreadf901_ DMDAVECRESTOREARRAYREADF901
16 #define dmdavecgetarrayreadf902_ DMDAVECGETARRAYREADF902
17 #define dmdavecrestorearrayreadf902_ DMDAVECRESTOREARRAYREADF902
18 #define dmdavecgetarrayreadf903_ DMDAVECGETARRAYREADF903
19 #define dmdavecrestorearrayreadf903_ DMDAVECRESTOREARRAYREADF903
20 #define dmdavecgetarrayreadf904_ DMDAVECGETARRAYREADF904
21 #define dmdavecrestorearrayreadf904_ DMDAVECRESTOREARRAYREADF904
22 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
23 #define dmdagetlocalinfof90_ dmdagetlocalinfof90
24 #define dmdavecgetarrayf901_ dmdavecgetarrayf901
25 #define dmdavecrestorearrayf901_ dmdavecrestorearrayf901
26 #define dmdavecgetarrayf902_ dmdavecgetarrayf902
27 #define dmdavecrestorearrayf902_ dmdavecrestorearrayf902
28 #define dmdavecgetarrayf903_ dmdavecgetarrayf903
29 #define dmdavecrestorearrayf903_ dmdavecrestorearrayf903
30 #define dmdavecgetarrayf904_ dmdavecgetarrayf904
31 #define dmdavecrestorearrayf904_ dmdavecrestorearrayf904
32 #define dmdavecgetarrayreadf901_ dmdavecgetarrayreadf901
33 #define dmdavecrestorearrayreadf901_ dmdavecrestorearrayreadf901
34 #define dmdavecgetarrayreadf902_ dmdavecgetarrayreadf902
35 #define dmdavecrestorearrayreadf902_ dmdavecrestorearrayreadf902
36 #define dmdavecgetarrayreadf903_ dmdavecgetarrayreadf903
37 #define dmdavecrestorearrayreadf903_ dmdavecrestorearrayreadf903
38 #define dmdavecgetarrayreadf904_ dmdavecgetarrayreadf904
39 #define dmdavecrestorearrayreadf904_ dmdavecrestorearrayreadf904
40 #endif
41
dmdagetlocalinfof90_(DM * da,DMDALocalInfo * info,PetscErrorCode * ierr)42 PETSC_EXTERN void dmdagetlocalinfof90_(DM *da,DMDALocalInfo *info,PetscErrorCode *ierr)
43 {
44 *ierr = DMDAGetLocalInfo(*da,info);
45 }
46
dmdavecgetarrayf901_(DM * da,Vec * v,F90Array1d * a,PetscErrorCode * ierr PETSC_F90_2PTR_PROTO (ptrd))47 PETSC_EXTERN void dmdavecgetarrayf901_(DM *da,Vec *v,F90Array1d *a,PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
48 {
49 PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
50 PetscScalar *aa;
51
52 *ierr = DMDAGetCorners(*da,&xs,&ys,&zs,&xm,&ym,&zm);if (*ierr) return;
53 *ierr = DMDAGetGhostCorners(*da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);if (*ierr) return;
54 *ierr = DMDAGetInfo(*da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);if (*ierr) return;
55
56 /* Handle case where user passes in global vector as opposed to local */
57 *ierr = VecGetLocalSize(*v,&N);if (*ierr) return;
58 if (N == xm*ym*zm*dof) {
59 gxm = xm;
60 gym = ym;
61 gzm = zm;
62 gxs = xs;
63 gys = ys;
64 gzs = zs;
65 } else if (N != gxm*gym*gzm*dof) {
66 *ierr = PETSC_ERR_ARG_INCOMP; return;
67 }
68 *ierr = VecGetArray(*v,&aa);if (*ierr) return;
69 *ierr = F90Array1dCreate(aa,MPIU_SCALAR,gxs,gxm,a PETSC_F90_2PTR_PARAM(ptrd));if (*ierr) return;
70 }
71
dmdavecrestorearrayf901_(DM * da,Vec * v,F90Array1d * a,PetscErrorCode * ierr PETSC_F90_2PTR_PROTO (ptrd))72 PETSC_EXTERN void dmdavecrestorearrayf901_(DM *da,Vec *v,F90Array1d *a,PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
73 {
74 PetscScalar *fa;
75 *ierr = F90Array1dAccess(a,MPIU_SCALAR,(void**)&fa PETSC_F90_2PTR_PARAM(ptrd));
76 *ierr = VecRestoreArray(*v,&fa);if (*ierr) return;
77 *ierr = F90Array1dDestroy(a,MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
78 }
79
dmdavecgetarrayf902_(DM * da,Vec * v,F90Array2d * a,PetscErrorCode * ierr PETSC_F90_2PTR_PROTO (ptrd))80 PETSC_EXTERN void dmdavecgetarrayf902_(DM *da,Vec *v,F90Array2d *a,PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
81 {
82 PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
83 PetscScalar *aa;
84
85 *ierr = DMDAGetCorners(*da,&xs,&ys,&zs,&xm,&ym,&zm);if (*ierr) return;
86 *ierr = DMDAGetGhostCorners(*da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);if (*ierr) return;
87 *ierr = DMDAGetInfo(*da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);if (*ierr) return;
88
89 /* Handle case where user passes in global vector as opposed to local */
90 *ierr = VecGetLocalSize(*v,&N);if (*ierr) return;
91 if (N == xm*ym*zm*dof) {
92 gxm = xm;
93 gym = ym;
94 gzm = zm;
95 gxs = xs;
96 gys = ys;
97 gzs = zs;
98 } else if (N != gxm*gym*gzm*dof) {
99 *ierr = PETSC_ERR_ARG_INCOMP; return;
100 }
101 if (dim == 1) {
102 gys = gxs;
103 gym = gxm;
104 gxs = 0;
105 gxm = dof;
106 }
107 *ierr = VecGetArray(*v,&aa);if (*ierr) return;
108 *ierr = F90Array2dCreate(aa,MPIU_SCALAR,gxs,gxm,gys,gym,a PETSC_F90_2PTR_PARAM(ptrd));if (*ierr) return;
109 }
110
dmdavecrestorearrayf902_(DM * da,Vec * v,F90Array2d * a,PetscErrorCode * ierr PETSC_F90_2PTR_PROTO (ptrd))111 PETSC_EXTERN void dmdavecrestorearrayf902_(DM *da,Vec *v,F90Array2d *a,PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
112 {
113 PetscScalar *fa;
114 *ierr = F90Array2dAccess(a,MPIU_SCALAR,(void**)&fa PETSC_F90_2PTR_PARAM(ptrd));
115 *ierr = VecRestoreArray(*v,&fa);if (*ierr) return;
116 *ierr = F90Array2dDestroy(a,MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
117 }
118
dmdavecgetarrayf903_(DM * da,Vec * v,F90Array3d * a,PetscErrorCode * ierr PETSC_F90_2PTR_PROTO (ptrd))119 PETSC_EXTERN void dmdavecgetarrayf903_(DM *da,Vec *v,F90Array3d *a,PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
120 {
121 PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
122 PetscScalar *aa;
123
124 *ierr = DMDAGetCorners(*da,&xs,&ys,&zs,&xm,&ym,&zm);if (*ierr) return;
125 *ierr = DMDAGetGhostCorners(*da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);if (*ierr) return;
126 *ierr = DMDAGetInfo(*da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);if (*ierr) return;
127
128 /* Handle case where user passes in global vector as opposed to local */
129 *ierr = VecGetLocalSize(*v,&N);if (*ierr) return;
130 if (N == xm*ym*zm*dof) {
131 gxm = xm;
132 gym = ym;
133 gzm = zm;
134 gxs = xs;
135 gys = ys;
136 gzs = zs;
137 } else if (N != gxm*gym*gzm*dof) {
138 *ierr = PETSC_ERR_ARG_INCOMP; return;
139 }
140 if (dim == 2) {
141 gzs = gys;
142 gzm = gym;
143 gys = gxs;
144 gym = gxm;
145 gxs = 0;
146 gxm = dof;
147 }
148 *ierr = VecGetArray(*v,&aa);if (*ierr) return;
149 *ierr = F90Array3dCreate(aa,MPIU_SCALAR,gxs,gxm,gys,gym,gzs,gzm,a PETSC_F90_2PTR_PARAM(ptrd));if (*ierr) return;
150 }
151
dmdavecrestorearrayf903_(DM * da,Vec * v,F90Array3d * a,PetscErrorCode * ierr PETSC_F90_2PTR_PROTO (ptrd))152 PETSC_EXTERN void dmdavecrestorearrayf903_(DM *da,Vec *v,F90Array3d *a,PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
153 {
154 PetscScalar *fa;
155 *ierr = F90Array3dAccess(a,MPIU_SCALAR,(void**)&fa PETSC_F90_2PTR_PARAM(ptrd));
156 *ierr = VecRestoreArray(*v,&fa);if (*ierr) return;
157 *ierr = F90Array3dDestroy(a,MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
158 }
159
dmdavecgetarrayf904_(DM * da,Vec * v,F90Array4d * a,PetscErrorCode * ierr PETSC_F90_2PTR_PROTO (ptrd))160 PETSC_EXTERN void dmdavecgetarrayf904_(DM *da,Vec *v,F90Array4d *a,PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
161 {
162 PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof,zero = 0;
163 PetscScalar *aa;
164
165 *ierr = DMDAGetCorners(*da,&xs,&ys,&zs,&xm,&ym,&zm);if (*ierr) return;
166 *ierr = DMDAGetGhostCorners(*da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);if (*ierr) return;
167 *ierr = DMDAGetInfo(*da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);if (*ierr) return;
168
169 /* Handle case where user passes in global vector as opposed to local */
170 *ierr = VecGetLocalSize(*v,&N);if (*ierr) return;
171 if (N == xm*ym*zm*dof) {
172 gxm = xm;
173 gym = ym;
174 gzm = zm;
175 gxs = xs;
176 gys = ys;
177 gzs = zs;
178 } else if (N != gxm*gym*gzm*dof) {
179 *ierr = PETSC_ERR_ARG_INCOMP; return;
180 }
181 *ierr = VecGetArray(*v,&aa);if (*ierr) return;
182 *ierr = F90Array4dCreate(aa,MPIU_SCALAR,zero,dof,gxs,gxm,gys,gym,gzs,gzm,a PETSC_F90_2PTR_PARAM(ptrd));if (*ierr) return;
183 }
184
dmdavecrestorearrayf904_(DM * da,Vec * v,F90Array4d * a,PetscErrorCode * ierr PETSC_F90_2PTR_PROTO (ptrd))185 PETSC_EXTERN void dmdavecrestorearrayf904_(DM *da,Vec *v,F90Array4d *a,PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
186 {
187 PetscScalar *fa;
188 /*
189 F90Array4dAccess is not implemented, so the following call would fail
190 */
191 *ierr = F90Array4dAccess(a,MPIU_SCALAR,(void**)&fa PETSC_F90_2PTR_PARAM(ptrd));
192 *ierr = VecRestoreArray(*v,&fa);if (*ierr) return;
193 *ierr = F90Array4dDestroy(a,MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
194 }
195
dmdavecgetarrayreadf901_(DM * da,Vec * v,F90Array1d * a,PetscErrorCode * ierr PETSC_F90_2PTR_PROTO (ptrd))196 PETSC_EXTERN void dmdavecgetarrayreadf901_(DM *da,Vec *v,F90Array1d *a,PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
197 {
198 PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
199 const PetscScalar *aa;
200
201 *ierr = DMDAGetCorners(*da,&xs,&ys,&zs,&xm,&ym,&zm);if (*ierr) return;
202 *ierr = DMDAGetGhostCorners(*da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);if (*ierr) return;
203 *ierr = DMDAGetInfo(*da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);if (*ierr) return;
204
205 /* Handle case where user passes in global vector as opposed to local */
206 *ierr = VecGetLocalSize(*v,&N);if (*ierr) return;
207 if (N == xm*ym*zm*dof) {
208 gxm = xm;
209 gym = ym;
210 gzm = zm;
211 gxs = xs;
212 gys = ys;
213 gzs = zs;
214 } else if (N != gxm*gym*gzm*dof) {
215 *ierr = PETSC_ERR_ARG_INCOMP; return;
216 }
217 *ierr = VecGetArrayRead(*v,&aa);if (*ierr) return;
218 *ierr = F90Array1dCreate((void*)aa,MPIU_SCALAR,gxs,gxm,a PETSC_F90_2PTR_PARAM(ptrd));if (*ierr) return;
219 }
220
dmdavecrestorearrayreadf901_(DM * da,Vec * v,F90Array1d * a,PetscErrorCode * ierr PETSC_F90_2PTR_PROTO (ptrd))221 PETSC_EXTERN void dmdavecrestorearrayreadf901_(DM *da,Vec *v,F90Array1d *a,PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
222 {
223 const PetscScalar *fa;
224 *ierr = F90Array1dAccess(a,MPIU_SCALAR,(void**)&fa PETSC_F90_2PTR_PARAM(ptrd));
225 *ierr = VecRestoreArrayRead(*v,&fa);if (*ierr) return;
226 *ierr = F90Array1dDestroy(a,MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
227 }
228
dmdavecgetarrayreadf902_(DM * da,Vec * v,F90Array2d * a,PetscErrorCode * ierr PETSC_F90_2PTR_PROTO (ptrd))229 PETSC_EXTERN void dmdavecgetarrayreadf902_(DM *da,Vec *v,F90Array2d *a,PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
230 {
231 PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
232 const PetscScalar *aa;
233
234 *ierr = DMDAGetCorners(*da,&xs,&ys,&zs,&xm,&ym,&zm);if (*ierr) return;
235 *ierr = DMDAGetGhostCorners(*da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);if (*ierr) return;
236 *ierr = DMDAGetInfo(*da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);if (*ierr) return;
237
238 /* Handle case where user passes in global vector as opposed to local */
239 *ierr = VecGetLocalSize(*v,&N);if (*ierr) return;
240 if (N == xm*ym*zm*dof) {
241 gxm = xm;
242 gym = ym;
243 gzm = zm;
244 gxs = xs;
245 gys = ys;
246 gzs = zs;
247 } else if (N != gxm*gym*gzm*dof) {
248 *ierr = PETSC_ERR_ARG_INCOMP; return;
249 }
250 if (dim == 1) {
251 gys = gxs;
252 gym = gxm;
253 gxs = 0;
254 gxm = dof;
255 }
256 *ierr = VecGetArrayRead(*v,&aa);if (*ierr) return;
257 *ierr = F90Array2dCreate((void*)aa,MPIU_SCALAR,gxs,gxm,gys,gym,a PETSC_F90_2PTR_PARAM(ptrd));if (*ierr) return;
258 }
259
dmdavecrestorearrayreadf902_(DM * da,Vec * v,F90Array2d * a,PetscErrorCode * ierr PETSC_F90_2PTR_PROTO (ptrd))260 PETSC_EXTERN void dmdavecrestorearrayreadf902_(DM *da,Vec *v,F90Array2d *a,PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
261 {
262 const PetscScalar *fa;
263 *ierr = F90Array2dAccess(a,MPIU_SCALAR,(void**)&fa PETSC_F90_2PTR_PARAM(ptrd));
264 *ierr = VecRestoreArrayRead(*v,&fa);if (*ierr) return;
265 *ierr = F90Array2dDestroy(a,MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
266 }
267
dmdavecgetarrayreadf903_(DM * da,Vec * v,F90Array3d * a,PetscErrorCode * ierr PETSC_F90_2PTR_PROTO (ptrd))268 PETSC_EXTERN void dmdavecgetarrayreadf903_(DM *da,Vec *v,F90Array3d *a,PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
269 {
270 PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
271 const PetscScalar *aa;
272
273 *ierr = DMDAGetCorners(*da,&xs,&ys,&zs,&xm,&ym,&zm);if (*ierr) return;
274 *ierr = DMDAGetGhostCorners(*da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);if (*ierr) return;
275 *ierr = DMDAGetInfo(*da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);if (*ierr) return;
276
277 /* Handle case where user passes in global vector as opposed to local */
278 *ierr = VecGetLocalSize(*v,&N);if (*ierr) return;
279 if (N == xm*ym*zm*dof) {
280 gxm = xm;
281 gym = ym;
282 gzm = zm;
283 gxs = xs;
284 gys = ys;
285 gzs = zs;
286 } else if (N != gxm*gym*gzm*dof) {
287 *ierr = PETSC_ERR_ARG_INCOMP; return;
288 }
289 if (dim == 2) {
290 gzs = gys;
291 gzm = gym;
292 gys = gxs;
293 gym = gxm;
294 gxs = 0;
295 gxm = dof;
296 }
297 *ierr = VecGetArrayRead(*v,&aa);if (*ierr) return;
298 *ierr = F90Array3dCreate((void*)aa,MPIU_SCALAR,gxs,gxm,gys,gym,gzs,gzm,a PETSC_F90_2PTR_PARAM(ptrd));if (*ierr) return;
299 }
300
dmdavecrestorearrayreadf903_(DM * da,Vec * v,F90Array3d * a,PetscErrorCode * ierr PETSC_F90_2PTR_PROTO (ptrd))301 PETSC_EXTERN void dmdavecrestorearrayreadf903_(DM *da,Vec *v,F90Array3d *a,PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
302 {
303 const PetscScalar *fa;
304 *ierr = F90Array3dAccess(a,MPIU_SCALAR,(void**)&fa PETSC_F90_2PTR_PARAM(ptrd));
305 *ierr = VecRestoreArrayRead(*v,&fa);if (*ierr) return;
306 *ierr = F90Array3dDestroy(a,MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
307 }
308
dmdavecgetarrayreadf904_(DM * da,Vec * v,F90Array4d * a,PetscErrorCode * ierr PETSC_F90_2PTR_PROTO (ptrd))309 PETSC_EXTERN void dmdavecgetarrayreadf904_(DM *da,Vec *v,F90Array4d *a,PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
310 {
311 PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof,zero = 0;
312 const PetscScalar *aa;
313
314 *ierr = DMDAGetCorners(*da,&xs,&ys,&zs,&xm,&ym,&zm);if (*ierr) return;
315 *ierr = DMDAGetGhostCorners(*da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);if (*ierr) return;
316 *ierr = DMDAGetInfo(*da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);if (*ierr) return;
317
318 /* Handle case where user passes in global vector as opposed to local */
319 *ierr = VecGetLocalSize(*v,&N);if (*ierr) return;
320 if (N == xm*ym*zm*dof) {
321 gxm = xm;
322 gym = ym;
323 gzm = zm;
324 gxs = xs;
325 gys = ys;
326 gzs = zs;
327 } else if (N != gxm*gym*gzm*dof) {
328 *ierr = PETSC_ERR_ARG_INCOMP; return;
329 }
330 *ierr = VecGetArrayRead(*v,&aa);if (*ierr) return;
331 *ierr = F90Array4dCreate((void*)aa,MPIU_SCALAR,zero,dof,gxs,gxm,gys,gym,gzs,gzm,a PETSC_F90_2PTR_PARAM(ptrd));if (*ierr) return;
332 }
333
dmdavecrestorearrayreadf904_(DM * da,Vec * v,F90Array4d * a,PetscErrorCode * ierr PETSC_F90_2PTR_PROTO (ptrd))334 PETSC_EXTERN void dmdavecrestorearrayreadf904_(DM *da,Vec *v,F90Array4d *a,PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
335 {
336 const PetscScalar *fa;
337 /*
338 F90Array4dAccess is not implemented, so the following call would fail
339 */
340 *ierr = F90Array4dAccess(a,MPIU_SCALAR,(void**)&fa PETSC_F90_2PTR_PARAM(ptrd));
341 *ierr = VecRestoreArrayRead(*v,&fa);if (*ierr) return;
342 *ierr = F90Array4dDestroy(a,MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
343 }
344