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: } DM_Shell;
15: /*@
16: DMGlobalToLocalBeginDefaultShell - Uses the GlobalToLocal VecScatter context set by the user to begin a global to local scatter
17: Collective
19: Input Parameters:
20: + dm - shell DM
21: . g - global vector
22: . mode - InsertMode
23: - l - local vector
25: Level: advanced
27: 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.
29: .seealso: DMGlobalToLocalEndDefaultShell()
30: @*/
31: PetscErrorCode DMGlobalToLocalBeginDefaultShell(DM dm,Vec g,InsertMode mode,Vec l)
32: {
33: DM_Shell *shell = (DM_Shell*)dm->data;
36: VecScatterBegin(shell->gtol,g,l,mode,SCATTER_FORWARD);
37: return 0;
38: }
40: /*@
41: DMGlobalToLocalEndDefaultShell - Uses the GlobalToLocal VecScatter context set by the user to end a global to local scatter
42: Collective
44: Input Parameters:
45: + dm - shell DM
46: . g - global vector
47: . mode - InsertMode
48: - l - local vector
50: Level: advanced
52: .seealso: DMGlobalToLocalBeginDefaultShell()
53: @*/
54: PetscErrorCode DMGlobalToLocalEndDefaultShell(DM dm,Vec g,InsertMode mode,Vec l)
55: {
56: DM_Shell *shell = (DM_Shell*)dm->data;
59: VecScatterEnd(shell->gtol,g,l,mode,SCATTER_FORWARD);
60: return 0;
61: }
63: /*@
64: DMLocalToGlobalBeginDefaultShell - Uses the LocalToGlobal VecScatter context set by the user to begin a local to global scatter
65: Collective
67: Input Parameters:
68: + dm - shell DM
69: . l - local vector
70: . mode - InsertMode
71: - g - global vector
73: Level: advanced
75: 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.
77: .seealso: DMLocalToGlobalEndDefaultShell()
78: @*/
79: PetscErrorCode DMLocalToGlobalBeginDefaultShell(DM dm,Vec l,InsertMode mode,Vec g)
80: {
81: DM_Shell *shell = (DM_Shell*)dm->data;
84: VecScatterBegin(shell->ltog,l,g,mode,SCATTER_FORWARD);
85: return 0;
86: }
88: /*@
89: DMLocalToGlobalEndDefaultShell - Uses the LocalToGlobal VecScatter context set by the user to end a local to global scatter
90: Collective
92: Input Parameters:
93: + dm - shell DM
94: . l - local vector
95: . mode - InsertMode
96: - g - global vector
98: Level: advanced
100: .seealso: DMLocalToGlobalBeginDefaultShell()
101: @*/
102: PetscErrorCode DMLocalToGlobalEndDefaultShell(DM dm,Vec l,InsertMode mode,Vec g)
103: {
104: DM_Shell *shell = (DM_Shell*)dm->data;
107: VecScatterEnd(shell->ltog,l,g,mode,SCATTER_FORWARD);
108: return 0;
109: }
111: /*@
112: DMLocalToLocalBeginDefaultShell - Uses the LocalToLocal VecScatter context set by the user to begin a local to local scatter
113: Collective
115: Input Parameters:
116: + dm - shell DM
117: . g - the original local vector
118: - mode - InsertMode
120: Output Parameter:
121: . l - the local vector with correct ghost values
123: Level: advanced
125: 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.
127: .seealso: DMLocalToLocalEndDefaultShell()
128: @*/
129: PetscErrorCode DMLocalToLocalBeginDefaultShell(DM dm,Vec g,InsertMode mode,Vec l)
130: {
131: DM_Shell *shell = (DM_Shell*)dm->data;
134: VecScatterBegin(shell->ltol,g,l,mode,SCATTER_FORWARD);
135: return 0;
136: }
138: /*@
139: DMLocalToLocalEndDefaultShell - Uses the LocalToLocal VecScatter context set by the user to end a local to local scatter
140: Collective
142: Input Parameters:
143: + dm - shell DM
144: . g - the original local vector
145: - mode - InsertMode
147: Output Parameter:
148: . l - the local vector with correct ghost values
150: Level: advanced
152: .seealso: DMLocalToLocalBeginDefaultShell()
153: @*/
154: PetscErrorCode DMLocalToLocalEndDefaultShell(DM dm,Vec g,InsertMode mode,Vec l)
155: {
156: DM_Shell *shell = (DM_Shell*)dm->data;
159: VecScatterEnd(shell->ltol,g,l,mode,SCATTER_FORWARD);
160: return 0;
161: }
163: static PetscErrorCode DMCreateMatrix_Shell(DM dm,Mat *J)
164: {
165: DM_Shell *shell = (DM_Shell*)dm->data;
166: Mat A;
170: if (!shell->A) {
171: if (shell->Xglobal) {
172: PetscInt m,M;
173: PetscInfo(dm,"Naively creating matrix using global vector distribution without preallocation\n");
174: VecGetSize(shell->Xglobal,&M);
175: VecGetLocalSize(shell->Xglobal,&m);
176: MatCreate(PetscObjectComm((PetscObject)dm),&shell->A);
177: MatSetSizes(shell->A,m,m,M,M);
178: MatSetType(shell->A,dm->mattype);
179: MatSetUp(shell->A);
180: } else SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_USER,"Must call DMShellSetMatrix(), DMShellSetCreateMatrix(), or provide a vector");
181: }
182: A = shell->A;
183: MatDuplicate(A,MAT_SHARE_NONZERO_PATTERN,J);
184: MatSetDM(*J,dm);
185: return 0;
186: }
188: PetscErrorCode DMCreateGlobalVector_Shell(DM dm,Vec *gvec)
189: {
190: DM_Shell *shell = (DM_Shell*)dm->data;
191: Vec X;
195: *gvec = NULL;
196: X = shell->Xglobal;
198: /* Need to create a copy in order to attach the DM to the vector */
199: VecDuplicate(X,gvec);
200: VecZeroEntries(*gvec);
201: VecSetDM(*gvec,dm);
202: return 0;
203: }
205: PetscErrorCode DMCreateLocalVector_Shell(DM dm,Vec *gvec)
206: {
207: DM_Shell *shell = (DM_Shell*)dm->data;
208: Vec X;
212: *gvec = NULL;
213: X = shell->Xlocal;
215: /* Need to create a copy in order to attach the DM to the vector */
216: VecDuplicate(X,gvec);
217: VecZeroEntries(*gvec);
218: VecSetDM(*gvec,dm);
219: return 0;
220: }
222: /*@
223: DMShellSetContext - set some data to be usable by this DM
225: Collective
227: Input Parameters:
228: + dm - shell DM
229: - ctx - the context
231: Level: advanced
233: .seealso: DMCreateMatrix(), DMShellGetContext()
234: @*/
235: PetscErrorCode DMShellSetContext(DM dm,void *ctx)
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->ctx = ctx;
244: return 0;
245: }
247: /*@
248: DMShellGetContext - Returns the user-provided context associated to the DM
250: Collective
252: Input Parameter:
253: . dm - shell DM
255: Output Parameter:
256: . ctx - the context
258: Level: advanced
260: .seealso: DMCreateMatrix(), DMShellSetContext()
261: @*/
262: PetscErrorCode DMShellGetContext(DM dm,void *ctx)
263: {
264: DM_Shell *shell = (DM_Shell*)dm->data;
265: PetscBool isshell;
268: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
270: *(void**)ctx = shell->ctx;
271: return 0;
272: }
274: /*@
275: DMShellSetMatrix - sets a template matrix associated with the DMShell
277: Collective
279: Input Parameters:
280: + dm - shell DM
281: - J - template matrix
283: Level: advanced
285: Developer Notes:
286: 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.
288: .seealso: DMCreateMatrix(), DMShellSetCreateMatrix(), DMShellSetContext(), DMShellGetContext()
289: @*/
290: PetscErrorCode DMShellSetMatrix(DM dm,Mat J)
291: {
292: DM_Shell *shell = (DM_Shell*)dm->data;
293: PetscBool isshell;
294: DM mdm;
298: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
299: if (!isshell) return 0;
300: if (J == shell->A) return 0;
301: MatGetDM(J,&mdm);
302: PetscObjectReference((PetscObject)J);
303: MatDestroy(&shell->A);
304: if (mdm == dm) {
305: MatDuplicate(J,MAT_SHARE_NONZERO_PATTERN,&shell->A);
306: MatSetDM(shell->A,NULL);
307: } else shell->A = J;
308: return 0;
309: }
311: /*@C
312: DMShellSetCreateMatrix - sets the routine to create a matrix associated with the shell DM
314: Logically Collective on dm
316: Input Parameters:
317: + dm - the shell DM
318: - func - the function to create a matrix
320: Level: advanced
322: .seealso: DMCreateMatrix(), DMShellSetMatrix(), DMShellSetContext(), DMShellGetContext()
323: @*/
324: PetscErrorCode DMShellSetCreateMatrix(DM dm,PetscErrorCode (*func)(DM,Mat*))
325: {
327: dm->ops->creatematrix = func;
328: return 0;
329: }
331: /*@
332: DMShellSetGlobalVector - sets a template global vector associated with the DMShell
334: Logically Collective on dm
336: Input Parameters:
337: + dm - shell DM
338: - X - template vector
340: Level: advanced
342: .seealso: DMCreateGlobalVector(), DMShellSetMatrix(), DMShellSetCreateGlobalVector()
343: @*/
344: PetscErrorCode DMShellSetGlobalVector(DM dm,Vec X)
345: {
346: DM_Shell *shell = (DM_Shell*)dm->data;
347: PetscBool isshell;
348: DM vdm;
352: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
353: if (!isshell) return 0;
354: VecGetDM(X,&vdm);
355: /*
356: if the vector proposed as the new base global vector for the DM is a DM vector associated
357: with the same DM then the current base global vector for the DM is ok and if we replace it with the new one
358: we get a circular dependency that prevents the DM from being destroy when it should be.
359: This occurs when SNESSet/GetNPC() is used with a SNES that does not have a user provided
360: DM attached to it since the inner SNES (which shares the DM with the outer SNES) tries
361: to set its input vector (which is associated with the DM) as the base global vector.
362: Thanks to Juan P. Mendez Granado Re: [petsc-maint] Nonlinear conjugate gradien
363: for pointing out the problem.
364: */
365: if (vdm == dm) return 0;
366: PetscObjectReference((PetscObject)X);
367: VecDestroy(&shell->Xglobal);
368: shell->Xglobal = X;
369: return 0;
370: }
372: /*@
373: DMShellGetGlobalVector - Returns the template global vector associated with the DMShell, or NULL if it was not set
375: Not collective
377: Input Parameters:
378: + dm - shell DM
379: - X - template vector
381: Level: advanced
383: .seealso: DMShellSetGlobalVector(), DMShellSetCreateGlobalVector(), DMCreateGlobalVector()
384: @*/
385: PetscErrorCode DMShellGetGlobalVector(DM dm, Vec *X)
386: {
387: DM_Shell *shell = (DM_Shell *) dm->data;
388: PetscBool isshell;
392: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
393: if (!isshell) return 0;
394: *X = shell->Xglobal;
395: return 0;
396: }
398: /*@C
399: DMShellSetCreateGlobalVector - sets the routine to create a global vector associated with the shell DM
401: Logically Collective
403: Input Parameters:
404: + dm - the shell DM
405: - func - the creation routine
407: Level: advanced
409: .seealso: DMShellSetGlobalVector(), DMShellSetCreateMatrix(), DMShellSetContext(), DMShellGetContext()
410: @*/
411: PetscErrorCode DMShellSetCreateGlobalVector(DM dm,PetscErrorCode (*func)(DM,Vec*))
412: {
414: dm->ops->createglobalvector = func;
415: return 0;
416: }
418: /*@
419: DMShellSetLocalVector - sets a template local vector associated with the DMShell
421: Logically Collective on dm
423: Input Parameters:
424: + dm - shell DM
425: - X - template vector
427: Level: advanced
429: .seealso: DMCreateLocalVector(), DMShellSetMatrix(), DMShellSetCreateLocalVector()
430: @*/
431: PetscErrorCode DMShellSetLocalVector(DM dm,Vec X)
432: {
433: DM_Shell *shell = (DM_Shell*)dm->data;
434: PetscBool isshell;
435: DM vdm;
439: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
440: if (!isshell) return 0;
441: VecGetDM(X,&vdm);
442: /*
443: if the vector proposed as the new base global vector for the DM is a DM vector associated
444: with the same DM then the current base global vector for the DM is ok and if we replace it with the new one
445: we get a circular dependency that prevents the DM from being destroy when it should be.
446: This occurs when SNESSet/GetNPC() is used with a SNES that does not have a user provided
447: DM attached to it since the inner SNES (which shares the DM with the outer SNES) tries
448: to set its input vector (which is associated with the DM) as the base global vector.
449: Thanks to Juan P. Mendez Granado Re: [petsc-maint] Nonlinear conjugate gradien
450: for pointing out the problem.
451: */
452: if (vdm == dm) return 0;
453: PetscObjectReference((PetscObject)X);
454: VecDestroy(&shell->Xlocal);
455: shell->Xlocal = X;
456: return 0;
457: }
459: /*@C
460: DMShellSetCreateLocalVector - sets the routine to create a local vector associated with the shell DM
462: Logically Collective
464: Input Parameters:
465: + dm - the shell DM
466: - func - the creation routine
468: Level: advanced
470: .seealso: DMShellSetLocalVector(), DMShellSetCreateMatrix(), DMShellSetContext(), DMShellGetContext()
471: @*/
472: PetscErrorCode DMShellSetCreateLocalVector(DM dm,PetscErrorCode (*func)(DM,Vec*))
473: {
475: dm->ops->createlocalvector = func;
476: return 0;
477: }
479: /*@C
480: DMShellSetGlobalToLocal - Sets the routines used to perform a global to local scatter
482: Logically Collective on dm
484: Input Parameters
485: + dm - the shell DM
486: . begin - the routine that begins the global to local scatter
487: - end - the routine that ends the global to local scatter
489: Notes:
490: If these functions are not provided but DMShellSetGlobalToLocalVecScatter() is called then
491: DMGlobalToLocalBeginDefaultShell()/DMGlobalToLocalEndDefaultShell() are used to to perform the transfers
493: Level: advanced
495: .seealso: DMShellSetLocalToGlobal(), DMGlobalToLocalBeginDefaultShell(), DMGlobalToLocalEndDefaultShell()
496: @*/
497: PetscErrorCode DMShellSetGlobalToLocal(DM dm,PetscErrorCode (*begin)(DM,Vec,InsertMode,Vec),PetscErrorCode (*end)(DM,Vec,InsertMode,Vec))
498: {
500: dm->ops->globaltolocalbegin = begin;
501: dm->ops->globaltolocalend = end;
502: return 0;
503: }
505: /*@C
506: DMShellSetLocalToGlobal - Sets the routines used to perform a local to global scatter
508: Logically Collective on dm
510: Input Parameters
511: + dm - the shell DM
512: . begin - the routine that begins the local to global scatter
513: - end - the routine that ends the local to global scatter
515: Notes:
516: If these functions are not provided but DMShellSetLocalToGlobalVecScatter() is called then
517: DMLocalToGlobalBeginDefaultShell()/DMLocalToGlobalEndDefaultShell() are used to to perform the transfers
519: Level: advanced
521: .seealso: DMShellSetGlobalToLocal()
522: @*/
523: PetscErrorCode DMShellSetLocalToGlobal(DM dm,PetscErrorCode (*begin)(DM,Vec,InsertMode,Vec),PetscErrorCode (*end)(DM,Vec,InsertMode,Vec))
524: {
526: dm->ops->localtoglobalbegin = begin;
527: dm->ops->localtoglobalend = end;
528: return 0;
529: }
531: /*@C
532: DMShellSetLocalToLocal - Sets the routines used to perform a local to local scatter
534: Logically Collective on dm
536: Input Parameters
537: + dm - the shell DM
538: . begin - the routine that begins the local to local scatter
539: - end - the routine that ends the local to local scatter
541: Notes:
542: If these functions are not provided but DMShellSetLocalToLocalVecScatter() is called then
543: DMLocalToLocalBeginDefaultShell()/DMLocalToLocalEndDefaultShell() are used to to perform the transfers
545: Level: advanced
547: .seealso: DMShellSetGlobalToLocal(), DMLocalToLocalBeginDefaultShell(), DMLocalToLocalEndDefaultShell()
548: @*/
549: PetscErrorCode DMShellSetLocalToLocal(DM dm,PetscErrorCode (*begin)(DM,Vec,InsertMode,Vec),PetscErrorCode (*end)(DM,Vec,InsertMode,Vec))
550: {
552: dm->ops->localtolocalbegin = begin;
553: dm->ops->localtolocalend = end;
554: return 0;
555: }
557: /*@
558: DMShellSetGlobalToLocalVecScatter - Sets a VecScatter context for global to local communication
560: Logically Collective on dm
562: Input Parameters
563: + dm - the shell DM
564: - gtol - the global to local VecScatter context
566: Level: advanced
568: .seealso: DMShellSetGlobalToLocal(), DMGlobalToLocalBeginDefaultShell(), DMGlobalToLocalEndDefaultShell()
569: @*/
570: PetscErrorCode DMShellSetGlobalToLocalVecScatter(DM dm, VecScatter gtol)
571: {
572: DM_Shell *shell = (DM_Shell*)dm->data;
576: PetscObjectReference((PetscObject)gtol);
577: VecScatterDestroy(&shell->gtol);
578: shell->gtol = gtol;
579: return 0;
580: }
582: /*@
583: DMShellSetLocalToGlobalVecScatter - Sets a VecScatter context for local to global communication
585: Logically Collective on dm
587: Input Parameters
588: + dm - the shell DM
589: - ltog - the local to global VecScatter context
591: Level: advanced
593: .seealso: DMShellSetLocalToGlobal(), DMLocalToGlobalBeginDefaultShell(), DMLocalToGlobalEndDefaultShell()
594: @*/
595: PetscErrorCode DMShellSetLocalToGlobalVecScatter(DM dm, VecScatter ltog)
596: {
597: DM_Shell *shell = (DM_Shell*)dm->data;
601: PetscObjectReference((PetscObject)ltog);
602: VecScatterDestroy(&shell->ltog);
603: shell->ltog = ltog;
604: return 0;
605: }
607: /*@
608: DMShellSetLocalToLocalVecScatter - Sets a VecScatter context for local to local communication
610: Logically Collective on dm
612: Input Parameters
613: + dm - the shell DM
614: - ltol - the local to local VecScatter context
616: Level: advanced
618: .seealso: DMShellSetLocalToLocal(), DMLocalToLocalBeginDefaultShell(), DMLocalToLocalEndDefaultShell()
619: @*/
620: PetscErrorCode DMShellSetLocalToLocalVecScatter(DM dm, VecScatter ltol)
621: {
622: DM_Shell *shell = (DM_Shell*)dm->data;
626: PetscObjectReference((PetscObject)ltol);
627: VecScatterDestroy(&shell->ltol);
628: shell->ltol = ltol;
629: return 0;
630: }
632: /*@C
633: DMShellSetCoarsen - Set the routine used to coarsen the shell DM
635: Logically Collective on dm
637: Input Parameters
638: + dm - the shell DM
639: - coarsen - the routine that coarsens the DM
641: Level: advanced
643: .seealso: DMShellSetRefine(), DMCoarsen(), DMShellGetCoarsen(), DMShellSetContext(), DMShellGetContext()
644: @*/
645: PetscErrorCode DMShellSetCoarsen(DM dm, PetscErrorCode (*coarsen)(DM,MPI_Comm,DM*))
646: {
647: PetscBool isshell;
650: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
651: if (!isshell) return 0;
652: dm->ops->coarsen = coarsen;
653: return 0;
654: }
656: /*@C
657: DMShellGetCoarsen - Get the routine used to coarsen the shell DM
659: Logically Collective on dm
661: Input Parameter:
662: . dm - the shell DM
664: Output Parameter:
665: . coarsen - the routine that coarsens the DM
667: Level: advanced
669: .seealso: DMShellSetCoarsen(), DMCoarsen(), DMShellSetRefine(), DMRefine()
670: @*/
671: PetscErrorCode DMShellGetCoarsen(DM dm, PetscErrorCode (**coarsen)(DM,MPI_Comm,DM*))
672: {
673: PetscBool isshell;
676: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
678: *coarsen = dm->ops->coarsen;
679: return 0;
680: }
682: /*@C
683: DMShellSetRefine - Set the routine used to refine the shell DM
685: Logically Collective on dm
687: Input Parameters
688: + dm - the shell DM
689: - refine - the routine that refines the DM
691: Level: advanced
693: .seealso: DMShellSetCoarsen(), DMRefine(), DMShellGetRefine(), DMShellSetContext(), DMShellGetContext()
694: @*/
695: PetscErrorCode DMShellSetRefine(DM dm, PetscErrorCode (*refine)(DM,MPI_Comm,DM*))
696: {
697: PetscBool isshell;
700: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
701: if (!isshell) return 0;
702: dm->ops->refine = refine;
703: return 0;
704: }
706: /*@C
707: DMShellGetRefine - Get the routine used to refine the shell DM
709: Logically Collective on dm
711: Input Parameter:
712: . dm - the shell DM
714: Output Parameter:
715: . refine - the routine that refines the DM
717: Level: advanced
719: .seealso: DMShellSetCoarsen(), DMCoarsen(), DMShellSetRefine(), DMRefine()
720: @*/
721: PetscErrorCode DMShellGetRefine(DM dm, PetscErrorCode (**refine)(DM,MPI_Comm,DM*))
722: {
723: PetscBool isshell;
726: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
728: *refine = dm->ops->refine;
729: return 0;
730: }
732: /*@C
733: DMShellSetCreateInterpolation - Set the routine used to create the interpolation operator
735: Logically Collective on dm
737: Input Parameters
738: + dm - the shell DM
739: - interp - the routine to create the interpolation
741: Level: advanced
743: .seealso: DMShellSetCreateInjection(), DMCreateInterpolation(), DMShellGetCreateInterpolation(), DMShellSetCreateRestriction(), DMShellSetContext(), DMShellGetContext()
744: @*/
745: PetscErrorCode DMShellSetCreateInterpolation(DM dm, PetscErrorCode (*interp)(DM,DM,Mat*,Vec*))
746: {
747: PetscBool isshell;
750: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
751: if (!isshell) return 0;
752: dm->ops->createinterpolation = interp;
753: return 0;
754: }
756: /*@C
757: DMShellGetCreateInterpolation - Get the routine used to create the interpolation operator
759: Logically Collective on dm
761: Input Parameter:
762: . dm - the shell DM
764: Output Parameter:
765: . interp - the routine to create the interpolation
767: Level: advanced
769: .seealso: DMShellGetCreateInjection(), DMCreateInterpolation(), DMShellGetCreateRestriction(), DMShellSetContext(), DMShellGetContext()
770: @*/
771: PetscErrorCode DMShellGetCreateInterpolation(DM dm, PetscErrorCode (**interp)(DM,DM,Mat*,Vec*))
772: {
773: PetscBool isshell;
776: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
778: *interp = dm->ops->createinterpolation;
779: return 0;
780: }
782: /*@C
783: DMShellSetCreateRestriction - Set the routine used to create the restriction operator
785: Logically Collective on dm
787: Input Parameters
788: + dm - the shell DM
789: - striction- the routine to create the restriction
791: Level: advanced
793: .seealso: DMShellSetCreateInjection(), DMCreateInterpolation(), DMShellGetCreateRestriction(), DMShellSetContext(), DMShellGetContext()
794: @*/
795: PetscErrorCode DMShellSetCreateRestriction(DM dm, PetscErrorCode (*restriction)(DM,DM,Mat*))
796: {
797: PetscBool isshell;
800: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
801: if (!isshell) return 0;
802: dm->ops->createrestriction = restriction;
803: return 0;
804: }
806: /*@C
807: DMShellGetCreateRestriction - Get the routine used to create the restriction operator
809: Logically Collective on dm
811: Input Parameter:
812: . dm - the shell DM
814: Output Parameter:
815: . restriction - the routine to create the restriction
817: Level: advanced
819: .seealso: DMShellSetCreateInjection(), DMCreateInterpolation(), DMShellSetContext(), DMShellGetContext()
820: @*/
821: PetscErrorCode DMShellGetCreateRestriction(DM dm, PetscErrorCode (**restriction)(DM,DM,Mat*))
822: {
823: PetscBool isshell;
826: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
828: *restriction = dm->ops->createrestriction;
829: return 0;
830: }
832: /*@C
833: DMShellSetCreateInjection - Set the routine used to create the injection operator
835: Logically Collective on dm
837: Input Parameters:
838: + dm - the shell DM
839: - inject - the routine to create the injection
841: Level: advanced
843: .seealso: DMShellSetCreateInterpolation(), DMCreateInjection(), DMShellGetCreateInjection(), DMShellSetContext(), DMShellGetContext()
844: @*/
845: PetscErrorCode DMShellSetCreateInjection(DM dm, PetscErrorCode (*inject)(DM,DM,Mat*))
846: {
847: PetscBool isshell;
850: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
851: if (!isshell) return 0;
852: dm->ops->createinjection = inject;
853: return 0;
854: }
856: /*@C
857: DMShellGetCreateInjection - Get the routine used to create the injection operator
859: Logically Collective on dm
861: Input Parameter:
862: . dm - the shell DM
864: Output Parameter:
865: . inject - the routine to create the injection
867: Level: advanced
869: .seealso: DMShellGetCreateInterpolation(), DMCreateInjection(), DMShellSetContext(), DMShellGetContext()
870: @*/
871: PetscErrorCode DMShellGetCreateInjection(DM dm, PetscErrorCode (**inject)(DM,DM,Mat*))
872: {
873: PetscBool isshell;
876: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
878: *inject = dm->ops->createinjection;
879: return 0;
880: }
882: /*@C
883: DMShellSetCreateFieldDecomposition - Set the routine used to create a decomposition of fields for the shell DM
885: Logically Collective on dm
887: Input Parameters:
888: + dm - the shell DM
889: - decomp - the routine to create the decomposition
891: Level: advanced
893: .seealso: DMCreateFieldDecomposition(), DMShellSetContext(), DMShellGetContext()
894: @*/
895: PetscErrorCode DMShellSetCreateFieldDecomposition(DM dm, PetscErrorCode (*decomp)(DM,PetscInt*,char***, IS**,DM**))
896: {
897: PetscBool isshell;
900: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
901: if (!isshell) return 0;
902: dm->ops->createfielddecomposition = decomp;
903: return 0;
904: }
906: /*@C
907: DMShellSetCreateDomainDecomposition - Set the routine used to create a domain decomposition for the shell DM
909: Logically Collective on dm
911: Input Parameters:
912: + dm - the shell DM
913: - decomp - the routine to create the decomposition
915: Level: advanced
917: .seealso: DMCreateDomainDecomposition(), DMShellSetContext(), DMShellGetContext()
918: @*/
919: PetscErrorCode DMShellSetCreateDomainDecomposition(DM dm, PetscErrorCode (*decomp)(DM,PetscInt*,char***, IS**,IS**,DM**))
920: {
921: PetscBool isshell;
924: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
925: if (!isshell) return 0;
926: dm->ops->createdomaindecomposition = decomp;
927: return 0;
928: }
930: /*@C
931: DMShellSetCreateDomainDecompositionScatters - Set the routine used to create the scatter contexts for domain decomposition with a shell DM
933: Logically Collective on dm
935: Input Parameters:
936: + dm - the shell DM
937: - scatter - the routine to create the scatters
939: Level: advanced
941: .seealso: DMCreateDomainDecompositionScatters(), DMShellSetContext(), DMShellGetContext()
942: @*/
943: PetscErrorCode DMShellSetCreateDomainDecompositionScatters(DM dm, PetscErrorCode (*scatter)(DM,PetscInt,DM*,VecScatter**,VecScatter**,VecScatter**))
944: {
945: PetscBool isshell;
948: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
949: if (!isshell) return 0;
950: dm->ops->createddscatters = scatter;
951: return 0;
952: }
954: /*@C
955: DMShellSetCreateSubDM - Set the routine used to create a sub DM from the shell DM
957: Logically Collective on dm
959: Input Parameters:
960: + dm - the shell DM
961: - subdm - the routine to create the decomposition
963: Level: advanced
965: .seealso: DMCreateSubDM(), DMShellGetCreateSubDM(), DMShellSetContext(), DMShellGetContext()
966: @*/
967: PetscErrorCode DMShellSetCreateSubDM(DM dm, PetscErrorCode (*subdm)(DM,PetscInt,const PetscInt[],IS*,DM*))
968: {
969: PetscBool isshell;
972: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
973: if (!isshell) return 0;
974: dm->ops->createsubdm = subdm;
975: return 0;
976: }
978: /*@C
979: DMShellGetCreateSubDM - Get the routine used to create a sub DM from the shell DM
981: Logically Collective on dm
983: Input Parameter:
984: . dm - the shell DM
986: Output Parameter:
987: . subdm - the routine to create the decomposition
989: Level: advanced
991: .seealso: DMCreateSubDM(), DMShellSetCreateSubDM(), DMShellSetContext(), DMShellGetContext()
992: @*/
993: PetscErrorCode DMShellGetCreateSubDM(DM dm, PetscErrorCode (**subdm)(DM,PetscInt,const PetscInt[],IS*,DM*))
994: {
995: PetscBool isshell;
998: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
1000: *subdm = dm->ops->createsubdm;
1001: return 0;
1002: }
1004: static PetscErrorCode DMDestroy_Shell(DM dm)
1005: {
1006: DM_Shell *shell = (DM_Shell*)dm->data;
1008: MatDestroy(&shell->A);
1009: VecDestroy(&shell->Xglobal);
1010: VecDestroy(&shell->Xlocal);
1011: VecScatterDestroy(&shell->gtol);
1012: VecScatterDestroy(&shell->ltog);
1013: VecScatterDestroy(&shell->ltol);
1014: /* This was originally freed in DMDestroy(), but that prevents reference counting of backend objects */
1015: PetscFree(shell);
1016: return 0;
1017: }
1019: static PetscErrorCode DMView_Shell(DM dm,PetscViewer v)
1020: {
1021: DM_Shell *shell = (DM_Shell*)dm->data;
1023: VecView(shell->Xglobal,v);
1024: return 0;
1025: }
1027: static PetscErrorCode DMLoad_Shell(DM dm,PetscViewer v)
1028: {
1029: DM_Shell *shell = (DM_Shell*)dm->data;
1031: VecCreate(PetscObjectComm((PetscObject)dm),&shell->Xglobal);
1032: VecLoad(shell->Xglobal,v);
1033: return 0;
1034: }
1036: PetscErrorCode DMCreateSubDM_Shell(DM dm, PetscInt numFields, const PetscInt fields[], IS *is, DM *subdm)
1037: {
1038: if (subdm) DMShellCreate(PetscObjectComm((PetscObject) dm), subdm);
1039: DMCreateSectionSubDM(dm, numFields, fields, is, subdm);
1040: return 0;
1041: }
1043: PETSC_EXTERN PetscErrorCode DMCreate_Shell(DM dm)
1044: {
1045: DM_Shell *shell;
1047: PetscNewLog(dm,&shell);
1048: dm->data = shell;
1050: dm->ops->destroy = DMDestroy_Shell;
1051: dm->ops->createglobalvector = DMCreateGlobalVector_Shell;
1052: dm->ops->createlocalvector = DMCreateLocalVector_Shell;
1053: dm->ops->creatematrix = DMCreateMatrix_Shell;
1054: dm->ops->view = DMView_Shell;
1055: dm->ops->load = DMLoad_Shell;
1056: dm->ops->globaltolocalbegin = DMGlobalToLocalBeginDefaultShell;
1057: dm->ops->globaltolocalend = DMGlobalToLocalEndDefaultShell;
1058: dm->ops->localtoglobalbegin = DMLocalToGlobalBeginDefaultShell;
1059: dm->ops->localtoglobalend = DMLocalToGlobalEndDefaultShell;
1060: dm->ops->localtolocalbegin = DMLocalToLocalBeginDefaultShell;
1061: dm->ops->localtolocalend = DMLocalToLocalEndDefaultShell;
1062: dm->ops->createsubdm = DMCreateSubDM_Shell;
1063: DMSetMatType(dm,MATDENSE);
1064: return 0;
1065: }
1067: /*@
1068: DMShellCreate - Creates a shell DM object, used to manage user-defined problem data
1070: Collective
1072: Input Parameter:
1073: . comm - the processors that will share the global vector
1075: Output Parameters:
1076: . shell - the shell DM
1078: Level: advanced
1080: .seealso DMDestroy(), DMCreateGlobalVector(), DMCreateLocalVector(), DMShellSetContext(), DMShellGetContext()
1081: @*/
1082: PetscErrorCode DMShellCreate(MPI_Comm comm,DM *dm)
1083: {
1085: DMCreate(comm,dm);
1086: DMSetType(*dm,DMSHELL);
1087: DMSetUp(*dm);
1088: return 0;
1089: }