Refinement: the twinned arm runs beside the untwinned one

DecideTwin's twinned refinement reads only the placed model and the data, never the
untwinned result, so RefineModel now starts it on a thread of its own as soon as the
untwinned arm has placed its engine and slots, on the next card where there is more
than one. DecideTwin touches the shared result only once the untwinned arm is done,
and decides nothing where that arm failed or was cancelled. The twinned arm logs to a
buffered logger that is replayed after the untwinned arm's lines, so the log reads as
before. The per-stage TIMING sums of the two arms are no longer separated.

The engines are made exactly as they would be alone: each arm's GPU/CPU decision uses
the same budget as before, and the arms overlap only where twice what the untwinned
arm holds of its busiest card fits that budget (OwnSlotPlaces now reports it);
otherwise the twinned arm waits for the untwinned one, as before. Trial outcomes do
not depend on the slots, so the result is unchanged.

Measured on one RTX 5080, -N 28, two runs each (Step 11 wall, s):
  cytochrome C (P3x21, twin law -h,-k,l):  6.0 / 6.1 -> 4.5 / 4.6
  thaumatin, lysozyme (no twin law):       unchanged
Every output (model, maps, MTZ, report, refinement log lines) is identical to the
rc184-port binary on all three sets. With two or more cards the twinned arm gets a
card of its own, so the step should drop to about the longer arm (~3.3 s here).

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01HVeDQsU5a8AJUdJPgWevAK
This commit is contained in:
leonarski_fandClaude Opus 5.5 committed 2026-11-08 21:49:09 +01:00
1 parent 12d6208126
commit 40e2462bde
1 file changed
+95 -18
@@ -10,6 +10,7 @@
#include <future>
#include <limits>
#include <map>
#include <mutex>
#include <optional>
#include <cmath>
#include <numeric>
@@ -1633,8 +1634,11 @@ struct SlotPlace {
// puts the long trials. A slot takes its engine's own part and a lane of its own where both fit, else its engine
// alone on a lane it shares (one evaluation at a time: lanes side by side keep the card busier, more slots overlap
// more host work); the first one on another card also the part the engines there share and a scale pool. On the
// CPU, the main one and those beside it are no more than the threads.
std::vector<SlotPlace> OwnSlotPlaces(const Context &main, int wanted, const RefinementSettings &settings) {
// CPU, the main one and those beside it are no more than the threads. `busiest`: what the main engine and the slots
// hold of the card that holds the most of them (0 on the CPU).
std::vector<SlotPlace> OwnSlotPlaces(const Context &main, int wanted, const RefinementSettings &settings,
size_t &busiest) {
busiest = 0;
#ifdef JFJOCH_USE_CUDA
if (main.engine->OnGPU()) {
const int n = get_gpu_count();
@@ -1656,6 +1660,7 @@ std::vector<SlotPlace> OwnSlotPlaces(const Context &main, int wanted, const Refi
places.push_back({d, false});
}
}
busiest = *std::max_element(used.begin(), used.end());
return places;
}
#endif
@@ -2237,10 +2242,20 @@ double RefinementDMin(double max_inv_d2) {
// One refinement of a model, with the solvent build in it where settings.solvent asks for it.
// `twin_laws`: against the data as twinned by them (one domain each beside the model's own), empty for untwinned;
// `null_first`: the first of them the null (TwinNull) replaces the mates of.
// What a refinement holds of the card that holds the most of it (OwnSlotPlaces), said once (RefineModel).
struct CardHeld {
std::promise<size_t> bytes;
std::once_flag said;
void Say(size_t b) {
std::call_once(said, [&] { bytes.set_value(b); });
}
};
RefinementResult RefineOnce(const gemmi::Structure &placed, const std::vector<MergedReflection> &merged,
const gemmi::UnitCell &cell, const gemmi::SpaceGroup &sg, double wavelength_A,
Logger &logger, const RefinementSettings &settings,
const std::vector<gemmi::Op> &twin_laws = {}, size_t null_first = 0) {
const std::vector<gemmi::Op> &twin_laws = {}, size_t null_first = 0,
CardHeld *held = nullptr) {
const auto start = std::chrono::steady_clock::now();
RefinementResult result;
Context ctx;
@@ -2485,8 +2500,10 @@ RefinementResult RefineOnce(const gemmi::Structure &placed, const std::vector<Me
settings.final_ladder.size(), size_t{3}}));
if (settings.max_trials_at_once > 0)
wanted = std::min(wanted, static_cast<int>(settings.max_trials_at_once));
const std::vector<SlotPlace> devices = settings.restrained ? OwnSlotPlaces(ctx, wanted - 1, settings)
: std::vector<SlotPlace>();
size_t busiest = 0;
const std::vector<SlotPlace> devices = OwnSlotPlaces(ctx, settings.restrained ? wanted - 1 : 0, settings, busiest);
if (held != nullptr)
held->Say(busiest);
const size_t slot_threads = std::max<size_t>(1, settings.nthreads / (devices.size() + 1));
ScalePools pools;
#ifdef JFJOCH_USE_CUDA
@@ -3297,15 +3314,31 @@ double FinalRWork(const RefinementResult &r) {
// standard errors, each added domain at a fraction of TWIN_DOMAIN_MIN at least.
void DecideTwin(const gemmi::Structure &placed, const std::vector<MergedReflection> &merged,
const gemmi::UnitCell &cell, const gemmi::SpaceGroup &sg, double wavelength_A, Logger &logger,
const RefinementSettings &settings, RefinementResult &result) {
const RefinementSettings &settings, RefinementResult &result, const std::shared_future<void> &untwinned) {
const RefinementStage stage("twin");
result.untwinned_r_work = FinalRWork(result);
result.untwinned_r_free = FinalRFree(result);
const std::vector<gemmi::Op> laws = RefinementTwinLaws(cell, sg);
const std::vector<gemmi::Op> domains = RefinementTwinDomains(cell, sg);
result.twin_laws = static_cast<int>(laws.size());
// The twinned arm reads only the placed model and the data, so it runs beside the untwinned refinement
// (RefineModel): `result` is read and written only once `untwinned` says that refinement is in it, and
// nothing is decided where it failed or was cut short.
std::string split_law;
int split_pairs = 0;
auto take_untwinned = [&]() {
untwinned.wait();
if (!result.ok || result.cancelled)
return false;
result.untwinned_r_work = FinalRWork(result);
result.untwinned_r_free = FinalRFree(result);
result.twin_laws = static_cast<int>(laws.size());
if (!split_law.empty()) {
result.twin_law = split_law;
result.twin_split_pairs = split_pairs;
}
return true;
};
if (laws.empty()) {
result.twin_decision = "NO_TWIN_LAW";
if (take_untwinned())
result.twin_decision = "NO_TWIN_LAW";
return;
}
// the laws whose twin pairs the free set keeps together
@@ -3318,15 +3351,16 @@ void DecideTwin(const gemmi::Structure &placed, const std::vector<MergedReflecti
else {
logger.Info("Refinement, twin law {}: not tested - {} of its twin pairs are split between the working and "
"the free set", law.as_hkl().triplet('h'), split);
result.twin_law = law.as_hkl().triplet('h');
result.twin_split_pairs = split;
split_law = law.as_hkl().triplet('h');
split_pairs = split;
}
}
if (usable.empty()) {
result.twin_decision = "FREE_SET_SPLIT";
if (take_untwinned())
result.twin_decision = "FREE_SET_SPLIT";
return;
}
result.twin_split_pairs = 0;
split_law.clear();
gemmi::Op law = usable[0];
if (usable.size() > 1) {
RefinementSettings score = settings;
@@ -3342,8 +3376,11 @@ void DecideTwin(const gemmi::Structure &placed, const std::vector<MergedReflecti
}
}
}
result.twin_law = law.as_hkl().triplet('h');
RefinementResult twinned = RefineOnce(placed, merged, cell, sg, wavelength_A, logger, settings, {law});
if (!take_untwinned())
return;
result.twin_split_pairs = 0;
result.twin_law = law.as_hkl().triplet('h');
if (!twinned.ok || twinned.cancelled) {
result.twin_decision = "UNTWINNED";
return;
@@ -3513,12 +3550,52 @@ RefinementResult RefineModel(const gemmi::Structure &placed, const std::vector<M
"reflections in one set - {} reflections changed set", moved);
}
const std::vector<MergedReflection> &data = moved > 0 ? closed : merged;
RefinementResult result = RefineOnce(placed, data, cell, sg, wavelength_A, logger, settings);
// Each arm, untwinned and twinned, is a whole refinement - the solvent built and the displacement model chosen
// against its own target - and the decision is taken on where they ended. A model only scored is scored both
// ways and decided by the same rule, so the R-free it reports is not an untwinned one by default.
if (settings.twin && result.ok && !result.cancelled)
DecideTwin(placed, data, cell, sg, wavelength_A, logger, settings, result);
// The arms run side by side, the twinned one on the next card where there is more than one, once the untwinned
// one has placed its engines and slots - and where twice what it holds of its busiest card does not fit the
// budget, after it, so that the twinned arm's engines are made as they would be alone. Its log is kept and
// follows the untwinned one's, as when they ran one after the other.
RefinementResult result;
if (!settings.twin)
result = RefineOnce(placed, data, cell, sg, wavelength_A, logger, settings);
else {
Logger twin_logger = Logger::Buffered();
CardHeld held;
std::future<size_t> held_bytes = held.bytes.get_future();
std::promise<void> untwinned_done;
const std::shared_future<void> untwinned = untwinned_done.get_future().share();
size_t budget = 0;
#ifdef JFJOCH_USE_CUDA
const int twin_card = get_gpu_count() > 1 ? (get_gpu() + 1) % get_gpu_count() : get_gpu();
if (get_gpu_count() > 0) {
const size_t usable = get_gpu_total_memory(twin_card) / 5 * 4;
budget = usable > settings.gpu_resident_bytes ? usable - settings.gpu_resident_bytes : 0;
}
#endif
std::future<void> twin_arm = std::async(std::launch::async, [&]() {
#ifdef JFJOCH_USE_CUDA
const ScopedGpuPin pin(twin_card);
#endif
if (2 * held_bytes.get() > budget)
untwinned.wait();
DecideTwin(placed, data, cell, sg, wavelength_A, twin_logger, settings, result, untwinned);
});
try {
result = RefineOnce(placed, data, cell, sg, wavelength_A, logger, settings, {}, 0, &held);
} catch (...) {
result = RefinementResult(); // not ok: the twinned arm decides nothing
held.Say(0);
untwinned_done.set_value();
twin_arm.wait();
throw;
}
held.Say(0); // where it stopped before placing its slots
untwinned_done.set_value();
twin_arm.get();
twin_logger.ReplayInto(logger);
}
result.twin_free_closed = moved;
result.map.reset();
return result;