Actual source code: ex2.c
1: static char help[] = "Demonstrates Conway's Game of Life using a 2d DMDA.\n\n";
3: /*
4: At each step in time, the following transitions occur:
6: Any live cell with fewer than two live neighbours dies, as if by underpopulation.
7: Any live cell with two or three live neighbours lives on to the next generation.
8: Any live cell with more than three live neighbours dies, as if by overpopulation.
9: Any dead cell with exactly three live neighbours becomes a live cell, as if by reproduction.
11: https://en.wikipedia.org/wiki/Conway%27s_Game_of_Life
12: */
14: #include <petscdm.h>
15: #include <petscdmda.h>
17: static const int GLIDER[3][3] = {
18: {0, 1, 0},
19: {0, 1, 1},
20: {1, 0, 1}
21: };
23: int main(int argc, char **argv)
24: {
25: DM da;
26: PetscViewer viewer;
27: Vec Xlocal, Xglobal;
28: PetscInt glider_loc[2] = {10, 20}, blinker_loc[2] = {20, 10}, two, steps = 100, viz_interval = 1;
29: PetscInt check_step_alive = -1, check_step_dead = -1;
30: PetscBool has_glider, has_blinker;
32: PetscFunctionBeginUser;
33: PetscCall(PetscInitialize(&argc, &argv, NULL, help));
34: PetscOptionsBegin(PETSC_COMM_WORLD, NULL, "Conway's Game of Life", "");
35: {
36: PetscCall(PetscOptionsIntArray("-glider", "Coordinate at which to center a glider", NULL, glider_loc, (two = 2, &two), &has_glider));
37: PetscCall(PetscOptionsIntArray("-blinker", "Coordinate at which to center a blinker", NULL, blinker_loc, (two = 2, &two), &has_blinker));
38: PetscCall(PetscOptionsInt("-steps", "Number of steps to take", NULL, steps, &steps, NULL));
39: PetscCall(PetscOptionsInt("-viz_interval", "Visualization interval", NULL, viz_interval, &viz_interval, NULL));
40: PetscCall(PetscOptionsInt("-check_step_alive", "Step on which to check that the simulation is alive", NULL, check_step_alive, &check_step_alive, NULL));
41: PetscCall(PetscOptionsInt("-check_step_dead", "Step on which to check that the simulation is dead", NULL, check_step_dead, &check_step_dead, NULL));
42: }
43: PetscOptionsEnd();
45: PetscCall(PetscViewerDrawOpen(PETSC_COMM_WORLD, NULL, "Life", PETSC_DECIDE, PETSC_DECIDE, 1000, 1000, &viewer));
47: /* Create distributed array and get vectors */
48: PetscCall(DMDACreate2d(PETSC_COMM_WORLD, DM_BOUNDARY_PERIODIC, DM_BOUNDARY_PERIODIC, DMDA_STENCIL_BOX, 30, 30, PETSC_DECIDE, PETSC_DECIDE, 1, 1, NULL, NULL, &da));
49: PetscCall(DMSetFromOptions(da));
50: PetscCall(DMSetUp(da));
51: PetscCall(DMCreateLocalVector(da, &Xlocal));
52: PetscCall(DMCreateGlobalVector(da, &Xglobal));
54: { /* Initialize */
55: DMDALocalInfo info;
56: PetscScalar **x;
57: PetscInt j;
59: PetscCall(DMDAGetLocalInfo(da, &info));
60: PetscCall(DMDAVecGetArray(da, Xlocal, &x));
61: for (j = info.ys; j < info.ys + info.ym; j++) {
62: for (PetscInt i = info.xs; i < info.xs + info.xm; i++) {
63: if (has_glider && i == glider_loc[0] && j == glider_loc[1]) {
64: PetscInt ii;
65: for (ii = -1; ii <= 1; ii++)
66: for (PetscInt jj = -1; jj <= 1; jj++) x[j + jj][i + ii] = GLIDER[1 - jj][ii + 1];
67: }
68: if (has_blinker && i == blinker_loc[0] && j == blinker_loc[1]) {
69: x[j - 1][i] = 1;
70: x[j][i] = 1;
71: x[j + 1][i] = 1;
72: }
73: }
74: }
75: PetscCall(DMDAVecRestoreArray(da, Xlocal, &x));
76: PetscCall(DMLocalToGlobal(da, Xlocal, ADD_VALUES, Xglobal));
77: }
79: /* View the initial condition */
80: PetscCall(VecView(Xglobal, viewer));
82: { /* Play */
83: for (PetscInt step = 0; step < steps; step++) {
84: const PetscScalar **x;
85: PetscScalar **y;
86: DMDALocalInfo info;
88: PetscCall(DMGlobalToLocal(da, Xglobal, INSERT_VALUES, Xlocal));
89: PetscCall(DMDAGetLocalInfo(da, &info));
90: PetscCall(DMDAVecGetArrayRead(da, Xlocal, (void *)&x));
91: PetscCall(DMDAVecGetArrayWrite(da, Xglobal, &y));
92: for (PetscInt j = info.ys; j < info.ys + info.ym; j++) {
93: for (PetscInt i = info.xs; i < info.xs + info.xm; i++) {
94: PetscInt live_neighbors = 0;
95: live_neighbors += PetscRealPart(x[j - 1][i - 1]) > 0;
96: live_neighbors += PetscRealPart(x[j - 1][i]) > 0;
97: live_neighbors += PetscRealPart(x[j - 1][i + 1]) > 0;
98: live_neighbors += PetscRealPart(x[j][i - 1]) > 0;
99: live_neighbors += PetscRealPart(x[j][i + 1]) > 0;
100: live_neighbors += PetscRealPart(x[j + 1][i - 1]) > 0;
101: live_neighbors += PetscRealPart(x[j + 1][i]) > 0;
102: live_neighbors += PetscRealPart(x[j + 1][i + 1]) > 0;
103: if (PetscRealPart(x[j][i]) > 0) { /* Live cell */
104: switch (live_neighbors) {
105: case 2:
106: case 3:
107: y[j][i] = 1; /* Survive */
108: break;
109: default:
110: y[j][i] = 0; /* Death */
111: }
112: } else { /* Dead cell */
113: if (live_neighbors == 3) y[j][i] = 1; /* Birth */
114: else y[j][i] = 0;
115: }
116: }
117: }
118: PetscCall(DMDAVecRestoreArrayRead(da, Xlocal, (void *)&x));
119: PetscCall(DMDAVecRestoreArrayWrite(da, Xglobal, &y));
120: if (step == check_step_alive || step == check_step_dead) {
121: PetscScalar sum;
122: PetscCall(VecSum(Xglobal, &sum));
123: if (PetscAbsScalar(sum) > 0.1) {
124: if (step == check_step_dead) PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Simulation alive at step %" PetscInt_FMT "\n", step));
125: } else if (step == check_step_alive) {
126: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Simulation dead at step %" PetscInt_FMT "\n", step));
127: }
128: }
129: if (step % viz_interval == 0) PetscCall(VecView(Xglobal, viewer));
130: }
131: }
133: PetscCall(PetscViewerDestroy(&viewer));
134: PetscCall(VecDestroy(&Xglobal));
135: PetscCall(VecDestroy(&Xlocal));
136: PetscCall(DMDestroy(&da));
137: PetscCall(PetscFinalize());
138: return 0;
139: }
141: /*TEST
143: test:
144: requires: x
145: nsize: 2
146: args: -glider 5,6 -blinker 12,12 -steps 35 -check_step_alive 31 -check_step_dead 32 -da_grid_x 20 -da_grid_y 20 -nox
147: output_file: output/empty.out
149: TEST*/