Twin-composite hypothesis: count the spots per class of hkl mod n, not per pair of forms
Whether a reflection lies on an index-n sub-lattice depends only on its indices mod n, so TwinCompositeAt histograms the validation spots over the n^3 classes once and counts each congruent pair of forms over the classes on the first form's sub-lattice. The counts are the same integers as before, so the decision and the log are unchanged; the index-2..10 loop drops from ~0.45 s to ~10 ms a cell, and it runs for the standing cell and every multiple in every first pass (on a median open set ~3 s of the first passes, single-threaded). Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01UqXYdL4qBGYRiSXSGXfCLu
This commit is contained in:
1 parent
6f8c77a6ef
commit
79a80f086c
1 file changed
+21
-4
@@ -352,6 +352,21 @@ std::optional<TwinCompositeEvidence> TwinCompositeAt(const CrystalLattice &prim,
|
||||
return true;
|
||||
return false;
|
||||
};
|
||||
// Whether a reflection is on a sub-lattice depends only on its indices mod n, so the spots are
|
||||
// counted once per class of (h, k, l) mod n, and the forms are asked of the classes.
|
||||
const int n_classes = n * n * n;
|
||||
const auto class_hkl = [n](int c) { return std::array<int, 3>{c / (n * n), c / n % n, c % n}; };
|
||||
std::vector<int> per_class(n_classes, 0);
|
||||
for (const auto &m : hkl)
|
||||
per_class[(mod(m[0]) * n + mod(m[1])) * n + mod(m[2])]++;
|
||||
std::vector<std::vector<int>> classes_on(forms.size()); // the classes on each form's sub-lattice
|
||||
std::vector<int> spots_on(forms.size(), 0);
|
||||
for (size_t i = 0; i < forms.size(); i++)
|
||||
for (int c = 0; c < n_classes; c++)
|
||||
if (on(forms[i], class_hkl(c))) {
|
||||
classes_on[i].push_back(c);
|
||||
spots_on[i] += per_class[c];
|
||||
}
|
||||
// The congruent pair whose union holds the most spots, the stronger domain first.
|
||||
std::optional<TwinCompositeEvidence> best;
|
||||
for (size_t i = 0; i < forms.size(); i++)
|
||||
@@ -359,10 +374,12 @@ std::optional<TwinCompositeEvidence> TwinCompositeAt(const CrystalLattice &prim,
|
||||
if (!congruent(cells[i], cells[j]))
|
||||
continue;
|
||||
TwinCompositeEvidence t{.f1 = forms[i], .f2 = forms[j]};
|
||||
for (const auto &m : hkl) {
|
||||
const bool o1 = on(forms[i], m), o2 = on(forms[j], m);
|
||||
(o1 && o2 ? t.both : o1 ? t.only1 : o2 ? t.only2 : t.neither)++;
|
||||
}
|
||||
for (const int c : classes_on[i])
|
||||
if (on(forms[j], class_hkl(c)))
|
||||
t.both += per_class[c];
|
||||
t.only1 = spots_on[i] - t.both;
|
||||
t.only2 = spots_on[j] - t.both;
|
||||
t.neither = static_cast<int>(hkl.size()) - spots_on[i] - spots_on[j] + t.both;
|
||||
if (!best || t.neither < best->neither)
|
||||
best = t;
|
||||
}
|
||||
|
||||
Reference in new issue
Block a user