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)) {
81: *ierr = MatRestoreRow(*B, *row, &n, NULL, NULL);
82: return;
83: }
84: if (!FORTRANNULLINTEGERPOINTER(ia)) {
85: *ierr = F90Array1dAccess(ia, MPIU_INT, (void **)&IA PETSC_F90_2PTR_PARAM(iad));
86: if (*ierr) return;
87: *ierr = F90Array1dDestroy(ia, MPIU_INT PETSC_F90_2PTR_PARAM(iad));
88: if (*ierr) return;
89: }
90: if (!FORTRANNULLSCALARPOINTER(a)) {
91: *ierr = F90Array1dAccess(a, MPIU_SCALAR, (void **)&A PETSC_F90_2PTR_PARAM(jad));
92: if (*ierr) return;
93: *ierr = F90Array1dDestroy(a, MPIU_INT PETSC_F90_2PTR_PARAM(jad));
94: if (*ierr) return;
95: }
96: if (FORTRANNULLINTEGERPOINTER(ia)) {
97: *ierr = MatRestoreRow(*B, *row, &n, NULL, &A);
98: } else if (FORTRANNULLSCALARPOINTER(a)) {
99: *ierr = MatRestoreRow(*B, *row, &n, &IA, NULL);
100: } else {
101: *ierr = MatRestoreRow(*B, *row, &n, &IA, &A);
102: }
103: }
104: PETSC_EXTERN void matgetghosts_(Mat *mat, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
105: {
106: const PetscInt *ghosts;
107: PetscInt N;
109: *ierr = MatGetGhosts(*mat, &N, &ghosts);
110: if (*ierr) return;
111: *ierr = F90Array1dCreate((PetscInt *)ghosts, MPIU_INT, 1, N, ptr PETSC_F90_2PTR_PARAM(ptrd));
112: }
113: PETSC_EXTERN void matdensegetarray2d_(Mat *mat, F90Array2d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
114: {
115: PetscScalar *fa;
116: PetscInt m, N, lda;
117: *ierr = MatDenseGetArray(*mat, &fa);
118: if (*ierr) return;
119: *ierr = MatGetLocalSize(*mat, &m, NULL);
120: if (*ierr) return;
121: *ierr = MatGetSize(*mat, NULL, &N);
122: if (*ierr) return;
123: *ierr = MatDenseGetLDA(*mat, &lda);
124: if (*ierr) return;
125: if (m != lda) { // TODO: add F90Array2dLDACreate()
126: *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);
127: return;
128: }
129: *ierr = F90Array2dCreate(fa, MPIU_SCALAR, 1, m, 1, N, ptr PETSC_F90_2PTR_PARAM(ptrd));
130: }
131: PETSC_EXTERN void matdenserestorearray2d_(Mat *mat, F90Array2d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
132: {
133: PetscScalar *fa;
134: *ierr = F90Array2dAccess(ptr, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd));
135: if (*ierr) return;
136: *ierr = F90Array2dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
137: if (*ierr) return;
138: *ierr = MatDenseRestoreArray(*mat, &fa);
139: }
140: PETSC_EXTERN void matdensegetarrayread2d_(Mat *mat, F90Array2d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
141: {
142: const PetscScalar *fa;
143: PetscInt m, N, lda;
144: *ierr = MatDenseGetArrayRead(*mat, &fa);
145: if (*ierr) return;
146: *ierr = MatGetLocalSize(*mat, &m, NULL);
147: if (*ierr) return;
148: *ierr = MatGetSize(*mat, NULL, &N);
149: if (*ierr) return;
150: *ierr = MatDenseGetLDA(*mat, &lda);
151: if (*ierr) return;
152: if (m != lda) { // TODO: add F90Array2dLDACreate()
153: *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);
154: return;
155: }
156: *ierr = F90Array2dCreate((void **)fa, MPIU_SCALAR, 1, m, 1, N, ptr PETSC_F90_2PTR_PARAM(ptrd));
157: }
158: PETSC_EXTERN void matdenserestorearrayread2d_(Mat *mat, F90Array2d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
159: {
160: const PetscScalar *fa;
161: *ierr = F90Array2dAccess(ptr, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd));
162: if (*ierr) return;
163: *ierr = F90Array2dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
164: if (*ierr) return;
165: *ierr = MatDenseRestoreArrayRead(*mat, &fa);
166: }
167: PETSC_EXTERN void matdensegetarraywrite2d_(Mat *mat, F90Array2d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
168: {
169: PetscScalar *fa;
170: PetscInt m, N, lda;
171: *ierr = MatDenseGetArrayWrite(*mat, &fa);
172: if (*ierr) return;
173: *ierr = MatGetLocalSize(*mat, &m, NULL);
174: if (*ierr) return;
175: *ierr = MatGetSize(*mat, NULL, &N);
176: if (*ierr) return;
177: *ierr = MatDenseGetLDA(*mat, &lda);
178: if (*ierr) return;
179: if (m != lda) { // TODO: add F90Array2dLDACreate()
180: *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);
181: return;
182: }
183: *ierr = F90Array2dCreate(fa, MPIU_SCALAR, 1, m, 1, N, ptr PETSC_F90_2PTR_PARAM(ptrd));
184: }
185: PETSC_EXTERN void matdenserestorearraywrite2d_(Mat *mat, F90Array2d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
186: {
187: PetscScalar *fa;
188: *ierr = F90Array2dAccess(ptr, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd));
189: if (*ierr) return;
190: *ierr = F90Array2dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
191: if (*ierr) return;
192: *ierr = MatDenseRestoreArrayWrite(*mat, &fa);
193: }
194: PETSC_EXTERN void matdensegetarray1d_(Mat *mat, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
195: {
196: PetscScalar *fa;
197: PetscInt m, N, lda;
198: *ierr = MatDenseGetArray(*mat, &fa);
199: if (*ierr) return;
200: *ierr = MatGetLocalSize(*mat, &m, NULL);
201: if (*ierr) return;
202: *ierr = MatGetSize(*mat, NULL, &N);
203: if (*ierr) return;
204: *ierr = MatDenseGetLDA(*mat, &lda);
205: if (*ierr) return;
206: if (m != lda) { // TODO: add F90Array1dLDACreate()
207: *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);
208: return;
209: }
210: *ierr = F90Array1dCreate(fa, MPIU_SCALAR, 1, m * N, ptr PETSC_F90_2PTR_PARAM(ptrd));
211: }
212: PETSC_EXTERN void matdenserestorearray1d_(Mat *mat, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
213: {
214: PetscScalar *fa;
215: *ierr = F90Array1dAccess(ptr, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd));
216: if (*ierr) return;
217: *ierr = F90Array1dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
218: if (*ierr) return;
219: *ierr = MatDenseRestoreArray(*mat, &fa);
220: }
221: PETSC_EXTERN void matdensegetarrayread1d_(Mat *mat, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
222: {
223: const PetscScalar *fa;
224: PetscInt m, N, lda;
225: *ierr = MatDenseGetArrayRead(*mat, &fa);
226: if (*ierr) return;
227: *ierr = MatGetLocalSize(*mat, &m, NULL);
228: if (*ierr) return;
229: *ierr = MatGetSize(*mat, NULL, &N);
230: if (*ierr) return;
231: *ierr = MatDenseGetLDA(*mat, &lda);
232: if (*ierr) return;
233: if (m != lda) { // TODO: add F90Array1dLDACreate()
234: *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);
235: return;
236: }
237: *ierr = F90Array1dCreate((void **)fa, MPIU_SCALAR, 1, m * N, ptr PETSC_F90_2PTR_PARAM(ptrd));
238: }
239: PETSC_EXTERN void matdenserestorearrayread1d_(Mat *mat, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
240: {
241: const PetscScalar *fa;
242: *ierr = F90Array1dAccess(ptr, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd));
243: if (*ierr) return;
244: *ierr = F90Array1dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
245: if (*ierr) return;
246: *ierr = MatDenseRestoreArrayRead(*mat, &fa);
247: }
248: PETSC_EXTERN void matdensegetarraywrite1d_(Mat *mat, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
249: {
250: PetscScalar *fa;
251: PetscInt m, N, lda;
252: *ierr = MatDenseGetArrayWrite(*mat, &fa);
253: if (*ierr) return;
254: *ierr = MatGetLocalSize(*mat, &m, NULL);
255: if (*ierr) return;
256: *ierr = MatGetSize(*mat, NULL, &N);
257: if (*ierr) return;
258: *ierr = MatDenseGetLDA(*mat, &lda);
259: if (*ierr) return;
260: if (m != lda) { // TODO: add F90Array1dLDACreate()
261: *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);
262: return;
263: }
264: *ierr = F90Array1dCreate(fa, MPIU_SCALAR, 1, m * N, ptr PETSC_F90_2PTR_PARAM(ptrd));
265: }
266: PETSC_EXTERN void matdenserestorearraywrite1d_(Mat *mat, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
267: {
268: PetscScalar *fa;
269: *ierr = F90Array1dAccess(ptr, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd));
270: if (*ierr) return;
271: *ierr = F90Array1dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
272: if (*ierr) return;
273: *ierr = MatDenseRestoreArrayWrite(*mat, &fa);
274: }
275: PETSC_EXTERN void matdensegetcolumn_(Mat *mat, PetscInt *col, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
276: {
277: PetscScalar *fa;
278: PetscInt m;
279: *ierr = MatDenseGetColumn(*mat, *col, &fa);
280: if (*ierr) return;
281: *ierr = MatGetLocalSize(*mat, &m, NULL);
282: if (*ierr) return;
283: *ierr = F90Array1dCreate(fa, MPIU_SCALAR, 1, m, ptr PETSC_F90_2PTR_PARAM(ptrd));
284: }
285: PETSC_EXTERN void matdenserestorecolumn_(Mat *mat, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
286: {
287: PetscScalar *fa;
288: *ierr = F90Array1dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
289: if (*ierr) return;
290: *ierr = MatDenseRestoreColumn(*mat, &fa);
291: }
292: PETSC_EXTERN void matseqaijgetarray_(Mat *mat, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
293: {
294: PetscScalar *fa;
295: PetscInt m, n;
296: *ierr = MatSeqAIJGetArray(*mat, &fa);
297: if (*ierr) return;
298: *ierr = MatGetLocalSize(*mat, &m, &n);
299: if (*ierr) return;
300: *ierr = F90Array1dCreate(fa, MPIU_SCALAR, 1, m * n, ptr PETSC_F90_2PTR_PARAM(ptrd));
301: }
302: PETSC_EXTERN void matseqaijrestorearray_(Mat *mat, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
303: {
304: PetscScalar *fa;
305: *ierr = F90Array1dAccess(ptr, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd));
306: if (*ierr) return;
307: *ierr = F90Array1dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
308: if (*ierr) return;
309: *ierr = MatSeqAIJRestoreArray(*mat, &fa);
310: }
311: 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))
312: {
313: const PetscInt *IA, *JA;
314: *ierr = MatGetRowIJ(*B, *shift, *sym, *blockcompressed, n, &IA, &JA, done);
315: if (*ierr) return;
316: if (!*done) return;
317: *ierr = F90Array1dCreate((PetscInt *)IA, MPIU_INT, 1, *n + 1, ia PETSC_F90_2PTR_PARAM(iad));
318: *ierr = F90Array1dCreate((PetscInt *)JA, MPIU_INT, 1, IA[*n], ja PETSC_F90_2PTR_PARAM(jad));
319: }
320: 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))
321: {
322: const PetscInt *IA, *JA;
323: *ierr = F90Array1dAccess(ia, MPIU_INT, (void **)&IA PETSC_F90_2PTR_PARAM(iad));
324: if (*ierr) return;
325: *ierr = F90Array1dDestroy(ia, MPIU_INT PETSC_F90_2PTR_PARAM(iad));
326: if (*ierr) return;
327: *ierr = F90Array1dAccess(ja, MPIU_INT, (void **)&JA PETSC_F90_2PTR_PARAM(jad));
328: if (*ierr) return;
329: *ierr = F90Array1dDestroy(ja, MPIU_INT PETSC_F90_2PTR_PARAM(jad));
330: if (*ierr) return;
331: *ierr = MatRestoreRowIJ(*B, *shift, *sym, *blockcompressed, n, &IA, &JA, done);
332: }
333: PETSC_EXTERN void matmpiaijgetseqaij_(Mat *mat, Mat *A, Mat *B, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
334: {
335: const PetscInt *fa;
336: PetscInt n;
337: *ierr = MatMPIAIJGetSeqAIJ(*mat, A, B, &fa);
338: if (*ierr) return;
339: *ierr = MatGetLocalSize(*B, NULL, &n);
340: if (*ierr) return;
341: *ierr = F90Array1dCreate((void *)fa, MPIU_INT, 1, n, ptr PETSC_F90_2PTR_PARAM(ptrd));
342: }
343: PETSC_EXTERN void matmpiaijrestoreseqaij_(Mat *mat, Mat *A, Mat *B, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
344: {
345: *ierr = F90Array1dDestroy(ptr, MPIU_INT PETSC_F90_2PTR_PARAM(ptrd));
346: if (*ierr) return;
347: }