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:   #define vecdestroy_             VECDESTROY
 12: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
 13:   #define vecgetarrayf90_         vecgetarrayf90
 14:   #define vecrestorearrayf90_     vecrestorearrayf90
 15:   #define vecgetarrayreadf90_     vecgetarrayreadf90
 16:   #define vecrestorearrayreadf90_ vecrestorearrayreadf90
 17:   #define vecduplicatevecsf90_    vecduplicatevecsf90
 18:   #define vecdestroyvecsf90_      vecdestroyvecsf90
 19:   #define vecdestroy_             vecdestroy
 20: #endif

 22: PETSC_EXTERN void vecgetarrayf90_(Vec *x, F90Array1d *ptr, int *__ierr PETSC_F90_2PTR_PROTO(ptrd))
 23: {
 24:   PetscScalar *fa;
 25:   PetscInt     len;
 26:   if (!ptr) {
 27:     *__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?");
 28:     return;
 29:   }
 30:   *__ierr = VecGetArray(*x, &fa);
 31:   if (*__ierr) return;
 32:   *__ierr = VecGetLocalSize(*x, &len);
 33:   if (*__ierr) return;
 34:   *__ierr = F90Array1dCreate(fa, MPIU_SCALAR, 1, len, ptr PETSC_F90_2PTR_PARAM(ptrd));
 35: }
 36: PETSC_EXTERN void vecrestorearrayf90_(Vec *x, F90Array1d *ptr, int *__ierr PETSC_F90_2PTR_PROTO(ptrd))
 37: {
 38:   PetscScalar *fa;
 39:   *__ierr = F90Array1dAccess(ptr, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd));
 40:   if (*__ierr) return;
 41:   *__ierr = F90Array1dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
 42:   if (*__ierr) return;
 43:   *__ierr = VecRestoreArray(*x, &fa);
 44: }

 46: PETSC_EXTERN void vecgetarrayreadf90_(Vec *x, F90Array1d *ptr, int *__ierr PETSC_F90_2PTR_PROTO(ptrd))
 47: {
 48:   const PetscScalar *fa;
 49:   PetscInt           len;
 50:   if (!ptr) {
 51:     *__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?");
 52:     return;
 53:   }
 54:   *__ierr = VecGetArrayRead(*x, &fa);
 55:   if (*__ierr) return;
 56:   *__ierr = VecGetLocalSize(*x, &len);
 57:   if (*__ierr) return;
 58:   *__ierr = F90Array1dCreate((PetscScalar *)fa, MPIU_SCALAR, 1, len, ptr PETSC_F90_2PTR_PARAM(ptrd));
 59: }
 60: PETSC_EXTERN void vecrestorearrayreadf90_(Vec *x, F90Array1d *ptr, int *__ierr PETSC_F90_2PTR_PROTO(ptrd))
 61: {
 62:   const PetscScalar *fa;
 63:   *__ierr = F90Array1dAccess(ptr, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd));
 64:   if (*__ierr) return;
 65:   *__ierr = F90Array1dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
 66:   if (*__ierr) return;
 67:   *__ierr = VecRestoreArrayRead(*x, &fa);
 68: }

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

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

 86: PETSC_EXTERN void vecdestroyvecsf90_(int *m, F90Array1d *ptr, int *__ierr PETSC_F90_2PTR_PROTO(ptrd))
 87: {
 88:   Vec *vecs;
 89:   int  i;

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

104: PETSC_EXTERN void vecdestroy_(Vec *x, int *ierr)
105: {
106:   PETSC_FORTRAN_OBJECT_F_DESTROYED_TO_C_NULL(x);
107:   *ierr = VecDestroy(x);
108:   if (*ierr) return;
109:   PETSC_FORTRAN_OBJECT_C_NULL_TO_F_DESTROYED(x);
110: }