Percolation/Bond percolation: Difference between revisions

Content added Content deleted
(Rewritten the D entry on the base of the C entry)
Line 161: Line 161:
</pre>
</pre>
=={{header|D}}==
=={{header|D}}==
{{trans|C}}
{{incorrect|D| Results from p=0.3 to p=0.6 are statistically implausible, please check code. }}
<lang d>import std.stdio, std.random, std.array, std.range, std.algorithm;

{{trans|Python}}
Compared to the Python entry, the <code>initialize</code> function is a performance optimization useful when nTries is high.
<lang d>import std.stdio, std.random, std.array, std.algorithm, std.range,
std.typecons;


struct Grid {
struct Grid {
// Not enforced by runtime and type system:
immutable int nr, nc;
// a Cell must contain only the flags bits.
bool[][] hWalls, vWalls, cells;
alias MaybeCol = Nullable!(int, int.min);
alias Cell = uint;

Xorshift rng;
enum : Cell { // Cell states (bit flags).
empty = 0,
filled = 1,
rightWall = 2,
bottomWall = 4
}

const size_t nc, nr;
Cell[] cells;


this(in uint nRows, in uint nCols, in uint seed) pure nothrow {
this(in size_t nRows, in size_t nCols) pure nothrow {
nr = nRows;
nr = nRows;
nc = nCols;
nc = nCols;
hWalls = new typeof(hWalls)(nr + 1, nc);
vWalls = new typeof(vWalls)(nr, nc + 1);
cells = new typeof(cells)(nr, nc);
rng.seed = seed;
}


// Allocate two addition rows to avoid checking bounds.
void initialize(in double probability) {
foreach (ref row; hWalls)
// Bottom row is also required by drippage.
foreach (ref x; row)
cells = new Cell[nc * (nr + 2)];
x = uniform(0.0, 1.0, rng) < probability;
foreach (ref row; vWalls)
foreach (immutable c, ref x; row)
x = (c == 0 || c == nc) ?
true :
(uniform(0.0, 1.0, rng) < probability);
foreach (ref row; cells)
row[] = false;
}
}


void show(in MaybeCol percolatedCol) const {
void initialize(in double prob, ref Xorshift rng) {
cells[0 .. nc] = bottomWall | rightWall; // First row.
// Horiz, vert, fill chars.
static immutable h = ["+ ", "+-"],
v = [" ", "|"],
f = [" ", "#", "X"];


foreach (immutable r; 0 .. nr) {
uint pos = nc;
foreach (immutable r; 1 .. nr + 1) {
writefln("%-(%s%)+", nc.iota.map!(c => h[hWalls[r][c]]));
writefln("%-(%s%)", iota(nc + 1)
foreach (immutable c; 1 .. nc)
.map!(c => v[vWalls[r][c]] ~
cells[pos++] =
f[c<nc ? cells[r][c] : 0]));
(uniform(0.0, 1.0, rng) < prob ?
bottomWall :
empty) |
(uniform(0.0, 1.0, rng) < prob ?
rightWall :
empty);
cells[pos++] = rightWall |
(uniform(0.0, 1.0, rng) < prob ?
bottomWall :
empty);
}
}

writefln("%-(%s%)+", nc.iota.map!(c => h[hWalls[nr][c]]));
cells[$ - nc .. $] = empty; // Last row.
if (!percolatedCol.isNull)
writeln(" ".replicate(percolatedCol), " ",f[2]);
}
}


MaybeCol pourOnTop() pure nothrow {
bool percolate() pure nothrow {
MaybeCol floodFill(in int r, in int c) pure nothrow {
bool fill(in size_t i) pure nothrow {
cells[r][c] = true; // Fill cell.
if (cells[i] & filled)
return false;


// Bottom.
cells[i] |= filled;
if (r < nr - 1 && !hWalls[r + 1][c] && !cells[r + 1][c])
return floodFill(r + 1, c);
else if (r == nr - 1 && !hWalls[r + 1][c]) // The bottom.
return MaybeCol(c);


// Left.
if (i >= cells.length - nc)
if (c && !vWalls[r][c] && !cells[r][c - 1])
return true; // Success: reached bottom row.
return floodFill(r, c - 1);


return (!(cells[i] & bottomWall) && fill(i + nc)) ||
// Right.
if (c < nc - 1 && !vWalls[r][c + 1] && !cells[r][c + 1])
(!(cells[i] & rightWall) && fill(i + 1)) ||
return floodFill(r, c + 1);
(!(cells[i - 1] & rightWall) && fill(i - 1)) ||
(!(cells[i - nc] & bottomWall) && fill(i - nc));
}


// Top.
return iota(nc, nc + nc).any!fill;
}
if (r && !hWalls[r][c] && !cells[r - 1][c])
return floodFill(r - 1, c);


void show() const {
return MaybeCol();
writeln("+-".replicate(nc), '+');
}


immutable r = 0; // From top.
foreach (immutable r; 1 .. nr + 2) {
foreach (immutable c; 0 .. nc)
write(r == nr + 1 ? ' ' : '|');
if (!hWalls[r][c]) {
foreach (immutable c; 0 .. nc) {
immutable percolatedCol = floodFill(r, c);
immutable cell = cells[r * nc + c];
if (!percolatedCol.isNull)
write((cell & filled) ? (r <= nr ? '#' : 'X') : ' ');
return percolatedCol;
write((cell & rightWall) ? '|' : ' ');
}
}
return MaybeCol();
writeln;

if (r == nr + 1)
return;

foreach (immutable c; 0 .. nc)
write((cells[r * nc + c] & bottomWall) ? "+-" : "+ ");
'+'.writeln;
}
}
}
}
}


