Actual source code: dmshell.c
petsc-3.12.2 2019-11-22
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: PetscObjectReference((PetscObject)A);
210: MatZeroEntries(A);
211: *J = A;
212: } else { /* Need to create a copy, could use MAT_SHARE_NONZERO_PATTERN in most cases */
213: MatDuplicate(A,MAT_DO_NOT_COPY_VALUES,J);
214: MatZeroEntries(*J);
215: }
216: return(0);
217: }
219: PetscErrorCode DMCreateGlobalVector_Shell(DM dm,Vec *gvec)
220: {
222: DM_Shell *shell = (DM_Shell*)dm->data;
223: Vec X;
228: *gvec = 0;
229: X = shell->Xglobal;
230: if (!X) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_USER,"Must call DMShellSetGlobalVector() or DMShellSetCreateGlobalVector()");
231: if (((PetscObject)X)->refct < 2) { /* We have an exclusive reference so we can give it out */
232: PetscObjectReference((PetscObject)X);
233: VecZeroEntries(X);
234: *gvec = X;
235: } else { /* Need to create a copy, could use MAT_SHARE_NONZERO_PATTERN in most cases */
236: VecDuplicate(X,gvec);
237: VecZeroEntries(*gvec);
238: }
239: VecSetDM(*gvec,dm);
240: return(0);
241: }
243: PetscErrorCode DMCreateLocalVector_Shell(DM dm,Vec *gvec)
244: {
246: DM_Shell *shell = (DM_Shell*)dm->data;
247: Vec X;
252: *gvec = 0;
253: X = shell->Xlocal;
254: if (!X) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_USER,"Must call DMShellSetLocalVector() or DMShellSetCreateLocalVector()");
255: if (((PetscObject)X)->refct < 2) { /* We have an exclusive reference so we can give it out */
256: PetscObjectReference((PetscObject)X);
257: VecZeroEntries(X);
258: *gvec = X;
259: } else { /* Need to create a copy, could use MAT_SHARE_NONZERO_PATTERN in most cases */
260: VecDuplicate(X,gvec);
261: VecZeroEntries(*gvec);
262: }
263: VecSetDM(*gvec,dm);
264: return(0);
265: }
267: /*@
268: DMShellSetContext - set some data to be usable by this DM
270: Collective
272: Input Arguments:
273: + dm - shell DM
274: - ctx - the context
276: Level: advanced
278: .seealso: DMCreateMatrix(), DMShellGetContext()
279: @*/
280: PetscErrorCode DMShellSetContext(DM dm,void *ctx)
281: {
282: DM_Shell *shell = (DM_Shell*)dm->data;
284: PetscBool isshell;
288: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
289: if (!isshell) return(0);
290: shell->ctx = ctx;
291: return(0);
292: }
294: /*@
295: DMShellGetContext - Returns the user-provided context associated to the DM
297: Collective
299: Input Argument:
300: . dm - shell DM
302: Output Argument:
303: . ctx - the context
305: Level: advanced
307: .seealso: DMCreateMatrix(), DMShellSetContext()
308: @*/
309: PetscErrorCode DMShellGetContext(DM dm,void **ctx)
310: {
311: DM_Shell *shell = (DM_Shell*)dm->data;
313: PetscBool isshell;
317: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
318: if (!isshell) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Can only use with DMSHELL type DMs");
319: *ctx = shell->ctx;
320: return(0);
321: }
323: /*@
324: DMShellSetMatrix - sets a template matrix associated with the DMShell
326: Collective
328: Input Arguments:
329: + dm - shell DM
330: - J - template matrix
332: Level: advanced
334: .seealso: DMCreateMatrix(), DMShellSetCreateMatrix(), DMShellSetContext(), DMShellGetContext()
335: @*/
336: PetscErrorCode DMShellSetMatrix(DM dm,Mat J)
337: {
338: DM_Shell *shell = (DM_Shell*)dm->data;
340: PetscBool isshell;
345: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
346: if (!isshell) return(0);
347: PetscObjectReference((PetscObject)J);
348: MatDestroy(&shell->A);
349: shell->A = J;
350: return(0);
351: }
353: /*@C
354: DMShellSetCreateMatrix - sets the routine to create a matrix associated with the shell DM
356: Logically Collective on dm
358: Input Arguments:
359: + dm - the shell DM
360: - func - the function to create a matrix
362: Level: advanced
364: .seealso: DMCreateMatrix(), DMShellSetMatrix(), DMShellSetContext(), DMShellGetContext()
365: @*/
366: PetscErrorCode DMShellSetCreateMatrix(DM dm,PetscErrorCode (*func)(DM,Mat*))
367: {
371: dm->ops->creatematrix = func;
372: return(0);
373: }
375: /*@
376: DMShellSetGlobalVector - sets a template global vector associated with the DMShell
378: Logically Collective on dm
380: Input Arguments:
381: + dm - shell DM
382: - X - template vector
384: Level: advanced
386: .seealso: DMCreateGlobalVector(), DMShellSetMatrix(), DMShellSetCreateGlobalVector()
387: @*/
388: PetscErrorCode DMShellSetGlobalVector(DM dm,Vec X)
389: {
390: DM_Shell *shell = (DM_Shell*)dm->data;
392: PetscBool isshell;
393: DM vdm;
398: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
399: if (!isshell) return(0);
400: VecGetDM(X,&vdm);
401: /*
402: if the vector proposed as the new base global vector for the DM is a DM vector associated
403: with the same DM then the current base global vector for the DM is ok and if we replace it with the new one
404: we get a circular dependency that prevents the DM from being destroy when it should be.
405: This occurs when SNESSet/GetNPC() is used with a SNES that does not have a user provided
406: DM attached to it since the inner SNES (which shares the DM with the outer SNES) tries
407: to set its input vector (which is associated with the DM) as the base global vector.
408: Thanks to Juan P. Mendez Granado Re: [petsc-maint] Nonlinear conjugate gradien
409: for pointing out the problem.
410: */
411: if (vdm == dm) return(0);
412: PetscObjectReference((PetscObject)X);
413: VecDestroy(&shell->Xglobal);
414: shell->Xglobal = X;
415: return(0);
416: }
418: /*@C
419: DMShellSetCreateGlobalVector - sets the routine to create a global vector associated with the shell DM
421: Logically Collective
423: Input Arguments:
424: + dm - the shell DM
425: - func - the creation routine
427: Level: advanced
429: .seealso: DMShellSetGlobalVector(), DMShellSetCreateMatrix(), DMShellSetContext(), DMShellGetContext()
430: @*/
431: PetscErrorCode DMShellSetCreateGlobalVector(DM dm,PetscErrorCode (*func)(DM,Vec*))
432: {
436: dm->ops->createglobalvector = func;
437: return(0);
438: }
440: /*@
441: DMShellSetLocalVector - sets a template local vector associated with the DMShell
443: Logically Collective on dm
445: Input Arguments:
446: + dm - shell DM
447: - X - template vector
449: Level: advanced
451: .seealso: DMCreateLocalVector(), DMShellSetMatrix(), DMShellSetCreateLocalVector()
452: @*/
453: PetscErrorCode DMShellSetLocalVector(DM dm,Vec X)
454: {
455: DM_Shell *shell = (DM_Shell*)dm->data;
457: PetscBool isshell;
458: DM vdm;
463: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
464: if (!isshell) return(0);
465: VecGetDM(X,&vdm);
466: /*
467: if the vector proposed as the new base global vector for the DM is a DM vector associated
468: with the same DM then the current base global vector for the DM is ok and if we replace it with the new one
469: we get a circular dependency that prevents the DM from being destroy when it should be.
470: This occurs when SNESSet/GetNPC() is used with a SNES that does not have a user provided
471: DM attached to it since the inner SNES (which shares the DM with the outer SNES) tries
472: to set its input vector (which is associated with the DM) as the base global vector.
473: Thanks to Juan P. Mendez Granado Re: [petsc-maint] Nonlinear conjugate gradien
474: for pointing out the problem.
475: */
476: if (vdm == dm) return(0);
477: PetscObjectReference((PetscObject)X);
478: VecDestroy(&shell->Xlocal);
479: shell->Xlocal = X;
480: return(0);
481: }
483: /*@C
484: DMShellSetCreateLocalVector - sets the routine to create a local vector associated with the shell DM
486: Logically Collective
488: Input Arguments:
489: + dm - the shell DM
490: - func - the creation routine
492: Level: advanced
494: .seealso: DMShellSetLocalVector(), DMShellSetCreateMatrix(), DMShellSetContext(), DMShellGetContext()
495: @*/
496: PetscErrorCode DMShellSetCreateLocalVector(DM dm,PetscErrorCode (*func)(DM,Vec*))
497: {
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 Arguments
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)) {
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 Arguments
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)) {
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 Arguments
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)) {
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 Arguments
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;
603: PetscObjectReference((PetscObject)gtol);
604: VecScatterDestroy(&shell->gtol);
605: shell->gtol = gtol;
606: return(0);
607: }
609: /*@
610: DMShellSetLocalToGlobalVecScatter - Sets a VecScatter context for local to global communication
612: Logically Collective on dm
614: Input Arguments
615: + dm - the shell DM
616: - ltog - the local to global VecScatter context
618: Level: advanced
620: .seealso: DMShellSetLocalToGlobal(), DMLocalToGlobalBeginDefaultShell(), DMLocalToGlobalEndDefaultShell()
621: @*/
622: PetscErrorCode DMShellSetLocalToGlobalVecScatter(DM dm, VecScatter ltog)
623: {
624: DM_Shell *shell = (DM_Shell*)dm->data;
630: PetscObjectReference((PetscObject)ltog);
631: VecScatterDestroy(&shell->ltog);
632: shell->ltog = ltog;
633: return(0);
634: }
636: /*@
637: DMShellSetLocalToLocalVecScatter - Sets a VecScatter context for local to local communication
639: Logically Collective on dm
641: Input Arguments
642: + dm - the shell DM
643: - ltol - the local to local VecScatter context
645: Level: advanced
647: .seealso: DMShellSetLocalToLocal(), DMLocalToLocalBeginDefaultShell(), DMLocalToLocalEndDefaultShell()
648: @*/
649: PetscErrorCode DMShellSetLocalToLocalVecScatter(DM dm, VecScatter ltol)
650: {
651: DM_Shell *shell = (DM_Shell*)dm->data;
657: PetscObjectReference((PetscObject)ltol);
658: VecScatterDestroy(&shell->ltol);
659: shell->ltol = ltol;
660: return(0);
661: }
663: /*@C
664: DMShellSetCoarsen - Set the routine used to coarsen the shell DM
666: Logically Collective on dm
668: Input Arguments
669: + dm - the shell DM
670: - coarsen - the routine that coarsens the DM
672: Level: advanced
674: .seealso: DMShellSetRefine(), DMCoarsen(), DMShellGetCoarsen(), DMShellSetContext(), DMShellGetContext()
675: @*/
676: PetscErrorCode DMShellSetCoarsen(DM dm, PetscErrorCode (*coarsen)(DM,MPI_Comm,DM*))
677: {
679: PetscBool isshell;
683: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
684: if (!isshell) return(0);
685: dm->ops->coarsen = coarsen;
686: return(0);
687: }
689: /*@C
690: DMShellGetCoarsen - Get the routine used to coarsen the shell DM
692: Logically Collective on dm
694: Input Argument:
695: . dm - the shell DM
697: Output Argument:
698: . coarsen - the routine that coarsens the DM
700: Level: advanced
702: .seealso: DMShellSetCoarsen(), DMCoarsen(), DMShellSetRefine(), DMRefine()
703: @*/
704: PetscErrorCode DMShellGetCoarsen(DM dm, PetscErrorCode (**coarsen)(DM,MPI_Comm,DM*))
705: {
707: PetscBool isshell;
711: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
712: if (!isshell) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Can only use with DMSHELL type DMs");
713: *coarsen = dm->ops->coarsen;
714: return(0);
715: }
717: /*@C
718: DMShellSetRefine - Set the routine used to refine the shell DM
720: Logically Collective on dm
722: Input Arguments
723: + dm - the shell DM
724: - refine - the routine that refines the DM
726: Level: advanced
728: .seealso: DMShellSetCoarsen(), DMRefine(), DMShellGetRefine(), DMShellSetContext(), DMShellGetContext()
729: @*/
730: PetscErrorCode DMShellSetRefine(DM dm, PetscErrorCode (*refine)(DM,MPI_Comm,DM*))
731: {
733: PetscBool isshell;
737: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
738: if (!isshell) return(0);
739: dm->ops->refine = refine;
740: return(0);
741: }
743: /*@C
744: DMShellGetRefine - Get the routine used to refine the shell DM
746: Logically Collective on dm
748: Input Argument:
749: . dm - the shell DM
751: Output Argument:
752: . refine - the routine that refines the DM
754: Level: advanced
756: .seealso: DMShellSetCoarsen(), DMCoarsen(), DMShellSetRefine(), DMRefine()
757: @*/
758: PetscErrorCode DMShellGetRefine(DM dm, PetscErrorCode (**refine)(DM,MPI_Comm,DM*))
759: {
761: PetscBool isshell;
765: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
766: if (!isshell) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Can only use with DMSHELL type DMs");
767: *refine = dm->ops->refine;
768: return(0);
769: }
771: /*@C
772: DMShellSetCreateInterpolation - Set the routine used to create the interpolation operator
774: Logically Collective on dm
776: Input Arguments
777: + dm - the shell DM
778: - interp - the routine to create the interpolation
780: Level: advanced
782: .seealso: DMShellSetCreateInjection(), DMCreateInterpolation(), DMShellGetCreateInterpolation(), DMShellSetCreateRestriction(), DMShellSetContext(), DMShellGetContext()
783: @*/
784: PetscErrorCode DMShellSetCreateInterpolation(DM dm, PetscErrorCode (*interp)(DM,DM,Mat*,Vec*))
785: {
787: PetscBool isshell;
791: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
792: if (!isshell) return(0);
793: dm->ops->createinterpolation = interp;
794: return(0);
795: }
797: /*@C
798: DMShellGetCreateInterpolation - Get the routine used to create the interpolation operator
800: Logically Collective on dm
802: Input Argument:
803: + dm - the shell DM
805: Output Argument:
806: - interp - the routine to create the interpolation
808: Level: advanced
810: .seealso: DMShellGetCreateInjection(), DMCreateInterpolation(), DMShellGetCreateRestriction(), DMShellSetContext(), DMShellGetContext()
811: @*/
812: PetscErrorCode DMShellGetCreateInterpolation(DM dm, PetscErrorCode (**interp)(DM,DM,Mat*,Vec*))
813: {
815: PetscBool isshell;
819: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
820: if (!isshell) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Can only use with DMSHELL type DMs");
821: *interp = dm->ops->createinterpolation;
822: return(0);
823: }
825: /*@C
826: DMShellSetCreateRestriction - Set the routine used to create the restriction operator
828: Logically Collective on dm
830: Input Arguments
831: + dm - the shell DM
832: - striction- the routine to create the restriction
834: Level: advanced
836: .seealso: DMShellSetCreateInjection(), DMCreateInterpolation(), DMShellGetCreateRestriction(), DMShellSetContext(), DMShellGetContext()
837: @*/
838: PetscErrorCode DMShellSetCreateRestriction(DM dm, PetscErrorCode (*restriction)(DM,DM,Mat*))
839: {
841: PetscBool isshell;
845: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
846: if (!isshell) return(0);
847: dm->ops->createrestriction = restriction;
848: return(0);
849: }
851: /*@C
852: DMShellGetCreateRestriction - Get the routine used to create the restriction operator
854: Logically Collective on dm
856: Input Argument:
857: + dm - the shell DM
859: Output Argument:
860: - restriction - the routine to create the restriction
862: Level: advanced
864: .seealso: DMShellSetCreateInjection(), DMCreateInterpolation(), DMShellSetContext(), DMShellGetContext()
865: @*/
866: PetscErrorCode DMShellGetCreateRestriction(DM dm, PetscErrorCode (**restriction)(DM,DM,Mat*))
867: {
869: PetscBool isshell;
873: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
874: if (!isshell) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Can only use with DMSHELL type DMs");
875: *restriction = dm->ops->createrestriction;
876: return(0);
877: }
879: /*@C
880: DMShellSetCreateInjection - Set the routine used to create the injection operator
882: Logically Collective on dm
884: Input Arguments
885: + dm - the shell DM
886: - inject - the routine to create the injection
888: Level: advanced
890: .seealso: DMShellSetCreateInterpolation(), DMCreateInjection(), DMShellGetCreateInjection(), DMShellSetContext(), DMShellGetContext()
891: @*/
892: PetscErrorCode DMShellSetCreateInjection(DM dm, PetscErrorCode (*inject)(DM,DM,Mat*))
893: {
895: PetscBool isshell;
899: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
900: if (!isshell) return(0);
901: dm->ops->createinjection = inject;
902: return(0);
903: }
905: /*@C
906: DMShellGetCreateInjection - Get the routine used to create the injection operator
908: Logically Collective on dm
910: Input Argument:
911: + dm - the shell DM
913: Output Argument:
914: - inject - the routine to create the injection
916: Level: advanced
918: .seealso: DMShellGetCreateInterpolation(), DMCreateInjection(), DMShellSetContext(), DMShellGetContext()
919: @*/
920: PetscErrorCode DMShellGetCreateInjection(DM dm, PetscErrorCode (**inject)(DM,DM,Mat*))
921: {
923: PetscBool isshell;
927: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
928: if (!isshell) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Can only use with DMSHELL type DMs");
929: *inject = dm->ops->createinjection;
930: return(0);
931: }
933: /*@C
934: DMShellSetCreateFieldDecomposition - Set the routine used to create a decomposition of fields for the shell DM
936: Logically Collective on dm
938: Input Arguments
939: + dm - the shell DM
940: - decomp - the routine to create the decomposition
942: Level: advanced
944: .seealso: DMCreateFieldDecomposition(), DMShellSetContext(), DMShellGetContext()
945: @*/
946: PetscErrorCode DMShellSetCreateFieldDecomposition(DM dm, PetscErrorCode (*decomp)(DM,PetscInt*,char***, IS**,DM**))
947: {
949: PetscBool isshell;
953: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
954: if (!isshell) return(0);
955: dm->ops->createfielddecomposition = decomp;
956: return(0);
957: }
959: /*@C
960: DMShellSetCreateDomainDecomposition - Set the routine used to create a domain decomposition for the shell DM
962: Logically Collective on dm
964: Input Arguments
965: + dm - the shell DM
966: - decomp - the routine to create the decomposition
968: Level: advanced
970: .seealso: DMCreateDomainDecomposition(), DMShellSetContext(), DMShellGetContext()
971: @*/
972: PetscErrorCode DMShellSetCreateDomainDecomposition(DM dm, PetscErrorCode (*decomp)(DM,PetscInt*,char***, IS**,IS**,DM**))
973: {
975: PetscBool isshell;
979: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
980: if (!isshell) return(0);
981: dm->ops->createdomaindecomposition = decomp;
982: return(0);
983: }
985: /*@C
986: DMShellSetCreateDomainDecompositionScatters - Set the routine used to create the scatter contexts for domain decomposition with a shell DM
988: Logically Collective on dm
990: Input Arguments
991: + dm - the shell DM
992: - scatter - the routine to create the scatters
994: Level: advanced
996: .seealso: DMCreateDomainDecompositionScatters(), DMShellSetContext(), DMShellGetContext()
997: @*/
998: PetscErrorCode DMShellSetCreateDomainDecompositionScatters(DM dm, PetscErrorCode (*scatter)(DM,PetscInt,DM*,VecScatter**,VecScatter**,VecScatter**))
999: {
1001: PetscBool isshell;
1005: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
1006: if (!isshell) return(0);
1007: dm->ops->createddscatters = scatter;
1008: return(0);
1009: }
1011: /*@C
1012: DMShellSetCreateSubDM - Set the routine used to create a sub DM from the shell DM
1014: Logically Collective on dm
1016: Input Arguments
1017: + dm - the shell DM
1018: - subdm - the routine to create the decomposition
1020: Level: advanced
1022: .seealso: DMCreateSubDM(), DMShellGetCreateSubDM(), DMShellSetContext(), DMShellGetContext()
1023: @*/
1024: PetscErrorCode DMShellSetCreateSubDM(DM dm, PetscErrorCode (*subdm)(DM,PetscInt,const PetscInt[],IS*,DM*))
1025: {
1027: PetscBool isshell;
1031: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
1032: if (!isshell) return(0);
1033: dm->ops->createsubdm = subdm;
1034: return(0);
1035: }
1037: /*@C
1038: DMShellGetCreateSubDM - Get the routine used to create a sub DM from the shell DM
1040: Logically Collective on dm
1042: Input Argument:
1043: + dm - the shell DM
1045: Output Argument:
1046: - subdm - the routine to create the decomposition
1048: Level: advanced
1050: .seealso: DMCreateSubDM(), DMShellSetCreateSubDM(), DMShellSetContext(), DMShellGetContext()
1051: @*/
1052: PetscErrorCode DMShellGetCreateSubDM(DM dm, PetscErrorCode (**subdm)(DM,PetscInt,const PetscInt[],IS*,DM*))
1053: {
1055: PetscBool isshell;
1059: PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
1060: if (!isshell) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Can only use with DMSHELL type DMs");
1061: *subdm = dm->ops->createsubdm;
1062: return(0);
1063: }
1065: static PetscErrorCode DMDestroy_Shell(DM dm)
1066: {
1068: DM_Shell *shell = (DM_Shell*)dm->data;
1071: MatDestroy(&shell->A);
1072: VecDestroy(&shell->Xglobal);
1073: VecDestroy(&shell->Xlocal);
1074: VecScatterDestroy(&shell->gtol);
1075: VecScatterDestroy(&shell->ltog);
1076: VecScatterDestroy(&shell->ltol);
1077: /* This was originally freed in DMDestroy(), but that prevents reference counting of backend objects */
1078: PetscFree(shell);
1079: return(0);
1080: }
1082: static PetscErrorCode DMView_Shell(DM dm,PetscViewer v)
1083: {
1085: DM_Shell *shell = (DM_Shell*)dm->data;
1088: VecView(shell->Xglobal,v);
1089: return(0);
1090: }
1092: static PetscErrorCode DMLoad_Shell(DM dm,PetscViewer v)
1093: {
1095: DM_Shell *shell = (DM_Shell*)dm->data;
1098: VecCreate(PetscObjectComm((PetscObject)dm),&shell->Xglobal);
1099: VecLoad(shell->Xglobal,v);
1100: return(0);
1101: }
1103: PetscErrorCode DMCreateSubDM_Shell(DM dm, PetscInt numFields, const PetscInt fields[], IS *is, DM *subdm)
1104: {
1108: if (subdm) {DMShellCreate(PetscObjectComm((PetscObject) dm), subdm);}
1109: DMCreateSectionSubDM(dm, numFields, fields, is, subdm);
1110: return(0);
1111: }
1113: PETSC_EXTERN PetscErrorCode DMCreate_Shell(DM dm)
1114: {
1116: DM_Shell *shell;
1119: PetscNewLog(dm,&shell);
1120: dm->data = shell;
1122: dm->ops->destroy = DMDestroy_Shell;
1123: dm->ops->createglobalvector = DMCreateGlobalVector_Shell;
1124: dm->ops->createlocalvector = DMCreateLocalVector_Shell;
1125: dm->ops->creatematrix = DMCreateMatrix_Shell;
1126: dm->ops->view = DMView_Shell;
1127: dm->ops->load = DMLoad_Shell;
1128: dm->ops->globaltolocalbegin = DMGlobalToLocalBeginDefaultShell;
1129: dm->ops->globaltolocalend = DMGlobalToLocalEndDefaultShell;
1130: dm->ops->localtoglobalbegin = DMLocalToGlobalBeginDefaultShell;
1131: dm->ops->localtoglobalend = DMLocalToGlobalEndDefaultShell;
1132: dm->ops->localtolocalbegin = DMLocalToLocalBeginDefaultShell;
1133: dm->ops->localtolocalend = DMLocalToLocalEndDefaultShell;
1134: dm->ops->createsubdm = DMCreateSubDM_Shell;
1135: return(0);
1136: }
1138: /*@
1139: DMShellCreate - Creates a shell DM object, used to manage user-defined problem data
1141: Collective
1143: Input Parameter:
1144: . comm - the processors that will share the global vector
1146: Output Parameters:
1147: . shell - the shell DM
1149: Level: advanced
1151: .seealso DMDestroy(), DMCreateGlobalVector(), DMCreateLocalVector(), DMShellSetContext(), DMShellGetContext()
1152: @*/
1153: PetscErrorCode DMShellCreate(MPI_Comm comm,DM *dm)
1154: {
1159: DMCreate(comm,dm);
1160: DMSetType(*dm,DMSHELL);
1161: DMSetUp(*dm);
1162: return(0);
1163: }