From 03169161df16ef5094a73e1436f4091c3d97a56e Mon Sep 17 00:00:00 2001 From: Benjamin Brock Date: Tue, 11 Aug 2026 23:37:23 +0000 Subject: [PATCH] Cut the peak memory of binsparse_to_ssmc_problem Reading a very large problem back needs several times the memory of writing it, because every intermediate in the conversion is an nnz-length array and the biggest matrices in the collection have nnz in the billions. For Sybrandt/AGATHA_2015 (n=183,964,077, nnz=11,588,725,964) the export runs but the read back fails with MATLAB:pmaxsize, which on Linux is an out-of-memory error in disguise. None of this changes what the converter produces; it changes how much scratch space it takes to produce it. - checked_indices validates on the stored integer type before widening to double. For an unsigned index array the finiteness and integrality tests are vacuous and the range tests collapse to min and max reductions, so the check no longer builds a second nnz-length temporary to hold fix(indices). - matrix_values keeps an ISO matrix in ISO form and returns a flag saying so. sparse() expands a scalar over the index vectors, so the single stored value never has to be repmat'ed to nnz length. expand_iso does the expansion in the three dense branches, which genuinely need an array. - expand_structure keeps the ISO form across mirroring whenever mirroring leaves the stored value alone, which covers every symmetric matrix along with the real hermitian and all-zero skew-symmetric ones. - convert_matrix takes a separate branch for ISO matrices, whose entries are either all zero or all nonzero, so it needs no mask and no masked copies of the index vectors. For the general case it tests values == 0 directly, as == already compares both parts of a complex value; the old form called imag() on real values and got back an nnz-length array of zeros. When nothing is an explicit zero it hands the index vectors to sparse() untouched rather than copying them through a mask. - require_ordered_pairs and require_strictly_increasing scan in blocks with a one-entry overlap. require_ordered_pairs also evaluated diff(first) twice. Tests cover the new ISO paths: a nonzero ISO COOR matrix, ISO with symmetric_lower staying ISO, ISO with skew_symmetric_lower and with hermitian_lower dropping it, and a dense ISO DVEC. The 14-matrix BSP test collection still round-trips byte-exact through ssread; five of those matrices, Pajek/IMDB among them, store ISO values and so exercise the new path on real data. Co-Authored-By: Claude Opus 5 (1M context) --- bindings/matlab/binsparse_to_ssmc_problem.m | 172 ++++++++++++++---- .../matlab/test_binsparse_to_ssmc_problem.m | 32 ++++ 2 files changed, 168 insertions(+), 36 deletions(-) diff --git a/bindings/matlab/binsparse_to_ssmc_problem.m b/bindings/matlab/binsparse_to_ssmc_problem.m index 8c06695..0555ee9 100644 --- a/bindings/matlab/binsparse_to_ssmc_problem.m +++ b/bindings/matlab/binsparse_to_ssmc_problem.m @@ -143,12 +143,12 @@ function validate_problem(bsp_problem) switch format case 'DVEC' - values = matrix_values(bsp, count); + [values, iso] = matrix_values(bsp, count); if n ~= 1 || count ~= m error('BinSparse:InvalidMatrix', ... 'DVEC dimensions do not match the stored value count'); end - matrix = reshape(values, [m, 1]); + matrix = reshape(expand_iso(values, iso, count), [m, 1]); Zeros = sparse(m, n); return; @@ -158,8 +158,8 @@ function validate_problem(bsp_problem) error('BinSparse:InvalidMatrix', ... 'DMATR dimensions do not match the stored value count'); end - values = matrix_values(bsp, count); - matrix = reshape(values, [n, m]).'; + [values, iso] = matrix_values(bsp, count); + matrix = reshape(expand_iso(values, iso, count), [n, m]).'; Zeros = sparse(m, n); return; @@ -169,24 +169,50 @@ function validate_problem(bsp_problem) error('BinSparse:InvalidMatrix', ... 'DMATC dimensions do not match the stored value count'); end - values = matrix_values(bsp, count); - matrix = reshape(values, [m, n]); + [values, iso] = matrix_values(bsp, count); + matrix = reshape(expand_iso(values, iso, count), [m, n]); Zeros = sparse(m, n); return; end [rows, cols] = sparse_indices(bsp, format, m, n, count); -values = matrix_values(bsp, count); -[rows, cols, values] = expand_structure(bsp, rows, cols, values, m, n); +[values, iso] = matrix_values(bsp, count); +[rows, cols, values, iso] = ... + expand_structure(bsp, rows, cols, values, iso, m, n); + +% An ISO matrix carries a single stored value, and sparse() expands a scalar +% over the index vectors, so its values are never materialized at nnz length. +% The entries are then either all zero or all nonzero, which also removes the +% need for a mask and for masked copies of the index vectors. +if iso + if split_zeros && values == 0 + matrix = sparse(m, n); + Zeros = sparse(rows, cols, 1, m, n); + else + matrix = sparse(rows, cols, values, m, n); + Zeros = sparse(m, n); + end + return; +end -is_zero = (real(values) == 0) & (imag(values) == 0); -if split_zeros - matrix = sparse(rows(~is_zero), cols(~is_zero), ... - values(~is_zero), m, n); - Zeros = sparse(rows(is_zero), cols(is_zero), 1, m, n); -else +if ~split_zeros matrix = sparse(rows, cols, values, m, n); Zeros = sparse(m, n); + return; +end + +% MATLAB's == compares both parts of a complex value, so this needs no +% separate imag() test; taking one would allocate an nnz-length array of +% zeros for the common case of real values. +is_zero = (values == 0); +if ~any(is_zero) + matrix = sparse(rows, cols, values, m, n); + Zeros = sparse(m, n); +else + keep = ~is_zero; + matrix = sparse(rows(keep), cols(keep), values(keep), m, n); + clear keep + Zeros = sparse(rows(is_zero), cols(is_zero), 1, m, n); end end @@ -247,8 +273,13 @@ function validate_problem(bsp_problem) end end -function values = matrix_values(bsp, count) +function [values, iso] = matrix_values(bsp, count) +% Returns the stored values and whether they are still in ISO form, that is a +% single value standing for all count entries. The caller expands it only +% where an nnz-length array is genuinely required; the sparse constructors +% take the scalar directly. values = bsp.values(:); +iso = false; if logical(bsp.is_iso) if count == 0 if numel(values) > 1 @@ -260,7 +291,7 @@ function validate_problem(bsp_problem) error('BinSparse:InvalidMatrix', ... 'An ISO matrix must contain exactly one value'); else - values = repmat(values, count, 1); + iso = true; end elseif numel(values) ~= count error('BinSparse:InvalidMatrix', ... @@ -269,7 +300,14 @@ function validate_problem(bsp_problem) values = exact_double(values); end -function [rows, cols, values] = expand_structure(bsp, rows, cols, values, m, n) +function values = expand_iso(values, iso, count) +if iso + values = repmat(values, count, 1); +end +end + +function [rows, cols, values, iso] = ... + expand_structure(bsp, rows, cols, values, iso, m, n) structure = 'general'; if isfield(bsp, 'structure') && ~isempty(bsp.structure) structure = lower(char(bsp.structure)); @@ -293,19 +331,33 @@ function validate_problem(bsp_problem) 'Stored entries do not match the declared matrix triangle'); end -off_diagonal = (rows ~= cols); -mirror_rows = cols(off_diagonal); -mirror_cols = rows(off_diagonal); -mirror_values = values(off_diagonal); if starts_with(structure, 'hermitian_') - mirror_values = conj(mirror_values); + mirror = @conj; elseif starts_with(structure, 'skew_symmetric_') - mirror_values = -mirror_values; -elseif ~starts_with(structure, 'symmetric_') + mirror = @uminus; +elseif starts_with(structure, 'symmetric_') + mirror = @(v) v; +else error('BinSparse:UnsupportedStructure', ... 'Unsupported Binsparse matrix structure: %s', structure); end +off_diagonal = (rows ~= cols); +mirror_rows = cols(off_diagonal); +mirror_cols = rows(off_diagonal); + +% Mirroring an ISO matrix leaves it ISO whenever it leaves the single stored +% value alone, which covers every symmetric matrix and the real hermitian and +% all-zero skew-symmetric ones. Only the rest have to be expanded here. +if iso && isequaln(mirror(values), values) + rows = [rows; mirror_rows]; + cols = [cols; mirror_cols]; + return; +end + +values = expand_iso(values, iso, numel(rows)); +iso = false; +mirror_values = mirror(values(off_diagonal)); rows = [rows; mirror_rows]; cols = [cols; mirror_cols]; values = [values; mirror_values]; @@ -359,14 +411,32 @@ function validate_matrix_header(bsp) if ~isnumeric(indices) || ~isreal(indices) error('BinSparse:InvalidMatrix', '%s must be real numeric data', name); end -indices = double(indices(:)); if ~isempty(expected) && numel(indices) ~= expected error('BinSparse:InvalidMatrix', ... '%s length does not match nnz', name); end -if any(~isfinite(indices)) || any(indices < 0) || ... - any(fix(indices) ~= indices) || any(indices >= limit) - error('BinSparse:InvalidMatrix', '%s contains an invalid index', name); + +if isinteger(indices) + % A Binsparse index array is normally stored as an integer type, where + % the finiteness and integrality tests are vacuous and the range tests + % reduce to min and max. Running them here, before the widening to + % double, keeps the check from building a second nnz-length temporary + % just to hold fix(indices). + if ~isempty(indices) + if double(max(indices(:))) >= limit || ... + (intmin(class(indices)) < 0 && double(min(indices(:))) < 0) + error('BinSparse:InvalidMatrix', ... + '%s contains an invalid index', name); + end + end + indices = double(indices(:)); +else + indices = double(indices(:)); + if any(~isfinite(indices)) || any(indices < 0) || ... + any(fix(indices) ~= indices) || any(indices >= limit) + error('BinSparse:InvalidMatrix', ... + '%s contains an invalid index', name); + end end end @@ -393,17 +463,35 @@ function validate_matrix_header(bsp) end function require_ordered_pairs(first, second, format) -if numel(first) < 2 +% Scanned in blocks, so the check costs a fixed few megabytes rather than +% several nnz-length temporaries. Consecutive blocks overlap by one entry so +% that no adjacent pair straddles a block boundary unchecked. +count = numel(first); +if count < 2 return; end -bad = diff(first) < 0 | ... - (diff(first) == 0 & diff(second) <= 0); -if any(bad) - error('BinSparse:InvalidMatrix', ... - '%s indices are not sorted and unique', format); +lo = 1; +while lo < count + hi = min(count, lo + scan_block_size()); + step = diff(first(lo:hi)); + bad = step < 0; + tied = (step == 0); + clear step + if any(tied) + bad = bad | (tied & diff(second(lo:hi)) <= 0); + end + if any(bad) + error('BinSparse:InvalidMatrix', ... + '%s indices are not sorted and unique', format); + end + lo = hi; end end +function n = scan_block_size() +n = 4194304; +end + function require_segment_order(indices, pointers, format) for k = 1:numel(pointers)-1 first = pointers(k) + 1; @@ -423,8 +511,20 @@ function require_nonempty_segments(pointers, format) end function require_strictly_increasing(values, label) -if numel(values) > 1 && any(diff(values) <= 0) - error('BinSparse:InvalidMatrix', '%s are not sorted and unique', label); +% Blocked for the same reason as require_ordered_pairs: a CVEC index array is +% nnz long, so diff() over the whole of it is not affordable at scale. +count = numel(values); +if count < 2 + return; +end +lo = 1; +while lo < count + hi = min(count, lo + scan_block_size()); + if any(diff(values(lo:hi)) <= 0) + error('BinSparse:InvalidMatrix', ... + '%s are not sorted and unique', label); + end + lo = hi; end end diff --git a/bindings/matlab/test_binsparse_to_ssmc_problem.m b/bindings/matlab/test_binsparse_to_ssmc_problem.m index 3f998bf..0cbcddd 100644 --- a/bindings/matlab/test_binsparse_to_ssmc_problem.m +++ b/bindings/matlab/test_binsparse_to_ssmc_problem.m @@ -33,6 +33,38 @@ assert(nnz(Problem.A) == 0); assert(isequal(Problem.Zeros, sparse(rows + 1, cols + 1, 1, 3, 3))); +% An ISO matrix keeps its single stored value all the way to sparse(), which +% expands it, so these check that the expansion still lands on every entry and +% that mirroring drops the ISO form exactly when it changes the value. +raw = struct('metadata', metadata, ... + 'A', make_matrix(4, rows, cols, [], 3, 3, 'COOR', true, 'general')); +Problem = binsparse_to_ssmc_problem(raw); +assert(isequal(Problem.A, sparse(rows + 1, cols + 1, 4, 3, 3))); +assert(~isfield(Problem, 'Zeros'), 'a nonzero ISO matrix has no Zeros'); + +iso_symmetric = make_matrix(7, [0; 1; 1], [0; 0; 1], ... + [], 2, 2, 'COO', true, 'symmetric_lower'); +raw = struct('metadata', metadata, 'A', iso_symmetric); +Problem = binsparse_to_ssmc_problem(raw); +assert(isequal(Problem.A, sparse(7 * ones(2, 2)))); + +iso_skew = make_matrix(5, 1, 0, [], 2, 2, 'COO', true, ... + 'skew_symmetric_lower'); +raw = struct('metadata', metadata, 'A', iso_skew); +Problem = binsparse_to_ssmc_problem(raw); +assert(isequal(Problem.A, sparse([0 -5; 5 0]))); + +iso_hermitian = make_matrix(2 + 3i, 1, 0, [], 2, 2, 'COO', true, ... + 'hermitian_lower'); +raw = struct('metadata', metadata, 'A', iso_hermitian); +Problem = binsparse_to_ssmc_problem(raw); +assert(isequal(Problem.A, sparse([0 2-3i; 2+3i 0]))); + +raw = struct('metadata', metadata, 'A', formats{1}, ... + 'b', make_matrix(9, [], [], [], 3, 1, 'DVEC', true, 'general')); +Problem = binsparse_to_ssmc_problem(raw); +assert(isequal(Problem.b, [9; 9; 9]), 'a dense ISO vector must expand'); + lower = make_matrix([1; 2; 0; 3], [0; 1; 1; 2], [0; 0; 1; 2], ... [], 3, 3, 'COO', false, 'symmetric_lower'); raw = struct('metadata', metadata, 'A', lower);