1 #include <petsc/private/ftnimpl.h>
2 #include <petscdmswarm.h>
3
4 #if defined(PETSC_HAVE_FORTRAN_CAPS)
5 #define dmswarmgetfield_ DMSWARMGETFIELD
6 #define dmswarmrestorefield_ DMSWARMRESTOREFIELD
7 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
8 #define dmswarmgetfield_ dmswarmgetfield
9 #define dmswarmrestorefield_ dmswarmrestorefield
10 #endif
11
12 /* Definitions of Fortran Wrapper routines */
13
dmswarmgetfield_(DM * dm,char * name,PetscInt * blocksize,PetscDataType * type,F90Array1d * ptr,int * ierr PETSC_F90_2PTR_PROTO (ptrd),PETSC_FORTRAN_CHARLEN_T lenN)14 PETSC_EXTERN void dmswarmgetfield_(DM *dm, char *name, PetscInt *blocksize, PetscDataType *type, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd), PETSC_FORTRAN_CHARLEN_T lenN)
15 {
16 PetscScalar *v;
17 PetscInt n;
18 char *fieldname;
19
20 FIXCHAR(name, lenN, fieldname);
21 *ierr = DMSwarmGetSize(*dm, &n);
22 if (*ierr) return;
23 *ierr = DMSwarmGetField(*dm, fieldname, blocksize, type, (void **)&v);
24 if (*ierr) return;
25 *ierr = F90Array1dCreate((void *)v, MPIU_SCALAR, 1, n, ptr PETSC_F90_2PTR_PARAM(ptrd));
26 FREECHAR(name, fieldname);
27 }
28
dmswarmrestorefield_(DM * dm,char * name,PetscInt * blocksize,PetscDataType * type,F90Array1d * ptr,int * ierr PETSC_F90_2PTR_PROTO (ptrd),PETSC_FORTRAN_CHARLEN_T lenN)29 PETSC_EXTERN void dmswarmrestorefield_(DM *dm, char *name, PetscInt *blocksize, PetscDataType *type, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd), PETSC_FORTRAN_CHARLEN_T lenN)
30 {
31 PetscScalar *v;
32 char *fieldname;
33
34 FIXCHAR(name, lenN, fieldname);
35 *ierr = F90Array1dAccess(ptr, MPIU_SCALAR, (void **)&v PETSC_F90_2PTR_PARAM(ptrd));
36 if (*ierr) return;
37 if (*ierr) return;
38 *ierr = DMSwarmRestoreField(*dm, fieldname, blocksize, type, (void **)&v);
39 if (*ierr) return;
40 *ierr = F90Array1dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
41 FREECHAR(name, fieldname);
42 }
43