1 #include <petsc/private/dmpatchimpl.h> /*I "petscdmpatch.h" I*/
2 #include <petscdmda.h>
3
DMSetFromOptions_Patch(DM dm,PetscOptionItems PetscOptionsObject)4 static PetscErrorCode DMSetFromOptions_Patch(DM dm, PetscOptionItems PetscOptionsObject)
5 {
6 /* DM_Patch *mesh = (DM_Patch*) dm->data; */
7
8 PetscFunctionBegin;
9 PetscOptionsHeadBegin(PetscOptionsObject, "DMPatch Options");
10 /* Handle associated vectors */
11 /* Handle viewing */
12 PetscOptionsHeadEnd();
13 PetscFunctionReturn(PETSC_SUCCESS);
14 }
15
16 /* External function declarations here */
17 extern PetscErrorCode DMSetUp_Patch(DM dm);
18 extern PetscErrorCode DMView_Patch(DM dm, PetscViewer viewer);
19 extern PetscErrorCode DMCreateGlobalVector_Patch(DM dm, Vec *g);
20 extern PetscErrorCode DMCreateLocalVector_Patch(DM dm, Vec *l);
21 extern PetscErrorCode DMDestroy_Patch(DM dm);
22 extern PetscErrorCode DMCreateSubDM_Patch(DM dm, PetscInt numFields, const PetscInt fields[], IS *is, DM *subdm);
23
DMInitialize_Patch(DM dm)24 static PetscErrorCode DMInitialize_Patch(DM dm)
25 {
26 PetscFunctionBegin;
27 dm->ops->view = DMView_Patch;
28 dm->ops->setfromoptions = DMSetFromOptions_Patch;
29 dm->ops->setup = DMSetUp_Patch;
30 dm->ops->createglobalvector = DMCreateGlobalVector_Patch;
31 dm->ops->createlocalvector = DMCreateLocalVector_Patch;
32 dm->ops->getlocaltoglobalmapping = NULL;
33 dm->ops->createfieldis = NULL;
34 dm->ops->getcoloring = NULL;
35 dm->ops->creatematrix = NULL;
36 dm->ops->createinterpolation = NULL;
37 dm->ops->createinjection = NULL;
38 dm->ops->refine = NULL;
39 dm->ops->coarsen = NULL;
40 dm->ops->refinehierarchy = NULL;
41 dm->ops->coarsenhierarchy = NULL;
42 dm->ops->globaltolocalbegin = NULL;
43 dm->ops->globaltolocalend = NULL;
44 dm->ops->localtoglobalbegin = NULL;
45 dm->ops->localtoglobalend = NULL;
46 dm->ops->destroy = DMDestroy_Patch;
47 dm->ops->createsubdm = DMCreateSubDM_Patch;
48 PetscFunctionReturn(PETSC_SUCCESS);
49 }
50
DMCreate_Patch(DM dm)51 PETSC_EXTERN PetscErrorCode DMCreate_Patch(DM dm)
52 {
53 DM_Patch *mesh;
54
55 PetscFunctionBegin;
56 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
57 PetscCall(PetscNew(&mesh));
58 dm->data = mesh;
59
60 mesh->refct = 1;
61 mesh->dmCoarse = NULL;
62 mesh->patchSize.i = 0;
63 mesh->patchSize.j = 0;
64 mesh->patchSize.k = 0;
65 mesh->patchSize.c = 0;
66
67 PetscCall(DMInitialize_Patch(dm));
68 PetscFunctionReturn(PETSC_SUCCESS);
69 }
70
71 /*@
72 DMPatchCreate - Creates a DMPatch object, which is a collections of DMs called patches.
73
74 Collective
75
76 Input Parameter:
77 . comm - The communicator for the DMPatch object
78
79 Output Parameter:
80 . mesh - The DMPatch object
81
82 Notes:
83
84 This code is incomplete and not used by other parts of PETSc.
85
86 Level: beginner
87
88 .seealso: `DMPatchZoom()`
89
90 @*/
DMPatchCreate(MPI_Comm comm,DM * mesh)91 PetscErrorCode DMPatchCreate(MPI_Comm comm, DM *mesh)
92 {
93 PetscFunctionBegin;
94 PetscAssertPointer(mesh, 2);
95 PetscCall(DMCreate(comm, mesh));
96 PetscCall(DMSetType(*mesh, DMPATCH));
97 PetscFunctionReturn(PETSC_SUCCESS);
98 }
99
DMPatchCreateGrid(MPI_Comm comm,PetscInt dim,MatStencil patchSize,MatStencil commSize,MatStencil gridSize,DM * dm)100 PetscErrorCode DMPatchCreateGrid(MPI_Comm comm, PetscInt dim, MatStencil patchSize, MatStencil commSize, MatStencil gridSize, DM *dm)
101 {
102 DM_Patch *mesh;
103 DM da;
104 PetscInt dof = 1, width = 1;
105
106 PetscFunctionBegin;
107 PetscCall(DMPatchCreate(comm, dm));
108 mesh = (DM_Patch *)(*dm)->data;
109 if (dim < 2) {
110 gridSize.j = 1;
111 patchSize.j = 1;
112 }
113 if (dim < 3) {
114 gridSize.k = 1;
115 patchSize.k = 1;
116 }
117 PetscCall(DMCreate(comm, &da));
118 PetscCall(DMSetType(da, DMDA));
119 PetscCall(DMSetDimension(da, dim));
120 PetscCall(DMDASetSizes(da, gridSize.i, gridSize.j, gridSize.k));
121 PetscCall(DMDASetBoundaryType(da, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE));
122 PetscCall(DMDASetDof(da, dof));
123 PetscCall(DMDASetStencilType(da, DMDA_STENCIL_BOX));
124 PetscCall(DMDASetStencilWidth(da, width));
125
126 mesh->dmCoarse = da;
127
128 PetscCall(DMPatchSetPatchSize(*dm, patchSize));
129 PetscCall(DMPatchSetCommSize(*dm, commSize));
130 PetscFunctionReturn(PETSC_SUCCESS);
131 }
132