xref: /petsc/src/vec/is/sf/interface/ftn-custom/zsf.c (revision 8e38371517a4c87f570b6f306e56fb899cb1fff8)
18af6ec1cSBarry Smith #include <petsc/private/f90impl.h>
2a6e9e4f7SMatthew G. Knepley #include <petsc/private/sfimpl.h>
3a6e9e4f7SMatthew G. Knepley 
4a6e9e4f7SMatthew G. Knepley #if defined(PETSC_HAVE_FORTRAN_CAPS)
55c87a30dSBarry Smith #define petscsfview_          PETSCSFVIEW
68d1b7334SBarry Smith #define petscsfgetgraph_      PETSCSFGETGRAPH
78af6ec1cSBarry Smith #define petscsfbcastbegin_    PETSCSFBCASTBEGIN
88af6ec1cSBarry Smith #define petscsfbcastend_      PETSCSFBCASTEND
97f139299SBarry Smith #define f90arraysfnodecreate_ F90ARRAYSFNODECREATE
10a6e9e4f7SMatthew G. Knepley #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
115c87a30dSBarry Smith #define petscsfgetgraph_      petscsfgetgraph
125c87a30dSBarry Smith #define petscsfview_          petscsfview
138af6ec1cSBarry Smith #define petscsfbcastbegin_    petscsfbcastbegin
148af6ec1cSBarry Smith #define petscsfbcastend_      petscsfbcastend
157f139299SBarry Smith #define f90arraysfnodecreate_ f90arraysfnodecreate
16a6e9e4f7SMatthew G. Knepley #endif
17a6e9e4f7SMatthew G. Knepley 
187f139299SBarry Smith PETSC_EXTERN void PETSC_STDCALL f90arraysfnodecreate_(const PetscInt *,PetscInt *,void * PETSC_F90_2PTR_PROTO_NOVAR);
197f139299SBarry Smith 
201943db53SBarry Smith PETSC_EXTERN void PETSC_STDCALL petscsfview_(PetscSF *sf, PetscViewer *vin, PetscErrorCode *ierr)
21a6e9e4f7SMatthew G. Knepley {
22a6e9e4f7SMatthew G. Knepley   PetscViewer v;
23a6e9e4f7SMatthew G. Knepley 
24a6e9e4f7SMatthew G. Knepley   PetscPatchDefaultViewers_Fortran(vin, v);
25a6e9e4f7SMatthew G. Knepley   *ierr = PetscSFView(*sf, v);
26a6e9e4f7SMatthew G. Knepley }
278af6ec1cSBarry Smith 
287f139299SBarry Smith 
295c87a30dSBarry Smith PETSC_EXTERN void PETSC_STDCALL  petscsfgetgraph_(PetscSF *sf,PetscInt *nroots,PetscInt *nleaves, F90Array1d  *ailocal, F90Array1d  *airemote, int *ierr PETSC_F90_2PTR_PROTO(pilocal) PETSC_F90_2PTR_PROTO(piremote))
305c87a30dSBarry Smith {
315c87a30dSBarry Smith   const PetscInt    *ilocal;
325c87a30dSBarry Smith   const PetscSFNode *iremote;
335c87a30dSBarry Smith 
345c87a30dSBarry Smith   *ierr = PetscSFGetGraph(*sf,nroots,nleaves,&ilocal,&iremote);if (*ierr) return;
356a43dea8SBarry Smith   *ierr = F90Array1dCreate((void*)ilocal,MPIU_INT,1,*nleaves, ailocal PETSC_F90_2PTR_PARAM(pilocal));
367f139299SBarry Smith   /* this creates a memory leak */
37*8e383715SSatish Balay   f90arraysfnodecreate_((PetscInt*)iremote,nleaves, airemote PETSC_F90_2PTR_PARAM(piremote));
385c87a30dSBarry Smith }
395c87a30dSBarry Smith 
408af6ec1cSBarry Smith PETSC_EXTERN void PETSC_STDCALL petscsfbcastbegin_(PetscSF *sf, MPI_Fint *unit,F90Array1d *rptr, F90Array1d *lptr, int *ierr PETSC_F90_2PTR_PROTO(rptrd) PETSC_F90_2PTR_PROTO(lptrd))
418af6ec1cSBarry Smith {
428af6ec1cSBarry Smith   MPI_Datatype dtype;
438af6ec1cSBarry Smith   const void   *rootdata;
448af6ec1cSBarry Smith   void         *leafdata;
458af6ec1cSBarry Smith 
465c87a30dSBarry Smith 
475c87a30dSBarry Smith   *ierr = PetscMPIFortranDatatypeToC(*unit,&dtype);if (*ierr) return;
488af6ec1cSBarry Smith   *ierr = F90Array1dAccess(rptr, dtype, (void**) &rootdata PETSC_F90_2PTR_PARAM(rptrd));if (*ierr) return;
498af6ec1cSBarry Smith   *ierr = F90Array1dAccess(lptr, dtype, (void**) &leafdata PETSC_F90_2PTR_PARAM(lptrd));if (*ierr) return;
508af6ec1cSBarry Smith   *ierr = PetscSFBcastBegin(*sf, dtype, rootdata, leafdata);
518af6ec1cSBarry Smith }
528af6ec1cSBarry Smith 
538af6ec1cSBarry Smith PETSC_EXTERN void PETSC_STDCALL petscsfbcastend_(PetscSF *sf, MPI_Fint *unit,F90Array1d *rptr, F90Array1d *lptr, int *ierr PETSC_F90_2PTR_PROTO(rptrd) PETSC_F90_2PTR_PROTO(lptrd))
548af6ec1cSBarry Smith {
558af6ec1cSBarry Smith   MPI_Datatype dtype;
568af6ec1cSBarry Smith   const void   *rootdata;
578af6ec1cSBarry Smith   void         *leafdata;
588af6ec1cSBarry Smith 
595c87a30dSBarry Smith   *ierr = PetscMPIFortranDatatypeToC(*unit,&dtype);if (*ierr) return;
608af6ec1cSBarry Smith   *ierr = F90Array1dAccess(rptr, dtype, (void**) &rootdata PETSC_F90_2PTR_PARAM(rptrd));if (*ierr) return;
618af6ec1cSBarry Smith   *ierr = F90Array1dAccess(lptr, dtype, (void**) &leafdata PETSC_F90_2PTR_PARAM(lptrd));if (*ierr) return;
628af6ec1cSBarry Smith   *ierr = PetscSFBcastEnd(*sf, dtype, rootdata, leafdata);
638af6ec1cSBarry Smith }
64