1 #include <petsc-private/dmpleximpl.h> /*I "petscdmplex.h" I*/ 2 3 #undef __FUNCT__ 4 #define __FUNCT__ "DMPlexGetPointLocal" 5 /*@ 6 DMPlexGetPointLocal - get location of point data in local Vec 7 8 Not Collective 9 10 Input Arguments: 11 + dm - DM defining the topological space 12 - point - topological point 13 14 Output Arguments: 15 + start - start of point data 16 - end - end of point data 17 18 Level: intermediate 19 20 .seealso: DMGetDefaultSection(), PetscSectionGetOffset(), PetscSectionGetDof(), DMPlexPointLocalRead(), DMPlexPointLocalRead(), DMPlexPointLocalRef() 21 @*/ 22 PetscErrorCode DMPlexGetPointLocal(DM dm,PetscInt point,PetscInt *start,PetscInt *end) 23 { 24 PetscErrorCode ierr; 25 PetscInt offset,dof; 26 27 PetscFunctionBegin; 28 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 29 ierr = PetscSectionGetOffset(dm->defaultSection,point,&offset);CHKERRQ(ierr); 30 ierr = PetscSectionGetDof(dm->defaultSection,point,&dof);CHKERRQ(ierr); 31 if (start) *start = offset; 32 if (end) *end = offset + dof; 33 PetscFunctionReturn(0); 34 } 35 36 #undef __FUNCT__ 37 #define __FUNCT__ "DMPlexPointLocalRead" 38 /*@ 39 DMPlexPointLocalRead - return read access to a point in local array 40 41 Not Collective 42 43 Input Arguments: 44 + dm - DM defining topological space 45 . point - topological point 46 - array - array to index into 47 48 Output Arguments: 49 . ptr - address of read reference to point data, type generic so user can place in structure 50 51 Level: intermediate 52 53 Note: 54 A common usage when data sizes are known statically: 55 56 $ const struct { PetscScalar foo,bar,baz; } *ptr; 57 $ DMPlexPointLocalRead(dm,point,array,&ptr); 58 $ x = 2*ptr->foo + 3*ptr->bar + 5*ptr->baz; 59 60 .seealso: DMGetDefaultSection(), PetscSectionGetOffset(), PetscSectionGetDof(), DMPlexGetPointLocal(), DMPlexPointGlobalRead() 61 @*/ 62 PetscErrorCode DMPlexPointLocalRead(DM dm,PetscInt point,const PetscScalar *array,const void *ptr) 63 { 64 PetscErrorCode ierr; 65 PetscInt start; 66 67 PetscFunctionBegin; 68 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 69 PetscValidScalarPointer(array,3); 70 PetscValidPointer(ptr,4); 71 ierr = DMPlexGetLocalOffset_Private(dm,point,&start);CHKERRQ(ierr); 72 *(const PetscScalar**)ptr = array + start; 73 PetscFunctionReturn(0); 74 } 75 76 #undef __FUNCT__ 77 #define __FUNCT__ "DMPlexPointLocalRef" 78 /*@ 79 DMPlexPointLocalRef - return read/write access to a point in local array 80 81 Not Collective 82 83 Input Arguments: 84 + dm - DM defining topological space 85 . point - topological point 86 - array - array to index into 87 88 Output Arguments: 89 . ptr - address of reference to point data, type generic so user can place in structure 90 91 Level: intermediate 92 93 Note: 94 A common usage when data sizes are known statically: 95 96 $ struct { PetscScalar foo,bar,baz; } *ptr; 97 $ DMPlexPointLocalRef(dm,point,array,&ptr); 98 $ ptr->foo = 2; ptr->bar = 3; ptr->baz = 5; 99 100 .seealso: DMGetDefaultSection(), PetscSectionGetOffset(), PetscSectionGetDof(), DMPlexGetPointLocal(), DMPlexPointGlobalRef() 101 @*/ 102 PetscErrorCode DMPlexPointLocalRef(DM dm,PetscInt point,PetscScalar *array,void *ptr) 103 { 104 PetscErrorCode ierr; 105 PetscInt start; 106 107 PetscFunctionBegin; 108 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 109 PetscValidScalarPointer(array,3); 110 PetscValidPointer(ptr,4); 111 ierr = DMPlexGetLocalOffset_Private(dm,point,&start);CHKERRQ(ierr); 112 *(PetscScalar**)ptr = array + start; 113 PetscFunctionReturn(0); 114 } 115 116 #undef __FUNCT__ 117 #define __FUNCT__ "DMPlexPointLocalFieldRead" 118 /*@ 119 DMPlexPointLocalFieldRead - return read access to a field on a point in local array 120 121 Not Collective 122 123 Input Arguments: 124 + dm - DM defining topological space 125 . point - topological point 126 . field - field number 127 - array - array to index into 128 129 Output Arguments: 130 . ptr - address of read reference to point data, type generic so user can place in structure 131 132 Level: intermediate 133 134 .seealso: DMGetDefaultSection(), PetscSectionGetOffset(), PetscSectionGetDof(), DMPlexGetPointLocal(), DMPlexPointGlobalRef() 135 @*/ 136 PetscErrorCode DMPlexPointLocalFieldRead(DM dm,PetscInt point,PetscInt field,const PetscScalar *array,const void *ptr) 137 { 138 PetscErrorCode ierr; 139 PetscInt start; 140 141 PetscFunctionBegin; 142 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 143 PetscValidScalarPointer(array,3); 144 PetscValidPointer(ptr,4); 145 ierr = DMPlexGetLocalFieldOffset_Private(dm,point,field,&start);CHKERRQ(ierr); 146 *(const PetscScalar**)ptr = array + start; 147 PetscFunctionReturn(0); 148 } 149 150 #undef __FUNCT__ 151 #define __FUNCT__ "DMPlexPointLocalFieldRef" 152 /*@ 153 DMPlexPointLocalFieldRef - return read/write access to a field on a point in local array 154 155 Not Collective 156 157 Input Arguments: 158 + dm - DM defining topological space 159 . point - topological point 160 . field - field number 161 - array - array to index into 162 163 Output Arguments: 164 . ptr - address of reference to point data, type generic so user can place in structure 165 166 Level: intermediate 167 168 .seealso: DMGetDefaultSection(), PetscSectionGetOffset(), PetscSectionGetDof(), DMPlexGetPointLocal(), DMPlexPointGlobalRef() 169 @*/ 170 PetscErrorCode DMPlexPointLocalFieldRef(DM dm,PetscInt point,PetscInt field,PetscScalar *array,void *ptr) 171 { 172 PetscErrorCode ierr; 173 PetscInt start; 174 175 PetscFunctionBegin; 176 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 177 PetscValidScalarPointer(array,3); 178 PetscValidPointer(ptr,4); 179 ierr = DMPlexGetLocalFieldOffset_Private(dm,point,field,&start);CHKERRQ(ierr); 180 *(PetscScalar**)ptr = array + start; 181 PetscFunctionReturn(0); 182 } 183 184 #undef __FUNCT__ 185 #define __FUNCT__ "DMPlexGetPointGlobal" 186 /*@ 187 DMPlexGetPointGlobal - get location of point data in global Vec 188 189 Not Collective 190 191 Input Arguments: 192 + dm - DM defining the topological space 193 - point - topological point 194 195 Output Arguments: 196 + start - start of point data; returns -(global_start+1) if point is not owned 197 - end - end of point data; returns -(global_end+1) if point is not owned 198 199 Level: intermediate 200 201 .seealso: DMGetDefaultSection(), PetscSectionGetOffset(), PetscSectionGetDof(), DMPlexPointGlobalRead(), DMPlexGetPointLocal(), DMPlexPointGlobalRead(), DMPlexPointGlobalRef() 202 @*/ 203 PetscErrorCode DMPlexGetPointGlobal(DM dm,PetscInt point,PetscInt *start,PetscInt *end) 204 { 205 PetscErrorCode ierr; 206 PetscInt offset,dof; 207 208 PetscFunctionBegin; 209 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 210 ierr = PetscSectionGetOffset(dm->defaultGlobalSection,point,&offset);CHKERRQ(ierr); 211 ierr = PetscSectionGetDof(dm->defaultGlobalSection,point,&dof);CHKERRQ(ierr); 212 if (start) *start = offset; 213 if (end) *end = offset + dof; 214 PetscFunctionReturn(0); 215 } 216 217 #undef __FUNCT__ 218 #define __FUNCT__ "DMPlexPointGlobalRead" 219 /*@ 220 DMPlexPointGlobalRead - return read access to a point in global array 221 222 Not Collective 223 224 Input Arguments: 225 + dm - DM defining topological space 226 . point - topological point 227 - array - array to index into 228 229 Output Arguments: 230 . ptr - address of read reference to point data, type generic so user can place in structure; returns NULL if global point is not owned 231 232 Level: intermediate 233 234 Note: 235 A common usage when data sizes are known statically: 236 237 $ const struct { PetscScalar foo,bar,baz; } *ptr; 238 $ DMPlexPointGlobalRead(dm,point,array,&ptr); 239 $ x = 2*ptr->foo + 3*ptr->bar + 5*ptr->baz; 240 241 .seealso: DMGetDefaultSection(), PetscSectionGetOffset(), PetscSectionGetDof(), DMPlexGetPointGlobal(), DMPlexPointLocalRead(), DMPlexPointGlobalRef() 242 @*/ 243 PetscErrorCode DMPlexPointGlobalRead(DM dm,PetscInt point,const PetscScalar *array,const void *ptr) 244 { 245 PetscInt start, end; 246 PetscErrorCode ierr; 247 248 PetscFunctionBegin; 249 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 250 PetscValidScalarPointer(array, 3); 251 PetscValidPointer(ptr, 4); 252 ierr = DMPlexGetGlobalOffset_Private(dm, point, &start, &end);CHKERRQ(ierr); 253 *(const PetscScalar**) ptr = (start < end) ? array + start - dm->map->rstart : NULL; 254 PetscFunctionReturn(0); 255 } 256 257 #undef __FUNCT__ 258 #define __FUNCT__ "DMPlexPointGlobalRef" 259 /*@ 260 DMPlexPointGlobalRef - return read/write access to a point in global array 261 262 Not Collective 263 264 Input Arguments: 265 + dm - DM defining topological space 266 . point - topological point 267 - array - array to index into 268 269 Output Arguments: 270 . ptr - address of reference to point data, type generic so user can place in structure; returns NULL if global point is not owned 271 272 Level: intermediate 273 274 Note: 275 A common usage when data sizes are known statically: 276 277 $ struct { PetscScalar foo,bar,baz; } *ptr; 278 $ DMPlexPointGlobalRef(dm,point,array,&ptr); 279 $ ptr->foo = 2; ptr->bar = 3; ptr->baz = 5; 280 281 .seealso: DMGetDefaultSection(), PetscSectionGetOffset(), PetscSectionGetDof(), DMPlexGetPointGlobal(), DMPlexPointLocalRef(), DMPlexPointGlobalRead() 282 @*/ 283 PetscErrorCode DMPlexPointGlobalRef(DM dm,PetscInt point,PetscScalar *array,void *ptr) 284 { 285 PetscInt start, end; 286 PetscErrorCode ierr; 287 288 PetscFunctionBegin; 289 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 290 PetscValidScalarPointer(array, 3); 291 PetscValidPointer(ptr, 4); 292 ierr = DMPlexGetGlobalOffset_Private(dm, point, &start, &end);CHKERRQ(ierr); 293 *(PetscScalar**) ptr = (start < end) ? array + start - dm->map->rstart : NULL; 294 PetscFunctionReturn(0); 295 } 296 297 #undef __FUNCT__ 298 #define __FUNCT__ "DMPlexPointGlobalFieldRead" 299 /*@ 300 DMPlexPointGlobalFieldRead - return read access to a field on a point in global array 301 302 Not Collective 303 304 Input Arguments: 305 + dm - DM defining topological space 306 . point - topological point 307 . field - field number 308 - array - array to index into 309 310 Output Arguments: 311 . ptr - address of read reference to point data, type generic so user can place in structure; returns NULL if global point is not owned 312 313 Level: intermediate 314 315 .seealso: DMGetDefaultSection(), PetscSectionGetOffset(), PetscSectionGetDof(), DMPlexGetPointGlobal(), DMPlexPointLocalRead(), DMPlexPointGlobalRef() 316 @*/ 317 PetscErrorCode DMPlexPointGlobalFieldRead(DM dm,PetscInt point,PetscInt field,const PetscScalar *array,const void *ptr) 318 { 319 PetscInt start, end; 320 PetscErrorCode ierr; 321 322 PetscFunctionBegin; 323 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 324 PetscValidScalarPointer(array, 3); 325 PetscValidPointer(ptr, 4); 326 ierr = DMPlexGetGlobalFieldOffset_Private(dm, point, field, &start, &end);CHKERRQ(ierr); 327 *(const PetscScalar**) ptr = (start < end) ? array + start - dm->map->rstart : NULL; 328 PetscFunctionReturn(0); 329 } 330 331 #undef __FUNCT__ 332 #define __FUNCT__ "DMPlexPointGlobalFieldRef" 333 /*@ 334 DMPlexPointGlobalFieldRef - return read/write access to a field on a point in global array 335 336 Not Collective 337 338 Input Arguments: 339 + dm - DM defining topological space 340 . point - topological point 341 . field - field number 342 - array - array to index into 343 344 Output Arguments: 345 . ptr - address of reference to point data, type generic so user can place in structure; returns NULL if global point is not owned 346 347 Level: intermediate 348 349 .seealso: DMGetDefaultSection(), PetscSectionGetOffset(), PetscSectionGetDof(), DMPlexGetPointGlobal(), DMPlexPointLocalRef(), DMPlexPointGlobalRead() 350 @*/ 351 PetscErrorCode DMPlexPointGlobalFieldRef(DM dm,PetscInt point,PetscInt field,PetscScalar *array,void *ptr) 352 { 353 PetscInt start, end; 354 PetscErrorCode ierr; 355 356 PetscFunctionBegin; 357 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 358 PetscValidScalarPointer(array, 3); 359 PetscValidPointer(ptr, 4); 360 ierr = DMPlexGetGlobalFieldOffset_Private(dm, point, field, &start, &end);CHKERRQ(ierr); 361 *(PetscScalar**) ptr = (start < end) ? array + start - dm->map->rstart : NULL; 362 PetscFunctionReturn(0); 363 } 364