Build Lab: Hirschberg's O(NM)-Time O(min(N,M))-Space Sequence Alignment & Diff Reconstructor

1. The Memory Bottleneck of Full Alignment Tables

When aligning two genomic sequences or source files of length $N = M = 100{,}000$, storing a 2D table of 32-bit parent pointers requires $10^{10} \times 4 \text{ bytes} = 40 \text{ GB}$ of RAM—causing an immediate Out-Of-Memory (OOM) crash. Yet a 1D rolling array of length $M + 1$ takes only $400 \text{ KB}$, at the cost of forgetting the path. How can we recover the full path using only $O(M)$ memory?

2. The Geometric Series Proof: Why Total Time Only Doubles ($2NM$)

At depth $0$ of the recursion, the forward and backward 1D passes do proportional to $\frac{n}{2} m + \frac{n}{2} m = n m$ work. At depth $1$, the two recursive calls process areas $\frac{n}{2} j^*$ and $\frac{n}{2} (m - j^*)$, whose sum is:

At depth $d$, the $2^d$ subproblems have total grid area $\frac{1}{2^d} n m$! Summing over all recursion levels gives the geometric series:

Meanwhile, by reusing preallocated scratch vectors (or freeing row buffers before recursing), peak memory is strictly $\Theta(\min(n, m) + \text{output size})$!

3. Complexity Matrix

4. Complete Production C++20 Build Lab Implementation

#include <algorithm>
#include <cstdint>
#include <string>
#include <string_view>
#include <vector>

// Hirschberg's Algorithm: Reconstruct exact Longest Common Subsequence in O(N*M) time and O(M) space
class HirschbergLCS {
    // 1D Rolling Forward LCS row for A[0..n-1] vs B[0..m-1]
    static std::vector<int> lcs_forward_row(std::string_view a, std::string_view b) {
        const int m = static_cast<int>(b.size());
        std::vector<int> dp(m + 1, 0);
        for (char ca : a) {
            int prev_diag = 0;
            for (int j = 1; j <= m; ++j) {
                int saved_up = dp[j];
                if (ca == b[j - 1]) dp[j] = prev_diag + 1;
                else dp[j] = std::max(dp[j], dp[j - 1]);
                prev_diag = saved_up;
            }
        }
        return dp;
    }

    // 1D Rolling Backward LCS row for reversed suffixes of A and B
    static std::vector<int> lcs_backward_row(std::string_view a, std::string_view b) {
        const int n = static_cast<int>(a.size());
        const int m = static_cast<int>(b.size());
        std::vector<int> dp(m + 1, 0);
        for (int i = n - 1; i >= 0; --i) {
            int prev_diag = 0;
            for (int j = m - 1; j >= 0; --j) {
                int saved_down = dp[j];
                if (a[i] == b[j]) dp[j] = prev_diag + 1;
                else dp[j] = std::max(dp[j], dp[j + 1]);
                prev_diag = saved_down;
            }
        }
        return dp;
    }

    static void solve_rec(std::string_view a, std::string_view b, std::string& out) {
        const int n = static_cast<int>(a.size());
        const int m = static_cast<int>(b.size());
        if (n == 0 || m == 0) return;

        // Base Case: Single character in a takes O(m) scan
        if (n == 1) {
            if (b.find(a[0]) != std::string_view::npos) {
                out.push_back(a[0]);
            }
            return;
        }

        const int mid = n / 2;
        int best_j = 0;
        {
            // Scope block ensures L and R are deallocated BEFORE recursing!
            std::vector<int> L = lcs_forward_row(a.substr(0, mid), b);
            std::vector<int> R = lcs_backward_row(a.substr(mid), b);
            int max_len = -1;
            for (int j = 0; j <= m; ++j) {
                if (L[j] + R[j] > max_len) {
                    max_len = L[j] + R[j];
                    best_j = j;
                }
            }
        }

        // Conquer Top-Left [0..mid) x [0..best_j) and Bottom-Right [mid..n) x [best_j..m)
        solve_rec(a.substr(0, mid), b.substr(0, best_j), out);
        solve_rec(a.substr(mid), b.substr(best_j), out);
    }

public:
    static std::string reconstruct_lcs(std::string_view a, std::string_view b) {
        std::string out;
        out.reserve(std::min(a.size(), b.size()));
        solve_rec(a, b, out);
        return out;
    }
};