Actual source code: ex7.c

  1: static char help[] = "Demonstrates calling a Fortran computational routine from C.\n\
  2: Also demonstrates passing  PETSc objects, MPI Communicators from C to Fortran\n\
  3: and from Fortran to C\n\n";

  5: #include <petscvec.h>
  6: /*
  7:   Ugly stuff to insure the function names match between Fortran
  8:   and C. This is out of our PETSc hands to cleanup.
  9: */
 10: #include <petsc/private/fortranimpl.h>
 11: #if defined(PETSC_HAVE_FORTRAN_CAPS)
 12:   #define ex7f_ EX7F
 13:   #define ex7c_ EX7C
 14: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
 15:   #define ex7f_ ex7f
 16:   #define ex7c_ ex7c
 17: #endif

 19: PETSC_INTERN void ex7f_(Vec *, int *);

 21: int main(int argc, char **args)
 22: {
 23:   PetscInt m = 10;
 24:   int      fcomm;
 25:   Vec      vec;

 27:   PetscFunctionBeginUser;
 28:   PetscCall(PetscInitialize(&argc, &args, (char *)0, help));
 29:   /* This function should be called to be able to use PETSc routines
 30:      from the FORTRAN subroutines needed by this program */

 32:   PetscCall(PetscInitializeFortran());

 34:   PetscCall(VecCreate(PETSC_COMM_WORLD, &vec));
 35:   PetscCall(VecSetSizes(vec, PETSC_DECIDE, m));
 36:   PetscCall(VecSetFromOptions(vec));

 38:   /*
 39:      Call Fortran routine - the use of MPI_Comm_c2f() allows
 40:      translation of the MPI_Comm from C so that it can be properly
 41:      interpreted from Fortran.
 42:   */
 43:   fcomm = MPI_Comm_c2f(PETSC_COMM_WORLD);

 45:   ex7f_(&vec, &fcomm);

 47:   PetscCall(VecView(vec, PETSC_VIEWER_STDOUT_WORLD));
 48:   PetscCall(VecDestroy(&vec));
 49:   PetscCall(PetscFinalize());
 50:   return 0;
 51: }

 53: PETSC_INTERN void ex7c_(Vec *fvec, int *fcomm, PetscErrorCode *ierr)
 54: {
 55:   MPI_Comm comm;
 56:   PetscInt vsize;

 58:   /*
 59:     Translate Fortran integer pointer back to C and
 60:     Fortran Communicator back to C communicator
 61:   */
 62:   comm = MPI_Comm_f2c(*fcomm);

 64:   /* Some PETSc/MPI operations on Vec/Communicator objects */
 65:   *ierr = VecGetSize(*fvec, &vsize);
 66:   if (*ierr) return;
 67:   if (MPI_Barrier(comm)) *ierr = PETSC_ERR_MPI;
 68: }

 70: /*TEST

 72:    build:
 73:      depends: ex7f.F
 74:      requires: fortran

 76:    test:
 77:       nsize: 3
 78:       filter: sort -b |grep -v " MPI process"
 79:       filter_output: sort -b

 81: TEST*/