1# --------------------------------------------------------------------
2
3class DMType(object):
4    DA        = S_(DMDA_type)
5    COMPOSITE = S_(DMCOMPOSITE)
6    SLICED    = S_(DMSLICED)
7    SHELL     = S_(DMSHELL)
8    PLEX      = S_(DMPLEX)
9    REDUNDANT = S_(DMREDUNDANT)
10    PATCH     = S_(DMPATCH)
11    MOAB      = S_(DMMOAB)
12    NETWORK   = S_(DMNETWORK)
13    FOREST    = S_(DMFOREST)
14    P4EST     = S_(DMP4EST)
15    P8EST     = S_(DMP8EST)
16    SWARM     = S_(DMSWARM)
17    PRODUCT   = S_(DMPRODUCT)
18    STAG      = S_(DMSTAG)
19
20class DMBoundaryType(object):
21    NONE     = DM_BOUNDARY_NONE
22    GHOSTED  = DM_BOUNDARY_GHOSTED
23    MIRROR   = DM_BOUNDARY_MIRROR
24    PERIODIC = DM_BOUNDARY_PERIODIC
25    TWIST    = DM_BOUNDARY_TWIST
26
27# --------------------------------------------------------------------
28
29cdef class DM(Object):
30
31    Type         = DMType
32    BoundaryType = DMBoundaryType
33
34    #
35
36    def __cinit__(self):
37        self.obj = <PetscObject*> &self.dm
38        self.dm  = NULL
39
40    def view(self, Viewer viewer=None):
41        cdef PetscViewer vwr = NULL
42        if viewer is not None: vwr = viewer.vwr
43        CHKERR( DMView(self.dm, vwr) )
44
45    def destroy(self):
46        CHKERR( DMDestroy(&self.dm) )
47        return self
48
49    def create(self, comm=None):
50        cdef MPI_Comm ccomm = def_Comm(comm, PETSC_COMM_DEFAULT)
51        cdef PetscDM newdm = NULL
52        CHKERR( DMCreate(ccomm, &newdm) )
53        PetscCLEAR(self.obj); self.dm = newdm
54        return self
55
56    def clone(self):
57        cdef DM dm = type(self)()
58        CHKERR( DMClone(self.dm, &dm.dm) )
59        return dm
60
61    def setType(self, dm_type):
62        cdef PetscDMType cval = NULL
63        dm_type = str2bytes(dm_type, &cval)
64        CHKERR( DMSetType(self.dm, cval) )
65
66    def getType(self):
67        cdef PetscDMType cval = NULL
68        CHKERR( DMGetType(self.dm, &cval) )
69        return bytes2str(cval)
70
71    def getDimension(self):
72        cdef PetscInt dim = 0
73        CHKERR( DMGetDimension(self.dm, &dim) )
74        return toInt(dim)
75
76    def setDimension(self, dim):
77        cdef PetscInt cdim = asInt(dim)
78        CHKERR( DMSetDimension(self.dm, cdim) )
79
80    def getCoordinateDim(self):
81        cdef PetscInt dim = 0
82        CHKERR( DMGetCoordinateDim(self.dm, &dim) )
83        return toInt(dim)
84
85    def setCoordinateDim(self, dim):
86        cdef PetscInt cdim = asInt(dim)
87        CHKERR( DMSetCoordinateDim(self.dm, cdim) )
88
89    def setOptionsPrefix(self, prefix):
90        cdef const char *cval = NULL
91        prefix = str2bytes(prefix, &cval)
92        CHKERR( DMSetOptionsPrefix(self.dm, cval) )
93
94    def setFromOptions(self):
95        CHKERR( DMSetFromOptions(self.dm) )
96
97    def setUp(self):
98        CHKERR( DMSetUp(self.dm) )
99        return self
100
101    # --- application context ---
102
103    def setAppCtx(self, appctx):
104        self.set_attr('__appctx__', appctx)
105
106    def getAppCtx(self):
107        return self.get_attr('__appctx__')
108
109    #
110
111    def setBasicAdjacency(self, useCone, useClosure):
112        cdef PetscBool uC  = useCone
113        cdef PetscBool uCl = useClosure
114        CHKERR( DMSetBasicAdjacency(self.dm, uC, uCl) )
115
116    def getBasicAdjacency(self):
117        cdef PetscBool uC  = PETSC_FALSE
118        cdef PetscBool uCl = PETSC_FALSE
119        CHKERR( DMGetBasicAdjacency(self.dm, &uC, &uCl) )
120        return toBool(uC), toBool(uCl)
121
122    def setFieldAdjacency(self, field, useCone, useClosure):
123        cdef PetscInt  f   = asInt(field)
124        cdef PetscBool uC  = useCone
125        cdef PetscBool uCl = useClosure
126        CHKERR( DMSetAdjacency(self.dm, f, uC, uCl) )
127
128    def getFieldAdjacency(self, field):
129        cdef PetscInt  f   = asInt(field)
130        cdef PetscBool uC  = PETSC_FALSE
131        cdef PetscBool uCl = PETSC_FALSE
132        CHKERR( DMGetAdjacency(self.dm, f, &uC, &uCl) )
133        return toBool(uC), toBool(uCl)
134
135    #
136
137    def setNumFields(self, numFields):
138        cdef PetscInt cnum = asInt(numFields)
139        CHKERR( DMSetNumFields(self.dm, cnum) )
140
141    def getNumFields(self):
142        cdef PetscInt cnum = 0
143        CHKERR( DMGetNumFields(self.dm, &cnum) )
144        return toInt(cnum)
145
146    def setField(self, index, Object field, label=None):
147        cdef PetscInt     cidx = asInt(index)
148        cdef PetscObject  cobj = field.obj[0]
149        cdef PetscDMLabel clbl = NULL
150        assert label is None
151        CHKERR( DMSetField(self.dm, cidx, clbl, cobj) )
152
153    def getField(self, index):
154        cdef PetscInt     cidx = asInt(index)
155        cdef PetscObject  cobj = NULL
156        cdef PetscDMLabel clbl = NULL
157        CHKERR( DMGetField(self.dm, cidx, &clbl, &cobj) )
158        assert clbl == NULL
159        cdef Object field = subtype_Object(cobj)()
160        field.obj[0] = cobj
161        PetscINCREF(field.obj)
162        return (field, None)
163
164    def addField(self, Object field, label=None):
165        cdef PetscObject  cobj = field.obj[0]
166        cdef PetscDMLabel clbl = NULL
167        assert label is None
168        CHKERR( DMAddField(self.dm, clbl, cobj) )
169
170    def copyFields(self, DM dm):
171        CHKERR( DMCopyFields(self.dm, dm.dm) )
172
173    def createDS(self):
174        CHKERR( DMCreateDS(self.dm) )
175
176    def clearDS(self):
177        CHKERR( DMClearDS(self.dm) )
178
179    def getDS(self):
180        cdef DS ds = DS()
181        CHKERR( DMGetDS(self.dm, &ds.ds) )
182        PetscINCREF(ds.obj)
183        return ds
184
185    def copyDS(self, DM dm):
186        CHKERR( DMCopyDS(self.dm, dm.dm) )
187
188    def copyDisc(self, DM dm):
189        CHKERR( DMCopyDisc(self.dm, dm.dm) )
190
191    #
192
193    def getBlockSize(self):
194        cdef PetscInt bs = 1
195        CHKERR( DMGetBlockSize(self.dm, &bs) )
196        return toInt(bs)
197
198    def setVecType(self, vec_type):
199        cdef PetscVecType vtype = NULL
200        vec_type = str2bytes(vec_type, &vtype)
201        CHKERR( DMSetVecType(self.dm, vtype) )
202
203    def createGlobalVec(self):
204        cdef Vec vg = Vec()
205        CHKERR( DMCreateGlobalVector(self.dm, &vg.vec) )
206        return vg
207
208    def createLocalVec(self):
209        cdef Vec vl = Vec()
210        CHKERR( DMCreateLocalVector(self.dm, &vl.vec) )
211        return vl
212
213    def getGlobalVec(self):
214        cdef Vec vg = Vec()
215        CHKERR( DMGetGlobalVector(self.dm, &vg.vec) )
216        PetscINCREF(vg.obj)
217        return vg
218
219    def restoreGlobalVec(self, Vec vg):
220        CHKERR( PetscObjectDereference(<PetscObject>vg.vec) )
221        CHKERR( DMRestoreGlobalVector(self.dm, &vg.vec) )
222
223    def getLocalVec(self):
224        cdef Vec vl = Vec()
225        CHKERR( DMGetLocalVector(self.dm, &vl.vec) )
226        PetscINCREF(vl.obj)
227        return vl
228
229    def restoreLocalVec(self, Vec vl):
230        CHKERR( PetscObjectDereference(<PetscObject>vl.vec) )
231        CHKERR( DMRestoreLocalVector(self.dm, &vl.vec) )
232
233    def globalToLocal(self, Vec vg, Vec vl, addv=None):
234        cdef PetscInsertMode im = insertmode(addv)
235        CHKERR( DMGlobalToLocalBegin(self.dm, vg.vec, im, vl.vec) )
236        CHKERR( DMGlobalToLocalEnd  (self.dm, vg.vec, im, vl.vec) )
237
238    def localToGlobal(self, Vec vl, Vec vg, addv=None):
239        cdef PetscInsertMode im = insertmode(addv)
240        CHKERR( DMLocalToGlobalBegin(self.dm, vl.vec, im, vg.vec) )
241        CHKERR( DMLocalToGlobalEnd(self.dm, vl.vec, im, vg.vec) )
242
243    def localToLocal(self, Vec vl, Vec vlg, addv=None):
244        cdef PetscInsertMode im = insertmode(addv)
245        CHKERR( DMLocalToLocalBegin(self.dm, vl.vec, im, vlg.vec) )
246        CHKERR( DMLocalToLocalEnd  (self.dm, vl.vec, im, vlg.vec) )
247
248    def getLGMap(self):
249        cdef LGMap lgm = LGMap()
250        CHKERR( DMGetLocalToGlobalMapping(self.dm, &lgm.lgm) )
251        PetscINCREF(lgm.obj)
252        return lgm
253
254    #
255
256    def getCoordinateDM(self):
257        cdef DM cdm = type(self)()
258        CHKERR( DMGetCoordinateDM(self.dm, &cdm.dm) )
259        PetscINCREF(cdm.obj)
260        return cdm
261
262    def getCoordinateSection(self):
263        cdef Section sec = Section()
264        CHKERR( DMGetCoordinateSection(self.dm, &sec.sec) )
265        PetscINCREF(sec.obj)
266        return sec
267
268    def setCoordinates(self, Vec c):
269        CHKERR( DMSetCoordinates(self.dm, c.vec) )
270
271    def getCoordinates(self):
272        cdef Vec c = Vec()
273        CHKERR( DMGetCoordinates(self.dm, &c.vec) )
274        PetscINCREF(c.obj)
275        return c
276
277    def setCoordinatesLocal(self, Vec c):
278        CHKERR( DMSetCoordinatesLocal(self.dm, c.vec) )
279
280    def getCoordinatesLocal(self):
281        cdef Vec c = Vec()
282        CHKERR( DMGetCoordinatesLocal(self.dm, &c.vec) )
283        PetscINCREF(c.obj)
284        return c
285
286    def getBoundingBox(self):
287        cdef PetscInt i,dim=0
288        CHKERR( DMGetCoordinateDim(self.dm, &dim) )
289        cdef PetscReal gmin[3], gmax[3]
290        CHKERR( DMGetBoundingBox(self.dm, gmin, gmax) )
291        return tuple([(toReal(gmin[i]), toReal(gmax[i]))
292                      for i from 0 <= i < dim])
293
294    def getLocalBoundingBox(self):
295        cdef PetscInt i,dim=0
296        CHKERR( DMGetCoordinateDim(self.dm, &dim) )
297        cdef PetscReal lmin[3], lmax[3]
298        CHKERR( DMGetLocalBoundingBox(self.dm, lmin, lmax) )
299        return tuple([(toReal(lmin[i]), toReal(lmax[i]))
300                      for i from 0 <= i < dim])
301
302    def localizeCoordinates(self):
303        CHKERR( DMLocalizeCoordinates(self.dm) )
304    #
305
306    def setMatType(self, mat_type):
307        """Set matrix type to be used by DM.createMat"""
308        cdef PetscMatType mtype = NULL
309        mat_type = str2bytes(mat_type, &mtype)
310        CHKERR( DMSetMatType(self.dm, mtype) )
311
312    def createMat(self):
313        cdef Mat mat = Mat()
314        CHKERR( DMCreateMatrix(self.dm, &mat.mat) )
315        return mat
316
317    def createInterpolation(self, DM dm):
318        cdef Mat A = Mat()
319        cdef Vec scale = Vec()
320        CHKERR( DMCreateInterpolation(self.dm, dm.dm,
321                                   &A.mat, &scale.vec))
322        return(A, scale)
323
324    def createInjection(self, DM dm):
325        cdef Mat inject = Mat()
326        CHKERR( DMCreateInjection(self.dm, dm.dm, &inject.mat) )
327        return inject
328
329    def createRestriction(self, DM dm):
330        cdef Mat mat = Mat()
331        CHKERR( DMCreateRestriction(self.dm, dm.dm, &mat.mat) )
332        return mat
333
334    def convert(self, dm_type):
335        cdef PetscDMType cval = NULL
336        dm_type = str2bytes(dm_type, &cval)
337        cdef PetscDM newdm = NULL
338        CHKERR( DMConvert(self.dm, cval, &newdm) )
339        cdef DM dm = <DM>subtype_DM(newdm)()
340        dm.dm = newdm
341        return dm
342
343    def refine(self, comm=None):
344        cdef MPI_Comm dmcomm = MPI_COMM_NULL
345        CHKERR( PetscObjectGetComm(<PetscObject>self.dm, &dmcomm) )
346        dmcomm = def_Comm(comm, dmcomm)
347        cdef PetscDM newdm = NULL
348        CHKERR( DMRefine(self.dm, dmcomm, &newdm) )
349        cdef DM dm = subtype_DM(newdm)()
350        dm.dm = newdm
351        return dm
352
353    def coarsen(self, comm=None):
354        cdef MPI_Comm dmcomm = MPI_COMM_NULL
355        CHKERR( PetscObjectGetComm(<PetscObject>self.dm, &dmcomm) )
356        dmcomm = def_Comm(comm, dmcomm)
357        cdef PetscDM newdm = NULL
358        CHKERR( DMCoarsen(self.dm, dmcomm, &newdm) )
359        cdef DM dm = subtype_DM(newdm)()
360        dm.dm = newdm
361        return dm
362
363    def refineHierarchy(self, nlevels):
364        cdef PetscInt i, n = asInt(nlevels)
365        cdef PetscDM *newdmf = NULL
366        cdef object tmp = oarray_p(empty_p(n), NULL, <void**>&newdmf)
367        CHKERR( DMRefineHierarchy(self.dm, n, newdmf) )
368        cdef DM dmf = None
369        cdef list hierarchy = []
370        for i from 0 <= i < n:
371            dmf = subtype_DM(newdmf[i])()
372            dmf.dm = newdmf[i]
373            hierarchy.append(dmf)
374        return hierarchy
375
376    def coarsenHierarchy(self, nlevels):
377        cdef PetscInt i, n = asInt(nlevels)
378        cdef PetscDM *newdmc = NULL
379        cdef object tmp = oarray_p(empty_p(n),NULL, <void**>&newdmc)
380        CHKERR( DMCoarsenHierarchy(self.dm, n, newdmc) )
381        cdef DM dmc = None
382        cdef list hierarchy = []
383        for i from 0 <= i < n:
384            dmc = subtype_DM(newdmc[i])()
385            dmc.dm = newdmc[i]
386            hierarchy.append(dmc)
387        return hierarchy
388
389    def getRefineLevel(self):
390        cdef PetscInt n = 0
391        CHKERR( DMGetRefineLevel(self.dm, &n) )
392        return toInt(n)
393
394    def setRefineLevel(self, level):
395        cdef PetscInt clevel = asInt(level)
396        CHKERR( DMSetRefineLevel(self.dm, clevel) )
397
398    def getCoarsenLevel(self):
399        cdef PetscInt n = 0
400        CHKERR( DMGetCoarsenLevel(self.dm, &n) )
401        return toInt(n)
402
403    #
404
405    def adaptLabel(self, label):
406        cdef const char *cval = NULL
407        cdef PetscDMLabel clbl = NULL
408        label = str2bytes(label, &cval)
409        CHKERR( DMGetLabel(self.dm, cval, &clbl) )
410        cdef DM newdm = DMPlex()
411        CHKERR( DMAdaptLabel(self.dm, clbl, &newdm.dm) )
412        return newdm
413
414    def adaptMetric(self, Vec metric, label=None):
415        cdef const char *cval = NULL
416        cdef PetscDMLabel clbl = NULL
417        label = str2bytes(label, &cval)
418        if cval == NULL: cval = b"" # XXX Should be fixed upstream
419        CHKERR( DMGetLabel(self.dm, cval, &clbl) )
420        cdef DM newdm = DMPlex()
421        CHKERR( DMAdaptMetric(self.dm, metric.vec, clbl, &newdm.dm) )
422        return newdm
423
424    def getLabel(self, name):
425        cdef const char *cname = NULL
426        cdef DMLabel dmlabel = DMLabel()
427        name = str2bytes(name, &cname)
428        CHKERR( DMGetLabel(self.dm, cname, &dmlabel.dmlabel) )
429        PetscINCREF(dmlabel.obj)
430        return dmlabel
431
432    #
433
434    def setSection(self, Section sec):
435        CHKERR( DMSetSection(self.dm, sec.sec) )
436
437    def getSection(self):
438        cdef Section sec = Section()
439        CHKERR( DMGetSection(self.dm, &sec.sec) )
440        PetscINCREF(sec.obj)
441        return sec
442
443    def setGlobalSection(self, Section sec):
444        CHKERR( DMSetGlobalSection(self.dm, sec.sec) )
445
446    def getGlobalSection(self):
447        cdef Section sec = Section()
448        CHKERR( DMGetGlobalSection(self.dm, &sec.sec) )
449        PetscINCREF(sec.obj)
450        return sec
451
452    setDefaultSection = setSection
453    getDefaultSection = getSection
454    setDefaultGlobalSection = setGlobalSection
455    getDefaultGlobalSection = getGlobalSection
456
457    def createSectionSF(self, Section localsec, Section globalsec):
458        CHKERR( DMCreateSectionSF(self.dm, localsec.sec, globalsec.sec) )
459
460    def getSectionSF(self):
461        cdef SF sf = SF()
462        CHKERR( DMGetSectionSF(self.dm, &sf.sf) )
463        PetscINCREF(sf.obj)
464        return sf
465
466    def setSectionSF(self, SF sf):
467        CHKERR( DMSetSectionSF(self.dm, sf.sf) )
468
469    createDefaultSF = createSectionSF
470    getDefaultSF = getSectionSF
471    setDefaultSF = setSectionSF
472
473    def getPointSF(self):
474        cdef SF sf = SF()
475        CHKERR( DMGetPointSF(self.dm, &sf.sf) )
476        PetscINCREF(sf.obj)
477        return sf
478
479    def setPointSF(self, SF sf):
480        CHKERR( DMSetPointSF(self.dm, sf.sf) )
481
482    def getNumLabels(self):
483        cdef PetscInt nLabels = 0
484        CHKERR( DMGetNumLabels(self.dm, &nLabels) )
485        return toInt(nLabels)
486
487    def getLabelName(self, index):
488        cdef PetscInt cindex = asInt(index)
489        cdef const char *cname = NULL
490        CHKERR( DMGetLabelName(self.dm, cindex, &cname) )
491        return bytes2str(cname)
492
493    def hasLabel(self, name):
494        cdef PetscBool flag = PETSC_FALSE
495        cdef const char *cname = NULL
496        name = str2bytes(name, &cname)
497        CHKERR( DMHasLabel(self.dm, cname, &flag) )
498        return toBool(flag)
499
500    def createLabel(self, name):
501        cdef const char *cname = NULL
502        name = str2bytes(name, &cname)
503        CHKERR( DMCreateLabel(self.dm, cname) )
504
505    def removeLabel(self, name):
506        cdef const char *cname = NULL
507        cdef PetscDMLabel clbl = NULL
508        name = str2bytes(name, &cname)
509        CHKERR( DMRemoveLabel(self.dm, cname, &clbl) )
510        # TODO: Once DMLabel is wrapped, this should return the label, like the C function.
511        CHKERR( DMLabelDestroy(&clbl) )
512
513    def getLabelValue(self, name, point):
514        cdef PetscInt cpoint = asInt(point), value = 0
515        cdef const char *cname = NULL
516        name = str2bytes(name, &cname)
517        CHKERR( DMGetLabelValue(self.dm, cname, cpoint, &value) )
518        return toInt(value)
519
520    def setLabelValue(self, name, point, value):
521        cdef PetscInt cpoint = asInt(point), cvalue = asInt(value)
522        cdef const char *cname = NULL
523        name = str2bytes(name, &cname)
524        CHKERR( DMSetLabelValue(self.dm, cname, cpoint, cvalue) )
525
526    def clearLabelValue(self, name, point, value):
527        cdef PetscInt cpoint = asInt(point), cvalue = asInt(value)
528        cdef const char *cname = NULL
529        name = str2bytes(name, &cname)
530        CHKERR( DMClearLabelValue(self.dm, cname, cpoint, cvalue) )
531
532    def getLabelSize(self, name):
533        cdef PetscInt size = 0
534        cdef const char *cname = NULL
535        name = str2bytes(name, &cname)
536        CHKERR( DMGetLabelSize(self.dm, cname, &size) )
537        return toInt(size)
538
539    def getLabelIdIS(self, name):
540        cdef const char *cname = NULL
541        name = str2bytes(name, &cname)
542        cdef IS lis = IS()
543        CHKERR( DMGetLabelIdIS(self.dm, cname, &lis.iset) )
544        return lis
545
546    def getStratumSize(self, name, value):
547        cdef PetscInt size = 0
548        cdef PetscInt cvalue = asInt(value)
549        cdef const char *cname = NULL
550        name = str2bytes(name, &cname)
551        CHKERR( DMGetStratumSize(self.dm, cname, cvalue, &size) )
552        return toInt(size)
553
554    def getStratumIS(self, name, value):
555        cdef PetscInt cvalue = asInt(value)
556        cdef const char *cname = NULL
557        name = str2bytes(name, &cname)
558        cdef IS sis = IS()
559        CHKERR( DMGetStratumIS(self.dm, cname, cvalue, &sis.iset) )
560        return sis
561
562    def clearLabelStratum(self, name, value):
563        cdef PetscInt cvalue = asInt(value)
564        cdef const char *cname = NULL
565        name = str2bytes(name, &cname)
566        CHKERR( DMClearLabelStratum(self.dm, cname, cvalue) )
567
568    def setLabelOutput(self, name, output):
569        cdef const char *cname = NULL
570        name = str2bytes(name, &cname)
571        cdef PetscBool coutput = output
572        CHKERR( DMSetLabelOutput(self.dm, cname, coutput) )
573
574    def getLabelOutput(self, name):
575        cdef const char *cname = NULL
576        name = str2bytes(name, &cname)
577        cdef PetscBool coutput = PETSC_FALSE
578        CHKERR( DMGetLabelOutput(self.dm, cname, &coutput) )
579        return coutput
580
581    # backward compatibility
582    createGlobalVector = createGlobalVec
583    createLocalVector = createLocalVec
584    getMatrix = createMatrix = createMat
585
586    def setKSPComputeOperators(self, operators, args=None, kargs=None):
587        if args  is None: args  = ()
588        if kargs is None: kargs = {}
589        context = (operators, args, kargs)
590        self.set_attr('__operators__', context)
591        CHKERR( DMKSPSetComputeOperators(self.dm, KSP_ComputeOps, <void*>context) )
592
593    def createFieldDecomposition(self):
594        cdef PetscInt clen = 0
595        cdef PetscIS *cis = NULL
596        cdef PetscDM *cdm = NULL
597        cdef char** cnamelist = NULL
598
599        CHKERR( DMCreateFieldDecomposition(self.dm, &clen, &cnamelist, &cis, &cdm) )
600
601        cdef list isets = [ref_IS(cis[i]) for i from 0 <= i < clen]
602        cdef list dms   = []
603        cdef list names = []
604        cdef DM dm = None
605
606        for i from 0 <= i < clen:
607            if cdm != NULL:
608                dm = subtype_DM(cdm[i])()
609                dm.dm = cdm[i]
610                PetscINCREF(dm.obj)
611                dms.append(dm)
612            else:
613                dms.append(None)
614
615            name = bytes2str(cnamelist[i])
616            names.append(name)
617            CHKERR( PetscFree(cnamelist[i]) )
618
619            CHKERR( ISDestroy(&cis[i]) )
620            CHKERR( DMDestroy(&cdm[i]) )
621
622        CHKERR( PetscFree(cis) )
623        CHKERR( PetscFree(cdm) )
624        CHKERR( PetscFree(cnamelist) )
625
626        return (names, isets, dms)
627
628    def setSNESFunction(self, function, args=None, kargs=None):
629        if function is not None:
630            if args  is None: args  = ()
631            if kargs is None: kargs = {}
632            context = (function, args, kargs)
633            self.set_attr('__function__', context)
634            CHKERR( DMSNESSetFunction(self.dm, SNES_Function, <void*>context) )
635        else:
636            CHKERR( DMSNESSetFunction(self.dm, NULL, NULL) )
637
638    def setSNESJacobian(self, jacobian, args=None, kargs=None):
639        if jacobian is not None:
640            if args  is None: args  = ()
641            if kargs is None: kargs = {}
642            context = (jacobian, args, kargs)
643            self.set_attr('__jacobian__', context)
644            CHKERR( DMSNESSetJacobian(self.dm, SNES_Jacobian, <void*>context) )
645        else:
646            CHKERR( DMSNESSetJacobian(self.dm, NULL, NULL) )
647
648    def addCoarsenHook(self, coarsenhook, restricthook, args=None, kargs=None):
649        if args  is None: args  = ()
650        if kargs is None: kargs = {}
651
652        if coarsenhook is not None:
653            coarsencontext = (coarsenhook, args, kargs)
654
655            coarsenhooks = self.get_attr('__coarsenhooks__')
656            if coarsenhooks is None:
657                coarsenhooks = [coarsencontext]
658                CHKERR( DMCoarsenHookAdd(self.dm, DM_PyCoarsenHook, NULL, <void*>NULL) )
659            else:
660                coarsenhooks.append(coarsencontext)
661            self.set_attr('__coarsenhooks__', coarsenhooks)
662
663        if restricthook is not None:
664            restrictcontext = (restricthook, args, kargs)
665
666            restricthooks = self.get_attr('__restricthooks__')
667            if restricthooks is None:
668                restricthooks = [restrictcontext]
669                CHKERR( DMCoarsenHookAdd(self.dm, NULL, DM_PyRestrictHook, <void*>NULL) )
670            else:
671                restricthooks.append(restrictcontext)
672            self.set_attr('__restricthooks__', restricthooks)
673
674    # --- application context ---
675
676    property appctx:
677        def __get__(self):
678            return self.getAppCtx()
679        def __set__(self, value):
680            self.setAppCtx(value)
681
682    # --- discretization space ---
683
684    property ds:
685        def __get__(self):
686            return self.getDS()
687        def __set__(self, value):
688            self.setDS(value)
689
690# --------------------------------------------------------------------
691
692del DMType
693del DMBoundaryType
694
695# --------------------------------------------------------------------
696