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