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: PetscCtx ctx;
13: PetscCtxDestroyFn *destroyctx;
14: } DM_Shell;
16: /*@
17: DMGlobalToLocalBeginDefaultShell - Uses the GlobalToLocal `VecScatter` context set by the user to begin a global to local scatter
19: Collective
21: Input Parameters:
22: + dm - `DMSHELL`
23: . g - global vector
24: . mode - `InsertMode`
25: - l - local vector
27: Level: advanced
29: Note:
30: 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.
32: .seealso: `DM`, `DMSHELL`, `DMGlobalToLocalEndDefaultShell()`
33: @*/
34: PetscErrorCode DMGlobalToLocalBeginDefaultShell(DM dm, Vec g, InsertMode mode, Vec l)
35: {
36: DM_Shell *shell = (DM_Shell *)dm->data;
38: PetscFunctionBegin;
39: PetscCheck(shell->gtol, ((PetscObject)dm)->comm, PETSC_ERR_ARG_WRONGSTATE, "Cannot be used without first setting the scatter context via DMShellSetGlobalToLocalVecScatter()");
40: PetscCall(VecScatterBegin(shell->gtol, g, l, mode, SCATTER_FORWARD));
41: PetscFunctionReturn(PETSC_SUCCESS);
42: }
44: /*@
45: DMGlobalToLocalEndDefaultShell - Uses the GlobalToLocal `VecScatter` context set by the user to end a global to local scatter
46: Collective
48: Input Parameters:
49: + dm - `DMSHELL`
50: . g - global vector
51: . mode - `InsertMode`
52: - l - local vector
54: Level: advanced
56: .seealso: `DM`, `DMSHELL`, `DMGlobalToLocalBeginDefaultShell()`
57: @*/
58: PetscErrorCode DMGlobalToLocalEndDefaultShell(DM dm, Vec g, InsertMode mode, Vec l)
59: {
60: DM_Shell *shell = (DM_Shell *)dm->data;
62: PetscFunctionBegin;
63: PetscCheck(shell->gtol, ((PetscObject)dm)->comm, PETSC_ERR_ARG_WRONGSTATE, "Cannot be used without first setting the scatter context via DMShellSetGlobalToLocalVecScatter()");
64: PetscCall(VecScatterEnd(shell->gtol, g, l, mode, SCATTER_FORWARD));
65: PetscFunctionReturn(PETSC_SUCCESS);
66: }
68: /*@
69: DMLocalToGlobalBeginDefaultShell - Uses the LocalToGlobal `VecScatter` context set by the user to begin a local to global scatter
70: Collective
72: Input Parameters:
73: + dm - `DMSHELL`
74: . l - local vector
75: . mode - `InsertMode`
76: - g - global vector
78: Level: advanced
80: Note:
81: 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.
83: .seealso: `DM`, `DMSHELL`, `DMLocalToGlobalEndDefaultShell()`
84: @*/
85: PetscErrorCode DMLocalToGlobalBeginDefaultShell(DM dm, Vec l, InsertMode mode, Vec g)
86: {
87: DM_Shell *shell = (DM_Shell *)dm->data;
89: PetscFunctionBegin;
90: PetscCheck(shell->ltog, ((PetscObject)dm)->comm, PETSC_ERR_ARG_WRONGSTATE, "Cannot be used without first setting the scatter context via DMShellSetLocalToGlobalVecScatter()");
91: PetscCall(VecScatterBegin(shell->ltog, l, g, mode, SCATTER_FORWARD));
92: PetscFunctionReturn(PETSC_SUCCESS);
93: }
95: /*@
96: DMLocalToGlobalEndDefaultShell - Uses the LocalToGlobal `VecScatter` context set by the user to end a local to global scatter
97: Collective
99: Input Parameters:
100: + dm - `DMSHELL`
101: . l - local vector
102: . mode - `InsertMode`
103: - g - global vector
105: Level: advanced
107: .seealso: `DM`, `DMSHELL`, `DMLocalToGlobalBeginDefaultShell()`
108: @*/
109: PetscErrorCode DMLocalToGlobalEndDefaultShell(DM dm, Vec l, InsertMode mode, Vec g)
110: {
111: DM_Shell *shell = (DM_Shell *)dm->data;
113: PetscFunctionBegin;
114: PetscCheck(shell->ltog, ((PetscObject)dm)->comm, PETSC_ERR_ARG_WRONGSTATE, "Cannot be used without first setting the scatter context via DMShellSetLocalToGlobalVecScatter()");
115: PetscCall(VecScatterEnd(shell->ltog, l, g, mode, SCATTER_FORWARD));
116: PetscFunctionReturn(PETSC_SUCCESS);
117: }
119: /*@
120: DMLocalToLocalBeginDefaultShell - Uses the LocalToLocal `VecScatter` context set by the user to begin a local to local scatter
121: Collective
123: Input Parameters:
124: + dm - `DMSHELL`
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:
134: 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.
136: .seealso: `DM`, `DMSHELL`, `DMLocalToLocalEndDefaultShell()`
137: @*/
138: PetscErrorCode DMLocalToLocalBeginDefaultShell(DM dm, Vec g, InsertMode mode, Vec l)
139: {
140: DM_Shell *shell = (DM_Shell *)dm->data;
142: PetscFunctionBegin;
143: PetscCheck(shell->ltol, ((PetscObject)dm)->comm, PETSC_ERR_ARG_WRONGSTATE, "Cannot be used without first setting the scatter context via DMShellSetLocalToLocalVecScatter()");
144: PetscCall(VecScatterBegin(shell->ltol, g, l, mode, SCATTER_FORWARD));
145: PetscFunctionReturn(PETSC_SUCCESS);
146: }
148: /*@
149: DMLocalToLocalEndDefaultShell - Uses the LocalToLocal `VecScatter` context set by the user to end a local to local scatter
150: Collective
152: Input Parameters:
153: + dm - `DMSHELL`
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: `DM`, `DMSHELL`, `DMLocalToLocalBeginDefaultShell()`
163: @*/
164: PetscErrorCode DMLocalToLocalEndDefaultShell(DM dm, Vec g, InsertMode mode, Vec l)
165: {
166: DM_Shell *shell = (DM_Shell *)dm->data;
168: PetscFunctionBegin;
169: PetscCheck(shell->ltol, ((PetscObject)dm)->comm, PETSC_ERR_ARG_WRONGSTATE, "Cannot be used without first setting the scatter context via DMShellSetGlobalToLocalVecScatter()");
170: PetscCall(VecScatterEnd(shell->ltol, g, l, mode, SCATTER_FORWARD));
171: PetscFunctionReturn(PETSC_SUCCESS);
172: }
174: static PetscErrorCode DMCreateMatrix_Shell(DM dm, Mat *J)
175: {
176: DM_Shell *shell = (DM_Shell *)dm->data;
177: Mat A;
179: PetscFunctionBegin;
181: PetscAssertPointer(J, 2);
182: if (!shell->A) {
183: PetscInt m, M;
185: PetscCheck(shell->Xglobal, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Must call DMShellSetMatrix(), DMShellSetCreateMatrix(), or provide a vector");
186: PetscCall(PetscInfo(dm, "Naively creating matrix using global vector distribution without preallocation\n"));
187: PetscCall(VecGetSize(shell->Xglobal, &M));
188: PetscCall(VecGetLocalSize(shell->Xglobal, &m));
189: PetscCall(MatCreate(PetscObjectComm((PetscObject)dm), &shell->A));
190: PetscCall(MatSetSizes(shell->A, m, m, M, M));
191: PetscCall(MatSetType(shell->A, dm->mattype));
192: PetscCall(MatSetUp(shell->A));
193: }
194: A = shell->A;
195: PetscCall(MatDuplicate(A, MAT_SHARE_NONZERO_PATTERN, J));
196: PetscCall(MatSetDM(*J, dm));
197: PetscFunctionReturn(PETSC_SUCCESS);
198: }
200: static PetscErrorCode DMCreateGlobalVector_Shell(DM dm, Vec *gvec)
201: {
202: DM_Shell *shell = (DM_Shell *)dm->data;
203: Vec X;
205: PetscFunctionBegin;
207: PetscAssertPointer(gvec, 2);
208: *gvec = NULL;
209: X = shell->Xglobal;
210: PetscCheck(X, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Must call DMShellSetGlobalVector() or DMShellSetCreateGlobalVector()");
211: /* Need to create a copy in order to attach the DM to the vector */
212: PetscCall(VecDuplicate(X, gvec));
213: PetscCall(VecZeroEntries(*gvec));
214: PetscCall(VecSetDM(*gvec, dm));
215: PetscFunctionReturn(PETSC_SUCCESS);
216: }
218: static PetscErrorCode DMCreateLocalVector_Shell(DM dm, Vec *gvec)
219: {
220: DM_Shell *shell = (DM_Shell *)dm->data;
221: Vec X;
223: PetscFunctionBegin;
225: PetscAssertPointer(gvec, 2);
226: *gvec = NULL;
227: X = shell->Xlocal;
228: PetscCheck(X, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Must call DMShellSetLocalVector() or DMShellSetCreateLocalVector()");
229: /* Need to create a copy in order to attach the DM to the vector */
230: PetscCall(VecDuplicate(X, gvec));
231: PetscCall(VecZeroEntries(*gvec));
232: PetscCall(VecSetDM(*gvec, dm));
233: PetscFunctionReturn(PETSC_SUCCESS);
234: }
236: /*@C
237: DMShellSetDestroyContext - set a function that destroys the context provided with `DMShellSetContext()`
239: Collective
241: Input Parameters:
242: + dm - the `DM` to attach the `destroyctx()` function to
243: - destroyctx - the function that destroys the context
245: Level: advanced
247: .seealso: `DM`, `DMSHELL`, `DMShellSetContext()`, `DMShellGetContext()`
248: @*/
249: PetscErrorCode DMShellSetDestroyContext(DM dm, PetscCtxDestroyFn *destroyctx)
250: {
251: DM_Shell *shell = (DM_Shell *)dm->data;
252: PetscBool isshell;
254: PetscFunctionBegin;
256: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell));
257: if (!isshell) PetscFunctionReturn(PETSC_SUCCESS);
258: shell->destroyctx = destroyctx;
259: PetscFunctionReturn(PETSC_SUCCESS);
260: }
262: /*@
263: DMShellSetContext - set some data to be usable by this `DMSHELL`
265: Collective
267: Input Parameters:
268: + dm - `DMSHELL`
269: - ctx - the context
271: Level: advanced
273: .seealso: `DM`, `DMSHELL`, `DMCreateMatrix()`, `DMShellGetContext()`
274: @*/
275: PetscErrorCode DMShellSetContext(DM dm, PetscCtx ctx)
276: {
277: DM_Shell *shell = (DM_Shell *)dm->data;
278: PetscBool isshell;
280: PetscFunctionBegin;
282: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell));
283: if (!isshell) PetscFunctionReturn(PETSC_SUCCESS);
284: shell->ctx = ctx;
285: PetscFunctionReturn(PETSC_SUCCESS);
286: }
288: /*@
289: DMShellGetContext - Returns the user-provided context associated to the `DMSHELL`
291: Collective
293: Input Parameter:
294: . dm - `DMSHELL`
296: Output Parameter:
297: . ctx - the context
299: Level: advanced
301: Fortran Notes:
302: This only works when the context is a Fortran derived type or a `PetscObject`. Declare `ctx` with
303: .vb
304: type(tUsertype), pointer :: ctx
305: .ve
307: .seealso: `DM`, `DMSHELL`, `DMCreateMatrix()`, `DMShellSetContext()`
308: @*/
309: PetscErrorCode DMShellGetContext(DM dm, PetscCtxRt ctx)
310: {
311: DM_Shell *shell = (DM_Shell *)dm->data;
312: PetscBool isshell;
314: PetscFunctionBegin;
316: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell));
317: PetscCheck(isshell, PetscObjectComm((PetscObject)dm), PETSC_ERR_SUP, "Can only use with DMSHELL type DMs");
318: *(void **)ctx = shell->ctx;
319: PetscFunctionReturn(PETSC_SUCCESS);
320: }
322: /*@
323: DMShellSetMatrix - sets a template matrix associated with the `DMSHELL`
325: Collective
327: Input Parameters:
328: + dm - `DMSHELL`
329: - J - template matrix
331: Level: advanced
333: Developer Notes:
334: 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.
336: .seealso: `DM`, `DMSHELL`, `DMCreateMatrix()`, `DMShellSetCreateMatrix()`, `DMShellSetContext()`, `DMShellGetContext()`
337: @*/
338: PetscErrorCode DMShellSetMatrix(DM dm, Mat J)
339: {
340: DM_Shell *shell = (DM_Shell *)dm->data;
341: PetscBool isshell;
342: DM mdm;
344: PetscFunctionBegin;
347: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell));
348: if (!isshell) PetscFunctionReturn(PETSC_SUCCESS);
349: if (J == shell->A) PetscFunctionReturn(PETSC_SUCCESS);
350: PetscCall(MatGetDM(J, &mdm));
351: PetscCall(PetscObjectReference((PetscObject)J));
352: PetscCall(MatDestroy(&shell->A));
353: if (mdm == dm) {
354: PetscCall(MatDuplicate(J, MAT_SHARE_NONZERO_PATTERN, &shell->A));
355: PetscCall(MatSetDM(shell->A, NULL));
356: } else shell->A = J;
357: PetscFunctionReturn(PETSC_SUCCESS);
358: }
360: /*@C
361: DMShellSetCreateMatrix - sets the routine to create a matrix associated with the `DMSHELL`
363: Logically Collective
365: Input Parameters:
366: + dm - the `DMSHELL`
367: - func - the function to create a matrix
369: Calling sequence of `func`:
370: + dm - the `DM`
371: - mat - the `Mat` to be created
373: Level: advanced
375: .seealso: `DM`, `DMSHELL`, `DMCreateMatrix()`, `DMShellSetMatrix()`, `DMShellSetContext()`, `DMShellGetContext()`
376: @*/
377: PetscErrorCode DMShellSetCreateMatrix(DM dm, PetscErrorCode (*func)(DM dm, Mat *mat))
378: {
379: PetscFunctionBegin;
381: dm->ops->creatematrix = func;
382: PetscFunctionReturn(PETSC_SUCCESS);
383: }
385: /*@
386: DMShellSetGlobalVector - sets a template global vector associated with the `DMSHELL`
388: Logically Collective
390: Input Parameters:
391: + dm - `DMSHELL`
392: - X - template vector
394: Level: advanced
396: .seealso: `DM`, `DMSHELL`, `DMCreateGlobalVector()`, `DMShellSetMatrix()`, `DMShellSetCreateGlobalVector()`
397: @*/
398: PetscErrorCode DMShellSetGlobalVector(DM dm, Vec X)
399: {
400: DM_Shell *shell = (DM_Shell *)dm->data;
401: PetscBool isshell;
402: DM vdm;
404: PetscFunctionBegin;
407: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell));
408: if (!isshell) PetscFunctionReturn(PETSC_SUCCESS);
409: PetscCall(VecGetDM(X, &vdm));
410: /*
411: if the vector proposed as the new base global vector for the DM is a DM vector associated
412: with the same DM then the current base global vector for the DM is ok and if we replace it with the new one
413: we get a circular dependency that prevents the DM from being destroy when it should be.
414: This occurs when SNESSet/GetNPC() is used with a SNES that does not have a user provided
415: DM attached to it since the inner SNES (which shares the DM with the outer SNES) tries
416: to set its input vector (which is associated with the DM) as the base global vector.
417: Thanks to Juan P. Mendez Granado Re: [petsc-maint] Nonlinear conjugate gradien
418: for pointing out the problem.
419: */
420: if (vdm == dm) PetscFunctionReturn(PETSC_SUCCESS);
421: PetscCall(PetscObjectReference((PetscObject)X));
422: PetscCall(VecDestroy(&shell->Xglobal));
423: shell->Xglobal = X;
424: PetscCall(DMClearGlobalVectors(dm));
425: PetscCall(DMClearNamedGlobalVectors(dm));
426: PetscFunctionReturn(PETSC_SUCCESS);
427: }
429: /*@
430: DMShellGetGlobalVector - Returns the template global vector associated with the `DMSHELL`, or `NULL` if it was not set
432: Not Collective
434: Input Parameters:
435: + dm - `DMSHELL`
436: - X - template vector
438: Level: advanced
440: .seealso: `DM`, `DMSHELL`, `DMShellSetGlobalVector()`, `DMShellSetCreateGlobalVector()`, `DMCreateGlobalVector()`
441: @*/
442: PetscErrorCode DMShellGetGlobalVector(DM dm, Vec *X)
443: {
444: DM_Shell *shell = (DM_Shell *)dm->data;
445: PetscBool isshell;
447: PetscFunctionBegin;
449: PetscAssertPointer(X, 2);
450: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell));
451: if (!isshell) PetscFunctionReturn(PETSC_SUCCESS);
452: *X = shell->Xglobal;
453: PetscFunctionReturn(PETSC_SUCCESS);
454: }
456: /*@C
457: DMShellSetCreateGlobalVector - sets the routine to create a global vector associated with the `DMSHELL`
459: Logically Collective
461: Input Parameters:
462: + dm - the `DMSHELL`
463: - func - the creation routine
465: Calling sequence of `func`:
466: + dm - the `DM`
467: - g - the global `Vec` to be created
469: Level: advanced
471: .seealso: `DM`, `DMSHELL`, `DMShellSetGlobalVector()`, `DMShellSetCreateMatrix()`, `DMShellSetContext()`, `DMShellGetContext()`
472: @*/
473: PetscErrorCode DMShellSetCreateGlobalVector(DM dm, PetscErrorCode (*func)(DM dm, Vec *g))
474: {
475: PetscFunctionBegin;
477: dm->ops->createglobalvector = func;
478: PetscFunctionReturn(PETSC_SUCCESS);
479: }
481: /*@
482: DMShellSetLocalVector - sets a template local vector associated with the `DMSHELL`
484: Logically Collective
486: Input Parameters:
487: + dm - `DMSHELL`
488: - X - template vector
490: Level: advanced
492: .seealso: `DM`, `DMSHELL`, `DMCreateLocalVector()`, `DMShellSetMatrix()`, `DMShellSetCreateLocalVector()`
493: @*/
494: PetscErrorCode DMShellSetLocalVector(DM dm, Vec X)
495: {
496: DM_Shell *shell = (DM_Shell *)dm->data;
497: PetscBool isshell;
498: DM vdm;
500: PetscFunctionBegin;
503: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell));
504: if (!isshell) PetscFunctionReturn(PETSC_SUCCESS);
505: PetscCall(VecGetDM(X, &vdm));
506: /*
507: if the vector proposed as the new base global vector for the DM is a DM vector associated
508: with the same DM then the current base global vector for the DM is ok and if we replace it with the new one
509: we get a circular dependency that prevents the DM from being destroy when it should be.
510: This occurs when SNESSet/GetNPC() is used with a SNES that does not have a user provided
511: DM attached to it since the inner SNES (which shares the DM with the outer SNES) tries
512: to set its input vector (which is associated with the DM) as the base global vector.
513: Thanks to Juan P. Mendez Granado Re: [petsc-maint] Nonlinear conjugate gradien
514: for pointing out the problem.
515: */
516: if (vdm == dm) PetscFunctionReturn(PETSC_SUCCESS);
517: PetscCall(PetscObjectReference((PetscObject)X));
518: PetscCall(VecDestroy(&shell->Xlocal));
519: shell->Xlocal = X;
520: PetscCall(DMClearLocalVectors(dm));
521: PetscCall(DMClearNamedLocalVectors(dm));
522: PetscFunctionReturn(PETSC_SUCCESS);
523: }
525: /*@C
526: DMShellSetCreateLocalVector - sets the routine to create a local vector associated with the `DMSHELL`
528: Logically Collective
530: Input Parameters:
531: + dm - the `DMSHELL`
532: - func - the creation routine
534: Calling sequence of `func`:
535: + dm - the `DM`
536: - l - the local `Vec` to be created
538: Level: advanced
540: .seealso: `DM`, `DMSHELL`, `DMShellSetLocalVector()`, `DMShellSetCreateMatrix()`, `DMShellSetContext()`, `DMShellGetContext()`
541: @*/
542: PetscErrorCode DMShellSetCreateLocalVector(DM dm, PetscErrorCode (*func)(DM dm, Vec *l))
543: {
544: PetscFunctionBegin;
546: dm->ops->createlocalvector = func;
547: PetscFunctionReturn(PETSC_SUCCESS);
548: }
550: /*@C
551: DMShellSetGlobalToLocal - Sets the routines used to perform a global to local scatter
553: Logically Collective
555: Input Parameters:
556: + dm - the `DMSHELL`
557: . begin - the routine that begins the global to local scatter
558: - end - the routine that ends the global to local scatter
560: Calling sequence of `begin`:
561: + dm - the `DM`
562: . global - the global `Vec` to be communicated
563: . mode - insert mode of the resulting vector
564: - local - the local `Vec` to receive the result
566: Calling sequence of `end`:
567: + dm - the `DM`
568: . global - the global `Vec` to be communicated
569: . mode - insert mode of the resulting vector
570: - local - the local `Vec` to receive the result
572: Level: advanced
574: Note:
575: If these functions are not provided but `DMShellSetGlobalToLocalVecScatter()` is called then
576: `DMGlobalToLocalBeginDefaultShell()`/`DMGlobalToLocalEndDefaultShell()` are used to perform the transfers
578: .seealso: `DM`, `DMSHELL`, `DMShellSetLocalToGlobal()`, `DMGlobalToLocalBeginDefaultShell()`, `DMGlobalToLocalEndDefaultShell()`
579: @*/
580: PetscErrorCode DMShellSetGlobalToLocal(DM dm, PetscErrorCode (*begin)(DM dm, Vec global, InsertMode mode, Vec local), PetscErrorCode (*end)(DM dm, Vec global, InsertMode mode, Vec local))
581: {
582: PetscFunctionBegin;
584: dm->ops->globaltolocalbegin = begin;
585: dm->ops->globaltolocalend = end;
586: PetscFunctionReturn(PETSC_SUCCESS);
587: }
589: /*@C
590: DMShellSetLocalToGlobal - Sets the routines used to perform a local to global scatter
592: Logically Collective
594: Input Parameters:
595: + dm - the `DMSHELL`
596: . begin - the routine that begins the local to global scatter
597: - end - the routine that ends the local to global scatter
599: Calling sequence of `begin`:
600: + dm - the `DM`
601: . local - the local `Vec` to be communicated
602: . mode - insert mode of the resulting vector
603: - global - the global `Vec` to receive the result
605: Calling sequence of `end`:
606: + dm - the `DM`
607: . local - the local `Vec` to be communicated
608: . mode - insert mode of the resulting vector
609: - global - the global `Vec` to receive the result
611: Level: advanced
613: Note:
614: If these functions are not provided but `DMShellSetLocalToGlobalVecScatter()` is called then
615: `DMLocalToGlobalBeginDefaultShell()`/`DMLocalToGlobalEndDefaultShell()` are used to perform the transfers
617: .seealso: `DM`, `DMSHELL`, `DMShellSetGlobalToLocal()`, `InsertMode`, `VecScatter`, `DMLocalToGlobal()`, `DMGlobalToLocal()`
618: @*/
619: PetscErrorCode DMShellSetLocalToGlobal(DM dm, PetscErrorCode (*begin)(DM dm, Vec local, InsertMode mode, Vec global), PetscErrorCode (*end)(DM dm, Vec local, InsertMode mode, Vec global))
620: {
621: PetscFunctionBegin;
623: dm->ops->localtoglobalbegin = begin;
624: dm->ops->localtoglobalend = end;
625: PetscFunctionReturn(PETSC_SUCCESS);
626: }
628: /*@C
629: DMShellSetLocalToLocal - Sets the routines used to perform a local to local scatter
631: Logically Collective
633: Input Parameters:
634: + dm - the `DMSHELL`
635: . begin - the routine that begins the local to local scatter
636: - end - the routine that ends the local to local scatter
638: Calling sequence of `begin`:
639: + dm - the `DM`
640: . local - the local `Vec` to be communicated
641: . mode - insert mode of the resulting vector
642: - nlocal - the local `Vec` to receive the result
644: Calling sequence of `end`:
645: + dm - the `DM`
646: . local - the local `Vec` to be communicated
647: . mode - insert mode of the resulting vector
648: - nlocal - the local `Vec` to receive the result
650: Level: advanced
652: Note:
653: If these functions are not provided but `DMShellSetLocalToLocalVecScatter()` is called then
654: `DMLocalToLocalBeginDefaultShell()`/`DMLocalToLocalEndDefaultShell()` are used to perform the transfers
656: .seealso: `DM`, `DMSHELL`, `DMShellSetGlobalToLocal()`, `DMLocalToLocalBeginDefaultShell()`, `DMLocalToLocalEndDefaultShell()`, `DMLocalToLocalBegin()`, `DMLocalToLocalEnd()`
657: @*/
658: PetscErrorCode DMShellSetLocalToLocal(DM dm, PetscErrorCode (*begin)(DM dm, Vec local, InsertMode mode, Vec nlocal), PetscErrorCode (*end)(DM dm, Vec local, InsertMode mode, Vec nlocal))
659: {
660: PetscFunctionBegin;
662: dm->ops->localtolocalbegin = begin;
663: dm->ops->localtolocalend = end;
664: PetscFunctionReturn(PETSC_SUCCESS);
665: }
667: /*@
668: DMShellSetGlobalToLocalVecScatter - Sets a `VecScatter` context for global to local communication
670: Logically Collective
672: Input Parameters:
673: + dm - the `DMSHELL`
674: - gtol - the global to local `VecScatter` context
676: Level: advanced
678: .seealso: `DM`, `DMSHELL`, `DMShellSetGlobalToLocal()`, `DMGlobalToLocalBeginDefaultShell()`, `DMGlobalToLocalEndDefaultShell()`
679: @*/
680: PetscErrorCode DMShellSetGlobalToLocalVecScatter(DM dm, VecScatter gtol)
681: {
682: DM_Shell *shell = (DM_Shell *)dm->data;
684: PetscFunctionBegin;
687: PetscCall(PetscObjectReference((PetscObject)gtol));
688: PetscCall(VecScatterDestroy(&shell->gtol));
689: shell->gtol = gtol;
690: PetscFunctionReturn(PETSC_SUCCESS);
691: }
693: /*@
694: DMShellSetLocalToGlobalVecScatter - Sets a` VecScatter` context for local to global communication
696: Logically Collective
698: Input Parameters:
699: + dm - the `DMSHELL`
700: - ltog - the local to global `VecScatter` context
702: Level: advanced
704: .seealso: `DM`, `DMSHELL`, `DMShellSetLocalToGlobal()`, `DMLocalToGlobalBeginDefaultShell()`, `DMLocalToGlobalEndDefaultShell()`
705: @*/
706: PetscErrorCode DMShellSetLocalToGlobalVecScatter(DM dm, VecScatter ltog)
707: {
708: DM_Shell *shell = (DM_Shell *)dm->data;
710: PetscFunctionBegin;
713: PetscCall(PetscObjectReference((PetscObject)ltog));
714: PetscCall(VecScatterDestroy(&shell->ltog));
715: shell->ltog = ltog;
716: PetscFunctionReturn(PETSC_SUCCESS);
717: }
719: /*@
720: DMShellSetLocalToLocalVecScatter - Sets a `VecScatter` context for local to local communication
722: Logically Collective
724: Input Parameters:
725: + dm - the `DMSHELL`
726: - ltol - the local to local `VecScatter` context
728: Level: advanced
730: .seealso: `DM`, `DMSHELL`, `DMShellSetLocalToLocal()`, `DMLocalToLocalBeginDefaultShell()`, `DMLocalToLocalEndDefaultShell()`
731: @*/
732: PetscErrorCode DMShellSetLocalToLocalVecScatter(DM dm, VecScatter ltol)
733: {
734: DM_Shell *shell = (DM_Shell *)dm->data;
736: PetscFunctionBegin;
739: PetscCall(PetscObjectReference((PetscObject)ltol));
740: PetscCall(VecScatterDestroy(&shell->ltol));
741: shell->ltol = ltol;
742: PetscFunctionReturn(PETSC_SUCCESS);
743: }
745: /*@C
746: DMShellSetCoarsen - Set the routine used to coarsen the `DMSHELL`
748: Logically Collective
750: Input Parameters:
751: + dm - the `DMSHELL`
752: - coarsen - the routine that coarsens the `DM`
754: Calling sequence of `coarsen`:
755: + fine - the `DM` to coarsen
756: . comm - the `MPI_Comm` to share the coarser `DM`
757: - coarse - the resulting coarse `DM`
759: Level: advanced
761: .seealso: `DM`, `DMSHELL`, `DMShellSetRefine()`, `DMCoarsen()`, `DMShellGetCoarsen()`, `DMShellSetContext()`, `DMShellGetContext()`
762: @*/
763: PetscErrorCode DMShellSetCoarsen(DM dm, PetscErrorCode (*coarsen)(DM fine, MPI_Comm comm, DM *coarse))
764: {
765: PetscBool isshell;
767: PetscFunctionBegin;
769: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell));
770: if (!isshell) PetscFunctionReturn(PETSC_SUCCESS);
771: dm->ops->coarsen = coarsen;
772: PetscFunctionReturn(PETSC_SUCCESS);
773: }
775: /*@C
776: DMShellGetCoarsen - Get the routine used to coarsen the `DMSHELL`
778: Logically Collective
780: Input Parameter:
781: . dm - the `DMSHELL`
783: Output Parameter:
784: . coarsen - the routine that coarsens the `DM`
786: Calling sequence of `coarsen`:
787: + fine - the `DM` to coarsen
788: . comm - the `MPI_Comm` to share the coarser `DM`
789: - coarse - the resulting coarse `DM`
791: Level: advanced
793: .seealso: `DM`, `DMSHELL`, `DMShellSetCoarsen()`, `DMCoarsen()`, `DMShellSetRefine()`, `DMRefine()`
794: @*/
795: PetscErrorCode DMShellGetCoarsen(DM dm, PetscErrorCode (**coarsen)(DM fine, MPI_Comm comm, DM *coarse))
796: {
797: PetscBool isshell;
799: PetscFunctionBegin;
801: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell));
802: PetscCheck(isshell, PetscObjectComm((PetscObject)dm), PETSC_ERR_SUP, "Can only use with DMSHELL type DMs");
803: *coarsen = dm->ops->coarsen;
804: PetscFunctionReturn(PETSC_SUCCESS);
805: }
807: /*@C
808: DMShellSetRefine - Set the routine used to refine the `DMSHELL`
810: Logically Collective
812: Input Parameters:
813: + dm - the `DMSHELL`
814: - refine - the routine that refines the `DM`
816: Calling sequence of `refine`:
817: + coarse - the `DM` to refine
818: . comm - the `MPI_Comm` to share the finer `DM`
819: - fine - the resulting fine `DM`
821: Level: advanced
823: .seealso: `DM`, `DMSHELL`, `DMShellSetCoarsen()`, `DMRefine()`, `DMShellGetRefine()`, `DMShellSetContext()`, `DMShellGetContext()`
824: @*/
825: PetscErrorCode DMShellSetRefine(DM dm, PetscErrorCode (*refine)(DM coarse, MPI_Comm comm, DM *fine))
826: {
827: PetscBool isshell;
829: PetscFunctionBegin;
831: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell));
832: if (!isshell) PetscFunctionReturn(PETSC_SUCCESS);
833: dm->ops->refine = refine;
834: PetscFunctionReturn(PETSC_SUCCESS);
835: }
837: /*@C
838: DMShellGetRefine - Get the routine used to refine the `DMSHELL`
840: Logically Collective
842: Input Parameter:
843: . dm - the `DMSHELL`
845: Output Parameter:
846: . refine - the routine that refines the `DM`
848: Calling sequence of `refine`:
849: + coarse - the `DM` to refine
850: . comm - the `MPI_Comm` to share the finer `DM`
851: - fine - the resulting fine `DM`
853: Level: advanced
855: .seealso: `DM`, `DMSHELL`, `DMShellSetCoarsen()`, `DMCoarsen()`, `DMShellSetRefine()`, `DMRefine()`
856: @*/
857: PetscErrorCode DMShellGetRefine(DM dm, PetscErrorCode (**refine)(DM coarse, MPI_Comm comm, DM *fine))
858: {
859: PetscBool isshell;
861: PetscFunctionBegin;
863: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell));
864: PetscCheck(isshell, PetscObjectComm((PetscObject)dm), PETSC_ERR_SUP, "Can only use with DMSHELL type DMs");
865: *refine = dm->ops->refine;
866: PetscFunctionReturn(PETSC_SUCCESS);
867: }
869: /*@C
870: DMShellSetCreateInterpolation - Set the routine used to create the interpolation operator
872: Logically Collective
874: Input Parameters:
875: + dm - the `DMSHELL`
876: - interp - the routine to create the interpolation
878: Calling sequence of `interp`:
879: + coarse - the `DM` to refine to
880: . fine - the fine `DM`
881: . interp - the output interpolation `Mat`
882: - rscale - an output scaling `Vec`, see `DMCreateInterpolationScale()`
884: Level: advanced
886: .seealso: `DM`, `DMSHELL`, `DMShellSetCreateInjection()`, `DMCreateInterpolation()`, `DMShellGetCreateInterpolation()`, `DMShellSetCreateRestriction()`, `DMShellSetContext()`, `DMShellGetContext()`
887: @*/
888: PetscErrorCode DMShellSetCreateInterpolation(DM dm, PetscErrorCode (*interp)(DM coarse, DM fine, Mat *interp, Vec *rscale))
889: {
890: PetscBool isshell;
892: PetscFunctionBegin;
894: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell));
895: if (!isshell) PetscFunctionReturn(PETSC_SUCCESS);
896: dm->ops->createinterpolation = interp;
897: PetscFunctionReturn(PETSC_SUCCESS);
898: }
900: /*@C
901: DMShellGetCreateInterpolation - Get the routine used to create the interpolation operator
903: Logically Collective
905: Input Parameter:
906: . dm - the `DMSHELL`
908: Output Parameter:
909: . interp - the routine to create the interpolation
911: Calling sequence of `interp`:
912: + coarse - the `DM` to refine to
913: . fine - the fine `DM`
914: . interp - the output interpolation `Mat`
915: - rscale - an output scaling `Vec`, see `DMCreateInterpolationScale()`
917: Level: advanced
919: .seealso: `DM`, `DMSHELL`, `DMShellGetCreateInjection()`, `DMCreateInterpolation()`, `DMShellGetCreateRestriction()`, `DMShellSetContext()`, `DMShellGetContext()`
920: @*/
921: PetscErrorCode DMShellGetCreateInterpolation(DM dm, PetscErrorCode (**interp)(DM coarse, DM fine, Mat *interp, Vec *rscale))
922: {
923: PetscBool isshell;
925: PetscFunctionBegin;
927: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell));
928: PetscCheck(isshell, PetscObjectComm((PetscObject)dm), PETSC_ERR_SUP, "Can only use with DMSHELL type DMs");
929: *interp = dm->ops->createinterpolation;
930: PetscFunctionReturn(PETSC_SUCCESS);
931: }
933: /*@C
934: DMShellSetCreateRestriction - Set the routine used to create the restriction operator
936: Logically Collective
938: Input Parameters:
939: + dm - the `DMSHELL`
940: - restriction - the routine to create the restriction
942: Calling sequence of `restriction`:
943: + fine - the fine `DM`
944: . coarse - the `DM` to restrict to
945: - restrct - the output restriction `Mat`
947: Level: advanced
949: .seealso: `DM`, `DMSHELL`, `DMShellSetCreateInjection()`, `DMCreateInterpolation()`, `DMShellGetCreateRestriction()`, `DMShellSetContext()`, `DMShellGetContext()`
950: @*/
951: PetscErrorCode DMShellSetCreateRestriction(DM dm, PetscErrorCode (*restriction)(DM fine, DM coarse, Mat *restrct))
952: {
953: PetscBool isshell;
955: PetscFunctionBegin;
957: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell));
958: if (!isshell) PetscFunctionReturn(PETSC_SUCCESS);
959: dm->ops->createrestriction = restriction;
960: PetscFunctionReturn(PETSC_SUCCESS);
961: }
963: /*@C
964: DMShellGetCreateRestriction - Get the routine used to create the restriction operator
966: Logically Collective
968: Input Parameter:
969: . dm - the `DMSHELL`
971: Output Parameter:
972: . restriction - the routine to create the restriction
974: Calling sequence of `restriction`:
975: + fine - the fine `DM`
976: . coarse - the `DM` to restrict to
977: - restrct - the output restriction `Mat`
979: Level: advanced
981: .seealso: `DM`, `DMSHELL`, `DMShellSetCreateInjection()`, `DMCreateInterpolation()`, `DMShellSetContext()`, `DMShellGetContext()`
982: @*/
983: PetscErrorCode DMShellGetCreateRestriction(DM dm, PetscErrorCode (**restriction)(DM fine, DM coarse, Mat *restrct))
984: {
985: PetscBool isshell;
987: PetscFunctionBegin;
989: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell));
990: PetscCheck(isshell, PetscObjectComm((PetscObject)dm), PETSC_ERR_SUP, "Can only use with DMSHELL type DMs");
991: *restriction = dm->ops->createrestriction;
992: PetscFunctionReturn(PETSC_SUCCESS);
993: }
995: /*@C
996: DMShellSetCreateInjection - Set the routine used to create the injection operator
998: Logically Collective
1000: Input Parameters:
1001: + dm - the `DMSHELL`
1002: - inject - the routine to create the injection
1004: Calling sequence of `inject`:
1005: + fine - the fine `DM`
1006: . coarse - the `DM` to inject to
1007: - inject - the output injection `Mat`
1009: Level: advanced
1011: .seealso: `DM`, `DMSHELL`, `DMShellSetCreateInterpolation()`, `DMCreateInjection()`, `DMShellGetCreateInjection()`, `DMShellSetContext()`, `DMShellGetContext()`
1012: @*/
1013: PetscErrorCode DMShellSetCreateInjection(DM dm, PetscErrorCode (*inject)(DM fine, DM coarse, Mat *inject))
1014: {
1015: PetscBool isshell;
1017: PetscFunctionBegin;
1019: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell));
1020: if (!isshell) PetscFunctionReturn(PETSC_SUCCESS);
1021: dm->ops->createinjection = inject;
1022: PetscFunctionReturn(PETSC_SUCCESS);
1023: }
1025: /*@C
1026: DMShellGetCreateInjection - Get the routine used to create the injection operator
1028: Logically Collective
1030: Input Parameter:
1031: . dm - the `DMSHELL`
1033: Output Parameter:
1034: . inject - the routine to create the injection
1036: Calling sequence of `inject`:
1037: + fine - the fine `DM`
1038: . coarse - the `DM` to inject to
1039: - inject - the output injection `Mat`
1041: Level: advanced
1043: .seealso: `DM`, `DMSHELL`, `DMShellGetCreateInterpolation()`, `DMCreateInjection()`, `DMShellSetContext()`, `DMShellGetContext()`
1044: @*/
1045: PetscErrorCode DMShellGetCreateInjection(DM dm, PetscErrorCode (**inject)(DM fine, DM coarse, Mat *inject))
1046: {
1047: PetscBool isshell;
1049: PetscFunctionBegin;
1051: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell));
1052: PetscCheck(isshell, PetscObjectComm((PetscObject)dm), PETSC_ERR_SUP, "Can only use with DMSHELL type DMs");
1053: *inject = dm->ops->createinjection;
1054: PetscFunctionReturn(PETSC_SUCCESS);
1055: }
1057: /*@C
1058: DMShellSetCreateFieldDecomposition - Set the routine used to create a decomposition of fields for the `DMSHELL`
1060: Logically Collective
1062: Input Parameters:
1063: + dm - the `DMSHELL`
1064: - decomp - the routine to create the decomposition
1066: Calling sequence of `decomp`:
1067: + dm - the `DM` to decompose into fields
1068: . len - output, the number of fields (or `NULL` if not requested)
1069: . namelist - output, the name for each field (or `NULL` if not requested)
1070: . islist - output, the global indices for each field (or `NULL` if not requested)
1071: - dmlist - output, the `DM`s for each field subproblem (or `NULL`, if not requested; if `NULL` is returned, no `DM`s are defined)
1073: Level: advanced
1075: .seealso: `DM`, `DMSHELL`, `DMCreateFieldDecomposition()`, `DMShellSetContext()`, `DMShellGetContext()`
1076: @*/
1077: PetscErrorCode DMShellSetCreateFieldDecomposition(DM dm, PetscErrorCode (*decomp)(DM dm, PetscInt *len, char **namelist[], IS *islist[], DM *dmlist[]))
1078: {
1079: PetscBool isshell;
1081: PetscFunctionBegin;
1083: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell));
1084: if (!isshell) PetscFunctionReturn(PETSC_SUCCESS);
1085: dm->ops->createfielddecomposition = decomp;
1086: PetscFunctionReturn(PETSC_SUCCESS);
1087: }
1089: /*@C
1090: DMShellSetCreateDomainDecomposition - Set the routine used to create a domain decomposition for the `DMSHELL`
1092: Logically Collective
1094: Input Parameters:
1095: + dm - the `DMSHELL`
1096: - decomp - the routine to create the decomposition
1098: Calling sequence of `decomp`:
1099: + dm - the `DM` to decompose into domains
1100: . len - output, the number of domains (or `NULL` if not requested)
1101: . namelist - output, the name for each domain (or `NULL` if not requested)
1102: . innerlist - output, the global indices for each domain's inner region (or `NULL` if not requested)
1103: . outerlist - output, the global indices for each domain's outer region (or `NULL` if not requested)
1104: - dmlist - output, the `DM`s for each field subproblem (or `NULL`, if not requested; if `NULL` is returned, no `DM`s are defined)
1106: Level: advanced
1108: .seealso: `DM`, `DMSHELL`, `DMCreateDomainDecomposition()`, `DMShellSetContext()`, `DMShellGetContext()`
1109: @*/
1110: PetscErrorCode DMShellSetCreateDomainDecomposition(DM dm, PetscErrorCode (*decomp)(DM dm, PetscInt *len, char **namelist[], IS *innerlist[], IS *outerlist[], DM *dmlist[]))
1111: {
1112: PetscBool isshell;
1114: PetscFunctionBegin;
1116: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell));
1117: if (!isshell) PetscFunctionReturn(PETSC_SUCCESS);
1118: dm->ops->createdomaindecomposition = decomp;
1119: PetscFunctionReturn(PETSC_SUCCESS);
1120: }
1122: /*@C
1123: DMShellSetCreateDomainDecompositionScatters - Set the routine used to create the scatter contexts for domain decomposition with a `DMSHELL`
1125: Logically Collective
1127: Input Parameters:
1128: + dm - the `DMSHELL`
1129: - scatter - the routine to create the scatters
1131: Calling sequence of `scatter`:
1132: + dm - the `DM` to decompose into domains
1133: . n - number of subdomains
1134: . subdms - the sub `DM`
1135: . iscat - output, the inner scatters for the subdomains
1136: . oscat - output, outer scatters for the subdomains
1137: - gscat - output, the global scatters for the subdomains
1139: Level: advanced
1141: .seealso: `DM`, `DMSHELL`, `DMCreateDomainDecompositionScatters()`, `DMShellSetContext()`, `DMShellGetContext()`
1142: @*/
1143: PetscErrorCode DMShellSetCreateDomainDecompositionScatters(DM dm, PetscErrorCode (*scatter)(DM dm, PetscInt n, DM subdms[], VecScatter *iscat[], VecScatter *oscat[], VecScatter *gscat[]))
1144: {
1145: PetscBool isshell;
1147: PetscFunctionBegin;
1149: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell));
1150: if (!isshell) PetscFunctionReturn(PETSC_SUCCESS);
1151: dm->ops->createddscatters = scatter;
1152: PetscFunctionReturn(PETSC_SUCCESS);
1153: }
1155: /*@C
1156: DMShellSetCreateSubDM - Set the routine used to create a sub `DM` from the `DMSHELL`
1158: Logically Collective
1160: Input Parameters:
1161: + dm - the `DMSHELL`
1162: - subdm - the routine to create the decomposition
1164: Calling sequence of `subdm`:
1165: + dm - the original `DM`
1166: . numFields - the number of fields to create
1167: . fields - the fields to create for
1168: . is - output, the `IS` defining the sub `DM`
1169: - subdm - the sub `DM`
1171: Level: advanced
1173: .seealso: `DM`, `DMSHELL`, `DMCreateSubDM()`, `DMShellGetCreateSubDM()`, `DMShellSetContext()`, `DMShellGetContext()`
1174: @*/
1175: PetscErrorCode DMShellSetCreateSubDM(DM dm, PetscErrorCode (*subdm)(DM dm, PetscInt numFields, const PetscInt fields[], IS *is, DM *subdm))
1176: {
1177: PetscBool isshell;
1179: PetscFunctionBegin;
1181: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell));
1182: if (!isshell) PetscFunctionReturn(PETSC_SUCCESS);
1183: dm->ops->createsubdm = subdm;
1184: PetscFunctionReturn(PETSC_SUCCESS);
1185: }
1187: /*@C
1188: DMShellGetCreateSubDM - Get the routine used to create a sub `DM` from the `DMSHELL`
1190: Logically Collective
1192: Input Parameter:
1193: . dm - the `DMSHELL`
1195: Output Parameter:
1196: . subdm - the routine to create the decomposition
1198: Calling sequence of `subdm`:
1199: + dm - the original `DM`
1200: . numFields - the number of fields to create
1201: . fields - the fields to create for
1202: . is - output, the `IS` defining the sub `DM`
1203: - subdm - the sub `DM`
1205: Level: advanced
1207: .seealso: `DM`, `DMSHELL`, `DMCreateSubDM()`, `DMShellSetCreateSubDM()`, `DMShellSetContext()`, `DMShellGetContext()`
1208: @*/
1209: PetscErrorCode DMShellGetCreateSubDM(DM dm, PetscErrorCode (**subdm)(DM dm, PetscInt numFields, const PetscInt fields[], IS *is, DM *subdm))
1210: {
1211: PetscBool isshell;
1213: PetscFunctionBegin;
1215: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMSHELL, &isshell));
1216: PetscCheck(isshell, PetscObjectComm((PetscObject)dm), PETSC_ERR_SUP, "Can only use with DMSHELL type DMs");
1217: *subdm = dm->ops->createsubdm;
1218: PetscFunctionReturn(PETSC_SUCCESS);
1219: }
1221: static PetscErrorCode DMDestroy_Shell(DM dm)
1222: {
1223: DM_Shell *shell = (DM_Shell *)dm->data;
1225: PetscFunctionBegin;
1226: if (shell->destroyctx) PetscCallBack("Destroy Context", (*shell->destroyctx)(&shell->ctx));
1227: PetscCall(MatDestroy(&shell->A));
1228: PetscCall(VecDestroy(&shell->Xglobal));
1229: PetscCall(VecDestroy(&shell->Xlocal));
1230: PetscCall(VecScatterDestroy(&shell->gtol));
1231: PetscCall(VecScatterDestroy(&shell->ltog));
1232: PetscCall(VecScatterDestroy(&shell->ltol));
1233: /* This was originally freed in DMDestroy(), but that prevents reference counting of backend objects */
1234: PetscCall(PetscFree(shell));
1235: PetscFunctionReturn(PETSC_SUCCESS);
1236: }
1238: static PetscErrorCode DMView_Shell(DM dm, PetscViewer v)
1239: {
1240: DM_Shell *shell = (DM_Shell *)dm->data;
1242: PetscFunctionBegin;
1243: if (shell->Xglobal) PetscCall(VecView(shell->Xglobal, v));
1244: PetscFunctionReturn(PETSC_SUCCESS);
1245: }
1247: static PetscErrorCode DMLoad_Shell(DM dm, PetscViewer v)
1248: {
1249: DM_Shell *shell = (DM_Shell *)dm->data;
1251: PetscFunctionBegin;
1252: PetscCall(VecCreate(PetscObjectComm((PetscObject)dm), &shell->Xglobal));
1253: PetscCall(VecLoad(shell->Xglobal, v));
1254: PetscFunctionReturn(PETSC_SUCCESS);
1255: }
1257: static PetscErrorCode DMCreateSubDM_Shell(DM dm, PetscInt numFields, const PetscInt fields[], IS *is, DM *subdm)
1258: {
1259: PetscFunctionBegin;
1260: if (subdm) PetscCall(DMShellCreate(PetscObjectComm((PetscObject)dm), subdm));
1261: PetscCall(DMCreateSectionSubDM(dm, numFields, fields, NULL, NULL, is, subdm));
1262: PetscFunctionReturn(PETSC_SUCCESS);
1263: }
1265: PETSC_EXTERN PetscErrorCode DMCreate_Shell(DM dm)
1266: {
1267: DM_Shell *shell;
1269: PetscFunctionBegin;
1270: PetscCall(PetscNew(&shell));
1271: dm->data = shell;
1273: dm->ops->destroy = DMDestroy_Shell;
1274: dm->ops->createglobalvector = DMCreateGlobalVector_Shell;
1275: dm->ops->createlocalvector = DMCreateLocalVector_Shell;
1276: dm->ops->creatematrix = DMCreateMatrix_Shell;
1277: dm->ops->view = DMView_Shell;
1278: dm->ops->load = DMLoad_Shell;
1279: dm->ops->globaltolocalbegin = DMGlobalToLocalBeginDefaultShell;
1280: dm->ops->globaltolocalend = DMGlobalToLocalEndDefaultShell;
1281: dm->ops->localtoglobalbegin = DMLocalToGlobalBeginDefaultShell;
1282: dm->ops->localtoglobalend = DMLocalToGlobalEndDefaultShell;
1283: dm->ops->localtolocalbegin = DMLocalToLocalBeginDefaultShell;
1284: dm->ops->localtolocalend = DMLocalToLocalEndDefaultShell;
1285: dm->ops->createsubdm = DMCreateSubDM_Shell;
1286: PetscCall(DMSetMatType(dm, MATDENSE));
1287: PetscFunctionReturn(PETSC_SUCCESS);
1288: }
1290: /*@
1291: DMShellCreate - Creates a `DMSHELL` object, used to manage user-defined problem data
1293: Collective
1295: Input Parameter:
1296: . comm - the processors that will share the global vector
1298: Output Parameter:
1299: . dm - the `DMSHELL`
1301: Level: advanced
1303: .seealso: `DMDestroy()`, `DMCreateGlobalVector()`, `DMCreateLocalVector()`, `DMShellSetContext()`, `DMShellGetContext()`
1304: @*/
1305: PetscErrorCode DMShellCreate(MPI_Comm comm, DM *dm)
1306: {
1307: PetscFunctionBegin;
1308: PetscAssertPointer(dm, 2);
1309: PetscCall(DMCreate(comm, dm));
1310: PetscCall(DMSetType(*dm, DMSHELL));
1311: PetscCall(DMSetUp(*dm));
1312: PetscFunctionReturn(PETSC_SUCCESS);
1313: }