1 #include <petsc/private/ftnimpl.h>
2 #include <petscdmshell.h> /*I "petscdmshell.h" I*/
3
4 #if defined(PETSC_HAVE_FORTRAN_CAPS)
5 #define dmshellsetcreatematrix_ DMSHELLSETCREATEMATRIX
6 #define dmshellsetcreateglobalvector_ DMSHELLSETCREATEGLOBALVECTOR
7 #define dmshellsetcreatelocalvector_ DMSHELLSETCREATELOCALVECTOR
8 #define dmshellsetglobaltolocal_ DMSHELLSETGLOBALTOLOCAL
9 #define dmshellsetlocaltoglobal_ DMSHELLSETLOCALTOGLOBAL
10 #define dmshellsetlocaltolocal_ DMSHELLSETLOCALTOLOCAL
11 #define dmshellsetcreatefielddecomposition_ DMSHELLSETCREATEFIELDDECOMPOSITION
12 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
13 #define dmshellsetcreatematrix_ dmshellsetcreatematrix
14 #define dmshellsetcreateglobalvector_ dmshellsetcreateglobalvector
15 #define dmshellsetcreatelocalvector_ dmshellsetcreatelocalvector
16 #define dmshellsetglobaltolocal_ dmshellsetglobaltolocal
17 #define dmshellsetlocaltoglobal_ dmshellsetlocaltoglobal
18 #define dmshellsetlocaltolocal_ dmshellsetlocaltolocal
19 #define dmshellsetcreatefielddecomposition_ dmshellsetcreatefielddecomposition
20 #endif
21
22 /*
23 * C routines are required for matrix and global vector creation. We define C routines here that call the corresponding
24 * Fortran routine (indexed by _cb) that was set by the user.
25 */
26
27 static struct {
28 PetscFortranCallbackId creatematrix;
29 PetscFortranCallbackId createglobalvector;
30 PetscFortranCallbackId createlocalvector;
31 PetscFortranCallbackId globaltolocalbegin;
32 PetscFortranCallbackId globaltolocalend;
33 PetscFortranCallbackId localtoglobalbegin;
34 PetscFortranCallbackId localtoglobalend;
35 PetscFortranCallbackId localtolocalbegin;
36 PetscFortranCallbackId localtolocalend;
37 PetscFortranCallbackId createfielddecomposition;
38 } _cb;
39
ourcreatematrix(DM dm,Mat * A)40 static PetscErrorCode ourcreatematrix(DM dm, Mat *A)
41 {
42 PetscObjectUseFortranCallbackSubType(dm, _cb.creatematrix, (DM *, Mat *, PetscErrorCode *), (&dm, A, &ierr));
43 }
44
ourcreateglobalvector(DM dm,Vec * v)45 static PetscErrorCode ourcreateglobalvector(DM dm, Vec *v)
46 {
47 PetscObjectUseFortranCallbackSubType(dm, _cb.createglobalvector, (DM *, Vec *, PetscErrorCode *), (&dm, v, &ierr));
48 }
49
ourcreatelocalvector(DM dm,Vec * v)50 static PetscErrorCode ourcreatelocalvector(DM dm, Vec *v)
51 {
52 PetscObjectUseFortranCallbackSubType(dm, _cb.createlocalvector, (DM *, Vec *, PetscErrorCode *), (&dm, v, &ierr));
53 }
54
ourglobaltolocalbegin(DM dm,Vec g,InsertMode mode,Vec l)55 static PetscErrorCode ourglobaltolocalbegin(DM dm, Vec g, InsertMode mode, Vec l)
56 {
57 PetscObjectUseFortranCallbackSubType(dm, _cb.globaltolocalbegin, (DM *, Vec *, InsertMode *, Vec *, PetscErrorCode *), (&dm, &g, &mode, &l, &ierr));
58 }
59
ourglobaltolocalend(DM dm,Vec g,InsertMode mode,Vec l)60 static PetscErrorCode ourglobaltolocalend(DM dm, Vec g, InsertMode mode, Vec l)
61 {
62 PetscObjectUseFortranCallbackSubType(dm, _cb.globaltolocalend, (DM *, Vec *, InsertMode *, Vec *, PetscErrorCode *), (&dm, &g, &mode, &l, &ierr));
63 }
64
ourlocaltoglobalbegin(DM dm,Vec l,InsertMode mode,Vec g)65 static PetscErrorCode ourlocaltoglobalbegin(DM dm, Vec l, InsertMode mode, Vec g)
66 {
67 PetscObjectUseFortranCallbackSubType(dm, _cb.localtoglobalbegin, (DM *, Vec *, InsertMode *, Vec *, PetscErrorCode *), (&dm, &l, &mode, &g, &ierr));
68 }
69
ourlocaltoglobalend(DM dm,Vec l,InsertMode mode,Vec g)70 static PetscErrorCode ourlocaltoglobalend(DM dm, Vec l, InsertMode mode, Vec g)
71 {
72 PetscObjectUseFortranCallbackSubType(dm, _cb.localtoglobalend, (DM *, Vec *, InsertMode *, Vec *, PetscErrorCode *), (&dm, &l, &mode, &g, &ierr));
73 }
74
ourlocaltolocalbegin(DM dm,Vec g,InsertMode mode,Vec l)75 static PetscErrorCode ourlocaltolocalbegin(DM dm, Vec g, InsertMode mode, Vec l)
76 {
77 PetscObjectUseFortranCallbackSubType(dm, _cb.localtolocalbegin, (DM *, Vec *, InsertMode *, Vec *, PetscErrorCode *), (&dm, &g, &mode, &l, &ierr));
78 }
79
ourlocaltolocalend(DM dm,Vec g,InsertMode mode,Vec l)80 static PetscErrorCode ourlocaltolocalend(DM dm, Vec g, InsertMode mode, Vec l)
81 {
82 PetscObjectUseFortranCallbackSubType(dm, _cb.localtolocalend, (DM *, Vec *, InsertMode *, Vec *, PetscErrorCode *), (&dm, &g, &mode, &l, &ierr));
83 }
84
ourcreatefielddecomposition(DM dm,PetscInt * nfields,char *** names,IS ** is,DM ** subdms)85 static PetscErrorCode ourcreatefielddecomposition(DM dm, PetscInt *nfields, char ***names, IS **is, DM **subdms)
86 {
87 PetscObjectUseFortranCallbackSubType(dm, _cb.createfielddecomposition, (DM *, PetscInt *, char ***, IS **, DM **, PetscErrorCode *), (&dm, nfields, names, is, subdms, &ierr));
88 }
89
dmshellsetcreatematrix_(DM * dm,void (* func)(DM *,Mat *,PetscErrorCode *,PETSC_FORTRAN_CHARLEN_T len),PetscErrorCode * ierr)90 PETSC_EXTERN void dmshellsetcreatematrix_(DM *dm, void (*func)(DM *, Mat *, PetscErrorCode *, PETSC_FORTRAN_CHARLEN_T len), PetscErrorCode *ierr)
91 {
92 *ierr = PetscObjectSetFortranCallback((PetscObject)*dm, PETSC_FORTRAN_CALLBACK_SUBTYPE, &_cb.creatematrix, (PetscVoidFn *)func, NULL);
93 if (*ierr) return;
94 *ierr = DMShellSetCreateMatrix(*dm, ourcreatematrix);
95 }
96
dmshellsetcreateglobalvector_(DM * dm,void (* func)(DM *,Vec *,PetscErrorCode *),PetscErrorCode * ierr)97 PETSC_EXTERN void dmshellsetcreateglobalvector_(DM *dm, void (*func)(DM *, Vec *, PetscErrorCode *), PetscErrorCode *ierr)
98 {
99 *ierr = PetscObjectSetFortranCallback((PetscObject)*dm, PETSC_FORTRAN_CALLBACK_SUBTYPE, &_cb.createglobalvector, (PetscVoidFn *)func, NULL);
100 if (*ierr) return;
101 *ierr = DMShellSetCreateGlobalVector(*dm, ourcreateglobalvector);
102 }
103
dmshellsetcreatelocalvector_(DM * dm,void (* func)(DM *,Vec *,PetscErrorCode *),PetscErrorCode * ierr)104 PETSC_EXTERN void dmshellsetcreatelocalvector_(DM *dm, void (*func)(DM *, Vec *, PetscErrorCode *), PetscErrorCode *ierr)
105 {
106 *ierr = PetscObjectSetFortranCallback((PetscObject)*dm, PETSC_FORTRAN_CALLBACK_SUBTYPE, &_cb.createlocalvector, (PetscVoidFn *)func, NULL);
107 if (*ierr) return;
108 *ierr = DMShellSetCreateLocalVector(*dm, ourcreatelocalvector);
109 }
110
dmshellsetglobaltolocal_(DM * dm,void (* begin)(DM *,Vec *,InsertMode *,Vec *,PetscErrorCode *),void (* end)(DM *,Vec *,InsertMode *,Vec *,PetscErrorCode *),PetscErrorCode * ierr)111 PETSC_EXTERN void dmshellsetglobaltolocal_(DM *dm, void (*begin)(DM *, Vec *, InsertMode *, Vec *, PetscErrorCode *), void (*end)(DM *, Vec *, InsertMode *, Vec *, PetscErrorCode *), PetscErrorCode *ierr)
112 {
113 *ierr = PetscObjectSetFortranCallback((PetscObject)*dm, PETSC_FORTRAN_CALLBACK_SUBTYPE, &_cb.globaltolocalbegin, (PetscVoidFn *)begin, NULL);
114 if (*ierr) return;
115 *ierr = PetscObjectSetFortranCallback((PetscObject)*dm, PETSC_FORTRAN_CALLBACK_SUBTYPE, &_cb.globaltolocalend, (PetscVoidFn *)end, NULL);
116 if (*ierr) return;
117 *ierr = DMShellSetGlobalToLocal(*dm, ourglobaltolocalbegin, ourglobaltolocalend);
118 }
119
dmshellsetlocaltoglobal_(DM * dm,void (* begin)(DM *,Vec *,InsertMode *,Vec *,PetscErrorCode *),void (* end)(DM *,Vec *,InsertMode *,Vec *,PetscErrorCode *),PetscErrorCode * ierr)120 PETSC_EXTERN void dmshellsetlocaltoglobal_(DM *dm, void (*begin)(DM *, Vec *, InsertMode *, Vec *, PetscErrorCode *), void (*end)(DM *, Vec *, InsertMode *, Vec *, PetscErrorCode *), PetscErrorCode *ierr)
121 {
122 *ierr = PetscObjectSetFortranCallback((PetscObject)*dm, PETSC_FORTRAN_CALLBACK_SUBTYPE, &_cb.localtoglobalbegin, (PetscVoidFn *)begin, NULL);
123 if (*ierr) return;
124 *ierr = PetscObjectSetFortranCallback((PetscObject)*dm, PETSC_FORTRAN_CALLBACK_SUBTYPE, &_cb.localtoglobalend, (PetscVoidFn *)end, NULL);
125 if (*ierr) return;
126 *ierr = DMShellSetLocalToGlobal(*dm, ourlocaltoglobalbegin, ourlocaltoglobalend);
127 }
128
dmshellsetlocaltolocal_(DM * dm,void (* begin)(DM *,Vec *,InsertMode *,Vec *,PetscErrorCode *),void (* end)(DM *,Vec *,InsertMode *,Vec *,PetscErrorCode *),PetscErrorCode * ierr)129 PETSC_EXTERN void dmshellsetlocaltolocal_(DM *dm, void (*begin)(DM *, Vec *, InsertMode *, Vec *, PetscErrorCode *), void (*end)(DM *, Vec *, InsertMode *, Vec *, PetscErrorCode *), PetscErrorCode *ierr)
130 {
131 *ierr = PetscObjectSetFortranCallback((PetscObject)*dm, PETSC_FORTRAN_CALLBACK_SUBTYPE, &_cb.localtolocalbegin, (PetscVoidFn *)begin, NULL);
132 if (*ierr) return;
133 *ierr = PetscObjectSetFortranCallback((PetscObject)*dm, PETSC_FORTRAN_CALLBACK_SUBTYPE, &_cb.localtolocalend, (PetscVoidFn *)end, NULL);
134 if (*ierr) return;
135 *ierr = DMShellSetLocalToLocal(*dm, ourlocaltolocalbegin, ourlocaltolocalend);
136 }
137
dmshellsetcreatefielddecomposition_(DM * dm,void (* func)(DM *,PetscInt *,char ***,IS **,DM **,PetscErrorCode *),PetscErrorCode * ierr)138 PETSC_EXTERN void dmshellsetcreatefielddecomposition_(DM *dm, void (*func)(DM *, PetscInt *, char ***, IS **, DM **, PetscErrorCode *), PetscErrorCode *ierr)
139 {
140 *ierr = PetscObjectSetFortranCallback((PetscObject)*dm, PETSC_FORTRAN_CALLBACK_SUBTYPE, &_cb.createfielddecomposition, (PetscVoidFn *)func, NULL);
141 if (*ierr) return;
142 *ierr = DMShellSetCreateFieldDecomposition(*dm, ourcreatefielddecomposition);
143 }
144