Actual source code: dmshell.c
petsc-3.12.5 2020-03-29
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 Arguments:
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: {
34: DM_Shell *shell = (DM_Shell*)dm->data;
37: if (!shell->gtol) SETERRQ(((PetscObject)dm)->comm,PETSC_ERR_ARG_WRONGSTATE, "Cannot be used without first setting the scatter context via DMShellSetGlobalToLocalVecScatter()");
38: VecScatterBegin(shell->gtol,g,l,mode,SCATTER_FORWARD);
39: return(0);
40: }
42: /*@
43: DMGlobalToLocalEndDefaultShell - Uses the GlobalToLocal VecScatter context set by the user to end a global to local scatter
44: Collective
46: Input Arguments:
47: + dm - shell DM
48: . g - global vector
49: . mode - InsertMode
50: - l - local vector
52: Level: advanced
54: .seealso: DMGlobalToLocalBeginDefaultShell()
55: @*/
56: PetscErrorCode DMGlobalToLocalEndDefaultShell(DM dm,Vec g,InsertMode mode,Vec l)
57: {
59: DM_Shell *shell = (DM_Shell*)dm->data;
62: if (!shell->gtol) SETERRQ(((PetscObject)dm)->comm,PETSC_ERR_ARG_WRONGSTATE, "Cannot be used without first setting the scatter context via DMShellSetGlobalToLocalVecScatter()");
63: VecScatterEnd(shell->gtol,g,l,mode,SCATTER_FORWARD);
64: return(0);
65: }
67: /*@
68: DMLocalToGlobalBeginDefaultShell - Uses the LocalToGlobal VecScatter context set by the user to begin a local to global scatter
69: Collective
71: Input Arguments:
72: + dm - shell DM
73: . l - local vector
74: . mode - InsertMode
75: - g - global vector
77: Level: advanced
79: 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.
81: .seealso: DMLocalToGlobalEndDefaultShell()
82: @*/
83: PetscErrorCode DMLocalToGlobalBeginDefaultShell(DM dm,Vec l,InsertMode mode,Vec g)
84: {
86: DM_Shell *shell = (DM_Shell*)dm->data;
89: if (!shell->ltog) SETERRQ(((PetscObject)dm)->comm,PETSC_ERR_ARG_WRONGSTATE, "Cannot be used without first setting the scatter context via DMShellSetLocalToGlobalVecScatter()");
90: VecScatterBegin(shell->ltog,l,g,mode,SCATTER_FORWARD);
91: return(0);
92: }
94: /*@
95: DMLocalToGlobalEndDefaultShell - Uses the LocalToGlobal VecScatter context set by the user to end a local to global scatter
96: Collective
98: Input Arguments:
99: + dm - shell DM
100: . l - local vector
101: . mode - InsertMode
102: - g - global vector
104: Level: advanced
106: .seealso: DMLocalToGlobalBeginDefaultShell()
107: @*/
108: PetscErrorCode DMLocalToGlobalEndDefaultShell(DM dm,Vec l,InsertMode mode,Vec g)
109: {
111: DM_Shell *shell = (DM_Shell*)dm->data;
114: if (!shell->ltog) SETERRQ(((PetscObject)dm)->comm,PETSC_ERR_ARG_WRONGSTATE, "Cannot be used without first setting the scatter context via DMShellSetLocalToGlobalVecScatter()");
115: VecScatterEnd(shell->ltog,l,g,mode,SCATTER_FORWARD);
116: return(0);
117: }
119: /*@
120: DMLocalToLocalBeginDefaultShell - Uses the LocalToLocal VecScatter context set by the user to begin a local to local scatter
121: Collective
123: Input Arguments:
124: + dm - shell DM
125: . g - the original local vector
126: - mode - InsertMode
128: Output Parameter:
129: . l - the local vector with correct ghost values
131: Level: advanced
133: 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.
135: .seealso: DMLocalToLocalEndDefaultShell()
136: @*/
137: PetscErrorCode DMLocalToLocalBeginDefaultShell(DM dm,Vec g,InsertMode mode,Vec l)
138: {
140: DM_Shell *shell = (DM_Shell*)dm->data;
143: if (!shell->ltol) SETERRQ(((PetscObject)dm)->comm,PETSC_ERR_ARG_WRONGSTATE, "Cannot be used without first setting the scatter context via DMShellSetLocalToLocalVecScatter()");
144: VecScatterBegin(shell->ltol,g,l,mode,SCATTER_FORWARD);
145: return(0);
146: }
148: /*@
149: DMLocalToLocalEndDefaultShell - Uses the LocalToLocal VecScatter context set by the user to end a local to local scatter
150: Collective
152: Input Arguments:
153: + dm - shell DM
154: . g - the original local vector
155: - mode - InsertMode
157: Output Parameter:
158: . l - the local vector with correct ghost values
160: Level: advanced
162: .seealso: DMLocalToLocalBeginDefaultShell()
163: @*/
164: PetscErrorCode DMLocalToLocalEndDefaultShell(DM dm,Vec g,InsertMode mode,Vec l)
165: {
167: DM_Shell *shell = (DM_Shell*)dm->data;
170: if (!shell->ltol) SETERRQ(((PetscObject)dm)->comm,PETSC_ERR_ARG_WRONGSTATE, "Cannot be used without first setting the scatter context via DMShellSetGlobalToLocalVecScatter()");
171: VecScatterEnd(shell->ltol,g,l,mode,SCATTER_FORWARD);
172: return(0);
173: }
175: static PetscErrorCode DMCreateMatrix_Shell(DM dm,Mat *J)
176: {
178: DM_Shell *shell = (DM_Shell*)dm->data;
179: Mat A;
184: if (!shell->A) {
185: if (shell->Xglobal) {
186: PetscInt m,M;
187: PetscInfo(dm,"Naively creating matrix using global vector distribution without preallocation\n");
188: VecGetSize(shell->Xglobal,&M);
189: VecGetLocalSize(shell->Xglobal,&m);
190: MatCreate(PetscObjectComm((PetscObject)dm),&shell->A);
191: MatSetSizes(shell->A,m,m,M,M);
192: MatSetType(shell->A,dm->mattype);
193: MatSetUp(shell->A);
194: } else SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_USER,"Must call DMShellSetMatrix(), DMShellSetCreateMatrix(), or provide a vector");
195: }
196: A = shell->A;
197: /* the check below is tacky and incomplete */
198: if (dm->mattype) {
199: PetscBool flg,aij,seqaij,mpiaij;
200: PetscObjectTypeCompare((PetscObject)A,dm->mattype,&flg);
201: PetscObjectTypeCompare((PetscObject)A,MATSEQAIJ,&seqaij);
202: PetscObjectTypeCompare((PetscObject)A,MATMPIAIJ,&mpiaij);
203: PetscStrcmp(dm->mattype,MATAIJ,&aij);
204: if (!flg) {
205: if (!(aij && (seqaij || mpiaij))) SETERRQ2(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_NOTSAMETYPE,"Requested matrix of type %s, but only %s available",dm->mattype,((PetscObject)A)->type_name);
206: }
207: }
208: if (((PetscObject)A)->refct < 2) { /* We have an exclusive reference so we can give it out */
209: PetscBool f;
211: PetscObjectReference((PetscObject)A);
212: /* MATSHELL does not implement MATOP_ZERO_ENTRIES */
213: MatHasOperation(A,MATOP_ZERO_ENTRIES,&f);
214: if (f) { MatZeroEntries(A); }
215: *J = A;
216: } else { /* Need to create a copy, could use MAT_SHARE_NONZERO_PATTERN in most cases */
217: MatDuplicate(A,MAT_DO_NOT_COPY_VALUES,J);
218: }
219: return(0);
220: }
222: PetscErrorCode DMCreateGlobalVector_Shell(DM dm,Vec *gvec)
223: {
225: DM_Shell *shell = (DM_Shell*)dm->data;
226: Vec X;
231: *gvec = 0;
232: X = shell->Xglobal;
233: if (!X) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_USER,"Must call DMShellSetGlobalVector() or DMShellSetCreateGlobalVector()");
234: if (((PetscObject)X)->refct < 2) { /* We have an exclusive reference so we can give it out */
235: PetscObjectReference((PetscObject)X);
236: VecZeroEntries(X);
237: *gvec = X;
238: } else { /* Need to create a copy */
239: VecDuplicate(X,gvec);
240: VecZeroEntries(*gvec);
241: }
242: VecSetDM(*gvec,dm);
243: return(0);
244: }
246: PetscErrorCode DMCreateLocalVector_Shell(DM dm,Vec *gvec)
247: {
249: DM_Shell *shell = (DM_Shell*)dm->data;
250: Vec X;
255: *gvec = 0;
256: X = shell->Xlocal;
257: if (!X) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_USER,"Must call DMShellSetLocalVector() or DMShellSetCreateLocalVector()");
258: if (((PetscObject)X)->refct < 2) { /* We have an exclusive reference so we can give it out */
259: PetscObjectReference((PetscObject)X);
260: VecZeroEntries(X);
261: *gvec = X;
262: } else { /* Need to create a copy, could use MAT_SHARE_NONZERO_PATTERN in most cases */
263: VecDuplicate(X,gvec);
264: VecZeroEntries(*gvec);
265: }
266: VecSetDM(*gvec,dm);
267: return(0);
268: }
270: /*@
271: DMShellSetContext - set some data to be usable by this DM
273: Collective
275: Input Arguments:
276: + dm - shell DM
277: - ctx - the context
279: Level: advanced
281: .seealso: DMCreateMatrix(), DMShellGetContext()
282: @*/
283: PetscErrorCode DMShellSetContext(DM dm,void *ctx)
284: {
285: DM_Shell *shell = (DM_Shell*)dm->data;
287: PetscBool isshell;
291: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
292: if (!isshell) return(0);
293: shell->ctx = ctx;
294: return(0);
295: }
297: /*@
298: DMShellGetContext - Returns the user-provided context associated to the DM
300: Collective
302: Input Argument:
303: . dm - shell DM
305: Output Argument:
306: . ctx - the context
308: Level: advanced
310: .seealso: DMCreateMatrix(), DMShellSetContext()
311: @*/
312: PetscErrorCode DMShellGetContext(DM dm,void **ctx)
313: {
314: DM_Shell *shell = (DM_Shell*)dm->data;
316: PetscBool isshell;
320: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
321: if (!isshell) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Can only use with DMSHELL type DMs");
322: *ctx = shell->ctx;
323: return(0);
324: }
326: /*@
327: DMShellSetMatrix - sets a template matrix associated with the DMShell
329: Collective
331: Input Arguments:
332: + dm - shell DM
333: - J - template matrix
335: Level: advanced
337: .seealso: DMCreateMatrix(), DMShellSetCreateMatrix(), DMShellSetContext(), DMShellGetContext()
338: @*/
339: PetscErrorCode DMShellSetMatrix(DM dm,Mat J)
340: {
341: DM_Shell *shell = (DM_Shell*)dm->data;
343: PetscBool isshell;
348: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
349: if (!isshell) return(0);
350: PetscObjectReference((PetscObject)J);
351: MatDestroy(&shell->A);
352: shell->A = J;
353: return(0);
354: }
356: /*@C
357: DMShellSetCreateMatrix - sets the routine to create a matrix associated with the shell DM
359: Logically Collective on dm
361: Input Arguments:
362: + dm - the shell DM
363: - func - the function to create a matrix
365: Level: advanced
367: .seealso: DMCreateMatrix(), DMShellSetMatrix(), DMShellSetContext(), DMShellGetContext()
368: @*/
369: PetscErrorCode DMShellSetCreateMatrix(DM dm,PetscErrorCode (*func)(DM,Mat*))
370: {
373: dm->ops->creatematrix = func;
374: return(0);
375: }
377: /*@
378: DMShellSetGlobalVector - sets a template global vector associated with the DMShell
380: Logically Collective on dm
382: Input Arguments:
383: + dm - shell DM
384: - X - template vector
386: Level: advanced
388: .seealso: DMCreateGlobalVector(), DMShellSetMatrix(), DMShellSetCreateGlobalVector()
389: @*/
390: PetscErrorCode DMShellSetGlobalVector(DM dm,Vec X)
391: {
392: DM_Shell *shell = (DM_Shell*)dm->data;
394: PetscBool isshell;
395: DM vdm;
400: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
401: if (!isshell) return(0);
402: VecGetDM(X,&vdm);
403: /*
404: if the vector proposed as the new base global vector for the DM is a DM vector associated
405: with the same DM then the current base global vector for the DM is ok and if we replace it with the new one
406: we get a circular dependency that prevents the DM from being destroy when it should be.
407: This occurs when SNESSet/GetNPC() is used with a SNES that does not have a user provided
408: DM attached to it since the inner SNES (which shares the DM with the outer SNES) tries
409: to set its input vector (which is associated with the DM) as the base global vector.
410: Thanks to Juan P. Mendez Granado Re: [petsc-maint] Nonlinear conjugate gradien
411: for pointing out the problem.
412: */
413: if (vdm == dm) return(0);
414: PetscObjectReference((PetscObject)X);
415: VecDestroy(&shell->Xglobal);
416: shell->Xglobal = X;
417: return(0);
418: }
420: /*@C
421: DMShellSetCreateGlobalVector - sets the routine to create a global vector associated with the shell DM
423: Logically Collective
425: Input Arguments:
426: + dm - the shell DM
427: - func - the creation routine
429: Level: advanced
431: .seealso: DMShellSetGlobalVector(), DMShellSetCreateMatrix(), DMShellSetContext(), DMShellGetContext()
432: @*/
433: PetscErrorCode DMShellSetCreateGlobalVector(DM dm,PetscErrorCode (*func)(DM,Vec*))
434: {
437: dm->ops->createglobalvector = func;
438: return(0);
439: }
441: /*@
442: DMShellSetLocalVector - sets a template local vector associated with the DMShell
444: Logically Collective on dm
446: Input Arguments:
447: + dm - shell DM
448: - X - template vector
450: Level: advanced
452: .seealso: DMCreateLocalVector(), DMShellSetMatrix(), DMShellSetCreateLocalVector()
453: @*/
454: PetscErrorCode DMShellSetLocalVector(DM dm,Vec X)
455: {
456: DM_Shell *shell = (DM_Shell*)dm->data;
458: PetscBool isshell;
459: 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 Arguments:
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: {
501: dm->ops->createlocalvector = func;
502: return(0);
503: }
505: /*@C
506: DMShellSetGlobalToLocal - Sets the routines used to perform a global to local scatter
508: Logically Collective on dm
510: Input Arguments
511: + dm - the shell DM
512: . begin - the routine that begins the global to local scatter
513: - end - the routine that ends the global to local scatter
515: Notes:
516: If these functions are not provided but DMShellSetGlobalToLocalVecScatter() is called then
517: DMGlobalToLocalBeginDefaultShell()/DMGlobalToLocalEndDefaultShell() are used to to perform the transfers
519: Level: advanced
521: .seealso: DMShellSetLocalToGlobal(), DMGlobalToLocalBeginDefaultShell(), DMGlobalToLocalEndDefaultShell()
522: @*/
523: PetscErrorCode DMShellSetGlobalToLocal(DM dm,PetscErrorCode (*begin)(DM,Vec,InsertMode,Vec),PetscErrorCode (*end)(DM,Vec,InsertMode,Vec)) {
526: dm->ops->globaltolocalbegin = begin;
527: dm->ops->globaltolocalend = end;
528: return(0);
529: }
531: /*@C
532: DMShellSetLocalToGlobal - Sets the routines used to perform a local to global scatter
534: Logically Collective on dm
536: Input Arguments
537: + dm - the shell DM
538: . begin - the routine that begins the local to global scatter
539: - end - the routine that ends the local to global scatter
541: Notes:
542: If these functions are not provided but DMShellSetLocalToGlobalVecScatter() is called then
543: DMLocalToGlobalBeginDefaultShell()/DMLocalToGlobalEndDefaultShell() are used to to perform the transfers
545: Level: advanced
547: .seealso: DMShellSetGlobalToLocal()
548: @*/
549: PetscErrorCode DMShellSetLocalToGlobal(DM dm,PetscErrorCode (*begin)(DM,Vec,InsertMode,Vec),PetscErrorCode (*end)(DM,Vec,InsertMode,Vec)) {
552: dm->ops->localtoglobalbegin = begin;
553: dm->ops->localtoglobalend = end;
554: return(0);
555: }
557: /*@C
558: DMShellSetLocalToLocal - Sets the routines used to perform a local to local scatter
560: Logically Collective on dm
562: Input Arguments
563: + dm - the shell DM
564: . begin - the routine that begins the local to local scatter
565: - end - the routine that ends the local to local scatter
567: Notes:
568: If these functions are not provided but DMShellSetLocalToLocalVecScatter() is called then
569: DMLocalToLocalBeginDefaultShell()/DMLocalToLocalEndDefaultShell() are used to to perform the transfers
571: Level: advanced
573: .seealso: DMShellSetGlobalToLocal(), DMLocalToLocalBeginDefaultShell(), DMLocalToLocalEndDefaultShell()
574: @*/
575: PetscErrorCode DMShellSetLocalToLocal(DM dm,PetscErrorCode (*begin)(DM,Vec,InsertMode,Vec),PetscErrorCode (*end)(DM,Vec,InsertMode,Vec)) {
578: dm->ops->localtolocalbegin = begin;
579: dm->ops->localtolocalend = end;
580: return(0);
581: }
583: /*@
584: DMShellSetGlobalToLocalVecScatter - Sets a VecScatter context for global to local communication
586: Logically Collective on dm
588: Input Arguments
589: + dm - the shell DM
590: - gtol - the global to local VecScatter context
592: Level: advanced
594: .seealso: DMShellSetGlobalToLocal(), DMGlobalToLocalBeginDefaultShell(), DMGlobalToLocalEndDefaultShell()
595: @*/
596: PetscErrorCode DMShellSetGlobalToLocalVecScatter(DM dm, VecScatter gtol)
597: {
598: DM_Shell *shell = (DM_Shell*)dm->data;
604: PetscObjectReference((PetscObject)gtol);
605: VecScatterDestroy(&shell->gtol);
606: shell->gtol = gtol;
607: return(0);
608: }
610: /*@
611: DMShellSetLocalToGlobalVecScatter - Sets a VecScatter context for local to global communication
613: Logically Collective on dm
615: Input Arguments
616: + dm - the shell DM
617: - ltog - the local to global VecScatter context
619: Level: advanced
621: .seealso: DMShellSetLocalToGlobal(), DMLocalToGlobalBeginDefaultShell(), DMLocalToGlobalEndDefaultShell()
622: @*/
623: PetscErrorCode DMShellSetLocalToGlobalVecScatter(DM dm, VecScatter ltog)
624: {
625: DM_Shell *shell = (DM_Shell*)dm->data;
631: PetscObjectReference((PetscObject)ltog);
632: VecScatterDestroy(&shell->ltog);
633: shell->ltog = ltog;
634: return(0);
635: }
637: /*@
638: DMShellSetLocalToLocalVecScatter - Sets a VecScatter context for local to local communication
640: Logically Collective on dm
642: Input Arguments
643: + dm - the shell DM
644: - ltol - the local to local VecScatter context
646: Level: advanced
648: .seealso: DMShellSetLocalToLocal(), DMLocalToLocalBeginDefaultShell(), DMLocalToLocalEndDefaultShell()
649: @*/
650: PetscErrorCode DMShellSetLocalToLocalVecScatter(DM dm, VecScatter ltol)
651: {
652: DM_Shell *shell = (DM_Shell*)dm->data;
658: PetscObjectReference((PetscObject)ltol);
659: VecScatterDestroy(&shell->ltol);
660: shell->ltol = ltol;
661: return(0);
662: }
664: /*@C
665: DMShellSetCoarsen - Set the routine used to coarsen the shell DM
667: Logically Collective on dm
669: Input Arguments
670: + dm - the shell DM
671: - coarsen - the routine that coarsens the DM
673: Level: advanced
675: .seealso: DMShellSetRefine(), DMCoarsen(), DMShellGetCoarsen(), DMShellSetContext(), DMShellGetContext()
676: @*/
677: PetscErrorCode DMShellSetCoarsen(DM dm, PetscErrorCode (*coarsen)(DM,MPI_Comm,DM*))
678: {
680: PetscBool isshell;
684: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
685: if (!isshell) return(0);
686: dm->ops->coarsen = coarsen;
687: return(0);
688: }
690: /*@C
691: DMShellGetCoarsen - Get the routine used to coarsen the shell DM
693: Logically Collective on dm
695: Input Argument:
696: . dm - the shell DM
698: Output Argument:
699: . coarsen - the routine that coarsens the DM
701: Level: advanced
703: .seealso: DMShellSetCoarsen(), DMCoarsen(), DMShellSetRefine(), DMRefine()
704: @*/
705: PetscErrorCode DMShellGetCoarsen(DM dm, PetscErrorCode (**coarsen)(DM,MPI_Comm,DM*))
706: {
708: PetscBool isshell;
712: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
713: if (!isshell) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Can only use with DMSHELL type DMs");
714: *coarsen = dm->ops->coarsen;
715: return(0);
716: }
718: /*@C
719: DMShellSetRefine - Set the routine used to refine the shell DM
721: Logically Collective on dm
723: Input Arguments
724: + dm - the shell DM
725: - refine - the routine that refines the DM
727: Level: advanced
729: .seealso: DMShellSetCoarsen(), DMRefine(), DMShellGetRefine(), DMShellSetContext(), DMShellGetContext()
730: @*/
731: PetscErrorCode DMShellSetRefine(DM dm, PetscErrorCode (*refine)(DM,MPI_Comm,DM*))
732: {
734: PetscBool isshell;
738: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
739: if (!isshell) return(0);
740: dm->ops->refine = refine;
741: return(0);
742: }
744: /*@C
745: DMShellGetRefine - Get the routine used to refine the shell DM
747: Logically Collective on dm
749: Input Argument:
750: . dm - the shell DM
752: Output Argument:
753: . refine - the routine that refines the DM
755: Level: advanced
757: .seealso: DMShellSetCoarsen(), DMCoarsen(), DMShellSetRefine(), DMRefine()
758: @*/
759: PetscErrorCode DMShellGetRefine(DM dm, PetscErrorCode (**refine)(DM,MPI_Comm,DM*))
760: {
762: PetscBool isshell;
766: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
767: if (!isshell) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Can only use with DMSHELL type DMs");
768: *refine = dm->ops->refine;
769: return(0);
770: }
772: /*@C
773: DMShellSetCreateInterpolation - Set the routine used to create the interpolation operator
775: Logically Collective on dm
777: Input Arguments
778: + dm - the shell DM
779: - interp - the routine to create the interpolation
781: Level: advanced
783: .seealso: DMShellSetCreateInjection(), DMCreateInterpolation(), DMShellGetCreateInterpolation(), DMShellSetCreateRestriction(), DMShellSetContext(), DMShellGetContext()
784: @*/
785: PetscErrorCode DMShellSetCreateInterpolation(DM dm, PetscErrorCode (*interp)(DM,DM,Mat*,Vec*))
786: {
788: PetscBool isshell;
792: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
793: if (!isshell) return(0);
794: dm->ops->createinterpolation = interp;
795: return(0);
796: }
798: /*@C
799: DMShellGetCreateInterpolation - Get the routine used to create the interpolation operator
801: Logically Collective on dm
803: Input Argument:
804: + dm - the shell DM
806: Output Argument:
807: - interp - the routine to create the interpolation
809: Level: advanced
811: .seealso: DMShellGetCreateInjection(), DMCreateInterpolation(), DMShellGetCreateRestriction(), DMShellSetContext(), DMShellGetContext()
812: @*/
813: PetscErrorCode DMShellGetCreateInterpolation(DM dm, PetscErrorCode (**interp)(DM,DM,Mat*,Vec*))
814: {
816: PetscBool isshell;
820: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
821: if (!isshell) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Can only use with DMSHELL type DMs");
822: *interp = dm->ops->createinterpolation;
823: return(0);
824: }
826: /*@C
827: DMShellSetCreateRestriction - Set the routine used to create the restriction operator
829: Logically Collective on dm
831: Input Arguments
832: + dm - the shell DM
833: - striction- the routine to create the restriction
835: Level: advanced
837: .seealso: DMShellSetCreateInjection(), DMCreateInterpolation(), DMShellGetCreateRestriction(), DMShellSetContext(), DMShellGetContext()
838: @*/
839: PetscErrorCode DMShellSetCreateRestriction(DM dm, PetscErrorCode (*restriction)(DM,DM,Mat*))
840: {
842: PetscBool isshell;
846: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
847: if (!isshell) return(0);
848: dm->ops->createrestriction = restriction;
849: return(0);
850: }
852: /*@C
853: DMShellGetCreateRestriction - Get the routine used to create the restriction operator
855: Logically Collective on dm
857: Input Argument:
858: + dm - the shell DM
860: Output Argument:
861: - restriction - the routine to create the restriction
863: Level: advanced
865: .seealso: DMShellSetCreateInjection(), DMCreateInterpolation(), DMShellSetContext(), DMShellGetContext()
866: @*/
867: PetscErrorCode DMShellGetCreateRestriction(DM dm, PetscErrorCode (**restriction)(DM,DM,Mat*))
868: {
870: PetscBool isshell;
874: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
875: if (!isshell) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Can only use with DMSHELL type DMs");
876: *restriction = dm->ops->createrestriction;
877: return(0);
878: }
880: /*@C
881: DMShellSetCreateInjection - Set the routine used to create the injection operator
883: Logically Collective on dm
885: Input Arguments
886: + dm - the shell DM
887: - inject - the routine to create the injection
889: Level: advanced
891: .seealso: DMShellSetCreateInterpolation(), DMCreateInjection(), DMShellGetCreateInjection(), DMShellSetContext(), DMShellGetContext()
892: @*/
893: PetscErrorCode DMShellSetCreateInjection(DM dm, PetscErrorCode (*inject)(DM,DM,Mat*))
894: {
896: PetscBool isshell;
900: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
901: if (!isshell) return(0);
902: dm->ops->createinjection = inject;
903: return(0);
904: }
906: /*@C
907: DMShellGetCreateInjection - Get the routine used to create the injection operator
909: Logically Collective on dm
911: Input Argument:
912: + dm - the shell DM
914: Output Argument:
915: - inject - the routine to create the injection
917: Level: advanced
919: .seealso: DMShellGetCreateInterpolation(), DMCreateInjection(), DMShellSetContext(), DMShellGetContext()
920: @*/
921: PetscErrorCode DMShellGetCreateInjection(DM dm, PetscErrorCode (**inject)(DM,DM,Mat*))
922: {
924: PetscBool isshell;
928: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
929: if (!isshell) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Can only use with DMSHELL type DMs");
930: *inject = dm->ops->createinjection;
931: return(0);
932: }
934: /*@C
935: DMShellSetCreateFieldDecomposition - Set the routine used to create a decomposition of fields for the shell DM
937: Logically Collective on dm
939: Input Arguments
940: + dm - the shell DM
941: - decomp - the routine to create the decomposition
943: Level: advanced
945: .seealso: DMCreateFieldDecomposition(), DMShellSetContext(), DMShellGetContext()
946: @*/
947: PetscErrorCode DMShellSetCreateFieldDecomposition(DM dm, PetscErrorCode (*decomp)(DM,PetscInt*,char***, IS**,DM**))
948: {
950: PetscBool isshell;
954: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
955: if (!isshell) return(0);
956: dm->ops->createfielddecomposition = decomp;
957: return(0);
958: }
960: /*@C
961: DMShellSetCreateDomainDecomposition - Set the routine used to create a domain decomposition for the shell DM
963: Logically Collective on dm
965: Input Arguments
966: + dm - the shell DM
967: - decomp - the routine to create the decomposition
969: Level: advanced
971: .seealso: DMCreateDomainDecomposition(), DMShellSetContext(), DMShellGetContext()
972: @*/
973: PetscErrorCode DMShellSetCreateDomainDecomposition(DM dm, PetscErrorCode (*decomp)(DM,PetscInt*,char***, IS**,IS**,DM**))
974: {
976: PetscBool isshell;
980: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
981: if (!isshell) return(0);
982: dm->ops->createdomaindecomposition = decomp;
983: return(0);
984: }
986: /*@C
987: DMShellSetCreateDomainDecompositionScatters - Set the routine used to create the scatter contexts for domain decomposition with a shell DM
989: Logically Collective on dm
991: Input Arguments
992: + dm - the shell DM
993: - scatter - the routine to create the scatters
995: Level: advanced
997: .seealso: DMCreateDomainDecompositionScatters(), DMShellSetContext(), DMShellGetContext()
998: @*/
999: PetscErrorCode DMShellSetCreateDomainDecompositionScatters(DM dm, PetscErrorCode (*scatter)(DM,PetscInt,DM*,VecScatter**,VecScatter**,VecScatter**))
1000: {
1002: PetscBool isshell;
1006: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
1007: if (!isshell) return(0);
1008: dm->ops->createddscatters = scatter;
1009: return(0);
1010: }
1012: /*@C
1013: DMShellSetCreateSubDM - Set the routine used to create a sub DM from the shell DM
1015: Logically Collective on dm
1017: Input Arguments
1018: + dm - the shell DM
1019: - subdm - the routine to create the decomposition
1021: Level: advanced
1023: .seealso: DMCreateSubDM(), DMShellGetCreateSubDM(), DMShellSetContext(), DMShellGetContext()
1024: @*/
1025: PetscErrorCode DMShellSetCreateSubDM(DM dm, PetscErrorCode (*subdm)(DM,PetscInt,const PetscInt[],IS*,DM*))
1026: {
1028: PetscBool isshell;
1032: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
1033: if (!isshell) return(0);
1034: dm->ops->createsubdm = subdm;
1035: return(0);
1036: }
1038: /*@C
1039: DMShellGetCreateSubDM - Get the routine used to create a sub DM from the shell DM
1041: Logically Collective on dm
1043: Input Argument:
1044: + dm - the shell DM
1046: Output Argument:
1047: - subdm - the routine to create the decomposition
1049: Level: advanced
1051: .seealso: DMCreateSubDM(), DMShellSetCreateSubDM(), DMShellSetContext(), DMShellGetContext()
1052: @*/
1053: PetscErrorCode DMShellGetCreateSubDM(DM dm, PetscErrorCode (**subdm)(DM,PetscInt,const PetscInt[],IS*,DM*))
1054: {
1056: PetscBool isshell;
1060: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
1061: if (!isshell) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Can only use with DMSHELL type DMs");
1062: *subdm = dm->ops->createsubdm;
1063: return(0);
1064: }
1066: static PetscErrorCode DMDestroy_Shell(DM dm)
1067: {
1069: DM_Shell *shell = (DM_Shell*)dm->data;
1072: MatDestroy(&shell->A);
1073: VecDestroy(&shell->Xglobal);
1074: VecDestroy(&shell->Xlocal);
1075: VecScatterDestroy(&shell->gtol);
1076: VecScatterDestroy(&shell->ltog);
1077: VecScatterDestroy(&shell->ltol);
1078: /* This was originally freed in DMDestroy(), but that prevents reference counting of backend objects */
1079: PetscFree(shell);
1080: return(0);
1081: }
1083: static PetscErrorCode DMView_Shell(DM dm,PetscViewer v)
1084: {
1086: DM_Shell *shell = (DM_Shell*)dm->data;
1089: VecView(shell->Xglobal,v);
1090: return(0);
1091: }
1093: static PetscErrorCode DMLoad_Shell(DM dm,PetscViewer v)
1094: {
1096: DM_Shell *shell = (DM_Shell*)dm->data;
1099: VecCreate(PetscObjectComm((PetscObject)dm),&shell->Xglobal);
1100: VecLoad(shell->Xglobal,v);
1101: return(0);
1102: }
1104: PetscErrorCode DMCreateSubDM_Shell(DM dm, PetscInt numFields, const PetscInt fields[], IS *is, DM *subdm)
1105: {
1109: if (subdm) {DMShellCreate(PetscObjectComm((PetscObject) dm), subdm);}
1110: DMCreateSectionSubDM(dm, numFields, fields, is, subdm);
1111: return(0);
1112: }
1114: PETSC_EXTERN PetscErrorCode DMCreate_Shell(DM dm)
1115: {
1117: DM_Shell *shell;
1120: PetscNewLog(dm,&shell);
1121: dm->data = shell;
1123: dm->ops->destroy = DMDestroy_Shell;
1124: dm->ops->createglobalvector = DMCreateGlobalVector_Shell;
1125: dm->ops->createlocalvector = DMCreateLocalVector_Shell;
1126: dm->ops->creatematrix = DMCreateMatrix_Shell;
1127: dm->ops->view = DMView_Shell;
1128: dm->ops->load = DMLoad_Shell;
1129: dm->ops->globaltolocalbegin = DMGlobalToLocalBeginDefaultShell;
1130: dm->ops->globaltolocalend = DMGlobalToLocalEndDefaultShell;
1131: dm->ops->localtoglobalbegin = DMLocalToGlobalBeginDefaultShell;
1132: dm->ops->localtoglobalend = DMLocalToGlobalEndDefaultShell;
1133: dm->ops->localtolocalbegin = DMLocalToLocalBeginDefaultShell;
1134: dm->ops->localtolocalend = DMLocalToLocalEndDefaultShell;
1135: dm->ops->createsubdm = DMCreateSubDM_Shell;
1136: return(0);
1137: }
1139: /*@
1140: DMShellCreate - Creates a shell DM object, used to manage user-defined problem data
1142: Collective
1144: Input Parameter:
1145: . comm - the processors that will share the global vector
1147: Output Parameters:
1148: . shell - the shell DM
1150: Level: advanced
1152: .seealso DMDestroy(), DMCreateGlobalVector(), DMCreateLocalVector(), DMShellSetContext(), DMShellGetContext()
1153: @*/
1154: PetscErrorCode DMShellCreate(MPI_Comm comm,DM *dm)
1155: {
1160: DMCreate(comm,dm);
1161: DMSetType(*dm,DMSHELL);
1162: DMSetUp(*dm);
1163: return(0);
1164: }