1 #include <petsc/private/dmpleximpl.h> /*I "petscdmplex.h" I*/ 2 3 #undef __FUNCT__ 4 #define __FUNCT__ "DMPlexCreateGlobalToNaturalSF" 5 /*@ 6 DMPlexCreateGlobalToNaturalSF - Creates the SF for mapping Global Vec to the Natural Vec 7 8 Input Parameters: 9 + dm - The DM 10 . section - The PetscSection before the mesh was distributed 11 - sfMigration - The PetscSF used to distribute the mesh 12 13 Output Parameters: 14 . sfNatural - PetscSF for mapping the Vec in PETSc ordering to the canonical ordering 15 16 Level: intermediate 17 18 .seealso: DMPlexDistribute(), DMPlexDistributeField() 19 @*/ 20 PetscErrorCode DMPlexCreateGlobalToNaturalSF(DM dm, PetscSection section, PetscSF sfMigration, PetscSF *sfNatural) 21 { 22 MPI_Comm comm; 23 Vec gv; 24 PetscSF sfEmbed, sfSeqToNatural, sfField, sfFieldInv; 25 PetscSection gSection, newSection; 26 PetscInt *spoints, *remoteOffsets; 27 PetscInt ssize, pStart, pEnd, p; 28 PetscErrorCode ierr; 29 30 PetscFunctionBegin; 31 ierr = PetscObjectGetComm((PetscObject) dm, &comm);CHKERRQ(ierr); 32 /* Get a pruned version of migration SF */ 33 ierr = DMGetDefaultGlobalSection(dm, &gSection);CHKERRQ(ierr); 34 ierr = PetscSectionGetChart(gSection, &pStart, &pEnd);CHKERRQ(ierr); 35 for (p = pStart, ssize = 0; p < pEnd; ++p) { 36 PetscInt dof, off; 37 38 ierr = PetscSectionGetDof(gSection, p, &dof);CHKERRQ(ierr); 39 ierr = PetscSectionGetOffset(gSection, p, &off);CHKERRQ(ierr); 40 if ((dof > 0) && (off >= 0)) ++ssize; 41 } 42 ierr = PetscMalloc1(ssize, &spoints);CHKERRQ(ierr); 43 for (p = pStart, ssize = 0; p < pEnd; ++p) { 44 PetscInt dof, off; 45 46 ierr = PetscSectionGetDof(gSection, p, &dof);CHKERRQ(ierr); 47 ierr = PetscSectionGetOffset(gSection, p, &off);CHKERRQ(ierr); 48 if ((dof > 0) && (off >= 0)) spoints[ssize++] = p; 49 } 50 ierr = PetscSFCreateEmbeddedLeafSF(sfMigration, ssize, spoints, &sfEmbed);CHKERRQ(ierr); 51 /* Create a new section from distributing the original section */ 52 ierr = PetscSectionCreate(comm, &newSection);CHKERRQ(ierr); 53 ierr = PetscSFDistributeSection(sfEmbed, section, &remoteOffsets, newSection);CHKERRQ(ierr); 54 ierr = PetscSFDestroy(&sfEmbed);CHKERRQ(ierr); 55 /* Create the SF for seq to natural */ 56 ierr = DMGetGlobalVector(dm, &gv);CHKERRQ(ierr); 57 ierr = PetscSFCreateFromZero(comm, gv, &sfSeqToNatural);CHKERRQ(ierr); 58 ierr = DMRestoreGlobalVector(dm, &gv);CHKERRQ(ierr); 59 /* Create the SF associated with this section */ 60 ierr = PetscSFCreateSectionSF(sfEmbed, section, remoteOffsets, newSection, &sfField);CHKERRQ(ierr); 61 ierr = PetscSFDestroy(&sfEmbed);CHKERRQ(ierr); 62 /* Invert the field SF so it's now from distributed to sequential */ 63 ierr = PetscSFCreateInverseSF(sfField, &sfFieldInv);CHKERRQ(ierr); 64 ierr = PetscSFDestroy(&sfField);CHKERRQ(ierr); 65 /* Multiply the sfFieldInv with the */ 66 ierr = PetscSFCompose(sfFieldInv, sfSeqToNatural, sfNatural);CHKERRQ(ierr); 67 /* Clean up */ 68 ierr = PetscSectionDestroy(&newSection);CHKERRQ(ierr); 69 ierr = PetscSFDestroy(&sfFieldInv);CHKERRQ(ierr); 70 ierr = PetscSFDestroy(&sfSeqToNatural);CHKERRQ(ierr); 71 PetscFunctionReturn(0); 72 } 73 74 #undef __FUNCT__ 75 #define __FUNCT__ "DMPlexGlobalToNaturalBegin" 76 /*@ 77 DMPlexGlobalToNaturalBegin - Rearranges a global Vector in the natural order. 78 79 Collective on dm 80 81 Input Parameters: 82 + dm - The distributed DMPlex 83 - gv - The global Vec 84 85 Output Parameters: 86 . nv - Vec in the canonical ordering distributed over all processors associated with gv 87 88 Note: The user must call DMPlexSetUseNaturalSF(dm, PETSC_TRUE) before DMPlexDistribute(). 89 90 Level: intermediate 91 92 .seealso: DMPlexDistribute(), DMPlexDistributeField(), DMPlexNaturalToGlobalBegin(), DMPlexGlobalToNaturalEnd() 93 @*/ 94 PetscErrorCode DMPlexGlobalToNaturalBegin(DM dm, Vec gv, Vec nv) 95 { 96 const PetscScalar *inarray; 97 PetscScalar *outarray; 98 PetscErrorCode ierr; 99 100 PetscFunctionBegin; 101 ierr = PetscLogEventBegin(DMPLEX_GlobalToNaturalBegin,dm,0,0,0);CHKERRQ(ierr); 102 if (dm->sfNatural) { 103 ierr = VecGetArray(nv, &outarray);CHKERRQ(ierr); 104 ierr = VecGetArrayRead(gv, &inarray);CHKERRQ(ierr); 105 ierr = PetscSFBcastBegin(dm->sfNatural, MPIU_SCALAR, (PetscScalar *) inarray, outarray);CHKERRQ(ierr); 106 ierr = VecRestoreArrayRead(gv, &inarray);CHKERRQ(ierr); 107 ierr = VecRestoreArray(nv, &outarray);CHKERRQ(ierr); 108 } else SETERRQ(PetscObjectComm((PetscObject) dm), PETSC_ERR_ARG_WRONGSTATE, "DM global to natural SF was not created.\nYou must call DMPlexSetUseNaturalSF() before DMPlexDistribute().\n"); 109 ierr = PetscLogEventEnd(DMPLEX_GlobalToNaturalBegin,dm,0,0,0);CHKERRQ(ierr); 110 PetscFunctionReturn(0); 111 } 112 113 #undef __FUNCT__ 114 #define __FUNCT__ "DMPlexGlobalToNaturalEnd" 115 /*@ 116 DMPlexGlobalToNaturalEnd - Rearranges a global Vector in the natural order. 117 118 Collective on dm 119 120 Input Parameters: 121 + dm - The distributed DMPlex 122 - gv - The global Vec 123 124 Output Parameters: 125 . nv - The natural Vec 126 127 Note: The user must call DMPlexSetUseNaturalSF(dm, PETSC_TRUE) before DMPlexDistribute(). 128 129 Level: intermediate 130 131 .seealso: DMPlexDistribute(), DMPlexDistributeField(), DMPlexNaturalToGlobalBegin(), DMPlexGlobalToNaturalBegin() 132 @*/ 133 PetscErrorCode DMPlexGlobalToNaturalEnd(DM dm, Vec gv, Vec nv) 134 { 135 const PetscScalar *inarray; 136 PetscScalar *outarray; 137 PetscErrorCode ierr; 138 139 PetscFunctionBegin; 140 ierr = PetscLogEventBegin(DMPLEX_GlobalToNaturalEnd,dm,0,0,0);CHKERRQ(ierr); 141 if (dm->sfNatural) { 142 ierr = VecGetArrayRead(gv, &inarray);CHKERRQ(ierr); 143 ierr = VecGetArray(nv, &outarray);CHKERRQ(ierr);CHKERRQ(ierr); 144 ierr = PetscSFBcastEnd(dm->sfNatural, MPIU_SCALAR, (PetscScalar *) inarray, outarray);CHKERRQ(ierr); 145 ierr = VecRestoreArrayRead(gv, &inarray);CHKERRQ(ierr); 146 ierr = VecRestoreArray(nv, &outarray);CHKERRQ(ierr); 147 } 148 ierr = PetscLogEventEnd(DMPLEX_GlobalToNaturalEnd,dm,0,0,0);CHKERRQ(ierr); 149 PetscFunctionReturn(0); 150 } 151 152 #undef __FUNCT__ 153 #define __FUNCT__ "DMPlexNaturalToGlobalBegin" 154 /*@ 155 DMPlexNaturalToGlobalBegin - Rearranges a Vector in the natural order to the Global order. 156 157 Collective on dm 158 159 Input Parameters: 160 + dm - The distributed DMPlex 161 - nv - The natural Vec 162 163 Output Parameters: 164 . gv - The global Vec 165 166 Note: The user must call DMPlexSetUseNaturalSF(dm, PETSC_TRUE) before DMPlexDistribute(). 167 168 Level: intermediate 169 170 .seealso: DMPlexDistribute(), DMPlexDistributeField(), DMPlexNaturalToGlobalBegin(),DMPlexGlobalToNaturalEnd() 171 @*/ 172 PetscErrorCode DMPlexNaturalToGlobalBegin(DM dm, Vec nv, Vec gv) 173 { 174 const PetscScalar *inarray; 175 PetscScalar *outarray; 176 PetscErrorCode ierr; 177 178 PetscFunctionBegin; 179 ierr = PetscLogEventBegin(DMPLEX_NaturalToGlobalBegin,dm,0,0,0);CHKERRQ(ierr); 180 if (dm->sfNatural) { 181 ierr = VecGetArray(gv, &outarray);CHKERRQ(ierr); 182 ierr = VecGetArrayRead(nv, &inarray);CHKERRQ(ierr); 183 ierr = PetscSFReduceBegin(dm->sfNatural, MPIU_SCALAR, (PetscScalar *) inarray, outarray, MPI_SUM);CHKERRQ(ierr); 184 ierr = VecRestoreArrayRead(nv, &inarray);CHKERRQ(ierr); 185 ierr = VecRestoreArray(gv, &outarray);CHKERRQ(ierr); 186 } else SETERRQ(PetscObjectComm((PetscObject) dm), PETSC_ERR_ARG_WRONGSTATE, "DM global to natural SF was not created.\nYou must call DMPlexSetUseNaturalSF() before DMPlexDistribute().\n"); 187 ierr = PetscLogEventEnd(DMPLEX_NaturalToGlobalBegin,dm,0,0,0);CHKERRQ(ierr); 188 PetscFunctionReturn(0); 189 } 190 191 #undef __FUNCT__ 192 #define __FUNCT__ "DMPlexNaturalToGlobalEnd" 193 /*@ 194 DMPlexNaturalToGlobalEnd - Rearranges a Vector in the natural order to the Global order. 195 196 Collective on dm 197 198 Input Parameters: 199 + dm - The distributed DMPlex 200 - nv - The natural Vec 201 202 Output Parameters: 203 . gv - The global Vec 204 205 Note: The user must call DMPlexSetUseNaturalSF(dm, PETSC_TRUE) before DMPlexDistribute(). 206 207 Level: intermediate 208 209 .seealso: DMPlexDistribute(), DMPlexDistributeField(), DMPlexNaturalToGlobalBegin(), DMPlexGlobalToNaturalBegin() 210 @*/ 211 PetscErrorCode DMPlexNaturalToGlobalEnd(DM dm, Vec nv, Vec gv) 212 { 213 const PetscScalar *inarray; 214 PetscScalar *outarray; 215 PetscErrorCode ierr; 216 217 PetscFunctionBegin; 218 ierr = PetscLogEventBegin(DMPLEX_NaturalToGlobalEnd,dm,0,0,0);CHKERRQ(ierr); 219 if (dm->sfNatural) { 220 ierr = VecGetArrayRead(nv, &inarray);CHKERRQ(ierr); 221 ierr = VecGetArray(gv, &outarray);CHKERRQ(ierr);CHKERRQ(ierr); 222 ierr = PetscSFReduceEnd(dm->sfNatural, MPIU_SCALAR, (PetscScalar *) inarray, outarray, MPI_SUM);CHKERRQ(ierr); 223 ierr = VecRestoreArrayRead(nv, &inarray);CHKERRQ(ierr); 224 ierr = VecRestoreArray(gv, &outarray);CHKERRQ(ierr); 225 } 226 ierr = PetscLogEventEnd(DMPLEX_NaturalToGlobalEnd,dm,0,0,0);CHKERRQ(ierr); 227 PetscFunctionReturn(0); 228 } 229