Skip to content

Keep gramcover accurate when squares leave range - #87

Merged
timholy merged 1 commit into
mainfrom
teh/gramcover-range
Oct 13, 2026
Merged

timholy merged 1 commit into
mainfrom
teh/gramcover-range

Conversation

@timholy

@timholy timholy commented Oct 13, 2026

Copy link
Copy Markdown
Member

gramcover(A) and gramcover(A, w) sum abs(w[i]) * A[i,j]^2 per column and take one square root. An entry below about 1e-162 squares to zero and one above about 1e154 squares to Inf, so a column whose 2-norm is an ordinary number could get a scale of 0 or Inf: gramcover([1e-170 1.0]) returned [0.0, 1.0], which does not cover A'*A.

The plain sum of squares is kept as the first pass. A column whose sum is below floatmin, infinite, or NaN is accumulated again as a running hypot over its terms sqrt(abs(w[i])) * abs(A[i,j]). The scale is then accurate to the documented O(n_j * eps) whenever the largest term in the column is at least floatmin, and is Inf only when the norm exceeds floatmax. The existing inflation factor bounds the rounding of both forms.

The weight now enters each term as sqrt(abs(w[i])) before squaring, so weighted results change in the last bits. Entries are converted to floating point before squaring, which keeps the squares of large integer entries from wrapping.

When no column needs the second pass the time is unchanged: 1.43 ms for a dense 2000x500 matrix and 0.14 ms for a sparse matrix with 100000 nonzeros. When every column needs it, the sparse case takes 1.3 ms. Using hypot for every entry instead takes 13.6 ms and 1.1 ms.

Assisted-by: Claude Opus 5.5 noreply@anthropic.com

`gramcover(A)` and `gramcover(A, w)` sum `abs(w[i]) * A[i,j]^2` per
column and take one square root. An entry below about 1e-162 squares to
zero and one above about 1e154 squares to `Inf`, so a column whose
2-norm is an ordinary number could get a scale of 0 or `Inf`:
`gramcover([1e-170 1.0])` returned `[0.0, 1.0]`, which does not cover
`A'*A`.

The plain sum of squares is kept as the first pass. A column whose sum
is below `floatmin`, infinite, or `NaN` is accumulated again as a
running `hypot` over its terms `sqrt(abs(w[i])) * abs(A[i,j])`. The
scale is then accurate to the documented `O(n_j * eps)` whenever the
largest term in the column is at least `floatmin`, and is `Inf` only
when the norm exceeds `floatmax`. The existing inflation factor bounds
the rounding of both forms.

The weight now enters each term as `sqrt(abs(w[i]))` before squaring,
so weighted results change in the last bits. Entries are converted to
floating point before squaring, which keeps the squares of large
integer entries from wrapping.

When no column needs the second pass the time is unchanged: 1.43 ms for
a dense 2000x500 matrix and 0.14 ms for a sparse matrix with 100000
nonzeros. When every column needs it, the sparse case takes 1.3 ms.
Using `hypot` for every entry instead takes 13.6 ms and 1.1 ms.

Assisted-by: Claude Opus 5.5 <noreply@anthropic.com>
@codecov

codecov Bot commented Oct 13, 2026

Copy link
Copy Markdown

Codecov Report

✅ All modified and coverable lines are covered by tests.
✅ Project coverage is 98.79%. Comparing base (b80567c) to head (7d716bf).

Additional details and impacted files
@@           Coverage Diff           @@
##             main      #87   +/-   ##
=======================================
  Coverage   98.79%   98.79%           
=======================================
  Files          23       23           
  Lines        6286     6302   +16     
=======================================
+ Hits         6210     6226   +16     
  Misses         76       76           

☔ View full report in Codecov by Harness.
📢 Have feedback on the report? Share it here.

🚀 New features to boost your workflow:
  • ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.

@timholy
timholy merged commit cb4bc5e into main Oct 13, 2026
5 checks passed
@timholy
timholy deleted the teh/gramcover-range branch October 13, 2026 12:38
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

1 participant