rugnux --model: evaluate the rigid-body Jacobian's columns in parallel

The rigid-body target is an Fcalc + solvent mask + scale re-fit per
evaluation, and its forward-difference Jacobian made six of them one after the
other after the central one - 7 of the ~8 evaluations per LM iteration.

The six shifted evaluations now run in parallel, column j on its own Evaluator
over its own copy of the model (an evaluation moves every atom, so two cannot
share one), set to the same zone and given the central evaluation's bulk-solvent
pair - the pair the serial loop's shifted evaluations used, since the zone's
solvent is fitted on the zone's first evaluation and every Evaluate starts with
the central one. Each column's arithmetic is the serial one and the Jacobian is
assembled in column order, so the refinement is bit-identical; the evaluation
count is kept as the serial loop kept it (up to and including a failing column).
The per-zone solvent fit (FitModelScale) gets the thread count as well.

Inside a null replicate (a pool worker) the columns run inline, as before.

To check: RIGID_BODY_* and MODEL_* keys and md5 of maps/.mtz/_model.cif identical
with and without this commit on the audit set; rigid-body seconds on a large
model.

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_013nW6FNRP1bBJJ8pfHiByAT
This commit is contained in:
2026-09-20 18:45:03 +02:00
co-authored by Claude Opus 5
parent 542a1a987d
commit 13de0e7a87
4 changed files with 95 additions and 14 deletions
+40
View File
@@ -622,3 +622,43 @@ TEST_CASE("ModelValidation_MapToFPhiMatchesGemmi", "[ModelValidation]") {
for (size_t i = 0; i < a.v.size(); i++)
CHECK(a.v[i].hkl == b.v[i].hkl);
}
// The Jacobian's six columns are evaluated in parallel; the placement must be the serial one, bit for
// bit.
TEST_CASE("ModelValidation_RigidBodySameOnAnyNumberOfThreads", "[ModelValidation]") {
Logger logger("ModelValidation_RigidBodySameOnAnyNumberOfThreads");
const auto path = WriteTemp("rigid_body_threads_test.pdb", ClusterPdb().c_str());
gemmi::Structure st = gemmi::read_structure_gz(path, gemmi::CoorFormat::Detect);
const gemmi::SpaceGroup *sg = st.find_spacegroup();
REQUIRE(sg != nullptr);
st.setup_cell_images();
const auto ref = ModelReferenceIntensities(path, {}, {}, 3.0, logger);
REQUIRE_FALSE(ref.empty());
gemmi::AsuData<gemmi::ValueSigma<float>> fobs;
fobs.unit_cell_ = st.cell;
fobs.spacegroup_ = sg;
for (const auto &r : ref)
fobs.v.push_back({{{r.h, r.k, r.l}}, {std::sqrt(r.I), 1.0f}});
fobs.ensure_sorted();
std::vector<gemmi::Position> displaced;
for (const gemmi::Position &p : ModelPositions(st.models[0]))
displaced.emplace_back(p.x + 0.40, p.y - 0.30, p.z + 0.20);
gemmi::Model serial = st.models[0], parallel = st.models[0];
SetModelPositions(serial, displaced);
SetModelPositions(parallel, displaced);
const RigidBodyRefineResult r1 = RefineRigidBody(serial, st.cell, *sg, fobs, 3.0, logger, 1);
const RigidBodyRefineResult r6 = RefineRigidBody(parallel, st.cell, *sg, fobs, 3.0, logger, 6);
CHECK(r1.evaluations == r6.evaluations);
CHECK(r1.angle_deg == r6.angle_deg);
CHECK(r1.shift_A == r6.shift_A);
const auto p1 = ModelPositions(serial), p6 = ModelPositions(parallel);
REQUIRE(p1.size() == p6.size());
for (size_t i = 0; i < p1.size(); i++)
CHECK((p1[i].x == p6[i].x && p1[i].y == p6[i].y && p1[i].z == p6[i].z));
std::filesystem::remove(path);
}