Actual source code: plexnatural.c
petsc-3.10.5 2019-03-28
1: #include <petsc/private/dmpleximpl.h>
3: /*@
4: DMPlexSetMigrationSF - Sets the SF for migrating from a parent DM into this DM
6: + dm - The DM
7: . naturalSF - The PetscSF
8: Level: intermediate
10: .seealso: DMPlexDistribute(), DMPlexDistributeField(), DMPlexCreateMigrationSF(), DMPlexGetMigrationSF()
11: @*/
12: PetscErrorCode DMPlexSetMigrationSF(DM dm, PetscSF migrationSF)
13: {
16: dm->sfMigration = migrationSF;
17: PetscObjectReference((PetscObject) migrationSF);
18: return(0);
19: }
21: /*@
22: DMPlexGetMigrationSF - Gets the SF for migrating from a parent DM into this DM
24: + dm - The DM
25: . *migrationSF - The PetscSF
26: Level: intermediate
28: .seealso: DMPlexDistribute(), DMPlexDistributeField(), DMPlexCreateMigrationSF(), DMPlexSetMigrationSF
29: @*/
30: PetscErrorCode DMPlexGetMigrationSF(DM dm, PetscSF *migrationSF)
31: {
33: *migrationSF = dm->sfMigration;
34: return(0);
35: }
37: /*@
38: DMPlexSetGlobalToNaturalSF - Sets the SF for mapping Global Vec to the Natural Vec
40: + dm - The DM
41: . naturalSF - The PetscSF
42: Level: intermediate
44: .seealso: DMPlexDistribute(), DMPlexDistributeField(), DMPlexCreateGlobalToNaturalSF(), DMPlexGetGlobaltoNaturalSF()
45: @*/
46: PetscErrorCode DMPlexSetGlobalToNaturalSF(DM dm, PetscSF naturalSF)
47: {
50: dm->sfNatural = naturalSF;
51: PetscObjectReference((PetscObject) naturalSF);
52: dm->useNatural = PETSC_TRUE;
53: return(0);
54: }
56: /*@
57: DMPlexGetGlobalToNaturalSF - Gets the SF for mapping Global Vec to the Natural Vec
59: + dm - The DM
60: . *naturalSF - The PetscSF
61: Level: intermediate
63: .seealso: DMPlexDistribute(), DMPlexDistributeField(), DMPlexCreateGlobalToNaturalSF(), DMPlexSetGlobaltoNaturalSF
64: @*/
65: PetscErrorCode DMPlexGetGlobalToNaturalSF(DM dm, PetscSF *naturalSF)
66: {
68: *naturalSF = dm->sfNatural;
69: return(0);
70: }
72: /*@
73: DMPlexCreateGlobalToNaturalSF - Creates the SF for mapping Global Vec to the Natural Vec
75: Input Parameters:
76: + dm - The DM
77: . section - The PetscSection before the mesh was distributed
78: - sfMigration - The PetscSF used to distribute the mesh
80: Output Parameters:
81: . sfNatural - PetscSF for mapping the Vec in PETSc ordering to the canonical ordering
83: Level: intermediate
85: .seealso: DMPlexDistribute(), DMPlexDistributeField()
86: @*/
87: PetscErrorCode DMPlexCreateGlobalToNaturalSF(DM dm, PetscSection section, PetscSF sfMigration, PetscSF *sfNatural)
88: {
89: MPI_Comm comm;
90: Vec gv;
91: PetscSF sf, sfEmbed, sfSeqToNatural, sfField, sfFieldInv;
92: PetscSection gSection, sectionDist, gLocSection;
93: PetscInt *spoints, *remoteOffsets;
94: PetscInt ssize, pStart, pEnd, p;
98: PetscObjectGetComm((PetscObject) dm, &comm);
99: /* PetscPrintf(comm, "Point migration SF\n");
100: PetscSFView(sfMigration, 0); */
101: /* Create a new section from distributing the original section */
102: PetscSectionCreate(comm, §ionDist);
103: PetscSFDistributeSection(sfMigration, section, &remoteOffsets, sectionDist);
104: /* PetscPrintf(comm, "Distributed Section\n");
105: PetscSectionView(sectionDist, PETSC_VIEWER_STDOUT_WORLD); */
106: DMSetSection(dm, sectionDist);
107: /* Get a pruned version of migration SF */
108: DMGetGlobalSection(dm, &gSection);
109: PetscSectionGetChart(gSection, &pStart, &pEnd);
110: for (p = pStart, ssize = 0; p < pEnd; ++p) {
111: PetscInt dof, off;
113: PetscSectionGetDof(gSection, p, &dof);
114: PetscSectionGetOffset(gSection, p, &off);
115: if ((dof > 0) && (off >= 0)) ++ssize;
116: }
117: PetscMalloc1(ssize, &spoints);
118: for (p = pStart, ssize = 0; p < pEnd; ++p) {
119: PetscInt dof, off;
121: PetscSectionGetDof(gSection, p, &dof);
122: PetscSectionGetOffset(gSection, p, &off);
123: if ((dof > 0) && (off >= 0)) spoints[ssize++] = p;
124: }
125: PetscSFCreateEmbeddedLeafSF(sfMigration, ssize, spoints, &sfEmbed);
126: PetscFree(spoints);
127: /* PetscPrintf(comm, "Embedded SF\n");
128: PetscSFView(sfEmbed, 0); */
129: /* Create the SF for seq to natural */
130: DMGetGlobalVector(dm, &gv);
131: PetscSFCreateFromZero(comm, gv, &sfSeqToNatural);
132: DMRestoreGlobalVector(dm, &gv);
133: /* PetscPrintf(comm, "Seq-to-Natural SF\n");
134: PetscSFView(sfSeqToNatural, 0); */
135: /* Create the SF associated with this section */
136: DMGetPointSF(dm, &sf);
137: PetscSectionCreateGlobalSection(sectionDist, sf, PETSC_FALSE, PETSC_TRUE, &gLocSection);
138: PetscSFCreateSectionSF(sfEmbed, section, remoteOffsets, gLocSection, &sfField);
139: PetscFree(remoteOffsets);
140: PetscSFDestroy(&sfEmbed);
141: PetscSectionDestroy(&gLocSection);
142: PetscSectionDestroy(§ionDist);
143: /* PetscPrintf(comm, "Field SF\n");
144: PetscSFView(sfField, 0); */
145: /* Invert the field SF so it's now from distributed to sequential */
146: PetscSFCreateInverseSF(sfField, &sfFieldInv);
147: PetscSFDestroy(&sfField);
148: /* PetscPrintf(comm, "Inverse Field SF\n");
149: PetscSFView(sfFieldInv, 0); */
150: /* Multiply the sfFieldInv with the */
151: PetscSFCompose(sfFieldInv, sfSeqToNatural, sfNatural);
152: PetscObjectViewFromOptions((PetscObject) *sfNatural, NULL, "-globaltonatural_sf_view");
153: /* Clean up */
154: PetscSFDestroy(&sfFieldInv);
155: PetscSFDestroy(&sfSeqToNatural);
156: return(0);
157: }
159: /*@
160: DMPlexGlobalToNaturalBegin - Rearranges a global Vector in the natural order.
162: Collective on dm
164: Input Parameters:
165: + dm - The distributed DMPlex
166: - gv - The global Vec
168: Output Parameters:
169: . nv - Vec in the canonical ordering distributed over all processors associated with gv
171: Note: The user must call DMPlexSetUseNaturalSF(dm, PETSC_TRUE) before DMPlexDistribute().
173: Level: intermediate
175: .seealso: DMPlexDistribute(), DMPlexDistributeField(), DMPlexNaturalToGlobalBegin(), DMPlexGlobalToNaturalEnd()
176: @*/
177: PetscErrorCode DMPlexGlobalToNaturalBegin(DM dm, Vec gv, Vec nv)
178: {
179: const PetscScalar *inarray;
180: PetscScalar *outarray;
181: PetscErrorCode ierr;
184: PetscLogEventBegin(DMPLEX_GlobalToNaturalBegin,dm,0,0,0);
185: if (dm->sfNatural) {
186: VecGetArray(nv, &outarray);
187: VecGetArrayRead(gv, &inarray);
188: PetscSFBcastBegin(dm->sfNatural, MPIU_SCALAR, (PetscScalar *) inarray, outarray);
189: VecRestoreArrayRead(gv, &inarray);
190: VecRestoreArray(nv, &outarray);
191: } else SETERRQ(PetscObjectComm((PetscObject) dm), PETSC_ERR_ARG_WRONGSTATE, "DM global to natural SF was not created.\nYou must call DMSetUseNatural() before DMPlexDistribute().\n");
192: PetscLogEventEnd(DMPLEX_GlobalToNaturalBegin,dm,0,0,0);
193: return(0);
194: }
196: /*@
197: DMPlexGlobalToNaturalEnd - Rearranges a global Vector in the natural order.
199: Collective on dm
201: Input Parameters:
202: + dm - The distributed DMPlex
203: - gv - The global Vec
205: Output Parameters:
206: . nv - The natural Vec
208: Note: The user must call DMPlexSetUseNaturalSF(dm, PETSC_TRUE) before DMPlexDistribute().
210: Level: intermediate
212: .seealso: DMPlexDistribute(), DMPlexDistributeField(), DMPlexNaturalToGlobalBegin(), DMPlexGlobalToNaturalBegin()
213: @*/
214: PetscErrorCode DMPlexGlobalToNaturalEnd(DM dm, Vec gv, Vec nv)
215: {
216: const PetscScalar *inarray;
217: PetscScalar *outarray;
218: PetscErrorCode ierr;
221: PetscLogEventBegin(DMPLEX_GlobalToNaturalEnd,dm,0,0,0);
222: if (dm->sfNatural) {
223: VecGetArrayRead(gv, &inarray);
224: VecGetArray(nv, &outarray);
225: PetscSFBcastEnd(dm->sfNatural, MPIU_SCALAR, (PetscScalar *) inarray, outarray);
226: VecRestoreArrayRead(gv, &inarray);
227: VecRestoreArray(nv, &outarray);
228: }
229: PetscLogEventEnd(DMPLEX_GlobalToNaturalEnd,dm,0,0,0);
230: return(0);
231: }
233: /*@
234: DMPlexNaturalToGlobalBegin - Rearranges a Vector in the natural order to the Global order.
236: Collective on dm
238: Input Parameters:
239: + dm - The distributed DMPlex
240: - nv - The natural Vec
242: Output Parameters:
243: . gv - The global Vec
245: Note: The user must call DMPlexSetUseNaturalSF(dm, PETSC_TRUE) before DMPlexDistribute().
247: Level: intermediate
249: .seealso: DMPlexDistribute(), DMPlexDistributeField(), DMPlexNaturalToGlobalBegin(),DMPlexGlobalToNaturalEnd()
250: @*/
251: PetscErrorCode DMPlexNaturalToGlobalBegin(DM dm, Vec nv, Vec gv)
252: {
253: const PetscScalar *inarray;
254: PetscScalar *outarray;
255: PetscErrorCode ierr;
258: PetscLogEventBegin(DMPLEX_NaturalToGlobalBegin,dm,0,0,0);
259: if (dm->sfNatural) {
260: /* We only have acces to the SF that goes from Global to Natural.
261: Instead of inverting dm->sfNatural, we can call PetscSFReduceBegin/End with MPI_Op MPI_SUM.
262: Here the SUM really does nothing since sfNatural is one to one, as long as gV is set to zero first. */
263: VecZeroEntries(gv);
264: VecGetArray(gv, &outarray);
265: VecGetArrayRead(nv, &inarray);
266: PetscSFReduceBegin(dm->sfNatural, MPIU_SCALAR, (PetscScalar *) inarray, outarray, MPI_SUM);
267: VecRestoreArrayRead(nv, &inarray);
268: VecRestoreArray(gv, &outarray);
269: } else SETERRQ(PetscObjectComm((PetscObject) dm), PETSC_ERR_ARG_WRONGSTATE, "DM global to natural SF was not created.\nYou must call DMPlexSetUseNaturalSF() before DMPlexDistribute().\n");
270: PetscLogEventEnd(DMPLEX_NaturalToGlobalBegin,dm,0,0,0);
271: return(0);
272: }
274: /*@
275: DMPlexNaturalToGlobalEnd - Rearranges a Vector in the natural order to the Global order.
277: Collective on dm
279: Input Parameters:
280: + dm - The distributed DMPlex
281: - nv - The natural Vec
283: Output Parameters:
284: . gv - The global Vec
286: Note: The user must call DMPlexSetUseNaturalSF(dm, PETSC_TRUE) before DMPlexDistribute().
288: Level: intermediate
290: .seealso: DMPlexDistribute(), DMPlexDistributeField(), DMPlexNaturalToGlobalBegin(), DMPlexGlobalToNaturalBegin()
291: @*/
292: PetscErrorCode DMPlexNaturalToGlobalEnd(DM dm, Vec nv, Vec gv)
293: {
294: const PetscScalar *inarray;
295: PetscScalar *outarray;
296: PetscErrorCode ierr;
299: PetscLogEventBegin(DMPLEX_NaturalToGlobalEnd,dm,0,0,0);
300: if (dm->sfNatural) {
301: VecGetArrayRead(nv, &inarray);
302: VecGetArray(gv, &outarray);
303: PetscSFReduceEnd(dm->sfNatural, MPIU_SCALAR, (PetscScalar *) inarray, outarray, MPI_SUM);
304: VecRestoreArrayRead(nv, &inarray);
305: VecRestoreArray(gv, &outarray);
306: }
307: PetscLogEventEnd(DMPLEX_NaturalToGlobalEnd,dm,0,0,0);
308: return(0);
309: }