Actual source code: zvectorf90.c

  1: #include <petscvec.h>
  2: #include <petsc/private/f90impl.h>

  4: #if defined(PETSC_HAVE_FORTRAN_CAPS)
  5:   #define vecgetarrayf90_         VECGETARRAYF90
  6:   #define vecrestorearrayf90_     VECRESTOREARRAYF90
  7:   #define vecgetarrayreadf90_     VECGETARRAYREADF90
  8:   #define vecrestorearrayreadf90_ VECRESTOREARRAYREADF90
  9:   #define vecduplicatevecsf90_    VECDUPLICATEVECSF90
 10:   #define vecdestroyvecsf90_      VECDESTROYVECSF90
 11: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
 12:   #define vecgetarrayf90_         vecgetarrayf90
 13:   #define vecrestorearrayf90_     vecrestorearrayf90
 14:   #define vecgetarrayreadf90_     vecgetarrayreadf90
 15:   #define vecrestorearrayreadf90_ vecrestorearrayreadf90
 16:   #define vecduplicatevecsf90_    vecduplicatevecsf90
 17:   #define vecdestroyvecsf90_      vecdestroyvecsf90
 18: #endif

 20: PETSC_EXTERN void vecgetarrayf90_(Vec *x, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
 21: {
 22:   PetscScalar *fa;
 23:   PetscInt     len;
 24:   if (!ptr) {
 25:     *ierr = PetscError(((PetscObject)*x)->comm, __LINE__, PETSC_FUNCTION_NAME, __FILE__, PETSC_ERR_ARG_BADPTR, PETSC_ERROR_INITIAL, "ptr==NULL, maybe #include <petsc/finclude/petscvec.h> is missing?");
 26:     return;
 27:   }
 28:   *ierr = VecGetArray(*x, &fa);
 29:   if (*ierr) return;
 30:   *ierr = VecGetLocalSize(*x, &len);
 31:   if (*ierr) return;
 32:   *ierr = F90Array1dCreate(fa, MPIU_SCALAR, 1, len, ptr PETSC_F90_2PTR_PARAM(ptrd));
 33: }
 34: PETSC_EXTERN void vecrestorearrayf90_(Vec *x, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
 35: {
 36:   PetscScalar *fa;
 37:   *ierr = F90Array1dAccess(ptr, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd));
 38:   if (*ierr) return;
 39:   *ierr = F90Array1dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
 40:   if (*ierr) return;
 41:   *ierr = VecRestoreArray(*x, &fa);
 42: }

 44: PETSC_EXTERN void vecgetarrayreadf90_(Vec *x, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
 45: {
 46:   const PetscScalar *fa;
 47:   PetscInt           len;
 48:   if (!ptr) {
 49:     *ierr = PetscError(((PetscObject)*x)->comm, __LINE__, PETSC_FUNCTION_NAME, __FILE__, PETSC_ERR_ARG_BADPTR, PETSC_ERROR_INITIAL, "ptr==NULL, maybe #include <petsc/finclude/petscvec.h> is missing?");
 50:     return;
 51:   }
 52:   *ierr = VecGetArrayRead(*x, &fa);
 53:   if (*ierr) return;
 54:   *ierr = VecGetLocalSize(*x, &len);
 55:   if (*ierr) return;
 56:   *ierr = F90Array1dCreate((PetscScalar *)fa, MPIU_SCALAR, 1, len, ptr PETSC_F90_2PTR_PARAM(ptrd));
 57: }
 58: PETSC_EXTERN void vecrestorearrayreadf90_(Vec *x, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
 59: {
 60:   const PetscScalar *fa;
 61:   *ierr = F90Array1dAccess(ptr, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd));
 62:   if (*ierr) return;
 63:   *ierr = F90Array1dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
 64:   if (*ierr) return;
 65:   *ierr = VecRestoreArrayRead(*x, &fa);
 66: }

 68: PETSC_EXTERN void vecduplicatevecsf90_(Vec *v, int *m, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
 69: {
 70:   Vec              *lV;
 71:   PetscFortranAddr *newvint;
 72:   int               i;
 73:   *ierr = VecDuplicateVecs(*v, *m, &lV);
 74:   if (*ierr) return;
 75:   *ierr = PetscMalloc1(*m, &newvint);
 76:   if (*ierr) return;

 78:   for (i = 0; i < *m; i++) newvint[i] = (PetscFortranAddr)lV[i];
 79:   *ierr = PetscFree(lV);
 80:   if (*ierr) return;
 81:   *ierr = F90Array1dCreate(newvint, MPIU_FORTRANADDR, 1, *m, ptr PETSC_F90_2PTR_PARAM(ptrd));
 82: }

 84: PETSC_EXTERN void vecdestroyvecsf90_(int *m, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
 85: {
 86:   Vec *vecs;
 87:   int  i;

 89:   *ierr = F90Array1dAccess(ptr, MPIU_FORTRANADDR, (void **)&vecs PETSC_F90_2PTR_PARAM(ptrd));
 90:   if (*ierr) return;
 91:   for (i = 0; i < *m; i++) {
 92:     PETSC_FORTRAN_OBJECT_F_DESTROYED_TO_C_NULL(&vecs[i]);
 93:     *ierr = VecDestroy(&vecs[i]);
 94:     if (*ierr) return;
 95:     PETSC_FORTRAN_OBJECT_C_NULL_TO_F_DESTROYED(&vecs[i]);
 96:   }
 97:   *ierr = F90Array1dDestroy(ptr, MPIU_FORTRANADDR PETSC_F90_2PTR_PARAM(ptrd));
 98:   if (*ierr) return;
 99:   *ierr = PetscFree(vecs);
100: }