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 */ 119371c9d4SSatish Balay PetscErrorCode ex2_1(void) { 12c4762a1bSJed Brown DM dms; 13c4762a1bSJed Brown Vec x; 14c4762a1bSJed Brown PetscMPIInt rank; 15c4762a1bSJed Brown PetscInt p, bs, nlocal; 16c4762a1bSJed Brown 17362febeeSStefano Zampini PetscFunctionBegin; 189566063dSJacob Faibussowitsch PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank)); 19c4762a1bSJed Brown 209566063dSJacob Faibussowitsch PetscCall(DMCreate(PETSC_COMM_WORLD, &dms)); 219566063dSJacob Faibussowitsch PetscCall(DMSetType(dms, DMSWARM)); 226a5217c0SMatthew G. Knepley PetscCall(PetscObjectSetName((PetscObject)dms, "Particles")); 239566063dSJacob Faibussowitsch PetscCall(DMSwarmInitializeFieldRegister(dms)); 249566063dSJacob Faibussowitsch PetscCall(DMSwarmRegisterPetscDatatypeField(dms, "viscosity", 1, PETSC_REAL)); 259566063dSJacob Faibussowitsch PetscCall(DMSwarmRegisterPetscDatatypeField(dms, "strain", 3, PETSC_REAL)); 269566063dSJacob Faibussowitsch PetscCall(DMSwarmFinalizeFieldRegister(dms)); 279566063dSJacob Faibussowitsch PetscCall(DMSwarmSetLocalSizes(dms, 5 + rank, 4)); 289566063dSJacob Faibussowitsch PetscCall(DMView(dms, PETSC_VIEWER_STDOUT_WORLD)); 299566063dSJacob Faibussowitsch PetscCall(DMSwarmGetLocalSize(dms, &nlocal)); 30c4762a1bSJed Brown 31c4762a1bSJed Brown { 32c4762a1bSJed Brown PetscReal *array; 339566063dSJacob Faibussowitsch PetscCall(DMSwarmGetField(dms, "viscosity", &bs, NULL, (void **)&array)); 34*ad540459SPierre Jolivet for (p = 0; p < nlocal; p++) array[p] = 11.1 + p * 0.1 + rank * 100.0; 359566063dSJacob Faibussowitsch PetscCall(DMSwarmRestoreField(dms, "viscosity", &bs, NULL, (void **)&array)); 36c4762a1bSJed Brown } 37c4762a1bSJed Brown 38c4762a1bSJed Brown { 39c4762a1bSJed Brown PetscReal *array; 409566063dSJacob Faibussowitsch PetscCall(DMSwarmGetField(dms, "strain", &bs, NULL, (void **)&array)); 41c4762a1bSJed Brown for (p = 0; p < nlocal; p++) { 42c4762a1bSJed Brown array[bs * p + 0] = 2.0e-2 + p * 0.001 + rank * 1.0; 43c4762a1bSJed Brown array[bs * p + 1] = 2.0e-2 + p * 0.002 + rank * 1.0; 44c4762a1bSJed Brown array[bs * p + 2] = 2.0e-2 + p * 0.003 + rank * 1.0; 45c4762a1bSJed Brown } 469566063dSJacob Faibussowitsch PetscCall(DMSwarmRestoreField(dms, "strain", &bs, NULL, (void **)&array)); 47c4762a1bSJed Brown } 48c4762a1bSJed Brown 499566063dSJacob Faibussowitsch PetscCall(DMSwarmCreateGlobalVectorFromField(dms, "viscosity", &x)); 509566063dSJacob Faibussowitsch PetscCall(VecView(x, PETSC_VIEWER_STDOUT_WORLD)); 519566063dSJacob Faibussowitsch PetscCall(DMSwarmDestroyGlobalVectorFromField(dms, "viscosity", &x)); 52c4762a1bSJed Brown 539566063dSJacob Faibussowitsch PetscCall(DMSwarmCreateGlobalVectorFromField(dms, "strain", &x)); 549566063dSJacob Faibussowitsch PetscCall(VecView(x, PETSC_VIEWER_STDOUT_WORLD)); 559566063dSJacob Faibussowitsch PetscCall(DMSwarmDestroyGlobalVectorFromField(dms, "strain", &x)); 56c4762a1bSJed Brown 579566063dSJacob Faibussowitsch PetscCall(DMSwarmVectorDefineField(dms, "strain")); 589566063dSJacob Faibussowitsch PetscCall(DMCreateGlobalVector(dms, &x)); 599566063dSJacob Faibussowitsch PetscCall(VecView(x, PETSC_VIEWER_STDOUT_WORLD)); 609566063dSJacob Faibussowitsch PetscCall(VecDestroy(&x)); 619566063dSJacob Faibussowitsch PetscCall(DMDestroy(&dms)); 62c4762a1bSJed Brown PetscFunctionReturn(0); 63c4762a1bSJed Brown } 64c4762a1bSJed Brown 659371c9d4SSatish Balay int main(int argc, char **argv) { 66327415f7SBarry Smith PetscFunctionBeginUser; 679566063dSJacob Faibussowitsch PetscCall(PetscInitialize(&argc, &argv, (char *)0, help)); 689566063dSJacob Faibussowitsch PetscCall(ex2_1()); 699566063dSJacob Faibussowitsch PetscCall(PetscFinalize()); 70b122ec5aSJacob Faibussowitsch return 0; 71c4762a1bSJed Brown } 72c4762a1bSJed Brown 73c4762a1bSJed Brown /*TEST 74c4762a1bSJed Brown 75c4762a1bSJed Brown test: 76c4762a1bSJed Brown requires: !complex double 77c4762a1bSJed Brown nsize: 3 78c4762a1bSJed Brown filter: grep -v atomic 79c4762a1bSJed Brown filter_output: grep -v atomic 80c4762a1bSJed Brown 81c4762a1bSJed Brown TEST*/ 82