1# --------------------------------------------------------------------
2
3cdef extern from * nogil:
4
5    ctypedef const char* PetscDMType "DMType"
6    PetscDMType DMDA_type "DMDA"
7    PetscDMType DMCOMPOSITE
8    PetscDMType DMSLICED
9    PetscDMType DMSHELL
10    PetscDMType DMPLEX
11    PetscDMType DMREDUNDANT
12    PetscDMType DMPATCH
13    PetscDMType DMMOAB
14    PetscDMType DMNETWORK
15    PetscDMType DMFOREST
16    PetscDMType DMP4EST
17    PetscDMType DMP8EST
18    PetscDMType DMSWARM
19    PetscDMType DMPRODUCT
20    PetscDMType DMSTAG
21
22    ctypedef enum PetscDMBoundaryType"DMBoundaryType":
23        DM_BOUNDARY_NONE
24        DM_BOUNDARY_GHOSTED
25        DM_BOUNDARY_MIRROR
26        DM_BOUNDARY_PERIODIC
27        DM_BOUNDARY_TWIST
28
29    ctypedef int (*PetscDMCoarsenHook)(PetscDM,
30                                       PetscDM,
31                                       void*) except PETSC_ERR_PYTHON
32    ctypedef int (*PetscDMRestrictHook)(PetscDM,
33                                        PetscMat,
34                                        PetscVec,
35                                        PetscMat,
36                                        PetscDM,
37                                        void*) except PETSC_ERR_PYTHON
38
39    int DMCreate(MPI_Comm,PetscDM*)
40    int DMClone(PetscDM,PetscDM*)
41    int DMDestroy(PetscDM*)
42    int DMView(PetscDM,PetscViewer)
43    int DMSetType(PetscDM,PetscDMType)
44    int DMGetType(PetscDM,PetscDMType*)
45    int DMGetDimension(PetscDM,PetscInt*)
46    int DMSetDimension(PetscDM,PetscInt)
47    int DMSetOptionsPrefix(PetscDM,char[])
48    int DMSetFromOptions(PetscDM)
49    int DMSetUp(PetscDM)
50
51    int DMGetAdjacency(PetscDM,PetscInt,PetscBool*,PetscBool*)
52    int DMSetAdjacency(PetscDM,PetscInt,PetscBool,PetscBool)
53    int DMGetBasicAdjacency(PetscDM,PetscBool*,PetscBool*)
54    int DMSetBasicAdjacency(PetscDM,PetscBool,PetscBool)
55
56    int DMSetNumFields(PetscDM,PetscInt)
57    int DMGetNumFields(PetscDM,PetscInt*)
58    int DMSetField(PetscDM,PetscInt,PetscDMLabel,PetscObject)
59    int DMAddField(PetscDM,PetscDMLabel,PetscObject)
60    int DMGetField(PetscDM,PetscInt,PetscDMLabel*,PetscObject*)
61    int DMCopyFields(PetscDM,PetscDM)
62    int DMCreateDS(PetscDM)
63    int DMClearDS(PetscDM)
64    int DMGetDS(PetscDM,PetscDS*)
65    int DMCopyDS(PetscDM,PetscDM)
66    int DMCopyDisc(PetscDM,PetscDM)
67
68    int DMGetBlockSize(PetscDM,PetscInt*)
69    int DMSetVecType(PetscDM,PetscVecType)
70    int DMCreateLocalVector(PetscDM,PetscVec*)
71    int DMCreateGlobalVector(PetscDM,PetscVec*)
72    int DMGetLocalVector(PetscDM,PetscVec*)
73    int DMRestoreLocalVector(PetscDM,PetscVec*)
74    int DMGetGlobalVector(PetscDM,PetscVec*)
75    int DMRestoreGlobalVector(PetscDM,PetscVec*)
76    int DMSetMatType(PetscDM,PetscMatType)
77    int DMCreateMatrix(PetscDM,PetscMat*)
78
79    int DMGetCoordinateDM(PetscDM,PetscDM*)
80    int DMGetCoordinateSection(PetscDM,PetscSection*)
81    int DMSetCoordinates(PetscDM,PetscVec)
82    int DMGetCoordinates(PetscDM,PetscVec*)
83    int DMSetCoordinatesLocal(PetscDM,PetscVec)
84    int DMGetCoordinatesLocal(PetscDM,PetscVec*)
85    int DMGetCoordinateDim(PetscDM,PetscInt*)
86    int DMSetCoordinateDim(PetscDM,PetscInt)
87    int DMLocalizeCoordinates(PetscDM)
88
89    int DMCreateInterpolation(PetscDM,PetscDM,PetscMat*,PetscVec*)
90    int DMCreateInjection(PetscDM,PetscDM,PetscMat*)
91    int DMCreateRestriction(PetscDM,PetscDM,PetscMat*)
92
93    int DMConvert(PetscDM,PetscDMType,PetscDM*)
94    int DMRefine(PetscDM,MPI_Comm,PetscDM*)
95    int DMCoarsen(PetscDM,MPI_Comm,PetscDM*)
96    int DMRefineHierarchy(PetscDM,PetscInt,PetscDM[])
97    int DMCoarsenHierarchy(PetscDM,PetscInt,PetscDM[])
98    int DMGetRefineLevel(PetscDM,PetscInt*)
99    int DMSetRefineLevel(PetscDM,PetscInt)
100    int DMGetCoarsenLevel(PetscDM,PetscInt*)
101
102    int DMAdaptLabel(PetscDM,PetscDMLabel,PetscDM*)
103    int DMAdaptMetric(PetscDM,PetscVec,PetscDMLabel,PetscDM*)
104
105    int DMGlobalToLocalBegin(PetscDM,PetscVec,PetscInsertMode,PetscVec)
106    int DMGlobalToLocalEnd(PetscDM,PetscVec,PetscInsertMode,PetscVec)
107    int DMLocalToGlobalBegin(PetscDM,PetscVec,PetscInsertMode,PetscVec)
108    int DMLocalToGlobalEnd(PetscDM,PetscVec,PetscInsertMode,PetscVec)
109    int DMLocalToLocalBegin(PetscDM,PetscVec,PetscInsertMode,PetscVec)
110    int DMLocalToLocalEnd(PetscDM,PetscVec,PetscInsertMode,PetscVec)
111
112    int DMGetLocalToGlobalMapping(PetscDM,PetscLGMap*)
113
114    int DMSetSection(PetscDM,PetscSection)
115    int DMGetSection(PetscDM,PetscSection*)
116    int DMSetGlobalSection(PetscDM,PetscSection)
117    int DMGetGlobalSection(PetscDM,PetscSection*)
118    int DMCreateSectionSF(PetscDM,PetscSection,PetscSection)
119    int DMGetSectionSF(PetscDM,PetscSF*)
120    int DMSetSectionSF(PetscDM,PetscSF)
121    int DMGetPointSF(PetscDM,PetscSF*)
122    int DMSetPointSF(PetscDM,PetscSF)
123
124    int DMCreateLabel(PetscDM,const char[])
125    int DMGetLabelValue(PetscDM,const char[],PetscInt,PetscInt*)
126    int DMSetLabelValue(PetscDM,const char[],PetscInt,PetscInt)
127    int DMHasLabel(PetscDM,const char[],PetscBool*)
128    int DMClearLabelValue(PetscDM,const char[],PetscInt,PetscInt)
129    int DMGetLabelSize(PetscDM,const char[],PetscInt*)
130    int DMGetLabelIdIS(PetscDM,const char[],PetscIS*)
131    int DMGetStratumSize(PetscDM,const char[],PetscInt,PetscInt*)
132    int DMGetStratumIS(PetscDM,const char[],PetscInt,PetscIS*)
133    int DMClearLabelStratum(PetscDM,const char[],PetscInt)
134    int DMSetLabelOutput(PetscDM,const char[],PetscBool)
135    int DMGetLabelOutput(PetscDM,const char[],PetscBool*)
136    int DMGetNumLabels(PetscDM,PetscInt*)
137    int DMGetLabelName(PetscDM,PetscInt,const char**)
138    int DMHasLabel(PetscDM,const char[],PetscBool*)
139    int DMGetLabel(PetscDM,const char*,PetscDMLabel*)
140    int DMAddLabel(PetscDM,PetscDMLabel)
141    int DMRemoveLabel(PetscDM,const char[],PetscDMLabel*)
142    int DMLabelDestroy(PetscDMLabel *)
143    #int DMCopyLabels(PetscDM,PetscDM)
144
145    int DMShellSetGlobalVector(PetscDM,PetscVec)
146    int DMShellSetLocalVector(PetscDM,PetscVec)
147
148    int DMKSPSetComputeOperators(PetscDM,PetscKSPComputeOpsFunction,void*)
149
150    int DMCreateFieldDecomposition(PetscDM,PetscInt*,char***,PetscIS**,PetscDM**)
151
152    int DMSNESSetFunction(PetscDM,PetscSNESFunctionFunction,void*)
153    int DMSNESSetJacobian(PetscDM,PetscSNESJacobianFunction,void*)
154
155    int DMCoarsenHookAdd(PetscDM,PetscDMCoarsenHook,PetscDMRestrictHook,void*)
156
157# --------------------------------------------------------------------
158
159cdef inline PetscDMBoundaryType asBoundaryType(object boundary) \
160    except <PetscDMBoundaryType>(-1):
161    if boundary is None:
162        return DM_BOUNDARY_NONE
163    if boundary is False:
164        return DM_BOUNDARY_NONE
165    if boundary is True:
166        return DM_BOUNDARY_PERIODIC
167    if isinstance(boundary, str):
168        if boundary == 'none':
169            return DM_BOUNDARY_NONE
170        elif boundary == 'ghosted':
171            return DM_BOUNDARY_GHOSTED
172        elif boundary == 'mirror':
173            return DM_BOUNDARY_MIRROR
174        elif boundary == 'periodic':
175            return DM_BOUNDARY_PERIODIC
176        elif boundary == 'twist':
177            return DM_BOUNDARY_TWIST
178        else:
179            raise ValueError("unknown boundary type: %s" % boundary)
180    return boundary
181
182cdef inline PetscInt asBoundary(object boundary,
183                                PetscDMBoundaryType *_x,
184                                PetscDMBoundaryType *_y,
185                                PetscDMBoundaryType *_z) except -1:
186    cdef PetscInt dim = 0
187    cdef object x=None, y=None, z=None
188    if (boundary is None or
189        isinstance(boundary, str) or
190        isinstance(boundary, int)):
191        _x[0] = _y[0] = _z[0] = asBoundaryType(boundary)
192    else:
193        _x[0] = _y[0] = _z[0] = DM_BOUNDARY_NONE
194        boundary = tuple(boundary)
195        dim = <PetscInt>len(boundary)
196        if   dim == 0: pass
197        elif dim == 1: (x,) = boundary
198        elif dim == 2: (x, y) = boundary
199        elif dim == 3: (x, y, z) = boundary
200        if dim >= 1: _x[0] = asBoundaryType(x)
201        if dim >= 2: _y[0] = asBoundaryType(y)
202        if dim >= 3: _z[0] = asBoundaryType(z)
203    return dim
204
205cdef inline object toBoundary(PetscInt dim,
206                              PetscDMBoundaryType x,
207                              PetscDMBoundaryType y,
208                              PetscDMBoundaryType z):
209    if   dim == 0: return ()
210    elif dim == 1: return (x,)
211    elif dim == 2: return (x, y)
212    elif dim == 3: return (x, y, z)
213
214# -----------------------------------------------------------------------------
215
216cdef inline DM ref_DM(PetscDM dm):
217    cdef DM ob = <DM> DM()
218    ob.dm = dm
219    PetscINCREF(ob.obj)
220    return ob
221
222# --------------------------------------------------------------------
223
224cdef int DM_PyCoarsenHook(
225    PetscDM fine,
226    PetscDM coarse,
227    void*   ctx,
228    ) except PETSC_ERR_PYTHON with gil:
229
230    cdef DM Fine = ref_DM(fine)
231    cdef DM Coarse = ref_DM(coarse)
232    cdef object hooks = Fine.get_attr('__coarsenhooks__')
233    assert hooks is not None and type(hooks) is list
234    for hook in hooks:
235        (hookop, args, kargs) = hook
236        hookop(Fine, Coarse, *args, **kargs)
237    return 0
238
239cdef int DM_PyRestrictHook(
240    PetscDM fine,
241    PetscMat mrestrict,
242    PetscVec rscale,
243    PetscMat inject,
244    PetscDM coarse,
245    void*   ctx,
246    ) except PETSC_ERR_PYTHON with gil:
247
248    cdef DM  Fine = ref_DM(fine)
249    cdef Mat Mrestrict = ref_Mat(mrestrict)
250    cdef Vec Rscale = ref_Vec(rscale)
251    cdef Mat Inject = ref_Mat(inject)
252    cdef DM  Coarse = ref_DM(coarse)
253    cdef object hooks = Fine.get_attr('__restricthooks__')
254    assert hooks is not None and type(hooks) is list
255    for hook in hooks:
256        (hookop, args, kargs) = hook
257        hookop(Fine, Mrestrict, Rscale, Inject, Coarse, *args, **kargs)
258    return 0
259# -----------------------------------------------------------------------------
260