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 i, 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 (i = info.xs; i < info.xs + info.xm; i++) {
63: if (has_glider && i == glider_loc[0] && j == glider_loc[1]) {
64: PetscInt ii, jj;
65: for (ii = -1; ii <= 1; ii++)
66: for (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: PetscInt step;
85: for (step = 0; step < steps; step++) {
86: const PetscScalar **x;
87: PetscScalar **y;
88: DMDALocalInfo info;
89: PetscInt i, j;
91: PetscCall(DMGlobalToLocal(da, Xglobal, INSERT_VALUES, Xlocal));
92: PetscCall(DMDAGetLocalInfo(da, &info));
93: PetscCall(DMDAVecGetArrayRead(da, Xlocal, (void *)&x));
94: PetscCall(DMDAVecGetArrayWrite(da, Xglobal, &y));
95: for (j = info.ys; j < info.ys + info.ym; j++) {
96: for (i = info.xs; i < info.xs + info.xm; i++) {
97: PetscInt live_neighbors = 0;
98: live_neighbors += PetscRealPart(x[j - 1][i - 1]) > 0;
99: live_neighbors += PetscRealPart(x[j - 1][i]) > 0;
100: live_neighbors += PetscRealPart(x[j - 1][i + 1]) > 0;
101: live_neighbors += PetscRealPart(x[j][i - 1]) > 0;
102: live_neighbors += PetscRealPart(x[j][i + 1]) > 0;
103: live_neighbors += PetscRealPart(x[j + 1][i - 1]) > 0;
104: live_neighbors += PetscRealPart(x[j + 1][i]) > 0;
105: live_neighbors += PetscRealPart(x[j + 1][i + 1]) > 0;
106: if (PetscRealPart(x[j][i]) > 0) { /* Live cell */
107: switch (live_neighbors) {
108: case 2:
109: case 3:
110: y[j][i] = 1; /* Survive */
111: break;
112: default:
113: y[j][i] = 0; /* Death */
114: }
115: } else { /* Dead cell */
116: if (live_neighbors == 3) y[j][i] = 1; /* Birth */
117: else y[j][i] = 0;
118: }
119: }
120: }
121: PetscCall(DMDAVecRestoreArrayRead(da, Xlocal, (void *)&x));
122: PetscCall(DMDAVecRestoreArrayWrite(da, Xglobal, &y));
123: if (step == check_step_alive || step == check_step_dead) {
124: PetscScalar sum;
125: PetscCall(VecSum(Xglobal, &sum));
126: if (PetscAbsScalar(sum) > 0.1) {
127: if (step == check_step_dead) PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Simulation alive at step %" PetscInt_FMT "\n", step));
128: } else if (step == check_step_alive) {
129: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Simulation dead at step %" PetscInt_FMT "\n", step));
130: }
131: }
132: if (step % viz_interval == 0) PetscCall(VecView(Xglobal, viewer));
133: }
134: }
136: PetscCall(PetscViewerDestroy(&viewer));
137: PetscCall(VecDestroy(&Xglobal));
138: PetscCall(VecDestroy(&Xlocal));
139: PetscCall(DMDestroy(&da));
140: PetscCall(PetscFinalize());
141: return 0;
142: }
144: /*TEST
146: test:
147: requires: x
148: nsize: 2
149: 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
151: TEST*/