Actual source code: zf90vec.c
2: #include "petscis.h"
3: #include "petscvec.h"
4: #include "petscf90.h"
6: #ifdef PETSC_HAVE_FORTRAN_CAPS
7: #define isgetindicesf90_ ISGETINDICESF90
8: #define isblockgetindicesf90_ ISBLOCKGETINDICESF90
9: #define isrestoreindicesf90_ ISRESTOREINDICESF90
10: #define isblockrestoreindicesf90_ ISBLOCKRESTOREINDICESF90
11: #define iscoloringgetisf90_ ISCOLORINGGETISF90
12: #define iscoloringrestoreisf90_ ISCOLORINGRESTOREF90
13: #define vecgetarrayf90_ VECGETARRAYF90
14: #define vecrestorearrayf90_ VECRESTOREARRAYF90
15: #define vecduplicatevecsf90_ VECDUPLICATEVECSF90
16: #define vecdestroyvecsf90_ VECDESTROYVECSF90
17: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
18: #define isgetindicesf90_ isgetindicesf90
19: #define isblockgetindicesf90_ isblockgetindicesf90
20: #define isrestoreindicesf90_ isrestoreindicesf90
21: #define isblockrestoreindicesf90_ isblockrestoreindicesf90
22: #define iscoloringgetisf90_ iscoloringgetisf90
23: #define iscoloringrestoreisf90_ iscoloringrestoreisf90
24: #define vecgetarrayf90_ vecgetarrayf90
25: #define vecrestorearrayf90_ vecrestorearrayf90
26: #define vecduplicatevecsf90_ vecduplicatevecsf90
27: #define vecdestroyvecsf90_ vecdestroyvecsf90
28: #endif
32: /* --------------------------------------------------------------- */
34: void PETSC_STDCALL isgetindicesf90_(IS *x,F90Array1d *ptr,int *__ierr)
35: {
36: int *fa;
37: int len;
39: *__ISGetIndices(*x,&fa); if (*__ierr) return;
40: *__ISGetLocalSize(*x,&len); if (*__ierr) return;
41: *__F90Array1dCreate(fa,PETSC_INT,1,len,ptr);
42: }
43: void PETSC_STDCALL isrestoreindicesf90_(IS *x,F90Array1d *ptr,int *__ierr)
44: {
45: int *fa;
46: *__F90Array1dAccess(ptr,(void**)&fa);if (*__ierr) return;
47: *__F90Array1dDestroy(ptr);if (*__ierr) return;
48: *__ISRestoreIndices(*x,&fa);
49: }
51: void PETSC_STDCALL isblockgetindicesf90_(IS *x,F90Array1d *ptr,int *__ierr)
52: {
53: int *fa;
54: int len;
55: *__ISBlockGetIndices(*x,&fa); if (*__ierr) return;
56: *__ISBlockGetSize(*x,&len); if (*__ierr) return;
57: *__F90Array1dCreate(fa,PETSC_INT,1,len,ptr);
58: }
59: void PETSC_STDCALL isblockrestoreindicesf90_(IS *x,F90Array1d *ptr,int *__ierr)
60: {
61: int *fa;
62: *__F90Array1dAccess(ptr,(void**)&fa);if (*__ierr) return;
63: *__F90Array1dDestroy(ptr);if (*__ierr) return;
64: *__ISBlockRestoreIndices(*x,&fa);
65: }
68: void PETSC_STDCALL iscoloringgetisf90_(ISColoring *iscoloring,int *n,F90Array1d *ptr,int *__ierr)
69: {
70: IS *lis;
71: PetscFortranAddr *newisint;
72: int i;
73: *__ISColoringGetIS(*iscoloring,n,&lis); if (*__ierr) return;
74: *__PetscMalloc((*n)*sizeof(PetscFortranAddr),&newisint); if (*__ierr) return;
75: for (i=0; i<*n; i++) {
76: newisint[i] = (PetscFortranAddr)lis[i];
77: }
78: *__F90Array1dCreate(newisint,PETSC_FORTRANADDR,1,*n,ptr);
79: }
81: void PETSC_STDCALL iscoloringrestoreisf90_(ISColoring *iscoloring,F90Array1d *ptr,int *__ierr)
82: {
83: PetscFortranAddr *is;
85: *__F90Array1dAccess(ptr,(void**)&is);if (*__ierr) return;
86: *__F90Array1dDestroy(ptr);if (*__ierr) return;
87: *__ISColoringRestoreIS(*iscoloring,(IS **)is);if (*__ierr) return;
88: *__PetscFree(is);
89: }
91: /* ---------------------------------------------------------------*/
93: void PETSC_STDCALL vecgetarrayf90_(Vec *x,F90Array1d *ptr,int *__ierr)
94: {
95: PetscScalar *fa;
96: int len;
97: *__VecGetArray(*x,&fa); if (*__ierr) return;
98: *__VecGetLocalSize(*x,&len); if (*__ierr) return;
99: *__F90Array1dCreate(fa,PETSC_SCALAR,1,len,ptr);
100: }
101: void PETSC_STDCALL vecrestorearrayf90_(Vec *x,F90Array1d *ptr,int *__ierr)
102: {
103: PetscScalar *fa;
104: *__F90Array1dAccess(ptr,(void**)&fa);if (*__ierr) return;
105: *__F90Array1dDestroy(ptr);if (*__ierr) return;
106: *__VecRestoreArray(*x,&fa);
107: }
109: void PETSC_STDCALL vecduplicatevecsf90_(Vec *v,int *m,F90Array1d *ptr,int *__ierr)
110: {
111: Vec *lV;
112: PetscFortranAddr *newvint;
113: int i;
114: *__VecDuplicateVecs(*v,*m,&lV); if (*__ierr) return;
115: *__PetscMalloc((*m)*sizeof(PetscFortranAddr),&newvint); if (*__ierr) return;
117: for (i=0; i<*m; i++) {
118: newvint[i] = (PetscFortranAddr)lV[i];
119: }
120: *__PetscFree(lV); if (*__ierr) return;
121: *__F90Array1dCreate(newvint,PETSC_FORTRANADDR,1,*m,ptr);
122: }
124: void PETSC_STDCALL vecdestroyvecsf90_(F90Array1d *ptr,int *m,int *__ierr)
125: {
126: PetscFortranAddr *vecs;
127: int i;
129: *__F90Array1dAccess(ptr,(void**)&vecs);if (*__ierr) return;
130: for (i=0; i<*m; i++) {
131: *__VecDestroy((Vec)vecs[i]);
132: if (*__ierr) return;
133: }
134: *__F90Array1dDestroy(ptr);if (*__ierr) return;
135: *__PetscFree(vecs);
136: }