RosettaCodeData/Task/Percolation-Mean-cluster-density/D/percolation-mean-cluster-density.d
2015-02-20 09:02:09 -05:00

90 lines
2.7 KiB
D

import std.stdio, std.algorithm, std.random, std.math, std.array,
std.range, std.ascii;
alias Cell = ubyte;
alias Grid = Cell[][];
enum Cell notClustered = 1; // Filled cell, but not in a cluster.
Grid initialize(Grid grid, in double prob, ref Xorshift rng) nothrow {
foreach (row; grid)
foreach (ref cell; row)
cell = Cell(rng.uniform01 < prob);
return grid;
}
void show(in Grid grid) {
immutable static cell2char = " #" ~ letters;
writeln('+', "-".replicate(grid.length), '+');
foreach (row; grid) {
write('|');
row.map!(c => c < cell2char.length ? cell2char[c] : '@').write;
writeln('|');
}
writeln('+', "-".replicate(grid.length), '+');
}
size_t countClusters(bool justCount=false)(Grid grid)
pure nothrow @safe @nogc {
immutable side = grid.length;
static if (justCount)
enum Cell clusterID = 2;
else
Cell clusterID = 1;
void walk(in size_t r, in size_t c) nothrow @safe @nogc {
grid[r][c] = clusterID; // Fill grid.
if (r < side - 1 && grid[r + 1][c] == notClustered) // Down.
walk(r + 1, c);
if (c < side - 1 && grid[r][c + 1] == notClustered) // Right.
walk(r, c + 1);
if (c > 0 && grid[r][c - 1] == notClustered) // Left.
walk(r, c - 1);
if (r > 0 && grid[r - 1][c] == notClustered) // Up.
walk(r - 1, c);
}
size_t nClusters = 0;
foreach (immutable r; 0 .. side)
foreach (immutable c; 0 .. side)
if (grid[r][c] == notClustered) {
static if (!justCount)
clusterID++;
nClusters++;
walk(r, c);
}
return nClusters;
}
double clusterDensity(Grid grid, in double prob, ref Xorshift rng) {
return grid.initialize(prob, rng).countClusters!true /
double(grid.length ^^ 2);
}
void showDemo(in size_t side, in double prob, ref Xorshift rng) {
auto grid = new Grid(side, side);
grid.initialize(prob, rng);
writefln("Found %d clusters in this %d by %d grid:\n",
grid.countClusters, side, side);
grid.show;
}
void main() {
immutable prob = 0.5;
immutable nIters = 5;
auto rng = Xorshift(unpredictableSeed);
showDemo(15, prob, rng);
writeln;
foreach (immutable i; iota(4, 14, 2)) {
immutable side = 2 ^^ i;
auto grid = new Grid(side, side);
immutable density = nIters
.iota
.map!(_ => grid.clusterDensity(prob, rng))
.sum / nIters;
writefln("n_iters=%3d, p=%4.2f, n=%5d, sim=%7.8f",
nIters, prob, side, density);
}
}