Actual source code: zmatrixf90.c
1: #include <petscmat.h>
2: #include <petsc/private/ftnimpl.h>
4: #if defined(PETSC_HAVE_FORTRAN_CAPS)
5: #define matgetrow_ MATGETROW
6: #define matrestorerow_ MATRESTOREROW
7: #define matmpiaijgetseqaij_ MATMPIAIJGETSEQAIJ
8: #define matmpiaijrestoreseqaij_ MATMPIAIJRESTORESEQAIJ
9: #define matdensegetarray1d_ MATDENSEGETARRAY1D
10: #define matdenserestorearray1d_ MATDENSERESTOREARRAY1D
11: #define matdensegetarrayread1d_ MATDENSEGETARRAYREAD1D
12: #define matdenserestorearrayread1d_ MATDENSERESTOREARRAYREAD1D
13: #define matdensegetarraywrite1d_ MATDENSEGETARRAYWRITE1D
14: #define matdenserestorearraywrite1d_ MATDENSERESTOREARRAYWRITE1D
15: #define matdensegetarray2d_ MATDENSEGETARRAY2D
16: #define matdenserestorearray2d_ MATDENSERESTOREARRAY2D
17: #define matdensegetarrayread2d_ MATDENSEGETARRAYREAD2D
18: #define matdenserestorearrayread2d_ MATDENSERESTOREARRAYREAD2D
19: #define matdensegetarraywrite2d_ MATDENSEGETARRAYWRITE2D
20: #define matdenserestorearraywrite2d_ MATDENSERESTOREARRAYWRITE2D
21: #define matdensegetcolumn_ MATDENSEGETCOLUMN
22: #define matdenserestorecolumn_ MATDENSERESTORECOLUMN
23: #define matseqaijgetarray_ MATSEQAIJGETARRAY
24: #define matseqaijrestorearray_ MATSEQAIJRESTOREARRAY
25: #define matgetghosts_ MATGETGHOSTS
26: #define matgetrowij_ MATGETROWIJ
27: #define matrestorerowij_ MATRESTOREROWIJ
28: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
29: #define matgetrow_ matgetrow
30: #define matrestorerow_ matrestorerow
31: #define matmpiaijgetseqaij_ matmpiaijgetseqaij
32: #define matmpiaijrestoreseqaij_ matmpiaijrestoreseqaij
33: #define matdensegetarray1d_ matdensegetarray1d
34: #define matdenserestorearray1d_ matdenserestorearray1d
35: #define matdensegetarrayread1d_ matdensegetarrayread1d
36: #define matdenserestorearrayread1d_ matdenserestorearrayread1d
37: #define matdensegetarraywrite1d_ matdensegetarraywrite1d
38: #define matdenserestorearraywrite1d_ matdenserestorearraywrite1d
39: #define matdensegetarray2d_ matdensegetarray2d
40: #define matdenserestorearray2d_ matdenserestorearray2d
41: #define matdensegetarrayread2d_ matdensegetarrayread2d
42: #define matdenserestorearrayread2d_ matdenserestorearrayread2d
43: #define matdensegetarraywrite2d_ matdensegetarraywrite2d
44: #define matdenserestorearraywrite2d_ matdenserestorearraywrite2d
45: #define matdensegetcolumn_ matdensegetcolumn
46: #define matdenserestorecolumn_ matdenserestorecolumn
47: #define matseqaijgetarray_ matseqaijgetarray
48: #define matseqaijrestorearray_ matseqaijrestorearray
49: #define matgetghosts_ matgetghosts
50: #define matgetrowij_ matgetrowij
51: #define matrestorerowij_ matrestorerowij
52: #endif
54: PETSC_EXTERN void matgetrow_(Mat *B, PetscInt *row, PetscInt *N, F90Array1d *ia, F90Array1d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(iad) PETSC_F90_2PTR_PROTO(jad))
55: {
56: PetscInt n;
57: const PetscInt *II = NULL;
58: const PetscScalar *A = NULL;
60: if (FORTRANNULLINTEGERPOINTER(ia) && FORTRANNULLSCALARPOINTER(a)) {
61: *ierr = MatGetRow(*B, *row, &n, NULL, NULL);
62: } else if (FORTRANNULLINTEGERPOINTER(ia)) {
63: *ierr = MatGetRow(*B, *row, &n, NULL, &A);
64: } else if (FORTRANNULLSCALARPOINTER(a)) {
65: *ierr = MatGetRow(*B, *row, &n, &II, NULL);
66: } else {
67: *ierr = MatGetRow(*B, *row, &n, &II, &A);
68: }
69: if (*ierr) return;
70: if (II) { *ierr = F90Array1dCreate((void *)II, MPIU_INT, 1, n, ia PETSC_F90_2PTR_PARAM(iad)); }
71: if (A) { *ierr = F90Array1dCreate((void *)A, MPIU_SCALAR, 1, n, a PETSC_F90_2PTR_PARAM(jad)); }
72: if (!FORTRANNULLINTEGER(N)) *N = n;
73: }
74: PETSC_EXTERN void matrestorerow_(Mat *B, PetscInt *row, PetscInt *N, F90Array1d *ia, F90Array1d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(iad) PETSC_F90_2PTR_PROTO(jad))
75: {
76: const PetscInt *IA = NULL;
77: const PetscScalar *A = NULL;
78: PetscInt n;
80: if (FORTRANNULLINTEGERPOINTER(ia) && FORTRANNULLSCALARPOINTER(a)) return;
81: if (!FORTRANNULLINTEGERPOINTER(ia)) {
82: *ierr = F90Array1dAccess(ia, MPIU_INT, (void **)&IA PETSC_F90_2PTR_PARAM(iad));
83: if (*ierr) return;
84: *ierr = F90Array1dDestroy(ia, MPIU_INT PETSC_F90_2PTR_PARAM(iad));
85: if (*ierr) return;
86: }
87: if (!FORTRANNULLSCALARPOINTER(a)) {
88: *ierr = F90Array1dAccess(a, MPIU_SCALAR, (void **)&A PETSC_F90_2PTR_PARAM(jad));
89: if (*ierr) return;
90: *ierr = F90Array1dDestroy(a, MPIU_INT PETSC_F90_2PTR_PARAM(jad));
91: if (*ierr) return;
92: }
93: if (FORTRANNULLINTEGERPOINTER(ia)) {
94: *ierr = MatRestoreRow(*B, *row, &n, NULL, &A);
95: } else if (FORTRANNULLSCALARPOINTER(a)) {
96: *ierr = MatRestoreRow(*B, *row, &n, &IA, NULL);
97: } else {
98: *ierr = MatRestoreRow(*B, *row, &n, &IA, &A);
99: }
100: }
101: PETSC_EXTERN void matgetghosts_(Mat *mat, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
102: {
103: const PetscInt *ghosts;
104: PetscInt N;
106: *ierr = MatGetGhosts(*mat, &N, &ghosts);
107: if (*ierr) return;
108: *ierr = F90Array1dCreate((PetscInt *)ghosts, MPIU_INT, 1, N, ptr PETSC_F90_2PTR_PARAM(ptrd));
109: }
110: PETSC_EXTERN void matdensegetarray2d_(Mat *mat, F90Array2d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
111: {
112: PetscScalar *fa;
113: PetscInt m, N, lda;
114: *ierr = MatDenseGetArray(*mat, &fa);
115: if (*ierr) return;
116: *ierr = MatGetLocalSize(*mat, &m, NULL);
117: if (*ierr) return;
118: *ierr = MatGetSize(*mat, NULL, &N);
119: if (*ierr) return;
120: *ierr = MatDenseGetLDA(*mat, &lda);
121: if (*ierr) return;
122: if (m != lda) { // TODO: add F90Array2dLDACreate()
123: *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);
124: return;
125: }
126: *ierr = F90Array2dCreate(fa, MPIU_SCALAR, 1, m, 1, N, ptr PETSC_F90_2PTR_PARAM(ptrd));
127: }
128: PETSC_EXTERN void matdenserestorearray2d_(Mat *mat, F90Array2d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
129: {
130: PetscScalar *fa;
131: *ierr = F90Array2dAccess(ptr, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd));
132: if (*ierr) return;
133: *ierr = F90Array2dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
134: if (*ierr) return;
135: *ierr = MatDenseRestoreArray(*mat, &fa);
136: }
137: PETSC_EXTERN void matdensegetarrayread2d_(Mat *mat, F90Array2d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
138: {
139: const PetscScalar *fa;
140: PetscInt m, N, lda;
141: *ierr = MatDenseGetArrayRead(*mat, &fa);
142: if (*ierr) return;
143: *ierr = MatGetLocalSize(*mat, &m, NULL);
144: if (*ierr) return;
145: *ierr = MatGetSize(*mat, NULL, &N);
146: if (*ierr) return;
147: *ierr = MatDenseGetLDA(*mat, &lda);
148: if (*ierr) return;
149: if (m != lda) { // TODO: add F90Array2dLDACreate()
150: *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);
151: return;
152: }
153: *ierr = F90Array2dCreate((void **)fa, MPIU_SCALAR, 1, m, 1, N, ptr PETSC_F90_2PTR_PARAM(ptrd));
154: }
155: PETSC_EXTERN void matdenserestorearrayread2d_(Mat *mat, F90Array2d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
156: {
157: const PetscScalar *fa;
158: *ierr = F90Array2dAccess(ptr, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd));
159: if (*ierr) return;
160: *ierr = F90Array2dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
161: if (*ierr) return;
162: *ierr = MatDenseRestoreArrayRead(*mat, &fa);
163: }
164: PETSC_EXTERN void matdensegetarraywrite2d_(Mat *mat, F90Array2d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
165: {
166: PetscScalar *fa;
167: PetscInt m, N, lda;
168: *ierr = MatDenseGetArrayWrite(*mat, &fa);
169: if (*ierr) return;
170: *ierr = MatGetLocalSize(*mat, &m, NULL);
171: if (*ierr) return;
172: *ierr = MatGetSize(*mat, NULL, &N);
173: if (*ierr) return;
174: *ierr = MatDenseGetLDA(*mat, &lda);
175: if (*ierr) return;
176: if (m != lda) { // TODO: add F90Array2dLDACreate()
177: *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);
178: return;
179: }
180: *ierr = F90Array2dCreate(fa, MPIU_SCALAR, 1, m, 1, N, ptr PETSC_F90_2PTR_PARAM(ptrd));
181: }
182: PETSC_EXTERN void matdenserestorearraywrite2d_(Mat *mat, F90Array2d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
183: {
184: PetscScalar *fa;
185: *ierr = F90Array2dAccess(ptr, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd));
186: if (*ierr) return;
187: *ierr = F90Array2dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
188: if (*ierr) return;
189: *ierr = MatDenseRestoreArrayWrite(*mat, &fa);
190: }
191: PETSC_EXTERN void matdensegetarray1d_(Mat *mat, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
192: {
193: PetscScalar *fa;
194: PetscInt m, N, lda;
195: *ierr = MatDenseGetArray(*mat, &fa);
196: if (*ierr) return;
197: *ierr = MatGetLocalSize(*mat, &m, NULL);
198: if (*ierr) return;
199: *ierr = MatGetSize(*mat, NULL, &N);
200: if (*ierr) return;
201: *ierr = MatDenseGetLDA(*mat, &lda);
202: if (*ierr) return;
203: if (m != lda) { // TODO: add F90Array1dLDACreate()
204: *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);
205: return;
206: }
207: *ierr = F90Array1dCreate(fa, MPIU_SCALAR, 1, m * N, ptr PETSC_F90_2PTR_PARAM(ptrd));
208: }
209: PETSC_EXTERN void matdenserestorearray1d_(Mat *mat, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
210: {
211: PetscScalar *fa;
212: *ierr = F90Array1dAccess(ptr, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd));
213: if (*ierr) return;
214: *ierr = F90Array1dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
215: if (*ierr) return;
216: *ierr = MatDenseRestoreArray(*mat, &fa);
217: }
218: PETSC_EXTERN void matdensegetarrayread1d_(Mat *mat, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
219: {
220: const PetscScalar *fa;
221: PetscInt m, N, lda;
222: *ierr = MatDenseGetArrayRead(*mat, &fa);
223: if (*ierr) return;
224: *ierr = MatGetLocalSize(*mat, &m, NULL);
225: if (*ierr) return;
226: *ierr = MatGetSize(*mat, NULL, &N);
227: if (*ierr) return;
228: *ierr = MatDenseGetLDA(*mat, &lda);
229: if (*ierr) return;
230: if (m != lda) { // TODO: add F90Array1dLDACreate()
231: *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);
232: return;
233: }
234: *ierr = F90Array1dCreate((void **)fa, MPIU_SCALAR, 1, m * N, ptr PETSC_F90_2PTR_PARAM(ptrd));
235: }
236: PETSC_EXTERN void matdenserestorearrayread1d_(Mat *mat, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
237: {
238: const PetscScalar *fa;
239: *ierr = F90Array1dAccess(ptr, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd));
240: if (*ierr) return;
241: *ierr = F90Array1dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
242: if (*ierr) return;
243: *ierr = MatDenseRestoreArrayRead(*mat, &fa);
244: }
245: PETSC_EXTERN void matdensegetarraywrite1d_(Mat *mat, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
246: {
247: PetscScalar *fa;
248: PetscInt m, N, lda;
249: *ierr = MatDenseGetArrayWrite(*mat, &fa);
250: if (*ierr) return;
251: *ierr = MatGetLocalSize(*mat, &m, NULL);
252: if (*ierr) return;
253: *ierr = MatGetSize(*mat, NULL, &N);
254: if (*ierr) return;
255: *ierr = MatDenseGetLDA(*mat, &lda);
256: if (*ierr) return;
257: if (m != lda) { // TODO: add F90Array1dLDACreate()
258: *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);
259: return;
260: }
261: *ierr = F90Array1dCreate(fa, MPIU_SCALAR, 1, m * N, ptr PETSC_F90_2PTR_PARAM(ptrd));
262: }
263: PETSC_EXTERN void matdenserestorearraywrite1d_(Mat *mat, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
264: {
265: PetscScalar *fa;
266: *ierr = F90Array1dAccess(ptr, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd));
267: if (*ierr) return;
268: *ierr = F90Array1dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
269: if (*ierr) return;
270: *ierr = MatDenseRestoreArrayWrite(*mat, &fa);
271: }
272: PETSC_EXTERN void matdensegetcolumn_(Mat *mat, PetscInt *col, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
273: {
274: PetscScalar *fa;
275: PetscInt m;
276: *ierr = MatDenseGetColumn(*mat, *col, &fa);
277: if (*ierr) return;
278: *ierr = MatGetLocalSize(*mat, &m, NULL);
279: if (*ierr) return;
280: *ierr = F90Array1dCreate(fa, MPIU_SCALAR, 1, m, ptr PETSC_F90_2PTR_PARAM(ptrd));
281: }
282: PETSC_EXTERN void matdenserestorecolumn_(Mat *mat, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
283: {
284: PetscScalar *fa;
285: *ierr = F90Array1dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
286: if (*ierr) return;
287: *ierr = MatDenseRestoreColumn(*mat, &fa);
288: }
289: PETSC_EXTERN void matseqaijgetarray_(Mat *mat, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
290: {
291: PetscScalar *fa;
292: PetscInt m, n;
293: *ierr = MatSeqAIJGetArray(*mat, &fa);
294: if (*ierr) return;
295: *ierr = MatGetLocalSize(*mat, &m, &n);
296: if (*ierr) return;
297: *ierr = F90Array1dCreate(fa, MPIU_SCALAR, 1, m * n, ptr PETSC_F90_2PTR_PARAM(ptrd));
298: }
299: PETSC_EXTERN void matseqaijrestorearray_(Mat *mat, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
300: {
301: PetscScalar *fa;
302: *ierr = F90Array1dAccess(ptr, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd));
303: if (*ierr) return;
304: *ierr = F90Array1dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
305: if (*ierr) return;
306: *ierr = MatSeqAIJRestoreArray(*mat, &fa);
307: }
308: PETSC_EXTERN void matgetrowij_(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))
309: {
310: const PetscInt *IA, *JA;
311: *ierr = MatGetRowIJ(*B, *shift, *sym, *blockcompressed, n, &IA, &JA, done);
312: if (*ierr) return;
313: if (!*done) return;
314: *ierr = F90Array1dCreate((PetscInt *)IA, MPIU_INT, 1, *n + 1, ia PETSC_F90_2PTR_PARAM(iad));
315: *ierr = F90Array1dCreate((PetscInt *)JA, MPIU_INT, 1, IA[*n], ja PETSC_F90_2PTR_PARAM(jad));
316: }
317: PETSC_EXTERN void matrestorerowij_(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))
318: {
319: const PetscInt *IA, *JA;
320: *ierr = F90Array1dAccess(ia, MPIU_INT, (void **)&IA PETSC_F90_2PTR_PARAM(iad));
321: if (*ierr) return;
322: *ierr = F90Array1dDestroy(ia, MPIU_INT PETSC_F90_2PTR_PARAM(iad));
323: if (*ierr) return;
324: *ierr = F90Array1dAccess(ja, MPIU_INT, (void **)&JA PETSC_F90_2PTR_PARAM(jad));
325: if (*ierr) return;
326: *ierr = F90Array1dDestroy(ja, MPIU_INT PETSC_F90_2PTR_PARAM(jad));
327: if (*ierr) return;
328: *ierr = MatRestoreRowIJ(*B, *shift, *sym, *blockcompressed, n, &IA, &JA, done);
329: }
330: PETSC_EXTERN void matmpiaijgetseqaij_(Mat *mat, Mat *A, Mat *B, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
331: {
332: const PetscInt *fa;
333: PetscInt n;
334: *ierr = MatMPIAIJGetSeqAIJ(*mat, A, B, &fa);
335: if (*ierr) return;
336: *ierr = MatGetLocalSize(*B, NULL, &n);
337: if (*ierr) return;
338: *ierr = F90Array1dCreate((void *)fa, MPIU_INT, 1, n, ptr PETSC_F90_2PTR_PARAM(ptrd));
339: }
340: PETSC_EXTERN void matmpiaijrestoreseqaij_(Mat *mat, Mat *A, Mat *B, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
341: {
342: *ierr = F90Array1dDestroy(ptr, MPIU_INT PETSC_F90_2PTR_PARAM(ptrd));
343: if (*ierr) return;
344: }