1 #include <petscvec.h>
2 #include <petsc/private/ftnimpl.h>
3
4 #if defined(PETSC_HAVE_FORTRAN_CAPS)
5 #define vecgetarraywrite_ VECGETARRAYWRITE
6 #define vecrestorearraywrite_ VECRESTOREARRAYWRITE
7 #define vecgetarray_ VECGETARRAY
8 #define vecrestorearray_ VECRESTOREARRAY
9 #define vecgetarrayread_ VECGETARRAYREAD
10 #define vecrestorearrayread_ VECRESTOREARRAYREAD
11 #define vecduplicatevecs_ VECDUPLICATEVECS
12 #define vecdestroyvecs_ VECDESTROYVECS
13 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
14 #define vecgetarraywrite_ vecgetarraywrite
15 #define vecrestorearraywrite_ vecrestorearraywrite
16 #define vecgetarray_ vecgetarray
17 #define vecrestorearray_ vecrestorearray
18 #define vecgetarrayread_ vecgetarrayread
19 #define vecrestorearrayread_ vecrestorearrayread
20 #define vecduplicatevecs_ vecduplicatevecs
21 #define vecdestroyvecs_ vecdestroyvecs
22 #endif
23
vecgetarraywrite_(Vec * x,F90Array1d * ptr,int * ierr PETSC_F90_2PTR_PROTO (ptrd))24 PETSC_EXTERN void vecgetarraywrite_(Vec *x, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
25 {
26 PetscScalar *fa;
27 PetscInt len;
28 if (!ptr) {
29 *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?");
30 return;
31 }
32 *ierr = VecGetArrayWrite(*x, &fa);
33 if (*ierr) return;
34 *ierr = VecGetLocalSize(*x, &len);
35 if (*ierr) return;
36 *ierr = F90Array1dCreate(fa, MPIU_SCALAR, 1, len, ptr PETSC_F90_2PTR_PARAM(ptrd));
37 }
38
vecrestorearraywrite_(Vec * x,F90Array1d * ptr,int * ierr PETSC_F90_2PTR_PROTO (ptrd))39 PETSC_EXTERN void vecrestorearraywrite_(Vec *x, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
40 {
41 PetscScalar *fa;
42 *ierr = F90Array1dAccess(ptr, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd));
43 if (*ierr) return;
44 *ierr = F90Array1dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
45 if (*ierr) return;
46 *ierr = VecRestoreArrayWrite(*x, &fa);
47 }
48
vecgetarray_(Vec * x,F90Array1d * ptr,int * ierr PETSC_F90_2PTR_PROTO (ptrd))49 PETSC_EXTERN void vecgetarray_(Vec *x, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
50 {
51 PetscScalar *fa;
52 PetscInt len;
53 if (!ptr) {
54 *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?");
55 return;
56 }
57 *ierr = VecGetArray(*x, &fa);
58 if (*ierr) return;
59 *ierr = VecGetLocalSize(*x, &len);
60 if (*ierr) return;
61 *ierr = F90Array1dCreate(fa, MPIU_SCALAR, 1, len, ptr PETSC_F90_2PTR_PARAM(ptrd));
62 }
63
vecrestorearray_(Vec * x,F90Array1d * ptr,int * ierr PETSC_F90_2PTR_PROTO (ptrd))64 PETSC_EXTERN void vecrestorearray_(Vec *x, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
65 {
66 PetscScalar *fa;
67 *ierr = F90Array1dAccess(ptr, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd));
68 if (*ierr) return;
69 *ierr = F90Array1dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
70 if (*ierr) return;
71 *ierr = VecRestoreArray(*x, &fa);
72 }
73
vecgetarrayread_(Vec * x,F90Array1d * ptr,int * ierr PETSC_F90_2PTR_PROTO (ptrd))74 PETSC_EXTERN void vecgetarrayread_(Vec *x, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
75 {
76 const PetscScalar *fa;
77 PetscInt len;
78 if (!ptr) {
79 *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?");
80 return;
81 }
82 *ierr = VecGetArrayRead(*x, &fa);
83 if (*ierr) return;
84 *ierr = VecGetLocalSize(*x, &len);
85 if (*ierr) return;
86 *ierr = F90Array1dCreate((PetscScalar *)fa, MPIU_SCALAR, 1, len, ptr PETSC_F90_2PTR_PARAM(ptrd));
87 }
88
vecrestorearrayread_(Vec * x,F90Array1d * ptr,int * ierr PETSC_F90_2PTR_PROTO (ptrd))89 PETSC_EXTERN void vecrestorearrayread_(Vec *x, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
90 {
91 const PetscScalar *fa;
92 *ierr = F90Array1dAccess(ptr, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd));
93 if (*ierr) return;
94 *ierr = F90Array1dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
95 if (*ierr) return;
96 *ierr = VecRestoreArrayRead(*x, &fa);
97 }
98
vecduplicatevecs_(Vec * v,int * m,F90Array1d * ptr,int * ierr PETSC_F90_2PTR_PROTO (ptrd))99 PETSC_EXTERN void vecduplicatevecs_(Vec *v, int *m, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
100 {
101 Vec *lV;
102 PetscFortranAddr *newvint;
103 int i;
104 *ierr = VecDuplicateVecs(*v, *m, &lV);
105 if (*ierr) return;
106 *ierr = PetscMalloc1(*m, &newvint);
107 if (*ierr) return;
108
109 for (i = 0; i < *m; i++) newvint[i] = (PetscFortranAddr)lV[i];
110 *ierr = PetscFree(lV);
111 if (*ierr) return;
112 *ierr = F90Array1dCreate(newvint, MPIU_FORTRANADDR, 1, *m, ptr PETSC_F90_2PTR_PARAM(ptrd));
113 }
114
vecdestroyvecs_(int * m,F90Array1d * ptr,int * ierr PETSC_F90_2PTR_PROTO (ptrd))115 PETSC_EXTERN void vecdestroyvecs_(int *m, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
116 {
117 Vec *vecs;
118 int i;
119
120 *ierr = F90Array1dAccess(ptr, MPIU_FORTRANADDR, (void **)&vecs PETSC_F90_2PTR_PARAM(ptrd));
121 if (*ierr) return;
122 for (i = 0; i < *m; i++) {
123 PETSC_FORTRAN_OBJECT_F_DESTROYED_TO_C_NULL(&vecs[i]);
124 *ierr = VecDestroy(&vecs[i]);
125 if (*ierr) return;
126 PETSC_FORTRAN_OBJECT_C_NULL_TO_F_DESTROYED(&vecs[i]);
127 }
128 *ierr = F90Array1dDestroy(ptr, MPIU_FORTRANADDR PETSC_F90_2PTR_PARAM(ptrd));
129 if (*ierr) return;
130 *ierr = PetscFree(vecs);
131 }
132