summaryrefslogtreecommitdiff
diff options
context:
space:
mode:
-rw-r--r--Makefile2
-rw-r--r--uniform.cpp115
2 files changed, 113 insertions, 4 deletions
diff --git a/Makefile b/Makefile
index 33c8b6c..31d3891 100644
--- a/Makefile
+++ b/Makefile
@@ -4,7 +4,7 @@ BLOSSOM_DIRS := $(BLOSSOM_DIR) $(BLOSSOM_DIR)/MinCost $(BLOSSOM_DIR)/GEOM
BLOSSOM_SOURCES := $(filter-out $(BLOSSOM_DIR)/example.cpp, $(foreach dir, $(BLOSSOM_DIRS), $(wildcard $(dir)/*.cpp)))
BLOSSOM_OBJS := $(patsubst %.cpp, %.o, $(BLOSSOM_SOURCES))
-CFLAGS := -flto -fopenmp -O3 -Os -march=native -D_NDEBUG -DPERFECT_MATCHING_DOUBLE
+CFLAGS := -flto -fopenmp -O3 -Os -march=native -mtune=native -D_NDEBUG -DPERFECT_MATCHING_DOUBLE
CXX = clang++
LD = ld.lld
LIBS := -lrt
diff --git a/uniform.cpp b/uniform.cpp
index 07d794e..beaf038 100644
--- a/uniform.cpp
+++ b/uniform.cpp
@@ -38,13 +38,13 @@ public:
edges[i].index = i;
}
for (unsigned i = 1; i <= n; i++) {
- faces[i - 1].reserve(pow(i, 2));
+ faces[n - i].reserve(pow(i, 2));
for (unsigned j = 0; j < pow(i, 2); j++) {
unsigned x = j % (i);
unsigned y = j / (i);
unsigned x0 = n - i;
unsigned y0 = n - i;
- faces[i - 1].push_back(Face(
+ faces[n - i].push_back(Face(
edges[2 * n * (y0 + 2 * y) + x0 + 2 * x],
edges[2 * n * (y0 + 2 * y) + x0 + 2 * x + 1],
edges[2 * n * (y0 + 2 * y + 1) + x0 + 2 * x],
@@ -53,6 +53,62 @@ public:
}
}
}
+
+ void computeWeights() {
+ for (std::vector<Face>& fs : faces) {
+ for (Face& f : fs) {
+ double w = f.edges[0].get().weights.back();
+ double x = f.edges[1].get().weights.back();
+ double y = f.edges[2].get().weights.back();
+ double z = f.edges[3].get().weights.back();
+
+ double cellFactor = w * z + x * y;
+
+ f.edges[0].get().weights.push_back(z / cellFactor);
+ f.edges[1].get().weights.push_back(y / cellFactor);
+ f.edges[2].get().weights.push_back(x / cellFactor);
+ f.edges[3].get().weights.push_back(w / cellFactor);
+ }
+ }
+
+ for (Edge& e : edges) {
+ e.weights.pop_back();
+ }
+ }
+
+ void computeProbabilities() {
+ for (auto it = faces.rbegin(); it != faces.rend(); it++) {
+ for (Face& f : *it) {
+ double p = f.edges[0].get().probability;
+ double q = f.edges[1].get().probability;
+ double r = f.edges[2].get().probability;
+ double s = f.edges[3].get().probability;
+
+ double w = f.edges[0].get().weights.back();
+ double x = f.edges[1].get().weights.back();
+ double y = f.edges[2].get().weights.back();
+ double z = f.edges[3].get().weights.back();
+
+ std::cout << w << " " << x << " " << y << " " << z << std::endl;
+
+ double cellFactor = w * z + x * y;
+ double deficit = 1 - p - q - r - s;
+
+ f.edges[0].get().probability = s + deficit * w * z / cellFactor;
+ f.edges[1].get().probability = r + deficit * x * y / cellFactor;
+ f.edges[2].get().probability = q + deficit * x * y / cellFactor;
+ f.edges[3].get().probability = p + deficit * w * z / cellFactor;
+
+ for (Edge& e : f.edges) {
+ if (e.weights.empty()) {
+ std::cout << "No weights left!" << std::endl;
+ } else {
+ e.weights.pop_back();
+ }
+ }
+ }
+ }
+ }
};
int main(int argc, char* argv[]) {
@@ -75,7 +131,7 @@ int main(int argc, char* argv[]) {
}
Rng r;
- AztecDiamond a(n);
+// AztecDiamond a(n);
@@ -93,5 +149,58 @@ int main(int argc, char* argv[]) {
}
*/
+ AztecDiamond a(3);
+ for (Edge& e : a.edges) {
+ e.weights.push_back(1);
+ }
+
+ a.edges[1].weights.front() = 0;
+ a.edges[4].weights.front() = 0;
+ a.edges[6].weights.front() = 0;
+ a.edges[11].weights.front() = 0;
+ a.edges[24].weights.front() = 0;
+ a.edges[29].weights.front() = 0;
+ a.edges[31].weights.front() = 0;
+ a.edges[34].weights.front() = 0;
+
+ a.computeWeights();
+
+ for (unsigned i = 0; i < 3; i++) {
+ for (unsigned j = 0; j < a.edges.size(); j++) {
+ if (j % 6 == 0) {
+ std::cout << std::endl;
+ }
+ if (!a.edges[j].weights.empty()) {
+ std::cout << a.edges[j].weights.front() << " ";
+ a.edges[j].weights.pop_front();
+ }
+ }
+ }
+ std::cout << std::endl;
+
+ for (Edge& e : a.edges) {
+ e.weights.push_back(1);
+ }
+
+ a.edges[1].weights.front() = 0;
+ a.edges[4].weights.front() = 0;
+ a.edges[6].weights.front() = 0;
+ a.edges[11].weights.front() = 0;
+ a.edges[24].weights.front() = 0;
+ a.edges[29].weights.front() = 0;
+ a.edges[31].weights.front() = 0;
+ a.edges[34].weights.front() = 0;
+
+ a.computeWeights();
+ a.computeProbabilities();
+
+ for (unsigned i = 0; i < a.edges.size(); i++) {
+ if (i%6 == 0) {
+ std::cout << std::endl;
+ }
+ std::cout << a.edges[i].probability << " ";
+ }
+ std::cout << std::endl;
+
return 0;
}