diff options
| -rw-r--r-- | percolation.cpp | 10 |
1 files changed, 7 insertions, 3 deletions
diff --git a/percolation.cpp b/percolation.cpp index a5c2a36..f012b07 100644 --- a/percolation.cpp +++ b/percolation.cpp @@ -26,6 +26,7 @@ class Vertex { public: std::vector<HalfEdge> bonds; unsigned index; + std::vector<unsigned> coordinates; }; class Cluster : public std::list<std::reference_wrapper<const Vertex>> {}; @@ -60,8 +61,8 @@ public: unsigned squaredDistance(unsigned vi, unsigned vj) const { unsigned sd = 0; for (unsigned i = 0; i < D; i++) { - int x1 = (vi / uPow(L, i)) % uPow(L, i + 1); - int x2 = (vj / uPow(L, i)) % uPow(L, i + 1); + int x1 = vertices[vi].coordinates[i]; + int x2 = vertices[vj].coordinates[i]; unsigned Δx = abs(x1 - x2); if (Δx > L / 2) { Δx = L - Δx; @@ -72,11 +73,14 @@ public: } /* Initialize a square lattice */ - Graph(unsigned D, unsigned L) : D(D), L(L), vertices(uPow(L, D)), multiplicities(D * L * L / 4) { + Graph(unsigned D, unsigned L) : D(D), L(L), vertices(uPow(L, D)), multiplicities((D * L * L) / 4) { for (unsigned i = 0; i < Nv(); i++) { + vertices[i].coordinates.resize(D); vertices[i].index = i; vertices[i].bonds.reserve(2 * D); for (unsigned j = 0; j < D; j++) { + vertices[i].coordinates[j] = (i / uPow(L, j)) % L; + unsigned n1 = uPow(L, j + 1) * (i / uPow(L, j + 1)) + (i + uPow(L, j)) % uPow(L, j + 1); unsigned n2 = uPow(L, j + 1) * (i / uPow(L, j + 1)) + (uPow(L, j + 1) + i - uPow(L, j)) % uPow(L, j + 1); |