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: }