1 2 #define PETSCDM_DLL 3 #include <petsc/private/dmswarmimpl.h> /*I "petscdmswarm.h" I*/ 4 #include "data_bucket.h" 5 6 PetscErrorCode DMSwarmMigrate_Push_Basic(DM dm,PetscBool remove_sent_points); 7 8 9 //typedef PetscErrorCode (*swarm_project)(DM,DM,Vec) DMSwarmProjectMethod; /* swarm, geometry, result */ 10 11 //typedef enum { PROJECT_DMDA_AQ1=0, PROJECT_DMDA_P0 } DMSwarmDMDAProjectionType; 12 13 #if 0 14 15 /* Defines what the local space will be */ 16 PetscErrorCode DMSwarmSetOverlap(void) 17 { 18 19 PetscFunctionReturn(0); 20 } 21 22 23 /* coordinates */ 24 /* 25 DMGetCoordinateDM returns self 26 DMGetCoordinates and DMGetCoordinatesLocal return same thing 27 Local view could be used to define overlapping information 28 */ 29 30 #endif 31 32 #undef __FUNCT__ 33 #define __FUNCT__ "DMSwarmVectorDefineField" 34 PETSC_EXTERN PetscErrorCode DMSwarmVectorDefineField(DM dm,const char fieldname[]) 35 { 36 DM_Swarm *swarm = (DM_Swarm*)dm->data; 37 PetscErrorCode ierr; 38 PetscInt bs,n; 39 PetscScalar *array; 40 PetscDataType type; 41 42 if (!swarm->issetup) { ierr = DMSetUp(dm);CHKERRQ(ierr); } 43 ierr = DataBucketGetSizes(swarm->db,&n,NULL,NULL);CHKERRQ(ierr); 44 ierr = DMSwarmGetField(dm,fieldname,&bs,&type,(void**)&array);CHKERRQ(ierr); 45 46 /* Check all fields are of type PETSC_REAL or PETSC_SCALAR */ 47 if (type != PETSC_REAL) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Only valid for PETSC_REAL"); 48 49 PetscSNPrintf(swarm->vec_field_name,PETSC_MAX_PATH_LEN-1,"%s",fieldname); 50 swarm->vec_field_set = PETSC_TRUE; 51 swarm->vec_field_bs = bs; 52 swarm->vec_field_nlocal = n; 53 ierr = DMSwarmRestoreField(dm,fieldname,&bs,&type,(void**)&array);CHKERRQ(ierr); 54 55 PetscFunctionReturn(0); 56 } 57 58 /* requires DMSwarmDefineFieldVector has been called */ 59 #undef __FUNCT__ 60 #define __FUNCT__ "DMCreateGlobalVector_Swarm" 61 PetscErrorCode DMCreateGlobalVector_Swarm(DM dm,Vec *vec) 62 { 63 DM_Swarm *swarm = (DM_Swarm*)dm->data; 64 PetscErrorCode ierr; 65 Vec x; 66 char name[PETSC_MAX_PATH_LEN]; 67 68 if (!swarm->issetup) { ierr = DMSetUp(dm);CHKERRQ(ierr); } 69 if (!swarm->vec_field_set) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_USER,"Must call DMSwarmVectorDefineField first"); 70 PetscSNPrintf(name,PETSC_MAX_PATH_LEN-1,"DMSwarmField_%s",swarm->vec_field_name); 71 ierr = VecCreate(PetscObjectComm((PetscObject)dm),&x);CHKERRQ(ierr); 72 ierr = PetscObjectSetName((PetscObject)x,name);CHKERRQ(ierr); 73 ierr = VecSetSizes(x,swarm->db->L*swarm->vec_field_bs,PETSC_DETERMINE);CHKERRQ(ierr); 74 ierr = VecSetBlockSize(x,swarm->vec_field_bs);CHKERRQ(ierr); 75 ierr = VecSetFromOptions(x);CHKERRQ(ierr); 76 *vec = x; 77 78 PetscFunctionReturn(0); 79 } 80 81 /* requires DMSwarmDefineFieldVector has been called */ 82 #undef __FUNCT__ 83 #define __FUNCT__ "DMCreateLocalVector_Swarm" 84 PetscErrorCode DMCreateLocalVector_Swarm(DM dm,Vec *vec) 85 { 86 DM_Swarm *swarm = (DM_Swarm*)dm->data; 87 PetscErrorCode ierr; 88 Vec x; 89 char name[PETSC_MAX_PATH_LEN]; 90 91 if (!swarm->issetup) { ierr = DMSetUp(dm);CHKERRQ(ierr); } 92 if (!swarm->vec_field_set) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_USER,"Must call DMSwarmVectorDefineField first"); 93 PetscSNPrintf(name,PETSC_MAX_PATH_LEN-1,"DMSwarmField_%s",swarm->vec_field_name); 94 ierr = VecCreate(PETSC_COMM_SELF,&x);CHKERRQ(ierr); 95 ierr = PetscObjectSetName((PetscObject)x,name);CHKERRQ(ierr); 96 ierr = VecSetSizes(x,swarm->db->L*swarm->vec_field_bs,swarm->db->L);CHKERRQ(ierr); 97 ierr = VecSetBlockSize(x,swarm->vec_field_bs);CHKERRQ(ierr); 98 ierr = VecSetFromOptions(x);CHKERRQ(ierr); 99 *vec = x; 100 101 PetscFunctionReturn(0); 102 } 103 104 #undef __FUNCT__ 105 #define __FUNCT__ "DMSwarmCreateGlobalVectorFromField" 106 PETSC_EXTERN PetscErrorCode DMSwarmCreateGlobalVectorFromField(DM dm,const char fieldname[],Vec *vec) 107 { 108 DM_Swarm *swarm = (DM_Swarm*)dm->data; 109 PetscErrorCode ierr; 110 PetscInt bs,n; 111 PetscScalar *array; 112 Vec x; 113 PetscDataType type; 114 char name[PETSC_MAX_PATH_LEN]; 115 PetscMPIInt commsize; 116 117 if (!swarm->issetup) { ierr = DMSetUp(dm);CHKERRQ(ierr); } 118 ierr = DataBucketGetSizes(swarm->db,&n,NULL,NULL);CHKERRQ(ierr); 119 ierr = DMSwarmGetField(dm,fieldname,&bs,&type,(void**)&array);CHKERRQ(ierr); 120 121 /* Check all fields are of type PETSC_REAL or PETSC_SCALAR */ 122 if (type != PETSC_REAL) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Only valid for PETSC_REAL"); 123 124 ierr = MPI_Comm_size(PetscObjectComm((PetscObject)dm),&commsize);CHKERRQ(ierr); 125 if (commsize == 1) { 126 ierr = VecCreateSeqWithArray(PetscObjectComm((PetscObject)dm),bs,n*bs,array,&x);CHKERRQ(ierr); 127 } else { 128 ierr = VecCreateMPIWithArray(PetscObjectComm((PetscObject)dm),bs,n*bs,PETSC_DETERMINE,array,&x);CHKERRQ(ierr); 129 } 130 PetscSNPrintf(name,PETSC_MAX_PATH_LEN-1,"DMSwarmSharedField_%s",fieldname); 131 ierr = PetscObjectSetName((PetscObject)x,name);CHKERRQ(ierr); 132 133 /* Set guard */ 134 PetscSNPrintf(name,PETSC_MAX_PATH_LEN-1,"DMSwarm_VecFieldInPlace_%s",fieldname); 135 ierr = PetscObjectComposeFunction((PetscObject)x,name,DMSwarmDestroyGlobalVectorFromField);CHKERRQ(ierr); 136 137 *vec = x; 138 PetscFunctionReturn(0); 139 } 140 141 #undef __FUNCT__ 142 #define __FUNCT__ "DMSwarmDestroyGlobalVectorFromField" 143 PETSC_EXTERN PetscErrorCode DMSwarmDestroyGlobalVectorFromField(DM dm,const char fieldname[],Vec *vec) 144 { 145 DM_Swarm *swarm = (DM_Swarm*)dm->data; 146 PetscErrorCode ierr; 147 DataField gfield; 148 char name[PETSC_MAX_PATH_LEN]; 149 void (*fptr)(void); 150 151 /* get data field */ 152 ierr = DataBucketGetDataFieldByName(swarm->db,fieldname,&gfield);CHKERRQ(ierr); 153 154 /* check vector is an inplace array */ 155 PetscSNPrintf(name,PETSC_MAX_PATH_LEN-1,"DMSwarm_VecFieldInPlace_%s",fieldname); 156 ierr = PetscObjectQueryFunction((PetscObject)(*vec),name,&fptr);CHKERRQ(ierr); 157 if (!fptr) SETERRQ1(PetscObjectComm((PetscObject)dm),PETSC_ERR_USER,"Vector being destroyed was not created from DMSwarm field(%s)",fieldname); 158 159 /* restore data field */ 160 ierr = DataFieldRestoreAccess(gfield);CHKERRQ(ierr); 161 162 ierr = VecDestroy(vec);CHKERRQ(ierr); 163 164 PetscFunctionReturn(0); 165 } 166 167 /* 168 PETSC_EXTERN PetscErrorCode DMSwarmCreateGlobalVectorFromFields(DM dm,const PetscInt nf,const char *fieldnames[],Vec *vec) 169 { 170 PetscFunctionReturn(0); 171 } 172 173 PETSC_EXTERN PetscErrorCode DMSwarmRestoreGlobalVectorFromFields(DM dm,Vec *vec) 174 { 175 PetscFunctionReturn(0); 176 } 177 */ 178 179 #undef __FUNCT__ 180 #define __FUNCT__ "DMSwarmInitializeFieldRegister" 181 PETSC_EXTERN PetscErrorCode DMSwarmInitializeFieldRegister(DM dm) 182 { 183 DM_Swarm *swarm = (DM_Swarm*)dm->data; 184 PetscErrorCode ierr; 185 186 swarm->field_registration_initialized = PETSC_TRUE; 187 188 ierr = DMSwarmRegisterPetscDatatypeField(dm,"DMSwarm_pid",1,PETSC_LONG);CHKERRQ(ierr); /* unique identifer */ 189 ierr = DMSwarmRegisterPetscDatatypeField(dm,"DMSwarm_rank",1,PETSC_INT);CHKERRQ(ierr); /* used for communication */ 190 191 PetscFunctionReturn(0); 192 } 193 194 #undef __FUNCT__ 195 #define __FUNCT__ "DMSwarmFinalizeFieldRegister" 196 PETSC_EXTERN PetscErrorCode DMSwarmFinalizeFieldRegister(DM dm) 197 { 198 DM_Swarm *swarm = (DM_Swarm*)dm->data; 199 PetscErrorCode ierr; 200 201 swarm->field_registration_finalized = PETSC_TRUE; 202 ierr = DataBucketFinalize(swarm->db);CHKERRQ(ierr); 203 PetscFunctionReturn(0); 204 } 205 206 #undef __FUNCT__ 207 #define __FUNCT__ "DMSwarmSetLocalSizes" 208 PETSC_EXTERN PetscErrorCode DMSwarmSetLocalSizes(DM dm,PetscInt nlocal,PetscInt buffer) 209 { 210 DM_Swarm *swarm = (DM_Swarm*)dm->data; 211 PetscErrorCode ierr; 212 213 ierr = DataBucketSetSizes(swarm->db,nlocal,buffer);CHKERRQ(ierr); 214 215 PetscFunctionReturn(0); 216 } 217 218 #undef __FUNCT__ 219 #define __FUNCT__ "DMSwarmGetLocalSize" 220 PETSC_EXTERN PetscErrorCode DMSwarmGetLocalSize(DM dm,PetscInt *nlocal) 221 { 222 DM_Swarm *swarm = (DM_Swarm*)dm->data; 223 PetscErrorCode ierr; 224 225 if (nlocal) { 226 ierr = DataBucketGetSizes(swarm->db,nlocal,NULL,NULL);CHKERRQ(ierr); 227 } 228 229 PetscFunctionReturn(0); 230 } 231 232 #undef __FUNCT__ 233 #define __FUNCT__ "DMSwarmGetSize" 234 PETSC_EXTERN PetscErrorCode DMSwarmGetSize(DM dm,PetscInt *n) 235 { 236 DM_Swarm *swarm = (DM_Swarm*)dm->data; 237 PetscErrorCode ierr; 238 PetscInt nlocal,ng; 239 240 ierr = DataBucketGetSizes(swarm->db,&nlocal,NULL,NULL);CHKERRQ(ierr); 241 ierr = MPI_Allreduce(&nlocal,&ng,1,MPIU_INT,MPI_SUM,PetscObjectComm((PetscObject)dm));CHKERRQ(ierr); 242 if (n) { *n = ng; } 243 PetscFunctionReturn(0); 244 } 245 246 #undef __FUNCT__ 247 #define __FUNCT__ "DMSwarmRegisterPetscDatatypeField" 248 PETSC_EXTERN PetscErrorCode DMSwarmRegisterPetscDatatypeField(DM dm,const char fieldname[],PetscInt blocksize,PetscDataType type) 249 { 250 PetscErrorCode ierr; 251 DM_Swarm *swarm = (DM_Swarm*)dm->data; 252 size_t size; 253 254 if (!swarm->field_registration_initialized) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_USER,"Must call DMSwarmInitializeFieldRegister() first"); 255 if (swarm->field_registration_finalized) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_USER,"Cannot register additional fields after calling DMSwarmFinalizeFieldRegister() first"); 256 257 if (type == PETSC_OBJECT) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Valid for {char,short,int,long,float,double}"); 258 if (type == PETSC_FUNCTION) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Valid for {char,short,int,long,float,double}"); 259 if (type == PETSC_STRING) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Valid for {char,short,int,long,float,double}"); 260 if (type == PETSC_STRUCT) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Valid for {char,short,int,long,float,double}"); 261 if (type == PETSC_DATATYPE_UNKNOWN) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Valid for {char,short,int,long,float,double}"); 262 263 switch (type) { 264 case PETSC_CHAR: 265 size = sizeof(PetscChar); 266 break; 267 case PETSC_SHORT: 268 size = sizeof(PetscShort); 269 break; 270 case PETSC_INT: 271 size = sizeof(PetscInt); 272 break; 273 case PETSC_LONG: 274 size = sizeof(Petsc64bitInt); 275 break; 276 case PETSC_FLOAT: 277 size = sizeof(PetscFloat); 278 break; 279 case PETSC_DOUBLE: 280 size = sizeof(PetscReal); 281 break; 282 283 default: 284 SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Valid for {char,short,int,long,float,double}"); 285 break; 286 } 287 288 /* Load a specific data type into data bucket, specifying textual name and its size in bytes */ 289 ierr = DataBucketRegisterField(swarm->db,"DMSwarmRegisterPetscDatatypeField",fieldname,blocksize*size,NULL);CHKERRQ(ierr); 290 { 291 DataField gfield; 292 293 ierr = DataBucketGetDataFieldByName(swarm->db,fieldname,&gfield);CHKERRQ(ierr); 294 ierr = DataFieldSetBlockSize(gfield,blocksize);CHKERRQ(ierr); 295 } 296 swarm->db->field[swarm->db->nfields-1]->petsc_type = type; 297 298 PetscFunctionReturn(0); 299 } 300 301 #undef __FUNCT__ 302 #define __FUNCT__ "DMSwarmRegisterUserStructField" 303 PETSC_EXTERN PetscErrorCode DMSwarmRegisterUserStructField(DM dm,const char fieldname[],size_t size) 304 { 305 PetscErrorCode ierr; 306 DM_Swarm *swarm = (DM_Swarm*)dm->data; 307 308 ierr = DataBucketRegisterField(swarm->db,"DMSwarmRegisterUserStructField",fieldname,size,NULL);CHKERRQ(ierr); 309 swarm->db->field[swarm->db->nfields-1]->petsc_type = PETSC_STRUCT ; 310 311 PetscFunctionReturn(0); 312 } 313 314 #undef __FUNCT__ 315 #define __FUNCT__ "DMSwarmRegisterUserDatatypeField" 316 PETSC_EXTERN PetscErrorCode DMSwarmRegisterUserDatatypeField(DM dm,const char fieldname[],size_t size,PetscInt blocksize) 317 { 318 DM_Swarm *swarm = (DM_Swarm*)dm->data; 319 PetscErrorCode ierr; 320 321 ierr = DataBucketRegisterField(swarm->db,"DMSwarmRegisterUserDatatypeField",fieldname,blocksize*size,NULL);CHKERRQ(ierr); 322 { 323 DataField gfield; 324 325 ierr = DataBucketGetDataFieldByName(swarm->db,fieldname,&gfield);CHKERRQ(ierr); 326 ierr = DataFieldSetBlockSize(gfield,blocksize);CHKERRQ(ierr); 327 } 328 swarm->db->field[swarm->db->nfields-1]->petsc_type = PETSC_DATATYPE_UNKNOWN; 329 330 PetscFunctionReturn(0); 331 } 332 333 #undef __FUNCT__ 334 #define __FUNCT__ "DMSwarmGetField" 335 PETSC_EXTERN PetscErrorCode DMSwarmGetField(DM dm,const char fieldname[],PetscInt *blocksize,PetscDataType *type,void **data) 336 { 337 DM_Swarm *swarm = (DM_Swarm*)dm->data; 338 DataField gfield; 339 PetscErrorCode ierr; 340 341 if (!swarm->issetup) { ierr = DMSetUp(dm);CHKERRQ(ierr); } 342 343 ierr = DataBucketGetDataFieldByName(swarm->db,fieldname,&gfield);CHKERRQ(ierr); 344 ierr = DataFieldGetAccess(gfield);CHKERRQ(ierr); 345 ierr = DataFieldGetEntries(gfield,data);CHKERRQ(ierr); 346 if (blocksize) {*blocksize = gfield->bs; } 347 if (type) { *type = gfield->petsc_type; } 348 349 PetscFunctionReturn(0); 350 } 351 352 #undef __FUNCT__ 353 #define __FUNCT__ "DMSwarmRestoreField" 354 PETSC_EXTERN PetscErrorCode DMSwarmRestoreField(DM dm,const char fieldname[],PetscInt *blocksize,PetscDataType *type,void **data) 355 { 356 DM_Swarm *swarm = (DM_Swarm*)dm->data; 357 DataField gfield; 358 PetscErrorCode ierr; 359 360 ierr = DataBucketGetDataFieldByName(swarm->db,fieldname,&gfield);CHKERRQ(ierr); 361 ierr = DataFieldRestoreAccess(gfield);CHKERRQ(ierr); 362 if (data) *data = NULL; 363 364 PetscFunctionReturn(0); 365 } 366 367 #undef __FUNCT__ 368 #define __FUNCT__ "DMSwarmAddPoint" 369 PETSC_EXTERN PetscErrorCode DMSwarmAddPoint(DM dm) 370 { 371 DM_Swarm *swarm = (DM_Swarm*)dm->data; 372 PetscErrorCode ierr; 373 374 if (!swarm->issetup) { ierr = DMSetUp(dm);CHKERRQ(ierr); } 375 ierr = DataBucketAddPoint(swarm->db);CHKERRQ(ierr); 376 PetscFunctionReturn(0); 377 } 378 379 #undef __FUNCT__ 380 #define __FUNCT__ "DMSwarmAddNPoints" 381 PETSC_EXTERN PetscErrorCode DMSwarmAddNPoints(DM dm,PetscInt npoints) 382 { 383 DM_Swarm *swarm = (DM_Swarm*)dm->data; 384 PetscErrorCode ierr; 385 PetscInt nlocal; 386 387 ierr = DataBucketGetSizes(swarm->db,&nlocal,NULL,NULL);CHKERRQ(ierr); 388 nlocal = nlocal + npoints; 389 ierr = DataBucketSetSizes(swarm->db,nlocal,-1);CHKERRQ(ierr); 390 PetscFunctionReturn(0); 391 } 392 393 #undef __FUNCT__ 394 #define __FUNCT__ "DMSwarmRemovePoint" 395 PETSC_EXTERN PetscErrorCode DMSwarmRemovePoint(DM dm) 396 { 397 DM_Swarm *swarm = (DM_Swarm*)dm->data; 398 PetscErrorCode ierr; 399 400 ierr = DataBucketRemovePoint(swarm->db);CHKERRQ(ierr); 401 PetscFunctionReturn(0); 402 } 403 404 #undef __FUNCT__ 405 #define __FUNCT__ "DMSwarmRemovePointAtIndex" 406 PETSC_EXTERN PetscErrorCode DMSwarmRemovePointAtIndex(DM dm,PetscInt idx) 407 { 408 DM_Swarm *swarm = (DM_Swarm*)dm->data; 409 PetscErrorCode ierr; 410 411 ierr = DataBucketRemovePointAtIndex(swarm->db,idx);CHKERRQ(ierr); 412 PetscFunctionReturn(0); 413 } 414 415 #undef __FUNCT__ 416 #define __FUNCT__ "DMSwarmMigrate_Basic" 417 PetscErrorCode DMSwarmMigrate_Basic(DM dm,PetscBool remove_sent_points) 418 { 419 PetscErrorCode ierr; 420 ierr = DMSwarmMigrate_Push_Basic(dm,remove_sent_points);CHKERRQ(ierr); 421 PetscFunctionReturn(0); 422 } 423 424 #undef __FUNCT__ 425 #define __FUNCT__ "DMSwarmMigrate" 426 PETSC_EXTERN PetscErrorCode DMSwarmMigrate(DM dm,PetscBool remove_sent_points) 427 { 428 PetscErrorCode ierr; 429 ierr = DMSwarmMigrate_Basic(dm,remove_sent_points);CHKERRQ(ierr); 430 PetscFunctionReturn(0); 431 } 432 433 #undef __FUNCT__ 434 #define __FUNCT__ "DMSwarmGlobalToLocalViewCreate" 435 PETSC_EXTERN PetscErrorCode DMSwarmGlobalToLocalViewCreate(DM dm,InsertMode mode) 436 { 437 PetscErrorCode ierr; 438 DM_Swarm *swarm = (DM_Swarm*)dm->data; 439 PetscInt ng; 440 441 if (mode != INSERT_VALUES) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Only mode INSERT_VALUES is supported"); 442 443 ierr = DMSwarmMigrate_GlobalToLocal_Basic(dm,&ng);CHKERRQ(ierr); 444 swarm->view_ng = ng; 445 446 PetscFunctionReturn(0); 447 } 448 449 #undef __FUNCT__ 450 #define __FUNCT__ "DMSwarmGlobalToLocalViewDestroy" 451 PETSC_EXTERN PetscErrorCode DMSwarmGlobalToLocalViewDestroy(DM dm,InsertMode mode) 452 { 453 PetscErrorCode ierr; 454 DM_Swarm *swarm = (DM_Swarm*)dm->data; 455 456 if (mode != INSERT_VALUES) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Only mode INSERT_VALUES is supported"); 457 458 ierr = DMSwarmSetLocalSizes(dm,swarm->view_ng,-1);CHKERRQ(ierr); 459 460 PetscFunctionReturn(0); 461 } 462 463 #undef __FUNCT__ 464 #define __FUNCT__ "DMSetup_Swarm" 465 PetscErrorCode DMSetup_Swarm(DM dm) 466 { 467 DM_Swarm *swarm = (DM_Swarm*)dm->data; 468 PetscErrorCode ierr; 469 PetscMPIInt rank; 470 PetscInt p,npoints,*rankval; 471 472 if (swarm->issetup) PetscFunctionReturn(0); 473 474 PetscPrintf(PETSC_COMM_SELF,"Swarm setup \n"); 475 swarm->issetup = PETSC_TRUE; 476 477 /* check some fields were registered */ 478 if (swarm->db->nfields <= 2) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_USER,"At least one field user must be registered via DMSwarmRegisterXXX()"); 479 480 /* check local sizes were set */ 481 if (swarm->db->L == -1) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_USER,"Local sizes must be set via DMSwarmSetLocalSizes()"); 482 483 /* initialize values in pid and rank placeholders */ 484 /* TODO: [pid - use MPI_Scan] */ 485 486 ierr = MPI_Comm_rank(PetscObjectComm((PetscObject)dm),&rank);CHKERRQ(ierr); 487 ierr = DataBucketGetSizes(swarm->db,&npoints,NULL,NULL);CHKERRQ(ierr); 488 ierr = DMSwarmGetField(dm,"DMSwarm_rank",NULL,NULL,(void**)&rankval);CHKERRQ(ierr); 489 for (p=0; p<npoints; p++) { 490 rankval[p] = (PetscInt)rank; 491 } 492 ierr = DMSwarmRestoreField(dm,"DMSwarm_rank",NULL,NULL,(void**)&rankval);CHKERRQ(ierr); 493 494 PetscFunctionReturn(0); 495 } 496 497 #undef __FUNCT__ 498 #define __FUNCT__ "DMDestroy_Swarm" 499 PetscErrorCode DMDestroy_Swarm(DM dm) 500 { 501 DM_Swarm *swarm = (DM_Swarm*)dm->data; 502 PetscErrorCode ierr; 503 504 PetscFunctionBegin; 505 ierr = DataBucketDestroy(&swarm->db);CHKERRQ(ierr); 506 ierr = PetscFree(swarm);CHKERRQ(ierr); 507 PetscFunctionReturn(0); 508 } 509 510 #undef __FUNCT__ 511 #define __FUNCT__ "DMView_Swarm" 512 PetscErrorCode DMView_Swarm(DM dm, PetscViewer viewer) 513 { 514 DM_Swarm *swarm = (DM_Swarm*)dm->data; 515 PetscBool iascii,ibinary,ishdf5,isvtk; 516 PetscErrorCode ierr; 517 518 PetscFunctionBegin; 519 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 520 PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,2); 521 ierr = PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii);CHKERRQ(ierr); 522 ierr = PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERBINARY,&ibinary);CHKERRQ(ierr); 523 ierr = PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERVTK, &isvtk);CHKERRQ(ierr); 524 ierr = PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERHDF5, &ishdf5);CHKERRQ(ierr); 525 if (iascii) { 526 ierr = DataBucketView(PetscObjectComm((PetscObject)dm),swarm->db,NULL,DATABUCKET_VIEW_STDOUT);CHKERRQ(ierr); 527 } else if (ibinary) { 528 SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"NO VTK support"); 529 } else if (ishdf5) { 530 #if defined(PETSC_HAVE_HDF5) 531 SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"NO HDF5 support"); 532 #else 533 SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"HDF5 not supported. Please reconfigure using --download-hdf5"); 534 #endif 535 } else if (isvtk) { 536 SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"NO VTK support"); 537 } 538 PetscFunctionReturn(0); 539 } 540 541 #undef __FUNCT__ 542 #define __FUNCT__ "DMCreate_Swarm" 543 PETSC_EXTERN PetscErrorCode DMCreate_Swarm(DM dm) 544 { 545 DM_Swarm *swarm; 546 PetscErrorCode ierr; 547 548 PetscFunctionBegin; 549 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 550 ierr = PetscNewLog(dm,&swarm);CHKERRQ(ierr); 551 552 ierr = DataBucketCreate(&swarm->db);CHKERRQ(ierr); 553 swarm->vec_field_set = PETSC_FALSE; 554 swarm->issetup = PETSC_FALSE; 555 556 dm->dim = 0; 557 dm->data = swarm; 558 559 dm->ops->view = DMView_Swarm; 560 dm->ops->load = NULL; 561 dm->ops->setfromoptions = NULL; 562 dm->ops->clone = NULL; 563 dm->ops->setup = DMSetup_Swarm; 564 dm->ops->createdefaultsection = NULL; 565 dm->ops->createdefaultconstraints = NULL; 566 dm->ops->createglobalvector = DMCreateGlobalVector_Swarm; 567 dm->ops->createlocalvector = DMCreateLocalVector_Swarm; 568 dm->ops->getlocaltoglobalmapping = NULL; 569 dm->ops->createfieldis = NULL; 570 dm->ops->createcoordinatedm = NULL; 571 dm->ops->getcoloring = NULL; 572 dm->ops->creatematrix = NULL; 573 dm->ops->createinterpolation = NULL; 574 dm->ops->getaggregates = NULL; 575 dm->ops->getinjection = NULL; 576 dm->ops->refine = NULL; 577 dm->ops->coarsen = NULL; 578 dm->ops->refinehierarchy = NULL; 579 dm->ops->coarsenhierarchy = NULL; 580 dm->ops->globaltolocalbegin = NULL; 581 dm->ops->globaltolocalend = NULL; 582 dm->ops->localtoglobalbegin = NULL; 583 dm->ops->localtoglobalend = NULL; 584 dm->ops->destroy = DMDestroy_Swarm; 585 dm->ops->createsubdm = NULL; 586 dm->ops->getdimpoints = NULL; 587 dm->ops->locatepoints = NULL; 588 589 PetscFunctionReturn(0); 590 }