Actual source code: swarm_ex2.c
2: static char help[] = "Tests DMSwarm\n\n";
4: #include <petscdm.h>
5: #include <petscdmda.h>
6: #include <petscdmswarm.h>
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;
18: PetscFunctionBegin;
19: PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
21: PetscCall(DMCreate(PETSC_COMM_WORLD, &dms));
22: PetscCall(DMSetType(dms, DMSWARM));
23: PetscCall(PetscObjectSetName((PetscObject)dms, "Particles"));
24: PetscCall(DMSwarmInitializeFieldRegister(dms));
25: PetscCall(DMSwarmRegisterPetscDatatypeField(dms, "viscosity", 1, PETSC_REAL));
26: PetscCall(DMSwarmRegisterPetscDatatypeField(dms, "strain", 3, PETSC_REAL));
27: PetscCall(DMSwarmFinalizeFieldRegister(dms));
28: PetscCall(DMSwarmSetLocalSizes(dms, 5 + rank, 4));
29: PetscCall(DMView(dms, PETSC_VIEWER_STDOUT_WORLD));
30: PetscCall(DMSwarmGetLocalSize(dms, &nlocal));
32: {
33: PetscReal *array;
34: PetscCall(DMSwarmGetField(dms, "viscosity", &bs, NULL, (void **)&array));
35: for (p = 0; p < nlocal; p++) array[p] = 11.1 + p * 0.1 + rank * 100.0;
36: PetscCall(DMSwarmRestoreField(dms, "viscosity", &bs, NULL, (void **)&array));
37: }
39: {
40: PetscReal *array;
41: PetscCall(DMSwarmGetField(dms, "strain", &bs, NULL, (void **)&array));
42: for (p = 0; p < nlocal; p++) {
43: array[bs * p + 0] = 2.0e-2 + p * 0.001 + rank * 1.0;
44: array[bs * p + 1] = 2.0e-2 + p * 0.002 + rank * 1.0;
45: array[bs * p + 2] = 2.0e-2 + p * 0.003 + rank * 1.0;
46: }
47: PetscCall(DMSwarmRestoreField(dms, "strain", &bs, NULL, (void **)&array));
48: }
50: PetscCall(DMSwarmCreateGlobalVectorFromField(dms, "viscosity", &x));
51: PetscCall(VecView(x, PETSC_VIEWER_STDOUT_WORLD));
52: PetscCall(DMSwarmDestroyGlobalVectorFromField(dms, "viscosity", &x));
54: PetscCall(DMSwarmCreateGlobalVectorFromField(dms, "strain", &x));
55: PetscCall(VecView(x, PETSC_VIEWER_STDOUT_WORLD));
56: PetscCall(DMSwarmDestroyGlobalVectorFromField(dms, "strain", &x));
58: PetscCall(DMSwarmVectorDefineField(dms, "strain"));
59: PetscCall(DMCreateGlobalVector(dms, &x));
60: PetscCall(VecView(x, PETSC_VIEWER_STDOUT_WORLD));
61: PetscCall(VecDestroy(&x));
62: PetscCall(DMDestroy(&dms));
63: PetscFunctionReturn(PETSC_SUCCESS);
64: }
66: int main(int argc, char **argv)
67: {
68: PetscFunctionBeginUser;
69: PetscCall(PetscInitialize(&argc, &argv, (char *)0, help));
70: PetscCall(ex2_1());
71: PetscCall(PetscFinalize());
72: return 0;
73: }
75: /*TEST
77: test:
78: requires: !complex double
79: nsize: 3
80: filter: grep -v atomic
81: filter_output: grep -v atomic
83: TEST*/