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