1 #include <petsc-private/dmmbimpl.h> /*I "petscdm.h" I*/ 2 3 #include <petscdmmoab.h> 4 5 #undef __FUNCT__ 6 #define __FUNCT__ "DMMoabSetFieldVector" 7 PetscErrorCode DMMoabSetFieldVector(DM dm, PetscInt ifield, Vec fvec) 8 { 9 DM_Moab *dmmoab; 10 moab::Tag vtag,ntag; 11 const PetscScalar *varray; 12 PetscScalar *farray; 13 moab::ErrorCode merr; 14 PetscErrorCode ierr; 15 std::string tag_name; 16 17 PetscFunctionBegin; 18 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 19 dmmoab = (DM_Moab*)(dm)->data; 20 21 /* Create a tag in MOAB mesh to index and keep track of number of Petsc vec tags */ 22 merr = dmmoab->mbiface->tag_get_handle(dmmoab->fields[ifield],1,moab::MB_TYPE_DOUBLE,ntag, 23 moab::MB_TAG_DENSE|moab::MB_TAG_CREAT);MBERRNM(merr); 24 25 ierr = DMMoabGetVecTag(fvec,&vtag);CHKERRQ(ierr); 26 27 merr = dmmoab->mbiface->tag_get_name(vtag, tag_name); 28 if (!tag_name.length() && merr !=moab::MB_SUCCESS) { 29 ierr = VecGetArrayRead(fvec,&varray);CHKERRQ(ierr); 30 /* use the entity handle and the Dof index to set the right value */ 31 merr = dmmoab->mbiface->tag_set_data(ntag, *dmmoab->vowned, (const void*)varray);MBERRNM(merr); 32 ierr = VecRestoreArrayRead(fvec,&varray);CHKERRQ(ierr); 33 } 34 else { 35 ierr = PetscMalloc(dmmoab->nloc*sizeof(PetscScalar),&farray);CHKERRQ(ierr); 36 /* we are using a MOAB Vec - directly copy the tag data to new one */ 37 merr = dmmoab->mbiface->tag_get_data(vtag, *dmmoab->vowned, (void*)farray);MBERRNM(merr); 38 merr = dmmoab->mbiface->tag_set_data(ntag, *dmmoab->vowned, (const void*)farray);MBERRNM(merr); 39 /* make sure the parallel exchange for ghosts are done appropriately */ 40 ierr = PetscFree(farray);CHKERRQ(ierr); 41 } 42 merr = dmmoab->pcomm->exchange_tags(ntag, *dmmoab->vowned);MBERRNM(merr); 43 PetscFunctionReturn(0); 44 } 45 46 47 #undef __FUNCT__ 48 #define __FUNCT__ "DMMoabSetGlobalFieldVector" 49 PetscErrorCode DMMoabSetGlobalFieldVector(DM dm, Vec fvec) 50 { 51 DM_Moab *dmmoab; 52 moab::Tag vtag,ntag; 53 const PetscScalar *varray; 54 PetscScalar *farray; 55 moab::ErrorCode merr; 56 PetscErrorCode ierr; 57 PetscInt i,ifield; 58 std::string tag_name; 59 moab::Range::iterator iter; 60 61 PetscFunctionBegin; 62 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 63 dmmoab = (DM_Moab*)(dm)->data; 64 65 /* get the Tag corresponding to the global vector - possible that there is no tag associated.. */ 66 ierr = DMMoabGetVecTag(fvec,&vtag);CHKERRQ(ierr); 67 merr = dmmoab->mbiface->tag_get_name(vtag, tag_name); 68 ierr = PetscMalloc(dmmoab->nloc*sizeof(PetscScalar),&farray);CHKERRQ(ierr); 69 if (!tag_name.length() && merr !=moab::MB_SUCCESS) { 70 /* not a MOAB vector - use VecGetSubVector to get the parts as needed */ 71 72 ierr = VecGetArrayRead(fvec,&varray);CHKERRQ(ierr); 73 for (ifield=0; ifield<dmmoab->nfields; ++ifield) { 74 75 /* Create a tag in MOAB mesh to index and keep track of number of Petsc vec tags */ 76 merr = dmmoab->mbiface->tag_get_handle(dmmoab->fields[ifield],1,moab::MB_TYPE_DOUBLE,ntag, 77 moab::MB_TAG_DENSE|moab::MB_TAG_CREAT);MBERRNM(merr); 78 79 for(i=0;i<dmmoab->nloc;i++) { 80 if (dmmoab->bs == 1) 81 farray[i]=varray[ifield*dmmoab->nloc+i]; 82 else 83 farray[i]=varray[i*dmmoab->nfields+ifield]; 84 } 85 86 /* use the entity handle and the Dof index to set the right value */ 87 merr = dmmoab->mbiface->tag_set_data(ntag, *dmmoab->vowned, (const void*)farray);MBERRNM(merr); 88 } 89 ierr = VecRestoreArrayRead(fvec,&varray);CHKERRQ(ierr); 90 } 91 else { 92 ierr = PetscMalloc(dmmoab->nloc*dmmoab->bs*sizeof(PetscScalar),&varray);CHKERRQ(ierr); 93 94 /* we are using a MOAB Vec - directly copy the tag data to new one */ 95 merr = dmmoab->mbiface->tag_get_data(vtag, *dmmoab->vowned, (void*)varray);MBERRNM(merr); 96 for (ifield=0; ifield<dmmoab->nfields; ++ifield) { 97 98 /* Create a tag in MOAB mesh to index and keep track of number of Petsc vec tags */ 99 merr = dmmoab->mbiface->tag_get_handle(dmmoab->fields[ifield],1,moab::MB_TYPE_DOUBLE,ntag, 100 moab::MB_TAG_DENSE|moab::MB_TAG_CREAT);MBERRNM(merr); 101 102 /* we are using a MOAB Vec - directly copy the tag data to new one */ 103 for(i=0; i < dmmoab->nloc; i++) { 104 farray[i] = varray[i*dmmoab->bs+ifield]; 105 } 106 107 merr = dmmoab->mbiface->tag_set_data(ntag, *dmmoab->vowned, (const void*)farray);MBERRNM(merr); 108 /* make sure the parallel exchange for ghosts are done appropriately */ 109 merr = dmmoab->pcomm->exchange_tags(ntag, *dmmoab->vlocal);MBERRNM(merr); 110 } 111 ierr = PetscFree(varray);CHKERRQ(ierr); 112 } 113 ierr = PetscFree(farray);CHKERRQ(ierr); 114 PetscFunctionReturn(0); 115 } 116 117 118 #undef __FUNCT__ 119 #define __FUNCT__ "DMMoabSetFields" 120 PetscErrorCode DMMoabSetFields(DM dm,PetscInt nfields,const char** fields) 121 { 122 DM_Moab *dmmoab; 123 124 PetscFunctionBegin; 125 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 126 dmmoab = (DM_Moab*)(dm)->data; 127 128 dmmoab->fields = fields; 129 dmmoab->nfields = nfields; 130 PetscFunctionReturn(0); 131 } 132 133 134 #undef __FUNCT__ 135 #define __FUNCT__ "DMMoabGetFieldDof" 136 PetscErrorCode DMMoabGetFieldDof(DM dm,moab::EntityHandle point,PetscInt field,PetscInt* dof) 137 { 138 DM_Moab *dmmoab; 139 140 PetscFunctionBegin; 141 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 142 dmmoab = (DM_Moab*)(dm)->data; 143 144 *dof=dmmoab->gidmap[(PetscInt)point]; 145 PetscFunctionReturn(0); 146 } 147 148 149 #undef __FUNCT__ 150 #define __FUNCT__ "DMMoabGetFieldDofs" 151 PetscErrorCode DMMoabGetFieldDofs(DM dm,PetscInt npoints,const moab::EntityHandle* points,PetscInt field,PetscInt* dof) 152 { 153 PetscInt i; 154 PetscErrorCode ierr; 155 DM_Moab *dmmoab; 156 157 PetscFunctionBegin; 158 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 159 PetscValidPointer(points,2); 160 dmmoab = (DM_Moab*)(dm)->data; 161 162 if (!dof) { 163 ierr = PetscMalloc(sizeof(PetscInt)*npoints, &dof);CHKERRQ(ierr); 164 } 165 166 /* first get the local indices */ 167 if (dmmoab->bs > 1) { 168 for (i=0; i<npoints; ++i) 169 dof[i] = dmmoab->gidmap[(PetscInt)points[i]]*dmmoab->nfields+field; 170 } 171 else { 172 /* assume all fields have equal distribution */ 173 for (i=0; i<npoints; ++i) 174 dof[i] = dmmoab->gidmap[(PetscInt)points[i]]+field*dmmoab->n; 175 } 176 PetscFunctionReturn(0); 177 } 178 179 180 #undef __FUNCT__ 181 #define __FUNCT__ "DMMoabGetFieldDofsLocal" 182 PetscErrorCode DMMoabGetFieldDofsLocal(DM dm,PetscInt npoints,const moab::EntityHandle* points,PetscInt field,PetscInt* dof) 183 { 184 PetscInt i,offset; 185 PetscErrorCode ierr; 186 DM_Moab *dmmoab; 187 188 PetscFunctionBegin; 189 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 190 PetscValidPointer(points,2); 191 dmmoab = (DM_Moab*)(dm)->data; 192 193 if (!dof) { 194 ierr = PetscMalloc(sizeof(PetscInt)*npoints, &dof);CHKERRQ(ierr); 195 } 196 197 if (dmmoab->bs > 1) { 198 for (i=0; i<npoints; ++i) 199 dof[i] = dmmoab->lidmap[(PetscInt)points[i]]*dmmoab->nfields+field; 200 } 201 else { 202 offset = field*dmmoab->n; /* assume all fields have equal distribution */ 203 for (i=0; i<npoints; ++i) 204 dof[i] = dmmoab->lidmap[(PetscInt)points[i]]+offset; 205 } 206 PetscFunctionReturn(0); 207 } 208 209 210 #undef __FUNCT__ 211 #define __FUNCT__ "DMMoabGetDofs" 212 PetscErrorCode DMMoabGetDofs(DM dm,PetscInt npoints,const moab::EntityHandle* points,PetscInt* dof) 213 { 214 PetscInt i,field,offset; 215 PetscErrorCode ierr; 216 DM_Moab *dmmoab; 217 218 PetscFunctionBegin; 219 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 220 PetscValidPointer(points,2); 221 dmmoab = (DM_Moab*)(dm)->data; 222 223 if (!dof) { 224 ierr = PetscMalloc(sizeof(PetscInt)*dmmoab->nfields*npoints, &dof);CHKERRQ(ierr); 225 } 226 227 if (dmmoab->bs > 1) { 228 for (field=0; field<dmmoab->nfields; ++field) { 229 for (i=0; i<npoints; ++i) 230 dof[i*dmmoab->nfields+field] = dmmoab->gidmap[(PetscInt)points[i]]*dmmoab->nfields+field; 231 } 232 } 233 else { 234 for (field=0; field<dmmoab->nfields; ++field) { 235 offset = field*dmmoab->n; /* assume all fields have equal distribution */ 236 for (i=0; i<npoints; ++i) 237 dof[i*dmmoab->nfields+field] = dmmoab->gidmap[(PetscInt)points[i]]+offset; 238 } 239 } 240 PetscFunctionReturn(0); 241 } 242 243 244 #undef __FUNCT__ 245 #define __FUNCT__ "DMMoabGetDofsLocal" 246 PetscErrorCode DMMoabGetDofsLocal(DM dm,PetscInt npoints,const moab::EntityHandle* points,PetscInt* dof) 247 { 248 PetscInt i,field,offset; 249 PetscErrorCode ierr; 250 DM_Moab *dmmoab; 251 252 PetscFunctionBegin; 253 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 254 PetscValidPointer(points,2); 255 dmmoab = (DM_Moab*)(dm)->data; 256 257 if (!dof) { 258 ierr = PetscMalloc(sizeof(PetscInt)*dmmoab->nfields*npoints, &dof);CHKERRQ(ierr); 259 } 260 261 if (dmmoab->bs > 1) { 262 for (field=0; field<dmmoab->nfields; ++field) { 263 for (i=0; i<npoints; ++i) 264 dof[i*dmmoab->nfields+field] = dmmoab->lidmap[(PetscInt)points[i]]*dmmoab->nfields+field; 265 } 266 } 267 else { 268 for (field=0; field<dmmoab->nfields; ++field) { 269 offset = field*dmmoab->n; /* assume all fields have equal distribution */ 270 for (i=0; i<npoints; ++i) 271 dof[i*dmmoab->nfields+field] = dmmoab->lidmap[(PetscInt)points[i]]+offset; 272 } 273 } 274 PetscFunctionReturn(0); 275 } 276 277 278 #undef __FUNCT__ 279 #define __FUNCT__ "DMMoabGetDofsBlocked" 280 PetscErrorCode DMMoabGetDofsBlocked(DM dm,PetscInt npoints,const moab::EntityHandle* points,PetscInt* dof) 281 { 282 PetscInt i; 283 DM_Moab *dmmoab; 284 PetscErrorCode ierr; 285 286 PetscFunctionBegin; 287 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 288 PetscValidPointer(points,2); 289 dmmoab = (DM_Moab*)(dm)->data; 290 291 if (!dof) { 292 ierr = PetscMalloc(sizeof(PetscInt)*npoints, &dof);CHKERRQ(ierr); 293 } 294 295 for (i=0; i<npoints; ++i) { 296 dof[i]=dmmoab->gidmap[(PetscInt)points[i]]; 297 } 298 PetscFunctionReturn(0); 299 } 300 301 302 #undef __FUNCT__ 303 #define __FUNCT__ "DMMoabGetDofsBlockedLocal" 304 PetscErrorCode DMMoabGetDofsBlockedLocal(DM dm,PetscInt npoints,const moab::EntityHandle* points,PetscInt* dof) 305 { 306 PetscInt i; 307 DM_Moab *dmmoab; 308 PetscErrorCode ierr; 309 310 PetscFunctionBegin; 311 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 312 PetscValidPointer(points,2); 313 dmmoab = (DM_Moab*)(dm)->data; 314 315 if (!dof) { 316 ierr = PetscMalloc(sizeof(PetscInt)*npoints, &dof);CHKERRQ(ierr); 317 } 318 319 for (i=0; i<npoints; ++i) 320 dof[i] = dmmoab->lidmap[(PetscInt)points[i]]; 321 PetscFunctionReturn(0); 322 } 323 324 325 #undef __FUNCT__ 326 #define __FUNCT__ "DMMoabGetVertexDofsBlocked" 327 PetscErrorCode DMMoabGetVertexDofsBlocked(DM dm,PetscInt** dof) 328 { 329 DM_Moab *dmmoab; 330 331 PetscFunctionBegin; 332 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 333 dmmoab = (DM_Moab*)(dm)->data; 334 335 *dof = dmmoab->gidmap; 336 PetscFunctionReturn(0); 337 } 338 339 340 #undef __FUNCT__ 341 #define __FUNCT__ "DMMoabGetVertexDofsBlockedLocal" 342 PetscErrorCode DMMoabGetVertexDofsBlockedLocal(DM dm,PetscInt** dof) 343 { 344 DM_Moab *dmmoab; 345 346 PetscFunctionBegin; 347 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 348 PetscValidPointer(dof,2); 349 dmmoab = (DM_Moab*)(dm)->data; 350 351 *dof = dmmoab->lidmap; 352 PetscFunctionReturn(0); 353 } 354 355 356 357