Actual source code: zshellf.c
1: #include <petsc/private/ftnimpl.h>
2: #include <petscmat.h>
4: #if defined(PETSC_HAVE_FORTRAN_CAPS)
5: #define matshellsetoperation_ MATSHELLSETOPERATION
6: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
7: #define matshellsetoperation_ matshellsetoperation
8: #endif
10: /**
11: * Subset of MatOperation that is supported by the Fortran wrappers.
12: */
13: enum FortranMatOperation {
14: FORTRAN_MATOP_MULT = 0,
15: FORTRAN_MATOP_MULT_ADD = 1,
16: FORTRAN_MATOP_MULT_TRANSPOSE = 2,
17: FORTRAN_MATOP_MULT_TRANSPOSE_ADD = 3,
18: FORTRAN_MATOP_SOR = 4,
19: FORTRAN_MATOP_TRANSPOSE = 5,
20: FORTRAN_MATOP_GET_DIAGONAL = 6,
21: FORTRAN_MATOP_DIAGONAL_SCALE = 7,
22: FORTRAN_MATOP_ZERO_ENTRIES = 8,
23: FORTRAN_MATOP_AXPY = 9,
24: FORTRAN_MATOP_SHIFT = 10,
25: FORTRAN_MATOP_DIAGONAL_SET = 11,
26: FORTRAN_MATOP_DESTROY = 12,
27: FORTRAN_MATOP_VIEW = 13,
28: FORTRAN_MATOP_CREATE_VECS = 14,
29: FORTRAN_MATOP_GET_DIAGONAL_BLOCK = 15,
30: FORTRAN_MATOP_COPY = 16,
31: FORTRAN_MATOP_SCALE = 17,
32: FORTRAN_MATOP_SET_RANDOM = 18,
33: FORTRAN_MATOP_ASSEMBLY_BEGIN = 19,
34: FORTRAN_MATOP_ASSEMBLY_END = 20,
35: FORTRAN_MATOP_DUPLICATE = 21,
36: FORTRAN_MATOP_MULT_HT = 22,
37: FORTRAN_MATOP_MULT_HT_ADD = 23,
38: FORTRAN_MATOP_SIZE = 24
39: };
41: /*
42: The MatShell Matrix Vector product requires a C routine.
43: This C routine then calls the corresponding Fortran routine that was
44: set by the user.
45: */
46: static PetscErrorCode ourmult(Mat mat, Vec x, Vec y)
47: {
48: PetscCallFortranVoidFunction((*(void (*)(Mat *, Vec *, Vec *, PetscErrorCode *))(((PetscObject)mat)->fortran_func_pointers[FORTRAN_MATOP_MULT]))(&mat, &x, &y, &ierr));
49: return PETSC_SUCCESS;
50: }
52: static PetscErrorCode ourmultadd(Mat mat, Vec x, Vec y, Vec z)
53: {
54: PetscCallFortranVoidFunction((*(void (*)(Mat *, Vec *, Vec *, Vec *, PetscErrorCode *))(((PetscObject)mat)->fortran_func_pointers[FORTRAN_MATOP_MULT_ADD]))(&mat, &x, &y, &z, &ierr));
55: return PETSC_SUCCESS;
56: }
58: static PetscErrorCode ourmulttranspose(Mat mat, Vec x, Vec y)
59: {
60: PetscCallFortranVoidFunction((*(void (*)(Mat *, Vec *, Vec *, PetscErrorCode *))(((PetscObject)mat)->fortran_func_pointers[FORTRAN_MATOP_MULT_TRANSPOSE]))(&mat, &x, &y, &ierr));
61: return PETSC_SUCCESS;
62: }
64: static PetscErrorCode ourmulthermitiantranspose(Mat mat, Vec x, Vec y)
65: {
66: PetscCallFortranVoidFunction((*(void (*)(Mat *, Vec *, Vec *, PetscErrorCode *))(((PetscObject)mat)->fortran_func_pointers[FORTRAN_MATOP_MULT_HT]))(&mat, &x, &y, &ierr));
67: return PETSC_SUCCESS;
68: }
70: static PetscErrorCode ourmulttransposeadd(Mat mat, Vec x, Vec y, Vec z)
71: {
72: PetscCallFortranVoidFunction((*(void (*)(Mat *, Vec *, Vec *, Vec *, PetscErrorCode *))(((PetscObject)mat)->fortran_func_pointers[FORTRAN_MATOP_MULT_TRANSPOSE_ADD]))(&mat, &x, &y, &z, &ierr));
73: return PETSC_SUCCESS;
74: }
76: static PetscErrorCode ourmulthermitiantransposeadd(Mat mat, Vec x, Vec y, Vec z)
77: {
78: PetscCallFortranVoidFunction((*(void (*)(Mat *, Vec *, Vec *, Vec *, PetscErrorCode *))(((PetscObject)mat)->fortran_func_pointers[FORTRAN_MATOP_MULT_HT_ADD]))(&mat, &x, &y, &z, &ierr));
79: return PETSC_SUCCESS;
80: }
82: static PetscErrorCode oursor(Mat mat, Vec b, PetscReal omega, MatSORType flg, PetscReal shift, PetscInt its, PetscInt lits, Vec x)
83: {
84: PetscErrorCode ierr = PETSC_SUCCESS;
86: (*(void (*)(Mat *, Vec *, PetscReal *, MatSORType *, PetscReal *, PetscInt *, PetscInt *, Vec *, PetscErrorCode *))(((PetscObject)mat)->fortran_func_pointers[FORTRAN_MATOP_SOR]))(&mat, &b, &omega, &flg, &shift, &its, &lits, &x, &ierr);
87: return ierr;
88: }
90: static PetscErrorCode ourtranspose(Mat mat, MatReuse reuse, Mat *B)
91: {
92: Mat bb = (Mat)-1;
93: Mat *b = (!B ? &bb : B);
95: PetscCallFortranVoidFunction((*(void (*)(Mat *, MatReuse *, Mat *, PetscErrorCode *))(((PetscObject)mat)->fortran_func_pointers[FORTRAN_MATOP_TRANSPOSE]))(&mat, &reuse, b, &ierr));
96: return PETSC_SUCCESS;
97: }
99: static PetscErrorCode ourgetdiagonal(Mat mat, Vec x)
100: {
101: PetscCallFortranVoidFunction((*(void (*)(Mat *, Vec *, PetscErrorCode *))(((PetscObject)mat)->fortran_func_pointers[FORTRAN_MATOP_GET_DIAGONAL]))(&mat, &x, &ierr));
102: return PETSC_SUCCESS;
103: }
105: static PetscErrorCode ourdiagonalscale(Mat mat, Vec l, Vec r)
106: {
107: Vec aa = (Vec)-1;
108: Vec *a = (!l ? &aa : &l);
109: Vec *b = (!r ? &aa : &r);
111: PetscCallFortranVoidFunction((*(void (*)(Mat *, Vec *, Vec *, PetscErrorCode *))(((PetscObject)mat)->fortran_func_pointers[FORTRAN_MATOP_DIAGONAL_SCALE]))(&mat, a, b, &ierr));
112: return PETSC_SUCCESS;
113: }
115: static PetscErrorCode ourzeroentries(Mat mat)
116: {
117: PetscCallFortranVoidFunction((*(void (*)(Mat *, PetscErrorCode *))(((PetscObject)mat)->fortran_func_pointers[FORTRAN_MATOP_ZERO_ENTRIES]))(&mat, &ierr));
118: return PETSC_SUCCESS;
119: }
121: static PetscErrorCode ouraxpy(Mat mat, PetscScalar a, Mat X, MatStructure str)
122: {
123: PetscCallFortranVoidFunction((*(void (*)(Mat *, PetscScalar *, Mat *, MatStructure *, PetscErrorCode *))(((PetscObject)mat)->fortran_func_pointers[FORTRAN_MATOP_AXPY]))(&mat, &a, &X, &str, &ierr));
124: return PETSC_SUCCESS;
125: }
127: static PetscErrorCode ourshift(Mat mat, PetscScalar a)
128: {
129: PetscCallFortranVoidFunction((*(void (*)(Mat *, PetscScalar *, PetscErrorCode *))(((PetscObject)mat)->fortran_func_pointers[FORTRAN_MATOP_SHIFT]))(&mat, &a, &ierr));
130: return PETSC_SUCCESS;
131: }
133: static PetscErrorCode ourdiagonalset(Mat mat, Vec x, InsertMode ins)
134: {
135: PetscCallFortranVoidFunction((*(void (*)(Mat *, Vec *, InsertMode *, PetscErrorCode *))(((PetscObject)mat)->fortran_func_pointers[FORTRAN_MATOP_DIAGONAL_SET]))(&mat, &x, &ins, &ierr));
136: return PETSC_SUCCESS;
137: }
139: static PetscErrorCode ourdestroy(Mat mat)
140: {
141: PetscCallFortranVoidFunction((*(void (*)(Mat *, PetscErrorCode *))(((PetscObject)mat)->fortran_func_pointers[FORTRAN_MATOP_DESTROY]))(&mat, &ierr));
142: return PETSC_SUCCESS;
143: }
145: static PetscErrorCode ourview(Mat mat, PetscViewer v)
146: {
147: PetscCallFortranVoidFunction((*(void (*)(Mat *, PetscViewer *, PetscErrorCode *))(((PetscObject)mat)->fortran_func_pointers[FORTRAN_MATOP_VIEW]))(&mat, &v, &ierr));
148: return PETSC_SUCCESS;
149: }
151: static PetscErrorCode ourgetvecs(Mat mat, Vec *l, Vec *r)
152: {
153: Vec aa = (Vec)-1;
154: Vec *a = (!l ? &aa : l);
155: Vec *b = (!r ? &aa : r);
157: PetscCallFortranVoidFunction((*(void (*)(Mat *, Vec *, Vec *, PetscErrorCode *))(((PetscObject)mat)->fortran_func_pointers[FORTRAN_MATOP_CREATE_VECS]))(&mat, a, b, &ierr));
158: return PETSC_SUCCESS;
159: }
161: static PetscErrorCode ourgetdiagonalblock(Mat mat, Mat *l)
162: {
163: PetscCallFortranVoidFunction((*(void (*)(Mat *, Mat *, PetscErrorCode *))(((PetscObject)mat)->fortran_func_pointers[FORTRAN_MATOP_GET_DIAGONAL_BLOCK]))(&mat, l, &ierr));
164: return PETSC_SUCCESS;
165: }
167: static PetscErrorCode ourcopy(Mat mat, Mat B, MatStructure str)
168: {
169: PetscCallFortranVoidFunction((*(void (*)(Mat *, Mat *, MatStructure *, PetscErrorCode *))(((PetscObject)mat)->fortran_func_pointers[FORTRAN_MATOP_COPY]))(&mat, &B, &str, &ierr));
170: return PETSC_SUCCESS;
171: }
173: static PetscErrorCode ourscale(Mat mat, PetscScalar a)
174: {
175: PetscCallFortranVoidFunction((*(void (*)(Mat *, PetscScalar *, PetscErrorCode *))(((PetscObject)mat)->fortran_func_pointers[FORTRAN_MATOP_SCALE]))(&mat, &a, &ierr));
176: return PETSC_SUCCESS;
177: }
179: static PetscErrorCode oursetrandom(Mat mat, PetscRandom ctx)
180: {
181: PetscCallFortranVoidFunction((*(void (*)(Mat *, PetscRandom *, PetscErrorCode *))(((PetscObject)mat)->fortran_func_pointers[FORTRAN_MATOP_SET_RANDOM]))(&mat, &ctx, &ierr));
182: return PETSC_SUCCESS;
183: }
185: static PetscErrorCode ourassemblybegin(Mat mat, MatAssemblyType type)
186: {
187: PetscCallFortranVoidFunction((*(void (*)(Mat *, MatAssemblyType *, PetscErrorCode *))(((PetscObject)mat)->fortran_func_pointers[FORTRAN_MATOP_ASSEMBLY_BEGIN]))(&mat, &type, &ierr));
188: return PETSC_SUCCESS;
189: }
191: static PetscErrorCode ourassemblyend(Mat mat, MatAssemblyType type)
192: {
193: PetscCallFortranVoidFunction((*(void (*)(Mat *, MatAssemblyType *, PetscErrorCode *))(((PetscObject)mat)->fortran_func_pointers[FORTRAN_MATOP_ASSEMBLY_END]))(&mat, &type, &ierr));
194: return PETSC_SUCCESS;
195: }
197: static PetscErrorCode ourduplicate(Mat mat, MatDuplicateOption op, Mat *M)
198: {
199: *((void **)(M)) = (void *)-2; // Initialize matrix since it will be passed to Fortran
200: PetscCallFortranVoidFunction((*(void (*)(Mat *, MatDuplicateOption *, Mat *, PetscErrorCode *))(((PetscObject)mat)->fortran_func_pointers[FORTRAN_MATOP_DUPLICATE]))(&mat, &op, M, &ierr));
201: return PETSC_SUCCESS;
202: }
204: PETSC_EXTERN void matshellsetoperation_(Mat *mat, MatOperation *op, PetscErrorCode (*f)(Mat *, Vec *, Vec *, PetscErrorCode *), PetscErrorCode *ierr)
205: {
206: MPI_Comm comm;
208: *ierr = PetscObjectGetComm((PetscObject)*mat, &comm);
209: if (*ierr) return;
210: PetscObjectAllocateFortranPointers(*mat, FORTRAN_MATOP_SIZE);
212: switch (*op) {
213: case MATOP_MULT:
214: *ierr = MatShellSetOperation(*mat, *op, (PetscVoidFn *)ourmult);
215: ((PetscObject)*mat)->fortran_func_pointers[FORTRAN_MATOP_MULT] = (PetscVoidFn *)f;
216: break;
217: case MATOP_MULT_ADD:
218: *ierr = MatShellSetOperation(*mat, *op, (PetscVoidFn *)ourmultadd);
219: ((PetscObject)*mat)->fortran_func_pointers[FORTRAN_MATOP_MULT_ADD] = (PetscVoidFn *)f;
220: break;
221: case MATOP_MULT_TRANSPOSE:
222: *ierr = MatShellSetOperation(*mat, *op, (PetscVoidFn *)ourmulttranspose);
223: ((PetscObject)*mat)->fortran_func_pointers[FORTRAN_MATOP_MULT_TRANSPOSE] = (PetscVoidFn *)f;
224: break;
225: case MATOP_MULT_HERMITIAN_TRANSPOSE:
226: *ierr = MatShellSetOperation(*mat, *op, (PetscVoidFn *)ourmulthermitiantranspose);
227: ((PetscObject)*mat)->fortran_func_pointers[FORTRAN_MATOP_MULT_HT] = (PetscVoidFn *)f;
228: break;
229: case MATOP_MULT_TRANSPOSE_ADD:
230: *ierr = MatShellSetOperation(*mat, *op, (PetscVoidFn *)ourmulttransposeadd);
231: ((PetscObject)*mat)->fortran_func_pointers[FORTRAN_MATOP_MULT_TRANSPOSE_ADD] = (PetscVoidFn *)f;
232: break;
233: case MATOP_MULT_HERMITIAN_TRANS_ADD:
234: *ierr = MatShellSetOperation(*mat, *op, (PetscVoidFn *)ourmulthermitiantransposeadd);
235: ((PetscObject)*mat)->fortran_func_pointers[FORTRAN_MATOP_MULT_HT_ADD] = (PetscVoidFn *)f;
236: break;
237: case MATOP_SOR:
238: *ierr = MatShellSetOperation(*mat, *op, (PetscVoidFn *)oursor);
239: ((PetscObject)*mat)->fortran_func_pointers[FORTRAN_MATOP_SOR] = (PetscVoidFn *)f;
240: break;
241: case MATOP_TRANSPOSE:
242: *ierr = MatShellSetOperation(*mat, *op, (PetscVoidFn *)ourtranspose);
243: ((PetscObject)*mat)->fortran_func_pointers[FORTRAN_MATOP_TRANSPOSE] = (PetscVoidFn *)f;
244: break;
245: case MATOP_GET_DIAGONAL:
246: *ierr = MatShellSetOperation(*mat, *op, (PetscVoidFn *)ourgetdiagonal);
247: ((PetscObject)*mat)->fortran_func_pointers[FORTRAN_MATOP_GET_DIAGONAL] = (PetscVoidFn *)f;
248: break;
249: case MATOP_DIAGONAL_SCALE:
250: *ierr = MatShellSetOperation(*mat, *op, (PetscVoidFn *)ourdiagonalscale);
251: ((PetscObject)*mat)->fortran_func_pointers[FORTRAN_MATOP_DIAGONAL_SCALE] = (PetscVoidFn *)f;
252: break;
253: case MATOP_ZERO_ENTRIES:
254: *ierr = MatShellSetOperation(*mat, *op, (PetscVoidFn *)ourzeroentries);
255: ((PetscObject)*mat)->fortran_func_pointers[FORTRAN_MATOP_ZERO_ENTRIES] = (PetscVoidFn *)f;
256: break;
257: case MATOP_AXPY:
258: *ierr = MatShellSetOperation(*mat, *op, (PetscVoidFn *)ouraxpy);
259: ((PetscObject)*mat)->fortran_func_pointers[FORTRAN_MATOP_AXPY] = (PetscVoidFn *)f;
260: break;
261: case MATOP_SHIFT:
262: *ierr = MatShellSetOperation(*mat, *op, (PetscVoidFn *)ourshift);
263: ((PetscObject)*mat)->fortran_func_pointers[FORTRAN_MATOP_SHIFT] = (PetscVoidFn *)f;
264: break;
265: case MATOP_DIAGONAL_SET:
266: *ierr = MatShellSetOperation(*mat, *op, (PetscVoidFn *)ourdiagonalset);
267: ((PetscObject)*mat)->fortran_func_pointers[FORTRAN_MATOP_DIAGONAL_SET] = (PetscVoidFn *)f;
268: break;
269: case MATOP_DESTROY:
270: *ierr = MatShellSetOperation(*mat, *op, (PetscVoidFn *)ourdestroy);
271: ((PetscObject)*mat)->fortran_func_pointers[FORTRAN_MATOP_DESTROY] = (PetscVoidFn *)f;
272: break;
273: case MATOP_VIEW:
274: *ierr = MatShellSetOperation(*mat, *op, (PetscVoidFn *)ourview);
275: ((PetscObject)*mat)->fortran_func_pointers[FORTRAN_MATOP_VIEW] = (PetscVoidFn *)f;
276: break;
277: case MATOP_CREATE_VECS:
278: *ierr = MatShellSetOperation(*mat, *op, (PetscVoidFn *)ourgetvecs);
279: ((PetscObject)*mat)->fortran_func_pointers[FORTRAN_MATOP_CREATE_VECS] = (PetscVoidFn *)f;
280: break;
281: case MATOP_GET_DIAGONAL_BLOCK:
282: *ierr = MatShellSetOperation(*mat, *op, (PetscVoidFn *)ourgetdiagonalblock);
283: ((PetscObject)*mat)->fortran_func_pointers[FORTRAN_MATOP_GET_DIAGONAL_BLOCK] = (PetscVoidFn *)f;
284: break;
285: case MATOP_COPY:
286: *ierr = MatShellSetOperation(*mat, *op, (PetscVoidFn *)ourcopy);
287: ((PetscObject)*mat)->fortran_func_pointers[FORTRAN_MATOP_COPY] = (PetscVoidFn *)f;
288: break;
289: case MATOP_SCALE:
290: *ierr = MatShellSetOperation(*mat, *op, (PetscVoidFn *)ourscale);
291: ((PetscObject)*mat)->fortran_func_pointers[FORTRAN_MATOP_SCALE] = (PetscVoidFn *)f;
292: break;
293: case MATOP_SET_RANDOM:
294: *ierr = MatShellSetOperation(*mat, *op, (PetscVoidFn *)oursetrandom);
295: ((PetscObject)*mat)->fortran_func_pointers[FORTRAN_MATOP_SET_RANDOM] = (PetscVoidFn *)f;
296: break;
297: case MATOP_ASSEMBLY_BEGIN:
298: *ierr = MatShellSetOperation(*mat, *op, (PetscVoidFn *)ourassemblybegin);
299: ((PetscObject)*mat)->fortran_func_pointers[FORTRAN_MATOP_ASSEMBLY_BEGIN] = (PetscVoidFn *)f;
300: break;
301: case MATOP_ASSEMBLY_END:
302: *ierr = MatShellSetOperation(*mat, *op, (PetscVoidFn *)ourassemblyend);
303: ((PetscObject)*mat)->fortran_func_pointers[FORTRAN_MATOP_ASSEMBLY_END] = (PetscVoidFn *)f;
304: break;
305: case MATOP_DUPLICATE:
306: *ierr = MatShellSetOperation(*mat, *op, (PetscVoidFn *)ourduplicate);
307: ((PetscObject)*mat)->fortran_func_pointers[FORTRAN_MATOP_DUPLICATE] = (PetscVoidFn *)f;
308: break;
309: default:
310: *ierr = PetscError(comm, __LINE__, "MatShellSetOperation_Fortran", __FILE__, PETSC_ERR_ARG_WRONG, PETSC_ERROR_INITIAL, "Cannot set that matrix operation");
311: *ierr = PETSC_ERR_ARG_WRONG;
312: }
313: }