1 #include <petsc-private/dmpleximpl.h> /*I "petscdmplex.h" I*/ 2 #include <petsc-private/isimpl.h> /* for inline access to atlasOff */ 3 4 #undef __FUNCT__ 5 #define __FUNCT__ "DMPlexGetLocalOffset_Private" 6 PETSC_STATIC_INLINE PetscErrorCode DMPlexGetLocalOffset_Private(DM dm,PetscInt point,PetscInt *offset) 7 { 8 PetscFunctionBegin; 9 #if defined(PETSC_USE_DEBUG) 10 { 11 PetscErrorCode ierr; 12 ierr = PetscSectionGetOffset(dm->defaultSection,point,offset);CHKERRQ(ierr); 13 } 14 #else 15 { 16 PetscSection s = dm->defaultSection; 17 *offset = s->atlasOff[point - s->atlasLayout.pStart]; 18 } 19 #endif 20 PetscFunctionReturn(0); 21 } 22 23 #undef __FUNCT__ 24 #define __FUNCT__ "DMPlexGetGlobalOffset_Private" 25 PETSC_STATIC_INLINE PetscErrorCode DMPlexGetGlobalOffset_Private(DM dm,PetscInt point,PetscInt *offset) 26 { 27 PetscFunctionBegin; 28 #if defined(PETSC_USE_DEBUG) 29 { 30 PetscErrorCode ierr; 31 PetscInt dof,cdof; 32 ierr = PetscSectionGetOffset(dm->defaultGlobalSection,point,offset);CHKERRQ(ierr); 33 ierr = PetscSectionGetDof(dm->defaultGlobalSection,point,&dof);CHKERRQ(ierr); 34 ierr = PetscSectionGetConstraintDof(dm->defaultGlobalSection,point,&cdof);CHKERRQ(ierr); 35 if (dof-cdof <= 0) *offset = -1; /* Indicates no data */ 36 } 37 #else 38 { 39 PetscSection s = dm->defaultGlobalSection; 40 PetscInt dof,cdof; 41 *offset = s->atlasOff[point - s->atlasLayout.pStart]; 42 dof = s->atlasDof[point - s->atlasLayout.pStart]; 43 cdof = s->bc ? s->bc->atlasDof[point - s->bc->atlasLayout.pStart] : 0; 44 if (dof-cdof <= 0) *offset = -1; 45 } 46 #endif 47 PetscFunctionReturn(0); 48 } 49 50 #undef __FUNCT__ 51 #define __FUNCT__ "DMPlexGetPointLocal" 52 /*@ 53 DMPlexGetPointLocal - get location of point data in local Vec 54 55 Not Collective 56 57 Input Arguments: 58 + dm - DM defining the topological space 59 - point - topological point 60 61 Output Arguments: 62 + start - start of point data 63 - end - end of point data 64 65 Level: intermediate 66 67 .seealso: DMGetDefaultSection(), PetscSectionGetOffset(), PetscSectionGetDof(), DMPlexPointLocalRead(), DMPlexPointLocalRead(), DMPlexPointLocalRef() 68 @*/ 69 PetscErrorCode DMPlexGetPointLocal(DM dm,PetscInt point,PetscInt *start,PetscInt *end) 70 { 71 PetscErrorCode ierr; 72 PetscInt offset,dof; 73 74 PetscFunctionBegin; 75 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 76 ierr = PetscSectionGetOffset(dm->defaultSection,point,&offset);CHKERRQ(ierr); 77 ierr = PetscSectionGetDof(dm->defaultSection,point,&dof);CHKERRQ(ierr); 78 if (start) *start = offset; 79 if (end) *end = offset + dof; 80 PetscFunctionReturn(0); 81 } 82 83 #undef __FUNCT__ 84 #define __FUNCT__ "DMPlexPointLocalRead" 85 /*@ 86 DMPlexPointLocalRead - return read access to a point in local array 87 88 Not Collective 89 90 Input Arguments: 91 + dm - DM defining topological space 92 . point - topological point 93 - array - array to index into 94 95 Output Arguments: 96 . ptr - address of read reference to point data, type generic so user can place in structure 97 98 Level: intermediate 99 100 Note: 101 A common usage when data sizes are known statically: 102 103 $ const struct { PetscScalar foo,bar,baz; } *ptr; 104 $ DMPlexPointLocalRead(dm,point,array,&ptr); 105 $ x = 2*ptr->foo + 3*ptr->bar + 5*ptr->baz; 106 107 .seealso: DMGetDefaultSection(), PetscSectionGetOffset(), PetscSectionGetDof(), DMPlexGetPointLocal(), DMPlexPointGlobalRead() 108 @*/ 109 PetscErrorCode DMPlexPointLocalRead(DM dm,PetscInt point,const PetscScalar *array,const void *ptr) 110 { 111 PetscErrorCode ierr; 112 PetscInt start; 113 114 PetscFunctionBegin; 115 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 116 PetscValidScalarPointer(array,3); 117 PetscValidPointer(ptr,4); 118 ierr = DMPlexGetLocalOffset_Private(dm,point,&start);CHKERRQ(ierr); 119 *(const PetscScalar**)ptr = array + start; 120 PetscFunctionReturn(0); 121 } 122 123 #undef __FUNCT__ 124 #define __FUNCT__ "DMPlexPointLocalRef" 125 /*@ 126 DMPlexPointLocalRef - return read/write access to a point in local array 127 128 Not Collective 129 130 Input Arguments: 131 + dm - DM defining topological space 132 . point - topological point 133 - array - array to index into 134 135 Output Arguments: 136 . ptr - address of reference to point data, type generic so user can place in structure 137 138 Level: intermediate 139 140 Note: 141 A common usage when data sizes are known statically: 142 143 $ struct { PetscScalar foo,bar,baz; } *ptr; 144 $ DMPlexPointLocalRef(dm,point,array,&ptr); 145 $ ptr->foo = 2; ptr->bar = 3; ptr->baz = 5; 146 147 .seealso: DMGetDefaultSection(), PetscSectionGetOffset(), PetscSectionGetDof(), DMPlexGetPointLocal(), DMPlexPointGlobalRef() 148 @*/ 149 PetscErrorCode DMPlexPointLocalRef(DM dm,PetscInt point,PetscScalar *array,void *ptr) 150 { 151 PetscErrorCode ierr; 152 PetscInt start; 153 154 PetscFunctionBegin; 155 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 156 PetscValidScalarPointer(array,3); 157 PetscValidPointer(ptr,4); 158 ierr = DMPlexGetLocalOffset_Private(dm,point,&start);CHKERRQ(ierr); 159 *(PetscScalar**)ptr = array + start; 160 PetscFunctionReturn(0); 161 } 162 163 #undef __FUNCT__ 164 #define __FUNCT__ "DMPlexGetPointGlobal" 165 /*@ 166 DMPlexGetPointGlobal - get location of point data in global Vec 167 168 Not Collective 169 170 Input Arguments: 171 + dm - DM defining the topological space 172 - point - topological point 173 174 Output Arguments: 175 + start - start of point data; returns -(global_start+1) if point is not owned 176 - end - end of point data; returns -(global_end+1) if point is not owned 177 178 Level: intermediate 179 180 .seealso: DMGetDefaultSection(), PetscSectionGetOffset(), PetscSectionGetDof(), DMPlexPointGlobalRead(), DMPlexGetPointLocal(), DMPlexPointGlobalRead(), DMPlexPointGlobalRef() 181 @*/ 182 PetscErrorCode DMPlexGetPointGlobal(DM dm,PetscInt point,PetscInt *start,PetscInt *end) 183 { 184 PetscErrorCode ierr; 185 PetscInt offset,dof; 186 187 PetscFunctionBegin; 188 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 189 ierr = PetscSectionGetOffset(dm->defaultGlobalSection,point,&offset);CHKERRQ(ierr); 190 ierr = PetscSectionGetDof(dm->defaultGlobalSection,point,&dof);CHKERRQ(ierr); 191 if (start) *start = offset; 192 if (end) *end = offset + dof; 193 PetscFunctionReturn(0); 194 } 195 196 #undef __FUNCT__ 197 #define __FUNCT__ "DMPlexPointGlobalRead" 198 /*@ 199 DMPlexPointGlobalRead - return read access to a point in global array 200 201 Not Collective 202 203 Input Arguments: 204 + dm - DM defining topological space 205 . point - topological point 206 - array - array to index into 207 208 Output Arguments: 209 . ptr - address of read reference to point data, type generic so user can place in structure; returns NULL if global point is not owned 210 211 Level: intermediate 212 213 Note: 214 A common usage when data sizes are known statically: 215 216 $ const struct { PetscScalar foo,bar,baz; } *ptr; 217 $ DMPlexPointGlobalRead(dm,point,array,&ptr); 218 $ x = 2*ptr->foo + 3*ptr->bar + 5*ptr->baz; 219 220 .seealso: DMGetDefaultSection(), PetscSectionGetOffset(), PetscSectionGetDof(), DMPlexGetPointGlobal(), DMPlexPointLocalRead(), DMPlexPointGlobalRef() 221 @*/ 222 PetscErrorCode DMPlexPointGlobalRead(DM dm,PetscInt point,const PetscScalar *array,const void *ptr) 223 { 224 PetscErrorCode ierr; 225 PetscInt start; 226 227 PetscFunctionBegin; 228 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 229 PetscValidScalarPointer(array,3); 230 PetscValidPointer(ptr,4); 231 ierr = DMPlexGetGlobalOffset_Private(dm,point,&start);CHKERRQ(ierr); 232 *(const PetscScalar**)ptr = (start >= 0) ? array + start - dm->map->rstart : NULL; 233 PetscFunctionReturn(0); 234 } 235 236 #undef __FUNCT__ 237 #define __FUNCT__ "DMPlexPointGlobalRef" 238 /*@ 239 DMPlexPointGlobalRef - return read/write access to a point in global array 240 241 Not Collective 242 243 Input Arguments: 244 + dm - DM defining topological space 245 . point - topological point 246 - array - array to index into 247 248 Output Arguments: 249 . ptr - address of reference to point data, type generic so user can place in structure; returns NULL if global point is not owned 250 251 Level: intermediate 252 253 Note: 254 A common usage when data sizes are known statically: 255 256 $ struct { PetscScalar foo,bar,baz; } *ptr; 257 $ DMPlexPointGlobalRef(dm,point,array,&ptr); 258 $ ptr->foo = 2; ptr->bar = 3; ptr->baz = 5; 259 260 .seealso: DMGetDefaultSection(), PetscSectionGetOffset(), PetscSectionGetDof(), DMPlexGetPointGlobal(), DMPlexPointLocalRef(), DMPlexPointGlobalRead() 261 @*/ 262 PetscErrorCode DMPlexPointGlobalRef(DM dm,PetscInt point,PetscScalar *array,void *ptr) 263 { 264 PetscErrorCode ierr; 265 PetscInt start; 266 267 PetscFunctionBegin; 268 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 269 PetscValidScalarPointer(array,3); 270 PetscValidPointer(ptr,4); 271 ierr = DMPlexGetGlobalOffset_Private(dm,point,&start);CHKERRQ(ierr); 272 *(PetscScalar**)ptr = (start >= 0) ? array + start - dm->map->rstart : NULL; 273 PetscFunctionReturn(0); 274 } 275