1 2 #include <petsc/private/dmswarmimpl.h> /*I "petscdmswarm.h" I*/ 3 #include "data_bucket.h" 4 #include "data_ex.h" 5 6 7 #undef __FUNCT__ 8 #define __FUNCT__ "DMSwarmMigrate_Push_Basic" 9 PetscErrorCode DMSwarmMigrate_Push_Basic(DM dm,PetscBool remove_sent_points) 10 { 11 DM_Swarm *swarm = (DM_Swarm*)dm->data; 12 PetscErrorCode ierr; 13 DataEx de; 14 PetscInt p,npoints,*rankval,n_points_recv; 15 PetscMPIInt rank,nrank; 16 void *point_buffer,*recv_points; 17 size_t sizeof_dmswarm_point; 18 19 ierr = MPI_Comm_rank(PetscObjectComm((PetscObject)dm),&rank);CHKERRQ(ierr); 20 21 ierr = DataBucketGetSizes(swarm->db,&npoints,NULL,NULL);CHKERRQ(ierr); 22 ierr = DMSwarmGetField(dm,"DMSwarm_rank",NULL,NULL,(void**)&rankval);CHKERRQ(ierr); 23 24 de = DataExCreate(PetscObjectComm((PetscObject)dm),0);CHKERRQ(ierr); 25 26 ierr = DataExTopologyInitialize(de);CHKERRQ(ierr); 27 for (p=0; p<npoints; p++) { 28 nrank = rankval[p]; 29 if (nrank != rank) { 30 ierr = DataExTopologyAddNeighbour(de,nrank);CHKERRQ(ierr); 31 } 32 } 33 ierr = DataExTopologyFinalize(de);CHKERRQ(ierr); 34 35 ierr = DataExInitializeSendCount(de);CHKERRQ(ierr); 36 for (p=0; p<npoints; p++) { 37 nrank = rankval[p]; 38 if (nrank != rank) { 39 ierr = DataExAddToSendCount(de,nrank,1);CHKERRQ(ierr); 40 } 41 } 42 ierr = DataExFinalizeSendCount(de);CHKERRQ(ierr); 43 44 ierr = DataBucketCreatePackedArray(swarm->db,&sizeof_dmswarm_point,&point_buffer);CHKERRQ(ierr); 45 ierr = DataExPackInitialize(de,sizeof_dmswarm_point);CHKERRQ(ierr); 46 for (p=0; p<npoints; p++) { 47 nrank = rankval[p]; 48 if (nrank != rank) { 49 /* copy point into buffer */ 50 ierr = DataBucketFillPackedArray(swarm->db,p,point_buffer);CHKERRQ(ierr); 51 /* insert point buffer into DataExchanger */ 52 ierr = DataExPackData(de,nrank,1,point_buffer);CHKERRQ(ierr); 53 } 54 } 55 ierr = DataExPackFinalize(de);CHKERRQ(ierr); 56 57 58 if (remove_sent_points) { 59 /* remove points which left processor */ 60 ierr = DataBucketGetSizes(swarm->db,&npoints,NULL,NULL);CHKERRQ(ierr); 61 for (p=0; p<npoints; p++) { 62 nrank = rankval[p]; 63 if (nrank != rank) { 64 /* kill point */ 65 ierr = DataBucketRemovePointAtIndex(swarm->db,p);CHKERRQ(ierr); 66 DataBucketGetSizes(swarm->db,&npoints,NULL,NULL);CHKERRQ(ierr); /* you need to update npoints as the list size decreases! */ 67 p--; /* check replacement point */ 68 } 69 } 70 } 71 ierr = DMSwarmRestoreField(dm,"DMSwarm_rank",NULL,NULL,(void**)&rankval);CHKERRQ(ierr); 72 73 ierr = DataExBegin(de);CHKERRQ(ierr); 74 ierr = DataExEnd(de);CHKERRQ(ierr); 75 76 ierr = DataExGetRecvData(de,&n_points_recv,(void**)&recv_points);CHKERRQ(ierr); 77 78 ierr = DataBucketGetSizes(swarm->db,&npoints,NULL,NULL);CHKERRQ(ierr); 79 ierr = DataBucketSetSizes(swarm->db,npoints + n_points_recv,-1);CHKERRQ(ierr); 80 for (p=0; p<n_points_recv; p++) { 81 void *data_p = (void*)( (char*)recv_points + p*sizeof_dmswarm_point ); 82 83 ierr = DataBucketInsertPackedArray(swarm->db,npoints+p,data_p);CHKERRQ(ierr); 84 } 85 ierr = DataExView(de);CHKERRQ(ierr); 86 ierr = DataBucketDestroyPackedArray(swarm->db,&point_buffer);CHKERRQ(ierr); 87 ierr = DataExDestroy(de);CHKERRQ(ierr); 88 89 PetscFunctionReturn(0); 90 } 91