From 3aa89ad067cb0f81b571f6d9b0d39eea38d408b9 Mon Sep 17 00:00:00 2001 From: Filip Leonarski Date: Thu, 15 Oct 2026 21:42:16 +0200 Subject: [PATCH 1/3] Rotation scaling: merge working arrays handed back before the space-group search, the all-observation arm's engine emptied once it has merged, output unchanged The merge's per-group and per-full working arrays (accumulators, R_meas sums and staging, ISa scatter, error-model pools, the per-full fields) are engine members kept between merges to save faulting them in again. Each is resized and refilled before it is read, so handing them back leaves the next merge exactly where an engine's first merge starts. - ReleaseMergeScratch(): the main engine's, before the space-group search. Its next merge is in the group the search picks, with far fewer groups than P1. - ReleaseAllButPerFrameScale(): the CloneIngest engine of the all-observation arm merges once; right after its Run() everything on its host but the per-frame state WriteBackPerFrameScale reads goes (fulls, staging, raw-hkl runs, ingest copies, merge scratch). Measured (heap profile + RSS sampled at 4 Hz, GPU path): 4nwv peak RSS 19.10 -> 18.17 GiB; the space-group search phase 19.1 -> 16.9 GiB, so the peak is now the two concurrent P1 merges. Outputs byte-identical on myob_x10sa, cytc_x10sa, 4nwv (timing benchmark, MR) and 8qaw, 4nwv (open arm): p.hkl, MTZ column data, p_maps.mtz, maps, placed model, every report field. Co-Authored-By: Claude Opus 5.5 Claude-Session: https://claude.ai/code/session_01PXc1zD1yxHqaR1DGdCTHL4 --- .../scale_merge/RotationScaleMerge.h | 46 +++++++++++++++++++ rugnux/RugnuxScaleMerge.cpp | 13 +++++- 2 files changed, 57 insertions(+), 2 deletions(-) diff --git a/image_analysis/scale_merge/RotationScaleMerge.h b/image_analysis/scale_merge/RotationScaleMerge.h index 514378909..477d9bcf4 100644 --- a/image_analysis/scale_merge/RotationScaleMerge.h +++ b/image_analysis/scale_merge/RotationScaleMerge.h @@ -161,6 +161,52 @@ public: // beside others with it off, at the point where it would have run in series. void WriteBackPerFrameScale() { FinalizePerFrameScale(last_cc, last_cc_n, last_frame_scaled); } + // Hand back the merge's working arrays - MergeAndStats' accumulators, R_meas sums, ISa scatter, + // error-model pools and per-full fields. Each is resized and refilled before it is read, so the + // next merge starts from them as an engine's first merge does, and gives the same answer. On a + // large cell in P1 they are gigabytes, kept between merges only to save faulting them in again; + // this is for a gap between merges where something else sits on the run's memory peak. + void ReleaseMergeScratch() { + merge_fields = {}; + merge_acc = {}; + merge_accum = {}; + merge_merged_I = {}; + merge_rejected = {}; + merge_half_set = {}; + merge_obs_key = {}; + merge_wobs = {}; + merge_rmeas = {}; + merge_rmeas_staging = {}; + asymptote_scatter = {}; + em_samples = {}; + em_cut_pool = {}; + em_fit_pool = {}; + em_chi2 = {}; + } + + // For an engine that will not merge again: everything on the host but the per-frame state + // WriteBackPerFrameScale reads (partials_out, g_partial, frame_flux, mos_smooth, last_cc, + // last_cc_n, last_frame_scaled). Its card memory goes with the engine. + void ReleaseAllButPerFrameScale() { + ReleaseMergeScratch(); + partials = {}; + fulls = {}; + fulls_staging = {}; + fulls_group_start = {}; + fulls_group_perm = {}; + perm = {}; + rawrun_start = {}; + rawrun_count = {}; + rawrun_h = {}; + rawrun_k = {}; + rawrun_l = {}; + rawrun_d = {}; + rawrun_group = {}; + finite_ok = {}; + corr_ingested = {}; + obs_position.reset(); + } + // A second engine on this one's ingest, on `device` and logging to `log`, for a merge to run beside // this engine's: the host half of the ingest is copied, and the copy's IngestDevice(this) is still // to be called. Call after this engine's IngestDevice(). No merge depends on the merges run diff --git a/rugnux/RugnuxScaleMerge.cpp b/rugnux/RugnuxScaleMerge.cpp index f072be960..1814d489b 100644 --- a/rugnux/RugnuxScaleMerge.cpp +++ b/rugnux/RugnuxScaleMerge.cpp @@ -1183,8 +1183,13 @@ bool Rugnux::ScaleMergeAndSymmetry(PipelineLocals &p) { all_arm_ahead = std::async(std::launch::async, [&rsm_all, &all_arm_log] { pin_gpu(rsm_all->GetGpuDevice()); // the card and its NUMA node, as the guard's merge has them TimingMark mark(all_arm_log, "merge P1, all observations"); - return rsm_all->Run(/*for_search=*/true, /*full_stats=*/true, - /*measure_cc_before_corrections=*/false); + auto r = rsm_all->Run(/*for_search=*/true, /*full_stats=*/true, + /*measure_cc_before_corrections=*/false); + // This engine merges once: what is left for it to do is the write-back of its + // per-frame scales, so its observations and its merge's working arrays go now, + // rather than through the search the others are about to run. + rsm_all->ReleaseAllButPerFrameScale(); + return r; }); } } @@ -1593,6 +1598,10 @@ bool Rugnux::ScaleMergeAndSymmetry(PipelineLocals &p) { postrefine_gathered_future.wait(); // the search re-indexes the reflections the gather reads TimingMark search_mark(logger, "space-group search"); if (search_space_group) { + // The engine's next merge is in the group the search picks; until then its last merge's + // working arrays only sit beside the searches' own (gigabytes on a large cell in P1). + if (rsm) + rsm->ReleaseMergeScratch(); SearchSpaceGroupOptions sg_opts; sg_opts.nthreads = static_cast(std::max(1, config_.nthreads)); sg_opts.merge_friedel = experiment_.GetScalingSettings().GetMergeFriedel(); From 1138f3a0bb13926d88e7f5594dbe9d3d7eb12ee7 Mon Sep 17 00:00:00 2001 From: Filip Leonarski Date: Thu, 15 Oct 2026 21:42:16 +0200 Subject: [PATCH 2/3] Post-refinement: the gathered partials handed back once the fit's observations are selected, output unchanged The joint geometry fit reads at most MAX_OBS (20000) spot positions, chosen by I/sigma, through pointers into the gathered partial array (28 B a partial). That array and a full-length pointer array (8 B a partial) stayed alive until the fit returned. Now the ranking is built straight from the gathered partials (no pointer array), the kept observations are copied out, and the gathered array is freed before the fit. Same comparisons in the same order, so the same selection. 8qaw (122.7 M partials, GPU path): peak RSS 23.33 -> 21.92 GiB with the previous commit; the post-refinement, which runs beside the P1 merge, holds 3.3 GB less for most of its fit. Outputs byte-identical on myob_x10sa, cytc_x10sa, 4nwv, 8qaw. [PostRefine] tests pass. Co-Authored-By: Claude Opus 5.5 Claude-Session: https://claude.ai/code/session_01PXc1zD1yxHqaR1DGdCTHL4 --- image_analysis/geom_refinement/PostRefine.cpp | 38 ++++++++++++------- 1 file changed, 24 insertions(+), 14 deletions(-) diff --git a/image_analysis/geom_refinement/PostRefine.cpp b/image_analysis/geom_refinement/PostRefine.cpp index 43a0921fb..6c2bce937 100644 --- a/image_analysis/geom_refinement/PostRefine.cpp +++ b/image_analysis/geom_refinement/PostRefine.cpp @@ -473,10 +473,9 @@ PostRefineResult PostRefineRotationGeometry(PostRefineObservations observations, // The observed spot positions, one entry per partial that produced a centroid. Count first, // then fill, exactly as the partial gather above does and for the same reason: this walks the - // same tens of millions of partials, and a pointer vector grown by push_back copies itself - // every time it doubles. new[] rather than a sized vector so the array is not zeroed on one - // thread before the parallel fill overwrites it. The fill lands in the order the serial loop - // produced, so the selection below sees the same sequence it always did. + // same tens of millions of partials, and a vector grown by push_back copies itself every + // time it doubles. The fill lands in the order the serial loop produced, so the selection + // below sees the same sequence it always did. const int n_obs_chunks = static_cast(std::clamp(nthreads, 1, std::max(1, n_pts))); const size_t obs_chunk = (n_pts + n_obs_chunks - 1) / n_obs_chunks; @@ -491,28 +490,39 @@ PostRefineResult PostRefineRotationGeometry(PostRefineObservations observations, }); for (int c = 0; c < n_obs_chunks; ++c) obs_offset[c + 1] += obs_offset[c]; size_t n_obs = obs_offset[n_obs_chunks]; - std::unique_ptr obs(new const Partial *[n_obs]); - ParallelChunks(static_cast(n_pts), nthreads, [&](int lo, int hi) { - size_t at = obs_offset[static_cast(lo) / obs_chunk]; - for (int i = lo; i < hi; ++i) if (keep_obs(i)) obs[at++] = &pts[i]; - }); + // The fit reads at most MAX_OBS of them, the strongest by I/sigma. Those are copied out, and + // the gathered partials - gigabytes on a large cell, and read by nothing below - go back now + // rather than when the fit returns. constexpr size_t MAX_OBS = 20000; + std::vector kept; if (n_obs > MAX_OBS) { // By I/sigma, with the ratio worked out once per observation beside it rather than at - // every comparison through the pointer: the same comparisons, so the same selection in - // the same order. + // every comparison: the same comparisons, so the same selection in the same order. struct Ranked { double i_over_sigma; const Partial *p; }; std::vector ranked(n_obs); - ParallelChunks(static_cast(n_obs), nthreads, [&](int lo, int hi) { + ParallelChunks(static_cast(n_pts), nthreads, [&](int lo, int hi) { + size_t at = obs_offset[static_cast(lo) / obs_chunk]; for (int i = lo; i < hi; ++i) - ranked[i] = {obs[i]->I / std::max(1e-9, static_cast(obs[i]->sigma)), obs[i]}; + if (keep_obs(i)) + ranked[at++] = {pts[i].I / std::max(1e-9, static_cast(pts[i].sigma)), &pts[i]}; }); std::nth_element(ranked.begin(), ranked.begin() + MAX_OBS, ranked.end(), [](const Ranked &a, const Ranked &b) { return a.i_over_sigma > b.i_over_sigma; }); + kept.resize(MAX_OBS); for (size_t i = 0; i < MAX_OBS; ++i) - obs[i] = ranked[i].p; + kept[i] = *ranked[i].p; n_obs = MAX_OBS; + } else { + kept.resize(n_obs); + ParallelChunks(static_cast(n_pts), nthreads, [&](int lo, int hi) { + size_t at = obs_offset[static_cast(lo) / obs_chunk]; + for (int i = lo; i < hi; ++i) if (keep_obs(i)) kept[at++] = pts[i]; + }); } + observations.pts.reset(); + std::unique_ptr obs(new const Partial *[n_obs]); + for (size_t i = 0; i < n_obs; ++i) + obs[i] = &kept[i]; result.obs_used = static_cast(n_obs); // ===== The joint fit ===== From 18835f93f802088e962e3e21cc0801becbcbbc27 Mon Sep 17 00:00:00 2001 From: Filip Leonarski Date: Thu, 15 Oct 2026 22:53:48 +0200 Subject: [PATCH 3/3] Rotation scaling: merge working arrays also handed back after the serial all-observation merge and before the analyses, output unchanged The two other points between merges where the engine's MergeAndStats arrays (~330 B per ASU group and ~55 B per full; ~3 GiB per P1 merge on a 400 A cell) sat beside other work: - the all-observation arm merged in series (the CPU path, or a card too small for the arm beside), whose arrays were held through that arm's SearchSpaceGroup; - after the answer's merge, through the analyses (a twin-law or ledger merge refills them). -N 10 --no-gpu (the macOS path's thread count): myob_x10sa 7.73 -> 7.28, cytc_x10sa 7.96 -> 7.56, 4nwv (timing benchmark) 20.34 -> 20.07 GiB peak RSS. Outputs byte-identical (p.hkl, MTZ column data, maps, placed model, report fields) at -N 10 on the CPU path and on the GPU path for myob_x10sa, cytc_x10sa, 4nwv, 8qaw. Co-Authored-By: Claude Opus 5.5 Claude-Session: https://claude.ai/code/session_01PXc1zD1yxHqaR1DGdCTHL4 --- rugnux/RugnuxScaleMerge.cpp | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/rugnux/RugnuxScaleMerge.cpp b/rugnux/RugnuxScaleMerge.cpp index 1814d489b..658fadb27 100644 --- a/rugnux/RugnuxScaleMerge.cpp +++ b/rugnux/RugnuxScaleMerge.cpp @@ -1699,6 +1699,7 @@ bool Rugnux::ScaleMergeAndSymmetry(PipelineLocals &p) { rsm->SetSearchMinZeta(0.0); sm_all = scale_and_merge("P1, all observations", true); rsm->SetSearchMinZeta(zeta); + rsm->ReleaseMergeScratch(); // as before the first arm's search sg_opts.merge_isa = result.error_model_isa; // this arm's own error model TimingMark all_arm_mark(logger, "SearchSpaceGroup (all-observation arm)"); alt = SearchSpaceGroup(sm_all.merged, sg_opts); @@ -3309,6 +3310,9 @@ bool Rugnux::ScaleMergeAndSymmetry(PipelineLocals &p) { } adoption_mark.Stop(); TimingMark analyses_mark(logger, "analyses"); + // The answer's merge is made and in `sm`; a merge below (twin law, ledger) refills the arrays. + if (rsm) + rsm->ReleaseMergeScratch(); // --finalist-ledger (see where it is read, below) is folded from a merge of its own: the same // observations in P1. It needs every node of the lattice integrated, which a fixed centred group