Actual source code: zmatrixf.c
1: #include <petsc/private/fortranimpl.h>
2: #include <petsc/private/f90impl.h>
3: #include <petscmat.h>
4: #include <petscviewer.h>
6: #if defined(PETSC_HAVE_FORTRAN_CAPS)
7: #define matdestroymatrices_ MATDESTROYMATRICES
8: #define matdestroysubmatrices_ MATDESTROYSUBMATRICES
9: #define matgetrowij_ MATGETROWIJ
10: #define matrestorerowij_ MATRESTOREROWIJ
11: #define matgetrow_ MATGETROW
12: #define matrestorerow_ MATRESTOREROW
13: #define matseqaijgetarray_ MATSEQAIJGETARRAY
14: #define matseqaijrestorearray_ MATSEQAIJRESTOREARRAY
15: #define matdensegetarray_ MATDENSEGETARRAY
16: #define matdensegetarrayread_ MATDENSEGETARRAYREAD
17: #define matdenserestorearray_ MATDENSERESTOREARRAY
18: #define matdenserestorearrayread_ MATDENSERESTOREARRAYREAD
19: #define matcreatesubmatrices_ MATCREATESUBMATRICES
20: #define matcreatesubmatricesmpi_ MATCREATESUBMATRICESMPI
21: #define matnullspacesetfunction_ MATNULLSPACESETFUNCTION
22: #define matfindnonzerorows_ MATFINDNONZEROROWS
23: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
24: #define matdestroymatrices_ matdestroymatrices
25: #define matdestroysubmatrices_ matdestroysubmatrices
26: #define matgetrowij_ matgetrowij
27: #define matrestorerowij_ matrestorerowij
28: #define matgetrow_ matgetrow
29: #define matrestorerow_ matrestorerow
30: #define matseqaijgetarray_ matseqaijgetarray
31: #define matseqaijrestorearray_ matseqaijrestorearray
32: #define matdensegetarray_ matdensegetarray
33: #define matdensegetarrayread_ matdensegetarrayread
34: #define matdenserestorearray_ matdenserestorearray
35: #define matdenserestorearrayread_ matdenserestorearrayread
36: #define matcreatesubmatrices_ matcreatesubmatrices
37: #define matcreatesubmatricesmpi_ matcreatesubmatricesmpi
38: #define matnullspacesetfunction_ matnullspacesetfunction
39: #define matfindnonzerorows_ matfindnonzerorows
40: #endif
42: static PetscErrorCode ournullfunction(MatNullSpace sp, Vec x, void *ctx)
43: {
44: PetscCallFortranVoidFunction((*(void (*)(MatNullSpace *, Vec *, void *, PetscErrorCode *))(((PetscObject)sp)->fortran_func_pointers[0]))(&sp, &x, ctx, &ierr));
45: return PETSC_SUCCESS;
46: }
48: PETSC_EXTERN void matnullspacesetfunction_(MatNullSpace *sp, PetscErrorCode (*rem)(MatNullSpace, Vec, void *), void *ctx, PetscErrorCode *ierr)
49: {
50: PetscObjectAllocateFortranPointers(*sp, 1);
51: ((PetscObject)*sp)->fortran_func_pointers[0] = (PetscVoidFn *)rem;
53: *ierr = MatNullSpaceSetFunction(*sp, ournullfunction, ctx);
54: }
56: PETSC_EXTERN void matgetrowij_(Mat *B, PetscInt *shift, PetscBool *sym, PetscBool *blockcompressed, PetscInt *n, PetscInt *ia, size_t *iia, PetscInt *ja, size_t *jja, PetscBool *done, PetscErrorCode *ierr)
57: {
58: const PetscInt *IA, *JA;
59: *ierr = MatGetRowIJ(*B, *shift, *sym, *blockcompressed, n, &IA, &JA, done);
60: if (*ierr) return;
61: *iia = PetscIntAddressToFortran(ia, (PetscInt *)IA);
62: *jja = PetscIntAddressToFortran(ja, (PetscInt *)JA);
63: }
65: PETSC_EXTERN void matrestorerowij_(Mat *B, PetscInt *shift, PetscBool *sym, PetscBool *blockcompressed, PetscInt *n, PetscInt *ia, size_t *iia, PetscInt *ja, size_t *jja, PetscBool *done, PetscErrorCode *ierr)
66: {
67: const PetscInt *IA = PetscIntAddressFromFortran(ia, *iia), *JA = PetscIntAddressFromFortran(ja, *jja);
68: *ierr = MatRestoreRowIJ(*B, *shift, *sym, *blockcompressed, n, &IA, &JA, done);
69: }
71: /*
72: This is a poor way of storing the column and value pointers
73: generated by MatGetRow() to be returned with MatRestoreRow()
74: but there is not natural,good place else to store them. Hence
75: Fortran programmers can only have one outstanding MatGetRows()
76: at a time.
77: */
78: static int matgetrowactive = 0;
79: static const PetscInt *my_ocols = NULL;
80: static const PetscScalar *my_ovals = NULL;
82: PETSC_EXTERN void matgetrow_(Mat *mat, PetscInt *row, PetscInt *ncols, PetscInt *cols, PetscScalar *vals, PetscErrorCode *ierr)
83: {
84: const PetscInt **oocols = &my_ocols;
85: const PetscScalar **oovals = &my_ovals;
87: if (matgetrowactive) {
88: *ierr = PetscError(PETSC_COMM_SELF, __LINE__, "MatGetRow_Fortran", __FILE__, PETSC_ERR_ARG_WRONGSTATE, PETSC_ERROR_INITIAL, "Cannot have two MatGetRow() active simultaneously\n\
89: call MatRestoreRow() before calling MatGetRow() a second time");
90: *ierr = PETSC_ERR_ARG_WRONGSTATE;
91: return;
92: }
94: CHKFORTRANNULLINTEGER(cols);
95: if (!cols) oocols = NULL;
96: CHKFORTRANNULLSCALAR(vals);
97: if (!vals) oovals = NULL;
99: *ierr = MatGetRow(*mat, *row, ncols, oocols, oovals);
100: if (*ierr) return;
102: if (oocols) {
103: *ierr = PetscArraycpy(cols, my_ocols, *ncols);
104: if (*ierr) return;
105: }
106: if (oovals) {
107: *ierr = PetscArraycpy(vals, my_ovals, *ncols);
108: if (*ierr) return;
109: }
110: matgetrowactive = 1;
111: }
113: PETSC_EXTERN void matrestorerow_(Mat *mat, PetscInt *row, PetscInt *ncols, PetscInt *cols, PetscScalar *vals, PetscErrorCode *ierr)
114: {
115: const PetscInt **oocols = &my_ocols;
116: const PetscScalar **oovals = &my_ovals;
118: if (!matgetrowactive) {
119: *ierr = PetscError(PETSC_COMM_SELF, __LINE__, "MatRestoreRow_Fortran", __FILE__, PETSC_ERR_ARG_WRONGSTATE, PETSC_ERROR_INITIAL, "Must call MatGetRow() first");
120: *ierr = PETSC_ERR_ARG_WRONGSTATE;
121: return;
122: }
123: CHKFORTRANNULLINTEGER(cols);
124: if (!cols) oocols = NULL;
125: CHKFORTRANNULLSCALAR(vals);
126: if (!vals) oovals = NULL;
128: *ierr = MatRestoreRow(*mat, *row, ncols, oocols, oovals);
129: matgetrowactive = 0;
130: }
132: PETSC_EXTERN void matseqaijgetarray_(Mat *mat, PetscScalar *fa, size_t *ia, PetscErrorCode *ierr)
133: {
134: PetscScalar *mm;
135: PetscInt m, n;
137: *ierr = MatSeqAIJGetArray(*mat, &mm);
138: if (*ierr) return;
139: *ierr = MatGetSize(*mat, &m, &n);
140: if (*ierr) return;
141: *ierr = PetscScalarAddressToFortran((PetscObject)*mat, 1, fa, mm, m * n, ia);
142: if (*ierr) return;
143: }
145: PETSC_EXTERN void matseqaijrestorearray_(Mat *mat, PetscScalar *fa, size_t *ia, PetscErrorCode *ierr)
146: {
147: PetscScalar *lx;
148: PetscInt m, n;
150: *ierr = MatGetSize(*mat, &m, &n);
151: if (*ierr) return;
152: *ierr = PetscScalarAddressFromFortran((PetscObject)*mat, fa, *ia, m * n, &lx);
153: if (*ierr) return;
154: *ierr = MatSeqAIJRestoreArray(*mat, &lx);
155: if (*ierr) return;
156: }
158: PETSC_EXTERN void matdensegetarray_(Mat *mat, PetscScalar *fa, size_t *ia, PetscErrorCode *ierr)
159: {
160: PetscScalar *mm;
161: PetscInt m, n;
163: *ierr = MatDenseGetArray(*mat, &mm);
164: if (*ierr) return;
165: *ierr = MatGetSize(*mat, &m, &n);
166: if (*ierr) return;
167: *ierr = PetscScalarAddressToFortran((PetscObject)*mat, 1, fa, mm, m * n, ia);
168: if (*ierr) return;
169: }
171: PETSC_EXTERN void matdenserestorearray_(Mat *mat, PetscScalar *fa, size_t *ia, PetscErrorCode *ierr)
172: {
173: PetscScalar *lx;
174: PetscInt m, n;
176: *ierr = MatGetSize(*mat, &m, &n);
177: if (*ierr) return;
178: *ierr = PetscScalarAddressFromFortran((PetscObject)*mat, fa, *ia, m * n, &lx);
179: if (*ierr) return;
180: *ierr = MatDenseRestoreArray(*mat, &lx);
181: if (*ierr) return;
182: }
184: PETSC_EXTERN void matdensegetarrayread_(Mat *mat, PetscScalar *fa, size_t *ia, PetscErrorCode *ierr)
185: {
186: const PetscScalar *mm;
187: PetscInt m, n;
189: *ierr = MatDenseGetArrayRead(*mat, &mm);
190: if (*ierr) return;
191: *ierr = MatGetSize(*mat, &m, &n);
192: if (*ierr) return;
193: *ierr = PetscScalarAddressToFortran((PetscObject)*mat, 1, fa, (PetscScalar *)mm, m * n, ia);
194: if (*ierr) return;
195: }
197: PETSC_EXTERN void matdenserestorearrayread_(Mat *mat, PetscScalar *fa, size_t *ia, PetscErrorCode *ierr)
198: {
199: const PetscScalar *lx;
200: PetscInt m, n;
202: *ierr = MatGetSize(*mat, &m, &n);
203: if (*ierr) return;
204: *ierr = PetscScalarAddressFromFortran((PetscObject)*mat, fa, *ia, m * n, (PetscScalar **)&lx);
205: if (*ierr) return;
206: *ierr = MatDenseRestoreArrayRead(*mat, &lx);
207: if (*ierr) return;
208: }
210: /*
211: MatCreateSubmatrices() is slightly different from C since the
212: Fortran provides the array to hold the submatrix objects,while in C that
213: array is allocated by the MatCreateSubmatrices()
214: */
215: PETSC_EXTERN void matcreatesubmatrices_(Mat *mat, PetscInt *n, IS *isrow, IS *iscol, MatReuse *scall, Mat *smat, PetscErrorCode *ierr)
216: {
217: Mat *lsmat;
218: PetscInt i;
220: if (*scall == MAT_INITIAL_MATRIX) {
221: *ierr = MatCreateSubMatrices(*mat, *n, isrow, iscol, *scall, &lsmat);
222: for (i = 0; i <= *n; i++) { /* lsmat[*n] might be a dummy matrix for saving data structure */
223: smat[i] = lsmat[i];
224: }
225: *ierr = PetscFree(lsmat);
226: } else {
227: *ierr = MatCreateSubMatrices(*mat, *n, isrow, iscol, *scall, &smat);
228: }
229: }
231: /*
232: MatCreateSubmatrices() is slightly different from C since the
233: Fortran provides the array to hold the submatrix objects,while in C that
234: array is allocated by the MatCreateSubmatrices()
235: */
236: PETSC_EXTERN void matcreatesubmatricesmpi_(Mat *mat, PetscInt *n, IS *isrow, IS *iscol, MatReuse *scall, Mat *smat, PetscErrorCode *ierr)
237: {
238: Mat *lsmat;
239: PetscInt i;
241: if (*scall == MAT_INITIAL_MATRIX) {
242: *ierr = MatCreateSubMatricesMPI(*mat, *n, isrow, iscol, *scall, &lsmat);
243: for (i = 0; i <= *n; i++) { /* lsmat[*n] might be a dummy matrix for saving data structure */
244: smat[i] = lsmat[i];
245: }
246: *ierr = PetscFree(lsmat);
247: } else {
248: *ierr = MatCreateSubMatricesMPI(*mat, *n, isrow, iscol, *scall, &smat);
249: }
250: }
252: /*
253: MatDestroyMatrices() is slightly different from C since the
254: Fortran does not free the array of matrix objects, while in C that
255: the array is freed
256: */
257: PETSC_EXTERN void matdestroymatrices_(PetscInt *n, Mat *smat, PetscErrorCode *ierr)
258: {
259: PetscInt i;
261: for (i = 0; i < *n; i++) {
262: PETSC_FORTRAN_OBJECT_F_DESTROYED_TO_C_NULL(&smat[i]);
263: *ierr = MatDestroy(&smat[i]);
264: if (*ierr) return;
265: PETSC_FORTRAN_OBJECT_C_NULL_TO_F_DESTROYED(&smat[i]);
266: }
267: }
269: /*
270: MatDestroySubMatrices() is slightly different from C since the
271: Fortran provides the array to hold the submatrix objects, while in C that
272: array is allocated by the MatCreateSubmatrices()
274: An extra matrix may be stored at the end of the array, hence the check see
275: MatDestroySubMatrices_Dummy()
276: */
277: PETSC_EXTERN void matdestroysubmatrices_(PetscInt *n, Mat *smat, PetscErrorCode *ierr)
278: {
279: Mat *lsmat;
280: PetscInt i;
282: if (*n == 0) return;
283: *ierr = PetscMalloc1(*n + 1, &lsmat);
284: if (*ierr) return;
285: for (i = 0; i <= *n; i++) { lsmat[i] = smat[i]; }
286: *ierr = MatDestroySubMatrices(*n, &lsmat);
287: if (*ierr) return;
288: for (i = 0; i <= *n; i++) { PETSC_FORTRAN_OBJECT_C_NULL_TO_F_DESTROYED(&smat[i]); }
289: }