Actual source code: zmatrixf90.c
1: #include <petscmat.h>
2: #include <petsc/private/f90impl.h>
4: #if defined(PETSC_HAVE_FORTRAN_CAPS)
5: #define matmpiaijgetseqaijf90_ MATMPIAIJGETSEQAIJF90
6: #define matmpiaijrestoreseqaijf90_ MATMPIAIJRESTORESEQAIJF90
7: #define matdensegetarrayf90_ MATDENSEGETARRAYF90
8: #define matdenserestorearrayf90_ MATDENSERESTOREARRAYF90
9: #define matdensegetarrayreadf90_ MATDENSEGETARRAYREADF90
10: #define matdenserestorearrayreadf90_ MATDENSERESTOREARRAYREADF90
11: #define matdensegetarraywritef90_ MATDENSEGETARRAYWRITEF90
12: #define matdenserestorearraywritef90_ MATDENSERESTOREARRAYWRITEF90
13: #define matdensegetcolumnf90_ MATDENSEGETCOLUMNF90
14: #define matdenserestorecolumnf90_ MATDENSERESTORECOLUMNF90
15: #define matseqaijgetarrayf90_ MATSEQAIJGETARRAYF90
16: #define matseqaijrestorearrayf90_ MATSEQAIJRESTOREARRAYF90
17: #define matgetghostsf90_ MATGETGHOSTSF90
18: #define matgetrowijf90_ MATGETROWIJF90
19: #define matrestorerowijf90_ MATRESTOREROWIJF90
20: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
21: #define matmpiaijgetseqaijf90_ matmpiaijgetseqaijf90
22: #define matmpiaijrestoreseqaijf90_ matmpiaijrestoreseqaijf90
23: #define matdensegetarrayf90_ matdensegetarrayf90
24: #define matdenserestorearrayf90_ matdenserestorearrayf90
25: #define matdensegetarrayreadf90_ matdensegetarrayreadf90
26: #define matdenserestorearrayreadf90_ matdenserestorearrayreadf90
27: #define matdensegetarraywritef90_ matdensegetarraywritef90
28: #define matdenserestorearraywritef90_ matdenserestorearraywritef90
29: #define matdensegetcolumnf90_ matdensegetcolumnf90
30: #define matdenserestorecolumnf90_ matdenserestorecolumnf90
31: #define matseqaijgetarrayf90_ matseqaijgetarrayf90
32: #define matseqaijrestorearrayf90_ matseqaijrestorearrayf90
33: #define matgetghostsf90_ matgetghostsf90
34: #define matgetrowijf90_ matgetrowijf90
35: #define matrestorerowijf90_ matrestorerowijf90
36: #endif
38: PETSC_EXTERN void matgetghostsf90_(Mat *mat, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
39: {
40: const PetscInt *ghosts;
41: PetscInt N;
43: *ierr = MatGetGhosts(*mat, &N, &ghosts);
44: if (*ierr) return;
45: *ierr = F90Array1dCreate((PetscInt *)ghosts, MPIU_INT, 1, N, ptr PETSC_F90_2PTR_PARAM(ptrd));
46: }
47: PETSC_EXTERN void matdensegetarrayf90_(Mat *mat, F90Array2d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
48: {
49: PetscScalar *fa;
50: PetscInt m, N, lda;
51: *ierr = MatDenseGetArray(*mat, &fa);
52: if (*ierr) return;
53: *ierr = MatGetLocalSize(*mat, &m, NULL);
54: if (*ierr) return;
55: *ierr = MatGetSize(*mat, NULL, &N);
56: if (*ierr) return;
57: *ierr = MatDenseGetLDA(*mat, &lda);
58: if (*ierr) return;
59: if (m != lda) { // TODO: add F90Array2dLDACreate()
60: *ierr = PetscError(((PetscObject)*mat)->comm, __LINE__, PETSC_FUNCTION_NAME, __FILE__, PETSC_ERR_ARG_BADPTR, PETSC_ERROR_INITIAL, "Array lda %" PetscInt_FMT " must match number of local rows %" PetscInt_FMT, lda, m);
61: return;
62: }
63: *ierr = F90Array2dCreate(fa, MPIU_SCALAR, 1, m, 1, N, ptr PETSC_F90_2PTR_PARAM(ptrd));
64: }
65: PETSC_EXTERN void matdenserestorearrayf90_(Mat *mat, F90Array2d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
66: {
67: PetscScalar *fa;
68: *ierr = F90Array2dAccess(ptr, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd));
69: if (*ierr) return;
70: *ierr = F90Array2dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
71: if (*ierr) return;
72: *ierr = MatDenseRestoreArray(*mat, &fa);
73: }
74: PETSC_EXTERN void matdensegetarrayreadf90_(Mat *mat, F90Array2d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
75: {
76: const PetscScalar *fa;
77: PetscInt m, N, lda;
78: *ierr = MatDenseGetArrayRead(*mat, &fa);
79: if (*ierr) return;
80: *ierr = MatGetLocalSize(*mat, &m, NULL);
81: if (*ierr) return;
82: *ierr = MatGetSize(*mat, NULL, &N);
83: if (*ierr) return;
84: *ierr = MatDenseGetLDA(*mat, &lda);
85: if (*ierr) return;
86: if (m != lda) { // TODO: add F90Array2dLDACreate()
87: *ierr = PetscError(((PetscObject)*mat)->comm, __LINE__, PETSC_FUNCTION_NAME, __FILE__, PETSC_ERR_ARG_BADPTR, PETSC_ERROR_INITIAL, "Array lda %" PetscInt_FMT " must match number of local rows %" PetscInt_FMT, lda, m);
88: return;
89: }
90: *ierr = F90Array2dCreate((void **)fa, MPIU_SCALAR, 1, m, 1, N, ptr PETSC_F90_2PTR_PARAM(ptrd));
91: }
92: PETSC_EXTERN void matdenserestorearrayreadf90_(Mat *mat, F90Array2d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
93: {
94: const PetscScalar *fa;
95: *ierr = F90Array2dAccess(ptr, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd));
96: if (*ierr) return;
97: *ierr = F90Array2dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
98: if (*ierr) return;
99: *ierr = MatDenseRestoreArrayRead(*mat, &fa);
100: }
101: PETSC_EXTERN void matdensegetarraywritef90_(Mat *mat, F90Array2d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
102: {
103: PetscScalar *fa;
104: PetscInt m, N, lda;
105: *ierr = MatDenseGetArrayWrite(*mat, &fa);
106: if (*ierr) return;
107: *ierr = MatGetLocalSize(*mat, &m, NULL);
108: if (*ierr) return;
109: *ierr = MatGetSize(*mat, NULL, &N);
110: if (*ierr) return;
111: *ierr = MatDenseGetLDA(*mat, &lda);
112: if (*ierr) return;
113: if (m != lda) { // TODO: add F90Array2dLDACreate()
114: *ierr = PetscError(((PetscObject)*mat)->comm, __LINE__, PETSC_FUNCTION_NAME, __FILE__, PETSC_ERR_ARG_BADPTR, PETSC_ERROR_INITIAL, "Array lda %" PetscInt_FMT " must match number of local rows %" PetscInt_FMT, lda, m);
115: return;
116: }
117: *ierr = F90Array2dCreate(fa, MPIU_SCALAR, 1, m, 1, N, ptr PETSC_F90_2PTR_PARAM(ptrd));
118: }
119: PETSC_EXTERN void matdenserestorearraywritef90_(Mat *mat, F90Array2d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
120: {
121: PetscScalar *fa;
122: *ierr = F90Array2dAccess(ptr, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd));
123: if (*ierr) return;
124: *ierr = F90Array2dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
125: if (*ierr) return;
126: *ierr = MatDenseRestoreArrayWrite(*mat, &fa);
127: }
128: PETSC_EXTERN void matdensegetcolumnf90_(Mat *mat, PetscInt *col, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
129: {
130: PetscScalar *fa;
131: PetscInt m;
132: *ierr = MatDenseGetColumn(*mat, *col, &fa);
133: if (*ierr) return;
134: *ierr = MatGetLocalSize(*mat, &m, NULL);
135: if (*ierr) return;
136: *ierr = F90Array1dCreate(fa, MPIU_SCALAR, 1, m, ptr PETSC_F90_2PTR_PARAM(ptrd));
137: }
138: PETSC_EXTERN void matdenserestorecolumnf90_(Mat *mat, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
139: {
140: PetscScalar *fa;
141: *ierr = F90Array1dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
142: if (*ierr) return;
143: *ierr = MatDenseRestoreColumn(*mat, &fa);
144: }
145: PETSC_EXTERN void matseqaijgetarrayf90_(Mat *mat, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
146: {
147: PetscScalar *fa;
148: PetscInt m, n;
149: *ierr = MatSeqAIJGetArray(*mat, &fa);
150: if (*ierr) return;
151: *ierr = MatGetLocalSize(*mat, &m, &n);
152: if (*ierr) return;
153: *ierr = F90Array1dCreate(fa, MPIU_SCALAR, 1, m * n, ptr PETSC_F90_2PTR_PARAM(ptrd));
154: }
155: PETSC_EXTERN void matseqaijrestorearrayf90_(Mat *mat, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
156: {
157: PetscScalar *fa;
158: *ierr = F90Array1dAccess(ptr, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd));
159: if (*ierr) return;
160: *ierr = F90Array1dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
161: if (*ierr) return;
162: *ierr = MatSeqAIJRestoreArray(*mat, &fa);
163: }
164: PETSC_EXTERN void matgetrowijf90_(Mat *B, PetscInt *shift, PetscBool *sym, PetscBool *blockcompressed, PetscInt *n, F90Array1d *ia, F90Array1d *ja, PetscBool *done, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(iad) PETSC_F90_2PTR_PROTO(jad))
165: {
166: const PetscInt *IA, *JA;
167: *ierr = MatGetRowIJ(*B, *shift, *sym, *blockcompressed, n, &IA, &JA, done);
168: if (*ierr) return;
169: if (!*done) return;
170: *ierr = F90Array1dCreate((PetscInt *)IA, MPIU_INT, 1, *n + 1, ia PETSC_F90_2PTR_PARAM(iad));
171: *ierr = F90Array1dCreate((PetscInt *)JA, MPIU_INT, 1, IA[*n], ja PETSC_F90_2PTR_PARAM(jad));
172: }
173: PETSC_EXTERN void matrestorerowijf90_(Mat *B, PetscInt *shift, PetscBool *sym, PetscBool *blockcompressed, PetscInt *n, F90Array1d *ia, F90Array1d *ja, PetscBool *done, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(iad) PETSC_F90_2PTR_PROTO(jad))
174: {
175: const PetscInt *IA, *JA;
176: *ierr = F90Array1dAccess(ia, MPIU_INT, (void **)&IA PETSC_F90_2PTR_PARAM(iad));
177: if (*ierr) return;
178: *ierr = F90Array1dDestroy(ia, MPIU_INT PETSC_F90_2PTR_PARAM(iad));
179: if (*ierr) return;
180: *ierr = F90Array1dAccess(ja, MPIU_INT, (void **)&JA PETSC_F90_2PTR_PARAM(jad));
181: if (*ierr) return;
182: *ierr = F90Array1dDestroy(ja, MPIU_INT PETSC_F90_2PTR_PARAM(jad));
183: if (*ierr) return;
184: *ierr = MatRestoreRowIJ(*B, *shift, *sym, *blockcompressed, n, &IA, &JA, done);
185: }
186: PETSC_EXTERN void matmpiaijgetseqaijf90_(Mat *mat, Mat *A, Mat *B, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
187: {
188: const PetscInt *fa;
189: PetscInt n;
190: *ierr = MatMPIAIJGetSeqAIJ(*mat, A, B, &fa);
191: if (*ierr) return;
192: *ierr = MatGetLocalSize(*B, NULL, &n);
193: if (*ierr) return;
194: *ierr = F90Array1dCreate((void *)fa, MPIU_INT, 1, n, ptr PETSC_F90_2PTR_PARAM(ptrd));
195: }
196: PETSC_EXTERN void matmpiaijrestoreseqaijf90_(Mat *mat, Mat *A, Mat *B, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
197: {
198: *ierr = F90Array1dDestroy(ptr, MPIU_INT PETSC_F90_2PTR_PARAM(ptrd));
199: if (*ierr) return;
200: }