1 2 static char help[] = "Tests DMSwarm\n\n"; 3 4 #include <petscdm.h> 5 #include <petscdmda.h> 6 #include <petscdmswarm.h> 7 8 /* 9 Checks for variable blocksize 10 */ 11 PetscErrorCode ex2_1(void) 12 { 13 DM dms; 14 Vec x; 15 PetscMPIInt rank; 16 PetscInt p,bs,nlocal; 17 18 PetscFunctionBegin; 19 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD,&rank)); 20 21 PetscCall(DMCreate(PETSC_COMM_WORLD,&dms)); 22 PetscCall(DMSetType(dms,DMSWARM)); 23 PetscCall(DMSwarmInitializeFieldRegister(dms)); 24 PetscCall(DMSwarmRegisterPetscDatatypeField(dms,"viscosity",1,PETSC_REAL)); 25 PetscCall(DMSwarmRegisterPetscDatatypeField(dms,"strain",3,PETSC_REAL)); 26 PetscCall(DMSwarmFinalizeFieldRegister(dms)); 27 PetscCall(DMSwarmSetLocalSizes(dms,5+rank,4)); 28 PetscCall(DMView(dms,PETSC_VIEWER_STDOUT_WORLD)); 29 PetscCall(DMSwarmGetLocalSize(dms,&nlocal)); 30 31 { 32 PetscReal *array; 33 PetscCall(DMSwarmGetField(dms,"viscosity",&bs,NULL,(void**)&array)); 34 for (p=0; p<nlocal; p++) { 35 array[p] = 11.1 + p*0.1 + rank*100.0; 36 } 37 PetscCall(DMSwarmRestoreField(dms,"viscosity",&bs,NULL,(void**)&array)); 38 } 39 40 { 41 PetscReal *array; 42 PetscCall(DMSwarmGetField(dms,"strain",&bs,NULL,(void**)&array)); 43 for (p=0; p<nlocal; p++) { 44 array[bs*p+0] = 2.0e-2 + p*0.001 + rank*1.0; 45 array[bs*p+1] = 2.0e-2 + p*0.002 + rank*1.0; 46 array[bs*p+2] = 2.0e-2 + p*0.003 + rank*1.0; 47 } 48 PetscCall(DMSwarmRestoreField(dms,"strain",&bs,NULL,(void**)&array)); 49 } 50 51 PetscCall(DMSwarmCreateGlobalVectorFromField(dms,"viscosity",&x)); 52 PetscCall(VecView(x,PETSC_VIEWER_STDOUT_WORLD)); 53 PetscCall(DMSwarmDestroyGlobalVectorFromField(dms,"viscosity",&x)); 54 55 PetscCall(DMSwarmCreateGlobalVectorFromField(dms,"strain",&x)); 56 PetscCall(VecView(x,PETSC_VIEWER_STDOUT_WORLD)); 57 PetscCall(DMSwarmDestroyGlobalVectorFromField(dms,"strain",&x)); 58 59 PetscCall(DMSwarmVectorDefineField(dms,"strain")); 60 PetscCall(DMCreateGlobalVector(dms,&x)); 61 PetscCall(VecView(x,PETSC_VIEWER_STDOUT_WORLD)); 62 PetscCall(VecDestroy(&x)); 63 PetscCall(DMDestroy(&dms)); 64 PetscFunctionReturn(0); 65 } 66 67 int main(int argc,char **argv) 68 { 69 70 PetscCall(PetscInitialize(&argc,&argv,(char*)0,help)); 71 PetscCall(ex2_1()); 72 PetscCall(PetscFinalize()); 73 return 0; 74 } 75 76 /*TEST 77 78 test: 79 requires: !complex double 80 nsize: 3 81 filter: grep -v atomic 82 filter_output: grep -v atomic 83 84 TEST*/ 85