Actual source code: dmshell.c
1: #include <petscdmshell.h>
2: #include <petscmat.h>
3: #include <petsc/private/dmimpl.h>
5: typedef struct {
6: Vec Xglobal;
7: Vec Xlocal;
8: Mat A;
9: VecScatter gtol;
10: VecScatter ltog;
11: VecScatter ltol;
12: void *ctx;
13: PetscErrorCode (*destroyctx)(void *);
14: } DM_Shell;
16: /*@
17: DMGlobalToLocalBeginDefaultShell - Uses the GlobalToLocal VecScatter context set by the user to begin a global to local scatter
18: Collective
20: Input Parameters:
21: + dm - shell DM
22: . g - global vector
23: . mode - InsertMode
24: - l - local vector
26: Level: advanced
28: Note: This is not normally called directly by user code, generally user code calls DMGlobalToLocalBegin() and DMGlobalToLocalEnd(). If the user provides their own custom routines to DMShellSetLocalToGlobal() then those routines might have reason to call this function.
30: .seealso: `DMGlobalToLocalEndDefaultShell()`
31: @*/
32: PetscErrorCode DMGlobalToLocalBeginDefaultShell(DM dm, Vec g, InsertMode mode, Vec l)
33: {
34: DM_Shell *shell = (DM_Shell *)dm->data;
37: VecScatterBegin(shell->gtol, g, l, mode, SCATTER_FORWARD);
38: return 0;
39: }
41: /*@
42: DMGlobalToLocalEndDefaultShell - Uses the GlobalToLocal VecScatter context set by the user to end a global to local scatter
43: Collective
45: Input Parameters:
46: + dm - shell DM
47: . g - global vector
48: . mode - InsertMode
49: - l - local vector
51: Level: advanced
53: .seealso: `DMGlobalToLocalBeginDefaultShell()`
54: @*/
55: PetscErrorCode DMGlobalToLocalEndDefaultShell(DM dm, Vec g, InsertMode mode, Vec l)
56: {
57: DM_Shell *shell = (DM_Shell *)dm->data;
60: VecScatterEnd(shell->gtol, g, l, mode, SCATTER_FORWARD);
61: return 0;
62: }
64: /*@
65: DMLocalToGlobalBeginDefaultShell - Uses the LocalToGlobal VecScatter context set by the user to begin a local to global scatter
66: Collective
68: Input Parameters:
69: + dm - shell DM
70: . l - local vector
71: . mode - InsertMode
72: - g - global vector
74: Level: advanced
76: Note: This is not normally called directly by user code, generally user code calls DMLocalToGlobalBegin() and DMLocalToGlobalEnd(). If the user provides their own custom routines to DMShellSetLocalToGlobal() then those routines might have reason to call this function.
78: .seealso: `DMLocalToGlobalEndDefaultShell()`
79: @*/
80: PetscErrorCode DMLocalToGlobalBeginDefaultShell(DM dm, Vec l, InsertMode mode, Vec g)
81: {
82: DM_Shell *shell = (DM_Shell *)dm->data;
85: VecScatterBegin(shell->ltog, l, g, mode, SCATTER_FORWARD);
86: return 0;
87: }
89: /*@
90: DMLocalToGlobalEndDefaultShell - Uses the LocalToGlobal VecScatter context set by the user to end a local to global scatter
91: Collective
93: Input Parameters:
94: + dm - shell DM
95: . l - local vector
96: . mode - InsertMode
97: - g - global vector
99: Level: advanced
101: .seealso: `DMLocalToGlobalBeginDefaultShell()`
102: @*/
103: PetscErrorCode DMLocalToGlobalEndDefaultShell(DM dm, Vec l, InsertMode mode, Vec g)
104: {
105: DM_Shell *shell = (DM_Shell *)dm->data;
108: VecScatterEnd(shell->ltog, l, g, mode, SCATTER_FORWARD);
109: return 0;
110: }
112: /*@
113: DMLocalToLocalBeginDefaultShell - Uses the LocalToLocal VecScatter context set by the user to begin a local to local scatter
114: Collective
116: Input Parameters:
117: + dm - shell DM
118: . g - the original local vector
119: - mode - InsertMode
121: Output Parameter:
122: . l - the local vector with correct ghost values
124: Level: advanced
126: Note: This is not normally called directly by user code, generally user code calls DMLocalToLocalBegin() and DMLocalToLocalEnd(). If the user provides their own custom routines to DMShellSetLocalToLocal() then those routines might have reason to call this function.
128: .seealso: `DMLocalToLocalEndDefaultShell()`
129: @*/
130: PetscErrorCode DMLocalToLocalBeginDefaultShell(DM dm, Vec g, InsertMode mode, Vec l)
131: {
132: DM_Shell *shell = (DM_Shell *)dm->data;
135: VecScatterBegin(shell->ltol, g, l, mode, SCATTER_FORWARD);
136: return 0;
137: }
139: /*@
140: DMLocalToLocalEndDefaultShell - Uses the LocalToLocal VecScatter context set by the user to end a local to local scatter
141: Collective
143: Input Parameters:
144: + dm - shell DM
145: . g - the original local vector
146: - mode - InsertMode
148: Output Parameter:
149: . l - the local vector with correct ghost values
151: Level: advanced
153: .seealso: `DMLocalToLocalBeginDefaultShell()`
154: @*/
155: PetscErrorCode DMLocalToLocalEndDefaultShell(DM dm, Vec g, InsertMode mode, Vec l)
156: {
157: DM_Shell *shell = (DM_Shell *)dm->data;
160: VecScatterEnd(shell->ltol, g, l, mode, SCATTER_FORWARD);
161: return 0;
162: }
164: static PetscErrorCode DMCreateMatrix_Shell(DM dm, Mat *J)
165: {
166: DM_Shell *shell = (DM_Shell *)dm->data;
167: Mat A;
171: if (!shell->A) {
172: if (shell->Xglobal) {
173: PetscInt m, M;
174: PetscInfo(dm, "Naively creating matrix using global vector distribution without preallocation\n");
175: VecGetSize(shell->Xglobal, &M);
176: VecGetLocalSize(shell->Xglobal, &m);
177: MatCreate(PetscObjectComm((PetscObject)dm), &shell->A);
178: MatSetSizes(shell->A, m, m, M, M);
179: MatSetType(shell->A, dm->mattype);
180: MatSetUp(shell->A);
181: } else SETERRQ(PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Must call DMShellSetMatrix(), DMShellSetCreateMatrix(), or provide a vector");
182: }
183: A = shell->A;
184: MatDuplicate(A, MAT_SHARE_NONZERO_PATTERN, J);
185: MatSetDM(*J, dm);
186: return 0;
187: }
189: PetscErrorCode DMCreateGlobalVector_Shell(DM dm, Vec *gvec)
190: {
191: DM_Shell *shell = (DM_Shell *)dm->data;
192: Vec X;
196: *gvec = NULL;
197: X = shell->Xglobal;
199: /* Need to create a copy in order to attach the DM to the vector */
200: VecDuplicate(X, gvec);
201: VecZeroEntries(*gvec);
202: VecSetDM(*gvec, dm);
203: return 0;
204: }
206: PetscErrorCode DMCreateLocalVector_Shell(DM dm, Vec *gvec)
207: {
208: DM_Shell *shell = (DM_Shell *)dm->data;
209: Vec X;
213: *gvec = NULL;
214: X = shell->Xlocal;
216: /* Need to create a copy in order to attach the DM to the vector */
217: VecDuplicate(X, gvec);
218: VecZeroEntries(*gvec);
219: VecSetDM(*gvec, dm);
220: return 0;
221: }
223: /*@
224: DMShellSetDestroyContext - set a function that destroys the context provided with `DMShellSetContext()`
226: Collective
228: Input Parameters:
229: . ctx - the context
231: Level: advanced
233: .seealso: `DMShellSetContext()`, `DMShellGetContext()`
234: @*/
235: PetscErrorCode DMShellSetDestroyContext(DM dm, PetscErrorCode (*destroyctx)(void *))
236: {
237: DM_Shell *shell = (DM_Shell *)dm->data;
238: PetscBool isshell;
241: PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell);
242: if (!isshell) return 0;
243: shell->destroyctx = destroyctx;
244: return 0;
245: }
247: /*@
248: DMShellSetContext - set some data to be usable by this DM
250: Collective
252: Input Parameters:
253: + dm - shell DM
254: - ctx - the context
256: Level: advanced
258: .seealso: `DMCreateMatrix()`, `DMShellGetContext()`
259: @*/
260: PetscErrorCode DMShellSetContext(DM dm, void *ctx)
261: {
262: DM_Shell *shell = (DM_Shell *)dm->data;
263: PetscBool isshell;
266: PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell);
267: if (!isshell) return 0;
268: shell->ctx = ctx;
269: return 0;
270: }
272: /*@
273: DMShellGetContext - Returns the user-provided context associated to the DM
275: Collective
277: Input Parameter:
278: . dm - shell DM
280: Output Parameter:
281: . ctx - the context
283: Level: advanced
285: .seealso: `DMCreateMatrix()`, `DMShellSetContext()`
286: @*/
287: PetscErrorCode DMShellGetContext(DM dm, void *ctx)
288: {
289: DM_Shell *shell = (DM_Shell *)dm->data;
290: PetscBool isshell;
293: PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell);
295: *(void **)ctx = shell->ctx;
296: return 0;
297: }
299: /*@
300: DMShellSetMatrix - sets a template matrix associated with the DMShell
302: Collective
304: Input Parameters:
305: + dm - shell DM
306: - J - template matrix
308: Level: advanced
310: Developer Notes:
311: To avoid circular references, if J is already associated to the same DM, then MatDuplicate(SHARE_NONZERO_PATTERN) is called, followed by removing the DM reference from the private template.
313: .seealso: `DMCreateMatrix()`, `DMShellSetCreateMatrix()`, `DMShellSetContext()`, `DMShellGetContext()`
314: @*/
315: PetscErrorCode DMShellSetMatrix(DM dm, Mat J)
316: {
317: DM_Shell *shell = (DM_Shell *)dm->data;
318: PetscBool isshell;
319: DM mdm;
323: PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell);
324: if (!isshell) return 0;
325: if (J == shell->A) return 0;
326: MatGetDM(J, &mdm);
327: PetscObjectReference((PetscObject)J);
328: MatDestroy(&shell->A);
329: if (mdm == dm) {
330: MatDuplicate(J, MAT_SHARE_NONZERO_PATTERN, &shell->A);
331: MatSetDM(shell->A, NULL);
332: } else shell->A = J;
333: return 0;
334: }
336: /*@C
337: DMShellSetCreateMatrix - sets the routine to create a matrix associated with the shell DM
339: Logically Collective on dm
341: Input Parameters:
342: + dm - the shell DM
343: - func - the function to create a matrix
345: Level: advanced
347: .seealso: `DMCreateMatrix()`, `DMShellSetMatrix()`, `DMShellSetContext()`, `DMShellGetContext()`
348: @*/
349: PetscErrorCode DMShellSetCreateMatrix(DM dm, PetscErrorCode (*func)(DM, Mat *))
350: {
352: dm->ops->creatematrix = func;
353: return 0;
354: }
356: /*@
357: DMShellSetGlobalVector - sets a template global vector associated with the DMShell
359: Logically Collective on dm
361: Input Parameters:
362: + dm - shell DM
363: - X - template vector
365: Level: advanced
367: .seealso: `DMCreateGlobalVector()`, `DMShellSetMatrix()`, `DMShellSetCreateGlobalVector()`
368: @*/
369: PetscErrorCode DMShellSetGlobalVector(DM dm, Vec X)
370: {
371: DM_Shell *shell = (DM_Shell *)dm->data;
372: PetscBool isshell;
373: DM vdm;
377: PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell);
378: if (!isshell) return 0;
379: VecGetDM(X, &vdm);
380: /*
381: if the vector proposed as the new base global vector for the DM is a DM vector associated
382: with the same DM then the current base global vector for the DM is ok and if we replace it with the new one
383: we get a circular dependency that prevents the DM from being destroy when it should be.
384: This occurs when SNESSet/GetNPC() is used with a SNES that does not have a user provided
385: DM attached to it since the inner SNES (which shares the DM with the outer SNES) tries
386: to set its input vector (which is associated with the DM) as the base global vector.
387: Thanks to Juan P. Mendez Granado Re: [petsc-maint] Nonlinear conjugate gradien
388: for pointing out the problem.
389: */
390: if (vdm == dm) return 0;
391: PetscObjectReference((PetscObject)X);
392: VecDestroy(&shell->Xglobal);
393: shell->Xglobal = X;
394: return 0;
395: }
397: /*@
398: DMShellGetGlobalVector - Returns the template global vector associated with the DMShell, or NULL if it was not set
400: Not collective
402: Input Parameters:
403: + dm - shell DM
404: - X - template vector
406: Level: advanced
408: .seealso: `DMShellSetGlobalVector()`, `DMShellSetCreateGlobalVector()`, `DMCreateGlobalVector()`
409: @*/
410: PetscErrorCode DMShellGetGlobalVector(DM dm, Vec *X)
411: {
412: DM_Shell *shell = (DM_Shell *)dm->data;
413: PetscBool isshell;
417: PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell);
418: if (!isshell) return 0;
419: *X = shell->Xglobal;
420: return 0;
421: }
423: /*@C
424: DMShellSetCreateGlobalVector - sets the routine to create a global vector associated with the shell DM
426: Logically Collective
428: Input Parameters:
429: + dm - the shell DM
430: - func - the creation routine
432: Level: advanced
434: .seealso: `DMShellSetGlobalVector()`, `DMShellSetCreateMatrix()`, `DMShellSetContext()`, `DMShellGetContext()`
435: @*/
436: PetscErrorCode DMShellSetCreateGlobalVector(DM dm, PetscErrorCode (*func)(DM, Vec *))
437: {
439: dm->ops->createglobalvector = func;
440: return 0;
441: }
443: /*@
444: DMShellSetLocalVector - sets a template local vector associated with the DMShell
446: Logically Collective on dm
448: Input Parameters:
449: + dm - shell DM
450: - X - template vector
452: Level: advanced
454: .seealso: `DMCreateLocalVector()`, `DMShellSetMatrix()`, `DMShellSetCreateLocalVector()`
455: @*/
456: PetscErrorCode DMShellSetLocalVector(DM dm, Vec X)
457: {
458: DM_Shell *shell = (DM_Shell *)dm->data;
459: PetscBool isshell;
460: DM vdm;
464: PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell);
465: if (!isshell) return 0;
466: VecGetDM(X, &vdm);
467: /*
468: if the vector proposed as the new base global vector for the DM is a DM vector associated
469: with the same DM then the current base global vector for the DM is ok and if we replace it with the new one
470: we get a circular dependency that prevents the DM from being destroy when it should be.
471: This occurs when SNESSet/GetNPC() is used with a SNES that does not have a user provided
472: DM attached to it since the inner SNES (which shares the DM with the outer SNES) tries
473: to set its input vector (which is associated with the DM) as the base global vector.
474: Thanks to Juan P. Mendez Granado Re: [petsc-maint] Nonlinear conjugate gradien
475: for pointing out the problem.
476: */
477: if (vdm == dm) return 0;
478: PetscObjectReference((PetscObject)X);
479: VecDestroy(&shell->Xlocal);
480: shell->Xlocal = X;
481: return 0;
482: }
484: /*@C
485: DMShellSetCreateLocalVector - sets the routine to create a local vector associated with the shell DM
487: Logically Collective
489: Input Parameters:
490: + dm - the shell DM
491: - func - the creation routine
493: Level: advanced
495: .seealso: `DMShellSetLocalVector()`, `DMShellSetCreateMatrix()`, `DMShellSetContext()`, `DMShellGetContext()`
496: @*/
497: PetscErrorCode DMShellSetCreateLocalVector(DM dm, PetscErrorCode (*func)(DM, Vec *))
498: {
500: dm->ops->createlocalvector = func;
501: return 0;
502: }
504: /*@C
505: DMShellSetGlobalToLocal - Sets the routines used to perform a global to local scatter
507: Logically Collective on dm
509: Input Parameters
510: + dm - the shell DM
511: . begin - the routine that begins the global to local scatter
512: - end - the routine that ends the global to local scatter
514: Notes:
515: If these functions are not provided but DMShellSetGlobalToLocalVecScatter() is called then
516: DMGlobalToLocalBeginDefaultShell()/DMGlobalToLocalEndDefaultShell() are used to to perform the transfers
518: Level: advanced
520: .seealso: `DMShellSetLocalToGlobal()`, `DMGlobalToLocalBeginDefaultShell()`, `DMGlobalToLocalEndDefaultShell()`
521: @*/
522: PetscErrorCode DMShellSetGlobalToLocal(DM dm, PetscErrorCode (*begin)(DM, Vec, InsertMode, Vec), PetscErrorCode (*end)(DM, Vec, InsertMode, Vec))
523: {
525: dm->ops->globaltolocalbegin = begin;
526: dm->ops->globaltolocalend = end;
527: return 0;
528: }
530: /*@C
531: DMShellSetLocalToGlobal - Sets the routines used to perform a local to global scatter
533: Logically Collective on dm
535: Input Parameters
536: + dm - the shell DM
537: . begin - the routine that begins the local to global scatter
538: - end - the routine that ends the local to global scatter
540: Notes:
541: If these functions are not provided but DMShellSetLocalToGlobalVecScatter() is called then
542: DMLocalToGlobalBeginDefaultShell()/DMLocalToGlobalEndDefaultShell() are used to to perform the transfers
544: Level: advanced
546: .seealso: `DMShellSetGlobalToLocal()`
547: @*/
548: PetscErrorCode DMShellSetLocalToGlobal(DM dm, PetscErrorCode (*begin)(DM, Vec, InsertMode, Vec), PetscErrorCode (*end)(DM, Vec, InsertMode, Vec))
549: {
551: dm->ops->localtoglobalbegin = begin;
552: dm->ops->localtoglobalend = end;
553: return 0;
554: }
556: /*@C
557: DMShellSetLocalToLocal - Sets the routines used to perform a local to local scatter
559: Logically Collective on dm
561: Input Parameters
562: + dm - the shell DM
563: . begin - the routine that begins the local to local scatter
564: - end - the routine that ends the local to local scatter
566: Notes:
567: If these functions are not provided but DMShellSetLocalToLocalVecScatter() is called then
568: DMLocalToLocalBeginDefaultShell()/DMLocalToLocalEndDefaultShell() are used to to perform the transfers
570: Level: advanced
572: .seealso: `DMShellSetGlobalToLocal()`, `DMLocalToLocalBeginDefaultShell()`, `DMLocalToLocalEndDefaultShell()`
573: @*/
574: PetscErrorCode DMShellSetLocalToLocal(DM dm, PetscErrorCode (*begin)(DM, Vec, InsertMode, Vec), PetscErrorCode (*end)(DM, Vec, InsertMode, Vec))
575: {
577: dm->ops->localtolocalbegin = begin;
578: dm->ops->localtolocalend = end;
579: return 0;
580: }
582: /*@
583: DMShellSetGlobalToLocalVecScatter - Sets a VecScatter context for global to local communication
585: Logically Collective on dm
587: Input Parameters
588: + dm - the shell DM
589: - gtol - the global to local VecScatter context
591: Level: advanced
593: .seealso: `DMShellSetGlobalToLocal()`, `DMGlobalToLocalBeginDefaultShell()`, `DMGlobalToLocalEndDefaultShell()`
594: @*/
595: PetscErrorCode DMShellSetGlobalToLocalVecScatter(DM dm, VecScatter gtol)
596: {
597: DM_Shell *shell = (DM_Shell *)dm->data;
601: PetscObjectReference((PetscObject)gtol);
602: VecScatterDestroy(&shell->gtol);
603: shell->gtol = gtol;
604: return 0;
605: }
607: /*@
608: DMShellSetLocalToGlobalVecScatter - Sets a VecScatter context for local to global communication
610: Logically Collective on dm
612: Input Parameters
613: + dm - the shell DM
614: - ltog - the local to global VecScatter context
616: Level: advanced
618: .seealso: `DMShellSetLocalToGlobal()`, `DMLocalToGlobalBeginDefaultShell()`, `DMLocalToGlobalEndDefaultShell()`
619: @*/
620: PetscErrorCode DMShellSetLocalToGlobalVecScatter(DM dm, VecScatter ltog)
621: {
622: DM_Shell *shell = (DM_Shell *)dm->data;
626: PetscObjectReference((PetscObject)ltog);
627: VecScatterDestroy(&shell->ltog);
628: shell->ltog = ltog;
629: return 0;
630: }
632: /*@
633: DMShellSetLocalToLocalVecScatter - Sets a VecScatter context for local to local communication
635: Logically Collective on dm
637: Input Parameters
638: + dm - the shell DM
639: - ltol - the local to local VecScatter context
641: Level: advanced
643: .seealso: `DMShellSetLocalToLocal()`, `DMLocalToLocalBeginDefaultShell()`, `DMLocalToLocalEndDefaultShell()`
644: @*/
645: PetscErrorCode DMShellSetLocalToLocalVecScatter(DM dm, VecScatter ltol)
646: {
647: DM_Shell *shell = (DM_Shell *)dm->data;
651: PetscObjectReference((PetscObject)ltol);
652: VecScatterDestroy(&shell->ltol);
653: shell->ltol = ltol;
654: return 0;
655: }
657: /*@C
658: DMShellSetCoarsen - Set the routine used to coarsen the shell DM
660: Logically Collective on dm
662: Input Parameters
663: + dm - the shell DM
664: - coarsen - the routine that coarsens the DM
666: Level: advanced
668: .seealso: `DMShellSetRefine()`, `DMCoarsen()`, `DMShellGetCoarsen()`, `DMShellSetContext()`, `DMShellGetContext()`
669: @*/
670: PetscErrorCode DMShellSetCoarsen(DM dm, PetscErrorCode (*coarsen)(DM, MPI_Comm, DM *))
671: {
672: PetscBool isshell;
675: PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell);
676: if (!isshell) return 0;
677: dm->ops->coarsen = coarsen;
678: return 0;
679: }
681: /*@C
682: DMShellGetCoarsen - Get the routine used to coarsen the shell DM
684: Logically Collective on dm
686: Input Parameter:
687: . dm - the shell DM
689: Output Parameter:
690: . coarsen - the routine that coarsens the DM
692: Level: advanced
694: .seealso: `DMShellSetCoarsen()`, `DMCoarsen()`, `DMShellSetRefine()`, `DMRefine()`
695: @*/
696: PetscErrorCode DMShellGetCoarsen(DM dm, PetscErrorCode (**coarsen)(DM, MPI_Comm, DM *))
697: {
698: PetscBool isshell;
701: PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell);
703: *coarsen = dm->ops->coarsen;
704: return 0;
705: }
707: /*@C
708: DMShellSetRefine - Set the routine used to refine the shell DM
710: Logically Collective on dm
712: Input Parameters
713: + dm - the shell DM
714: - refine - the routine that refines the DM
716: Level: advanced
718: .seealso: `DMShellSetCoarsen()`, `DMRefine()`, `DMShellGetRefine()`, `DMShellSetContext()`, `DMShellGetContext()`
719: @*/
720: PetscErrorCode DMShellSetRefine(DM dm, PetscErrorCode (*refine)(DM, MPI_Comm, DM *))
721: {
722: PetscBool isshell;
725: PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell);
726: if (!isshell) return 0;
727: dm->ops->refine = refine;
728: return 0;
729: }
731: /*@C
732: DMShellGetRefine - Get the routine used to refine the shell DM
734: Logically Collective on dm
736: Input Parameter:
737: . dm - the shell DM
739: Output Parameter:
740: . refine - the routine that refines the DM
742: Level: advanced
744: .seealso: `DMShellSetCoarsen()`, `DMCoarsen()`, `DMShellSetRefine()`, `DMRefine()`
745: @*/
746: PetscErrorCode DMShellGetRefine(DM dm, PetscErrorCode (**refine)(DM, MPI_Comm, DM *))
747: {
748: PetscBool isshell;
751: PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell);
753: *refine = dm->ops->refine;
754: return 0;
755: }
757: /*@C
758: DMShellSetCreateInterpolation - Set the routine used to create the interpolation operator
760: Logically Collective on dm
762: Input Parameters
763: + dm - the shell DM
764: - interp - the routine to create the interpolation
766: Level: advanced
768: .seealso: `DMShellSetCreateInjection()`, `DMCreateInterpolation()`, `DMShellGetCreateInterpolation()`, `DMShellSetCreateRestriction()`, `DMShellSetContext()`, `DMShellGetContext()`
769: @*/
770: PetscErrorCode DMShellSetCreateInterpolation(DM dm, PetscErrorCode (*interp)(DM, DM, Mat *, Vec *))
771: {
772: PetscBool isshell;
775: PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell);
776: if (!isshell) return 0;
777: dm->ops->createinterpolation = interp;
778: return 0;
779: }
781: /*@C
782: DMShellGetCreateInterpolation - Get the routine used to create the interpolation operator
784: Logically Collective on dm
786: Input Parameter:
787: . dm - the shell DM
789: Output Parameter:
790: . interp - the routine to create the interpolation
792: Level: advanced
794: .seealso: `DMShellGetCreateInjection()`, `DMCreateInterpolation()`, `DMShellGetCreateRestriction()`, `DMShellSetContext()`, `DMShellGetContext()`
795: @*/
796: PetscErrorCode DMShellGetCreateInterpolation(DM dm, PetscErrorCode (**interp)(DM, DM, Mat *, Vec *))
797: {
798: PetscBool isshell;
801: PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell);
803: *interp = dm->ops->createinterpolation;
804: return 0;
805: }
807: /*@C
808: DMShellSetCreateRestriction - Set the routine used to create the restriction operator
810: Logically Collective on dm
812: Input Parameters
813: + dm - the shell DM
814: - striction- the routine to create the restriction
816: Level: advanced
818: .seealso: `DMShellSetCreateInjection()`, `DMCreateInterpolation()`, `DMShellGetCreateRestriction()`, `DMShellSetContext()`, `DMShellGetContext()`
819: @*/
820: PetscErrorCode DMShellSetCreateRestriction(DM dm, PetscErrorCode (*restriction)(DM, DM, Mat *))
821: {
822: PetscBool isshell;
825: PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell);
826: if (!isshell) return 0;
827: dm->ops->createrestriction = restriction;
828: return 0;
829: }
831: /*@C
832: DMShellGetCreateRestriction - Get the routine used to create the restriction operator
834: Logically Collective on dm
836: Input Parameter:
837: . dm - the shell DM
839: Output Parameter:
840: . restriction - the routine to create the restriction
842: Level: advanced
844: .seealso: `DMShellSetCreateInjection()`, `DMCreateInterpolation()`, `DMShellSetContext()`, `DMShellGetContext()`
845: @*/
846: PetscErrorCode DMShellGetCreateRestriction(DM dm, PetscErrorCode (**restriction)(DM, DM, Mat *))
847: {
848: PetscBool isshell;
851: PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell);
853: *restriction = dm->ops->createrestriction;
854: return 0;
855: }
857: /*@C
858: DMShellSetCreateInjection - Set the routine used to create the injection operator
860: Logically Collective on dm
862: Input Parameters:
863: + dm - the shell DM
864: - inject - the routine to create the injection
866: Level: advanced
868: .seealso: `DMShellSetCreateInterpolation()`, `DMCreateInjection()`, `DMShellGetCreateInjection()`, `DMShellSetContext()`, `DMShellGetContext()`
869: @*/
870: PetscErrorCode DMShellSetCreateInjection(DM dm, PetscErrorCode (*inject)(DM, DM, Mat *))
871: {
872: PetscBool isshell;
875: PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell);
876: if (!isshell) return 0;
877: dm->ops->createinjection = inject;
878: return 0;
879: }
881: /*@C
882: DMShellGetCreateInjection - Get the routine used to create the injection operator
884: Logically Collective on dm
886: Input Parameter:
887: . dm - the shell DM
889: Output Parameter:
890: . inject - the routine to create the injection
892: Level: advanced
894: .seealso: `DMShellGetCreateInterpolation()`, `DMCreateInjection()`, `DMShellSetContext()`, `DMShellGetContext()`
895: @*/
896: PetscErrorCode DMShellGetCreateInjection(DM dm, PetscErrorCode (**inject)(DM, DM, Mat *))
897: {
898: PetscBool isshell;
901: PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell);
903: *inject = dm->ops->createinjection;
904: return 0;
905: }
907: /*@C
908: DMShellSetCreateFieldDecomposition - Set the routine used to create a decomposition of fields for the shell DM
910: Logically Collective on dm
912: Input Parameters:
913: + dm - the shell DM
914: - decomp - the routine to create the decomposition
916: Level: advanced
918: .seealso: `DMCreateFieldDecomposition()`, `DMShellSetContext()`, `DMShellGetContext()`
919: @*/
920: PetscErrorCode DMShellSetCreateFieldDecomposition(DM dm, PetscErrorCode (*decomp)(DM, PetscInt *, char ***, IS **, DM **))
921: {
922: PetscBool isshell;
925: PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell);
926: if (!isshell) return 0;
927: dm->ops->createfielddecomposition = decomp;
928: return 0;
929: }
931: /*@C
932: DMShellSetCreateDomainDecomposition - Set the routine used to create a domain decomposition for the shell DM
934: Logically Collective on dm
936: Input Parameters:
937: + dm - the shell DM
938: - decomp - the routine to create the decomposition
940: Level: advanced
942: .seealso: `DMCreateDomainDecomposition()`, `DMShellSetContext()`, `DMShellGetContext()`
943: @*/
944: PetscErrorCode DMShellSetCreateDomainDecomposition(DM dm, PetscErrorCode (*decomp)(DM, PetscInt *, char ***, IS **, IS **, DM **))
945: {
946: PetscBool isshell;
949: PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell);
950: if (!isshell) return 0;
951: dm->ops->createdomaindecomposition = decomp;
952: return 0;
953: }
955: /*@C
956: DMShellSetCreateDomainDecompositionScatters - Set the routine used to create the scatter contexts for domain decomposition with a shell DM
958: Logically Collective on dm
960: Input Parameters:
961: + dm - the shell DM
962: - scatter - the routine to create the scatters
964: Level: advanced
966: .seealso: `DMCreateDomainDecompositionScatters()`, `DMShellSetContext()`, `DMShellGetContext()`
967: @*/
968: PetscErrorCode DMShellSetCreateDomainDecompositionScatters(DM dm, PetscErrorCode (*scatter)(DM, PetscInt, DM *, VecScatter **, VecScatter **, VecScatter **))
969: {
970: PetscBool isshell;
973: PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell);
974: if (!isshell) return 0;
975: dm->ops->createddscatters = scatter;
976: return 0;
977: }
979: /*@C
980: DMShellSetCreateSubDM - Set the routine used to create a sub DM from the shell DM
982: Logically Collective on dm
984: Input Parameters:
985: + dm - the shell DM
986: - subdm - the routine to create the decomposition
988: Level: advanced
990: .seealso: `DMCreateSubDM()`, `DMShellGetCreateSubDM()`, `DMShellSetContext()`, `DMShellGetContext()`
991: @*/
992: PetscErrorCode DMShellSetCreateSubDM(DM dm, PetscErrorCode (*subdm)(DM, PetscInt, const PetscInt[], IS *, DM *))
993: {
994: PetscBool isshell;
997: PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell);
998: if (!isshell) return 0;
999: dm->ops->createsubdm = subdm;
1000: return 0;
1001: }
1003: /*@C
1004: DMShellGetCreateSubDM - Get the routine used to create a sub DM from the shell DM
1006: Logically Collective on dm
1008: Input Parameter:
1009: . dm - the shell DM
1011: Output Parameter:
1012: . subdm - the routine to create the decomposition
1014: Level: advanced
1016: .seealso: `DMCreateSubDM()`, `DMShellSetCreateSubDM()`, `DMShellSetContext()`, `DMShellGetContext()`
1017: @*/
1018: PetscErrorCode DMShellGetCreateSubDM(DM dm, PetscErrorCode (**subdm)(DM, PetscInt, const PetscInt[], IS *, DM *))
1019: {
1020: PetscBool isshell;
1023: PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell);
1025: *subdm = dm->ops->createsubdm;
1026: return 0;
1027: }
1029: static PetscErrorCode DMDestroy_Shell(DM dm)
1030: {
1031: DM_Shell *shell = (DM_Shell *)dm->data;
1033: if (shell->destroyctx) PetscCallBack("Destroy Context", (*shell->destroyctx)(shell->ctx));
1034: MatDestroy(&shell->A);
1035: VecDestroy(&shell->Xglobal);
1036: VecDestroy(&shell->Xlocal);
1037: VecScatterDestroy(&shell->gtol);
1038: VecScatterDestroy(&shell->ltog);
1039: VecScatterDestroy(&shell->ltol);
1040: /* This was originally freed in DMDestroy(), but that prevents reference counting of backend objects */
1041: PetscFree(shell);
1042: return 0;
1043: }
1045: static PetscErrorCode DMView_Shell(DM dm, PetscViewer v)
1046: {
1047: DM_Shell *shell = (DM_Shell *)dm->data;
1049: if (shell->Xglobal) VecView(shell->Xglobal, v);
1050: return 0;
1051: }
1053: static PetscErrorCode DMLoad_Shell(DM dm, PetscViewer v)
1054: {
1055: DM_Shell *shell = (DM_Shell *)dm->data;
1057: VecCreate(PetscObjectComm((PetscObject)dm), &shell->Xglobal);
1058: VecLoad(shell->Xglobal, v);
1059: return 0;
1060: }
1062: PetscErrorCode DMCreateSubDM_Shell(DM dm, PetscInt numFields, const PetscInt fields[], IS *is, DM *subdm)
1063: {
1064: if (subdm) DMShellCreate(PetscObjectComm((PetscObject)dm), subdm);
1065: DMCreateSectionSubDM(dm, numFields, fields, is, subdm);
1066: return 0;
1067: }
1069: PETSC_EXTERN PetscErrorCode DMCreate_Shell(DM dm)
1070: {
1071: DM_Shell *shell;
1073: PetscNew(&shell);
1074: dm->data = shell;
1076: dm->ops->destroy = DMDestroy_Shell;
1077: dm->ops->createglobalvector = DMCreateGlobalVector_Shell;
1078: dm->ops->createlocalvector = DMCreateLocalVector_Shell;
1079: dm->ops->creatematrix = DMCreateMatrix_Shell;
1080: dm->ops->view = DMView_Shell;
1081: dm->ops->load = DMLoad_Shell;
1082: dm->ops->globaltolocalbegin = DMGlobalToLocalBeginDefaultShell;
1083: dm->ops->globaltolocalend = DMGlobalToLocalEndDefaultShell;
1084: dm->ops->localtoglobalbegin = DMLocalToGlobalBeginDefaultShell;
1085: dm->ops->localtoglobalend = DMLocalToGlobalEndDefaultShell;
1086: dm->ops->localtolocalbegin = DMLocalToLocalBeginDefaultShell;
1087: dm->ops->localtolocalend = DMLocalToLocalEndDefaultShell;
1088: dm->ops->createsubdm = DMCreateSubDM_Shell;
1089: DMSetMatType(dm, MATDENSE);
1090: return 0;
1091: }
1093: /*@
1094: DMShellCreate - Creates a shell DM object, used to manage user-defined problem data
1096: Collective
1098: Input Parameter:
1099: . comm - the processors that will share the global vector
1101: Output Parameters:
1102: . shell - the shell DM
1104: Level: advanced
1106: .seealso `DMDestroy()`, `DMCreateGlobalVector()`, `DMCreateLocalVector()`, `DMShellSetContext()`, `DMShellGetContext()`
1107: @*/
1108: PetscErrorCode DMShellCreate(MPI_Comm comm, DM *dm)
1109: {
1111: DMCreate(comm, dm);
1112: DMSetType(*dm, DMSHELL);
1113: DMSetUp(*dm);
1114: return 0;
1115: }