1 #include <petsc/private/ftnimpl.h> 2 #include <petscdmda.h> 3 4 #if defined(PETSC_HAVE_FORTRAN_CAPS) 5 #define dmdavecgetarrayf901_ DMDAVECGETARRAYF901 6 #define dmdavecrestorearrayf901_ DMDAVECRESTOREARRAYF901 7 #define dmdavecgetarrayf902_ DMDAVECGETARRAYF902 8 #define dmdavecrestorearrayf902_ DMDAVECRESTOREARRAYF902 9 #define dmdavecgetarrayf903_ DMDAVECGETARRAYF903 10 #define dmdavecrestorearrayf903_ DMDAVECRESTOREARRAYF903 11 #define dmdavecgetarrayf904_ DMDAVECGETARRAYF904 12 #define dmdavecrestorearrayf904_ DMDAVECRESTOREARRAYF904 13 #define dmdavecgetarrayreadf901_ DMDAVECGETARRAYREADF901 14 #define dmdavecrestorearrayreadf901_ DMDAVECRESTOREARRAYREADF901 15 #define dmdavecgetarrayreadf902_ DMDAVECGETARRAYREADF902 16 #define dmdavecrestorearrayreadf902_ DMDAVECRESTOREARRAYREADF902 17 #define dmdavecgetarrayreadf903_ DMDAVECGETARRAYREADF903 18 #define dmdavecrestorearrayreadf903_ DMDAVECRESTOREARRAYREADF903 19 #define dmdavecgetarrayreadf904_ DMDAVECGETARRAYREADF904 20 #define dmdavecrestorearrayreadf904_ DMDAVECRESTOREARRAYREADF904 21 #define dmdagetelements_ DMDAGETELEMENTS 22 #define dmdarestoreelements_ DMDARESTOREELEMENTS 23 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE) 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 #define dmdagetelements_ dmdagetelements 41 #define dmdarestoreelements_ dmdarestoreelements 42 #endif 43 44 PETSC_EXTERN void dmdagetelements_(DM *dm, PetscInt *nel, PetscInt *nen, F90Array1d *e, int *ierr PETSC_F90_2PTR_PROTO(ptrd)) 45 { 46 const PetscInt *fa; 47 48 if (!e) { 49 *ierr = PetscError(((PetscObject)e)->comm, __LINE__, PETSC_FUNCTION_NAME, __FILE__, PETSC_ERR_ARG_BADPTR, PETSC_ERROR_INITIAL, "e==NULL, maybe #include <petsc/finclude/petscvec.h> is missing?"); 50 return; 51 } 52 *ierr = DMDAGetElements(*dm, nel, nen, &fa); 53 if (*ierr) return; 54 *ierr = F90Array1dCreate((PetscInt *)fa, MPIU_INT, 1, (*nel) * (*nen), e PETSC_F90_2PTR_PARAM(ptrd)); 55 } 56 57 PETSC_EXTERN void dmdarestoreelements_(DM *dm, PetscInt *nel, PetscInt *nen, F90Array1d *e, int *ierr PETSC_F90_2PTR_PROTO(ptrd)) 58 { 59 if (!e) { 60 *ierr = PetscError(((PetscObject)e)->comm, __LINE__, PETSC_FUNCTION_NAME, __FILE__, PETSC_ERR_ARG_BADPTR, PETSC_ERROR_INITIAL, "e==NULL, maybe #include <petsc/finclude/petscvec.h> is missing?"); 61 return; 62 } 63 *ierr = F90Array1dDestroy(e, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd)); 64 } 65 66 PETSC_EXTERN void dmdavecgetarrayf901_(DM *da, Vec *v, F90Array1d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 67 { 68 PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof; 69 PetscScalar *aa; 70 71 *ierr = DMDAGetCorners(*da, &xs, &ys, &zs, &xm, &ym, &zm); 72 if (*ierr) return; 73 *ierr = DMDAGetGhostCorners(*da, &gxs, &gys, &gzs, &gxm, &gym, &gzm); 74 if (*ierr) return; 75 *ierr = DMDAGetInfo(*da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL); 76 if (*ierr) return; 77 78 /* Handle case where user passes in global vector as opposed to local */ 79 *ierr = VecGetLocalSize(*v, &N); 80 if (*ierr) return; 81 if (N == xm * ym * zm * dof) { 82 gxm = xm; 83 gym = ym; 84 gzm = zm; 85 gxs = xs; 86 gys = ys; 87 gzs = zs; 88 } else if (N != gxm * gym * gzm * dof) { 89 *ierr = PETSC_ERR_ARG_INCOMP; 90 return; 91 } 92 *ierr = VecGetArray(*v, &aa); 93 if (*ierr) return; 94 *ierr = F90Array1dCreate(aa, MPIU_SCALAR, gxs, gxm, a PETSC_F90_2PTR_PARAM(ptrd)); 95 if (*ierr) return; 96 } 97 98 PETSC_EXTERN void dmdavecrestorearrayf901_(DM *da, Vec *v, F90Array1d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 99 { 100 PetscScalar *fa; 101 *ierr = F90Array1dAccess(a, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd)); 102 *ierr = VecRestoreArray(*v, &fa); 103 if (*ierr) return; 104 *ierr = F90Array1dDestroy(a, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd)); 105 } 106 107 PETSC_EXTERN void dmdavecgetarrayf902_(DM *da, Vec *v, F90Array2d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 108 { 109 PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof; 110 PetscScalar *aa; 111 112 *ierr = DMDAGetCorners(*da, &xs, &ys, &zs, &xm, &ym, &zm); 113 if (*ierr) return; 114 *ierr = DMDAGetGhostCorners(*da, &gxs, &gys, &gzs, &gxm, &gym, &gzm); 115 if (*ierr) return; 116 *ierr = DMDAGetInfo(*da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL); 117 if (*ierr) return; 118 119 /* Handle case where user passes in global vector as opposed to local */ 120 *ierr = VecGetLocalSize(*v, &N); 121 if (*ierr) return; 122 if (N == xm * ym * zm * dof) { 123 gxm = xm; 124 gym = ym; 125 gzm = zm; 126 gxs = xs; 127 gys = ys; 128 gzs = zs; 129 } else if (N != gxm * gym * gzm * dof) { 130 *ierr = PETSC_ERR_ARG_INCOMP; 131 return; 132 } 133 if (dim == 1) { 134 gys = gxs; 135 gym = gxm; 136 gxs = 0; 137 gxm = dof; 138 } 139 *ierr = VecGetArray(*v, &aa); 140 if (*ierr) return; 141 *ierr = F90Array2dCreate(aa, MPIU_SCALAR, gxs, gxm, gys, gym, a PETSC_F90_2PTR_PARAM(ptrd)); 142 if (*ierr) return; 143 } 144 145 PETSC_EXTERN void dmdavecrestorearrayf902_(DM *da, Vec *v, F90Array2d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 146 { 147 PetscScalar *fa; 148 *ierr = F90Array2dAccess(a, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd)); 149 *ierr = VecRestoreArray(*v, &fa); 150 if (*ierr) return; 151 *ierr = F90Array2dDestroy(a, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd)); 152 } 153 154 PETSC_EXTERN void dmdavecgetarrayf903_(DM *da, Vec *v, F90Array3d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 155 { 156 PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof; 157 PetscScalar *aa; 158 159 *ierr = DMDAGetCorners(*da, &xs, &ys, &zs, &xm, &ym, &zm); 160 if (*ierr) return; 161 *ierr = DMDAGetGhostCorners(*da, &gxs, &gys, &gzs, &gxm, &gym, &gzm); 162 if (*ierr) return; 163 *ierr = DMDAGetInfo(*da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL); 164 if (*ierr) return; 165 166 /* Handle case where user passes in global vector as opposed to local */ 167 *ierr = VecGetLocalSize(*v, &N); 168 if (*ierr) return; 169 if (N == xm * ym * zm * dof) { 170 gxm = xm; 171 gym = ym; 172 gzm = zm; 173 gxs = xs; 174 gys = ys; 175 gzs = zs; 176 } else if (N != gxm * gym * gzm * dof) { 177 *ierr = PETSC_ERR_ARG_INCOMP; 178 return; 179 } 180 if (dim == 2) { 181 gzs = gys; 182 gzm = gym; 183 gys = gxs; 184 gym = gxm; 185 gxs = 0; 186 gxm = dof; 187 } 188 *ierr = VecGetArray(*v, &aa); 189 if (*ierr) return; 190 *ierr = F90Array3dCreate(aa, MPIU_SCALAR, gxs, gxm, gys, gym, gzs, gzm, a PETSC_F90_2PTR_PARAM(ptrd)); 191 if (*ierr) return; 192 } 193 194 PETSC_EXTERN void dmdavecrestorearrayf903_(DM *da, Vec *v, F90Array3d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 195 { 196 PetscScalar *fa; 197 *ierr = F90Array3dAccess(a, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd)); 198 *ierr = VecRestoreArray(*v, &fa); 199 if (*ierr) return; 200 *ierr = F90Array3dDestroy(a, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd)); 201 } 202 203 PETSC_EXTERN void dmdavecgetarrayf904_(DM *da, Vec *v, F90Array4d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 204 { 205 PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof, zero = 0; 206 PetscScalar *aa; 207 208 *ierr = DMDAGetCorners(*da, &xs, &ys, &zs, &xm, &ym, &zm); 209 if (*ierr) return; 210 *ierr = DMDAGetGhostCorners(*da, &gxs, &gys, &gzs, &gxm, &gym, &gzm); 211 if (*ierr) return; 212 *ierr = DMDAGetInfo(*da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL); 213 if (*ierr) return; 214 215 /* Handle case where user passes in global vector as opposed to local */ 216 *ierr = VecGetLocalSize(*v, &N); 217 if (*ierr) return; 218 if (N == xm * ym * zm * dof) { 219 gxm = xm; 220 gym = ym; 221 gzm = zm; 222 gxs = xs; 223 gys = ys; 224 gzs = zs; 225 } else if (N != gxm * gym * gzm * dof) { 226 *ierr = PETSC_ERR_ARG_INCOMP; 227 return; 228 } 229 *ierr = VecGetArray(*v, &aa); 230 if (*ierr) return; 231 *ierr = F90Array4dCreate(aa, MPIU_SCALAR, zero, dof, gxs, gxm, gys, gym, gzs, gzm, a PETSC_F90_2PTR_PARAM(ptrd)); 232 if (*ierr) return; 233 } 234 235 PETSC_EXTERN void dmdavecrestorearrayf904_(DM *da, Vec *v, F90Array4d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 236 { 237 PetscScalar *fa; 238 /* 239 F90Array4dAccess is not implemented, so the following call would fail 240 */ 241 *ierr = F90Array4dAccess(a, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd)); 242 *ierr = VecRestoreArray(*v, &fa); 243 if (*ierr) return; 244 *ierr = F90Array4dDestroy(a, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd)); 245 } 246 247 PETSC_EXTERN void dmdavecgetarrayreadf901_(DM *da, Vec *v, F90Array1d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 248 { 249 PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof; 250 const PetscScalar *aa; 251 252 *ierr = DMDAGetCorners(*da, &xs, &ys, &zs, &xm, &ym, &zm); 253 if (*ierr) return; 254 *ierr = DMDAGetGhostCorners(*da, &gxs, &gys, &gzs, &gxm, &gym, &gzm); 255 if (*ierr) return; 256 *ierr = DMDAGetInfo(*da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL); 257 if (*ierr) return; 258 259 /* Handle case where user passes in global vector as opposed to local */ 260 *ierr = VecGetLocalSize(*v, &N); 261 if (*ierr) return; 262 if (N == xm * ym * zm * dof) { 263 gxm = xm; 264 gym = ym; 265 gzm = zm; 266 gxs = xs; 267 gys = ys; 268 gzs = zs; 269 } else if (N != gxm * gym * gzm * dof) { 270 *ierr = PETSC_ERR_ARG_INCOMP; 271 return; 272 } 273 *ierr = VecGetArrayRead(*v, &aa); 274 if (*ierr) return; 275 *ierr = F90Array1dCreate((void *)aa, MPIU_SCALAR, gxs, gxm, a PETSC_F90_2PTR_PARAM(ptrd)); 276 if (*ierr) return; 277 } 278 279 PETSC_EXTERN void dmdavecrestorearrayreadf901_(DM *da, Vec *v, F90Array1d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 280 { 281 const PetscScalar *fa; 282 *ierr = F90Array1dAccess(a, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd)); 283 *ierr = VecRestoreArrayRead(*v, &fa); 284 if (*ierr) return; 285 *ierr = F90Array1dDestroy(a, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd)); 286 } 287 288 PETSC_EXTERN void dmdavecgetarrayreadf902_(DM *da, Vec *v, F90Array2d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 289 { 290 PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof; 291 const PetscScalar *aa; 292 293 *ierr = DMDAGetCorners(*da, &xs, &ys, &zs, &xm, &ym, &zm); 294 if (*ierr) return; 295 *ierr = DMDAGetGhostCorners(*da, &gxs, &gys, &gzs, &gxm, &gym, &gzm); 296 if (*ierr) return; 297 *ierr = DMDAGetInfo(*da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL); 298 if (*ierr) return; 299 300 /* Handle case where user passes in global vector as opposed to local */ 301 *ierr = VecGetLocalSize(*v, &N); 302 if (*ierr) return; 303 if (N == xm * ym * zm * dof) { 304 gxm = xm; 305 gym = ym; 306 gzm = zm; 307 gxs = xs; 308 gys = ys; 309 gzs = zs; 310 } else if (N != gxm * gym * gzm * dof) { 311 *ierr = PETSC_ERR_ARG_INCOMP; 312 return; 313 } 314 if (dim == 1) { 315 gys = gxs; 316 gym = gxm; 317 gxs = 0; 318 gxm = dof; 319 } 320 *ierr = VecGetArrayRead(*v, &aa); 321 if (*ierr) return; 322 *ierr = F90Array2dCreate((void *)aa, MPIU_SCALAR, gxs, gxm, gys, gym, a PETSC_F90_2PTR_PARAM(ptrd)); 323 if (*ierr) return; 324 } 325 326 PETSC_EXTERN void dmdavecrestorearrayreadf902_(DM *da, Vec *v, F90Array2d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 327 { 328 const PetscScalar *fa; 329 *ierr = F90Array2dAccess(a, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd)); 330 *ierr = VecRestoreArrayRead(*v, &fa); 331 if (*ierr) return; 332 *ierr = F90Array2dDestroy(a, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd)); 333 } 334 335 PETSC_EXTERN void dmdavecgetarrayreadf903_(DM *da, Vec *v, F90Array3d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 336 { 337 PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof; 338 const PetscScalar *aa; 339 340 *ierr = DMDAGetCorners(*da, &xs, &ys, &zs, &xm, &ym, &zm); 341 if (*ierr) return; 342 *ierr = DMDAGetGhostCorners(*da, &gxs, &gys, &gzs, &gxm, &gym, &gzm); 343 if (*ierr) return; 344 *ierr = DMDAGetInfo(*da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL); 345 if (*ierr) return; 346 347 /* Handle case where user passes in global vector as opposed to local */ 348 *ierr = VecGetLocalSize(*v, &N); 349 if (*ierr) return; 350 if (N == xm * ym * zm * dof) { 351 gxm = xm; 352 gym = ym; 353 gzm = zm; 354 gxs = xs; 355 gys = ys; 356 gzs = zs; 357 } else if (N != gxm * gym * gzm * dof) { 358 *ierr = PETSC_ERR_ARG_INCOMP; 359 return; 360 } 361 if (dim == 2) { 362 gzs = gys; 363 gzm = gym; 364 gys = gxs; 365 gym = gxm; 366 gxs = 0; 367 gxm = dof; 368 } 369 *ierr = VecGetArrayRead(*v, &aa); 370 if (*ierr) return; 371 *ierr = F90Array3dCreate((void *)aa, MPIU_SCALAR, gxs, gxm, gys, gym, gzs, gzm, a PETSC_F90_2PTR_PARAM(ptrd)); 372 if (*ierr) return; 373 } 374 375 PETSC_EXTERN void dmdavecrestorearrayreadf903_(DM *da, Vec *v, F90Array3d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 376 { 377 const PetscScalar *fa; 378 *ierr = F90Array3dAccess(a, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd)); 379 *ierr = VecRestoreArrayRead(*v, &fa); 380 if (*ierr) return; 381 *ierr = F90Array3dDestroy(a, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd)); 382 } 383 384 PETSC_EXTERN void dmdavecgetarrayreadf904_(DM *da, Vec *v, F90Array4d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 385 { 386 PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof, zero = 0; 387 const PetscScalar *aa; 388 389 *ierr = DMDAGetCorners(*da, &xs, &ys, &zs, &xm, &ym, &zm); 390 if (*ierr) return; 391 *ierr = DMDAGetGhostCorners(*da, &gxs, &gys, &gzs, &gxm, &gym, &gzm); 392 if (*ierr) return; 393 *ierr = DMDAGetInfo(*da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL); 394 if (*ierr) return; 395 396 /* Handle case where user passes in global vector as opposed to local */ 397 *ierr = VecGetLocalSize(*v, &N); 398 if (*ierr) return; 399 if (N == xm * ym * zm * dof) { 400 gxm = xm; 401 gym = ym; 402 gzm = zm; 403 gxs = xs; 404 gys = ys; 405 gzs = zs; 406 } else if (N != gxm * gym * gzm * dof) { 407 *ierr = PETSC_ERR_ARG_INCOMP; 408 return; 409 } 410 *ierr = VecGetArrayRead(*v, &aa); 411 if (*ierr) return; 412 *ierr = F90Array4dCreate((void *)aa, MPIU_SCALAR, zero, dof, gxs, gxm, gys, gym, gzs, gzm, a PETSC_F90_2PTR_PARAM(ptrd)); 413 if (*ierr) return; 414 } 415 416 PETSC_EXTERN void dmdavecrestorearrayreadf904_(DM *da, Vec *v, F90Array4d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 417 { 418 const PetscScalar *fa; 419 /* 420 F90Array4dAccess is not implemented, so the following call would fail 421 */ 422 *ierr = F90Array4dAccess(a, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd)); 423 *ierr = VecRestoreArrayRead(*v, &fa); 424 if (*ierr) return; 425 *ierr = F90Array4dDestroy(a, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd)); 426 } 427