1 #include <petsc/private/ftnimpl.h>
2 #include <petscdt.h>
3
4 #if defined(PETSC_HAVE_FORTRAN_CAPS)
5 #define petscquadraturegetdata_ PETSCQUADRATUREGETDATA
6 #define petscquadraturerestoredata_ PETSCQUADRATURERESTOREDATA
7 #define petscquadraturesetdata_ PETSCQUADRATURESETDATA
8 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
9 #define petscquadraturegetdata_ petscquadraturegetdata
10 #define petscquadraturerestoredata_ petscquadraturerestoredata
11 #define petscquadraturesetdata_ petscquadraturesetdata
12 #endif
13
petscquadraturegetdata_(PetscQuadrature * q,PetscInt * dim,PetscInt * Nc,PetscInt * npoints,F90Array1d * ptrP,F90Array1d * ptrW,PetscErrorCode * ierr PETSC_F90_2PTR_PROTO (ptrp)PETSC_F90_2PTR_PROTO (ptrw))14 PETSC_EXTERN void petscquadraturegetdata_(PetscQuadrature *q, PetscInt *dim, PetscInt *Nc, PetscInt *npoints, F90Array1d *ptrP, F90Array1d *ptrW, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrp) PETSC_F90_2PTR_PROTO(ptrw))
15 {
16 const PetscReal *points, *weights;
17
18 *ierr = PetscQuadratureGetData(*q, dim, Nc, npoints, &points, &weights);
19 if (*ierr) return;
20 *ierr = F90Array1dCreate((void *)points, MPIU_REAL, 1, (*npoints) * (*dim), ptrP PETSC_F90_2PTR_PARAM(ptrp));
21 if (*ierr) return;
22 *ierr = F90Array1dCreate((void *)weights, MPIU_REAL, 1, (*npoints) * (*Nc), ptrW PETSC_F90_2PTR_PARAM(ptrw));
23 }
24
petscquadraturerestoredata_(PetscQuadrature * q,PetscInt * dim,PetscInt * Nc,PetscInt * npoints,F90Array1d * ptrP,F90Array1d * ptrW,PetscErrorCode * ierr PETSC_F90_2PTR_PROTO (ptrp)PETSC_F90_2PTR_PROTO (ptrw))25 PETSC_EXTERN void petscquadraturerestoredata_(PetscQuadrature *q, PetscInt *dim, PetscInt *Nc, PetscInt *npoints, F90Array1d *ptrP, F90Array1d *ptrW, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrp) PETSC_F90_2PTR_PROTO(ptrw))
26 {
27 *ierr = F90Array1dDestroy(ptrP, MPIU_REAL PETSC_F90_2PTR_PARAM(ptrp));
28 if (*ierr) return;
29 *ierr = F90Array1dDestroy(ptrW, MPIU_REAL PETSC_F90_2PTR_PARAM(ptrw));
30 }
31
petscquadraturesetdata_(PetscQuadrature * q,PetscInt * dim,PetscInt * Nc,PetscInt * npoints,F90Array1d * ptrP,F90Array1d * ptrW,PetscErrorCode * ierr PETSC_F90_2PTR_PROTO (ptrp)PETSC_F90_2PTR_PROTO (ptrw))32 PETSC_EXTERN void petscquadraturesetdata_(PetscQuadrature *q, PetscInt *dim, PetscInt *Nc, PetscInt *npoints, F90Array1d *ptrP, F90Array1d *ptrW, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrp) PETSC_F90_2PTR_PROTO(ptrw))
33 {
34 PetscReal *points, *weights;
35
36 *ierr = F90Array1dAccess(ptrP, MPIU_REAL, (void **)&points PETSC_F90_2PTR_PARAM(ptrp));
37 if (*ierr) return;
38 *ierr = F90Array1dAccess(ptrW, MPIU_REAL, (void **)&weights PETSC_F90_2PTR_PARAM(ptrw));
39 if (*ierr) return;
40 *ierr = PetscQuadratureSetData(*q, *dim, *Nc, *npoints, points, weights);
41 }
42