void main() {
void main() {
enum int nr = 10, nc = 10; // N. rows and columns of the grid.
enum uint nr = 10, nc = 10; // N. rows and columns of the grid.
enum int nTries = 1_000; // N. simulations for each probability.
enum uint nTries = 10_000; // N. simulations for each probability.
enum int nStepsProb = 10; // N. steps of probability.
enum uint nStepsProb = 10; // N. steps of probability.


bool sampleShown = false;
auto rng = Xorshift(2);
alias Pair = Tuple!(double,"prob", uint,"co");
auto g = Grid(nr, nc);
g.initialize(0.5, rng);
Pair[] results;
g.percolate;
g.show;


writefln("\nRunning %dx%d grids %d times for each p:",
foreach (immutable i; 0 .. nStepsProb + 1) {
// Count down so sample print is interesting.
nr, nc, nTries);
immutable probability = (nStepsProb - i) /
foreach (immutable p; 0 .. nStepsProb) {
cast(double)nStepsProb;
immutable probability = p / cast(double)nStepsProb;
auto grid = Grid(nr, nc, unpredictableSeed);
uint nPercolated = 0;
uint lastCount = 0;
foreach (immutable i; 0 .. nTries) {
g.initialize(probability, rng);

foreach (immutable _; 0 .. nTries) {
nPercolated += g.percolate;
grid.initialize(probability);
immutable percolatedCol = grid.pourOnTop;
if (!percolatedCol.isNull) {
lastCount++;
if (!sampleShown) {
writefln("Percolating sample %dx%d grid," ~
" probability=%1.1f:", nc,nr,probability);
grid.show(percolatedCol);
sampleShown = true;
}
}
}
}
results ~= Pair(probability, lastCount);
writefln("p = %0.2f: %.4f",
probability, nPercolated / cast(double)nTries);
}
}

writefln("\nFraction of %d tries that percolate," ~
" varying probability p:", nTries);
foreach (const r; results)
writefln("p=%1.2f: %1.3f", r.prob, r.co / cast(double)nTries);
}</lang>
}</lang>
{{out}}
{{out}}
<pre>+-+-+-+-+-+-+-+-+-+-+
<pre>Percolating sample 10x10 grid, probability=0.6:
|#|#|#|#| | | |
+ + + + + +-+-+-+-+-+
+ +-+-+ +-+-+-+ +-+-+
|#|# # # #| | | |
+ + + + + + + + +-+-+
|#| | # | | | | |
+ +-+-+ + +-+-+ + +-+
|# #|# #|# #| | | | |
|#|# #|#| | | |
+ + + +-+-+ +-+ +-+ +
|# #|#| | |#| |
+ +-+ + +-+ + + +-+ +
|#|# #|#| | | | |
+-+ +-+ +-+ +-+ + +-+
+-+ + + + +-+-+-+-+-+
| |#| # # # | | |
|# # # # # | | | |
+-+-+ + +-+-+-+ + + +
| | | |# #| | | | |
+ + + + + + + +-+ +-+
+-+ + +-+ +-+ + +-+-+
|#|# # #|# # # | |
| | |# #| | | | |
+-+ + + +-+-+ + + + +
| |#|# #| | |# |
+-+-+-+ +-+ + +-+ + +
+-+-+-+-+ +-+ +-+-+-+
| | | |# | | |
| | | # #| |
+-+-+ + +-+-+-+-+-+-+
+-+-+-+ +-+ +-+-+-+ +
| |# | | | |
+ +-+-+ + +-+ + + + +
| | | # # # |
+ + +-+ +-+-+-+ +-+ +
| #| | | | | |
+-+ + + +-+ + +-+ + +
| | | |# |
| | #| | | |
+ +-+ + + + +-+ + + +
X
+ + + + +-+-+-+ +-+ +
X


Running 10x10 grids 10000 times for each p:
Fraction of 1000 tries that percolate, varying probability p:
p=1.00: 0.000
p = 0.00: 1.0000
p=0.90: 0.000
p = 0.10: 1.0000
p=0.80: 0.000
p = 0.20: 1.0000
p=0.70: 0.000
p = 0.30: 0.9973
p=0.60: 0.007
p = 0.40: 0.9177
p=0.50: 0.142
p = 0.50: 0.5050
p=0.40: 0.573
p = 0.60: 0.0880
p=0.30: 0.929
p = 0.70: 0.0035
p=0.20: 0.998
p = 0.80: 0.0001
p=0.10: 1.000
p = 0.90: 0.0000</pre>
With LDC2 compiler this code is almost three times slower than the C entry.
p=0.00: 1.000</pre>


=={{header|Python}}==
=={{header|Python}}==