Skip to content

Diagnose invalid angular grid sizes instead of returning zeros - #101

Open
susilehtola wants to merge 1 commit into
wavefunction91:masterfrom
susilehtola:fix/angular-invalid-npts
Open

Diagnose invalid angular grid sizes instead of returning zeros#101
susilehtola wants to merge 1 commit into
wavefunction91:masterfrom
susilehtola:fix/angular-invalid-npts

Conversation

@susilehtola

Copy link
Copy Markdown
Collaborator

The generate() dispatch for all four S2 quadrature families was a chain
of if/else-if with no terminal else. Requesting a point count that is not
tabulated left the value-initialised containers untouched, so the caller
received a grid whose weights were all zero and which therefore
integrates every function to zero. No exception, assertion or warning was
raised, and the README claims the opposite ("will fail if the grid order
is incompatible").

Add a terminal else to each family that reports the requested size and
points at npts_by_algebraic_order().

Two data problems surfaced once the check was in place:

  • lebedev_laikov_5294.hpp is shipped and included by the grids header but
    had no dispatch branch and no entry in the algebraic order tables, so
    the 5294-point Lebedev-Laikov grid was unreachable. Wire it up (order
    125, between 4802/119 and 5810/131).

  • The tabulated 552-point Ahrens-Beylkin grid is not a valid order-39
    rule. Its weights sum to 0.9632 * 4pi, and after renormalisation it
    still integrates degree-8 monomials to only 5.7e-02 relative accuracy,
    while every neighbouring size is exact to machine precision. Withdraw
    it from the dispatch and the order tables until the table can be
    regenerated from the original Ahrens-Beylkin data; next_algebraic_order
    now rounds order 39 up to 41.

The new test checks both directions: unsupported sizes must throw, and
every size advertised by the order tables must be dispatched and must
integrate even monomials exactly up to its claimed algebraic order. That
second check is what found the Ahrens-Beylkin defect; a sum-of-weights
check alone would not have, since the grid is exact at low order once
renormalised.

Co-Authored-By: Claude Opus 5 (1M context) noreply@anthropic.com
Claude-Session: https://claude.ai/code/session_01FDTFYJMQ76iujDFNHzZyXF

The generate() dispatch for all four S2 quadrature families was a chain
of if/else-if with no terminal else. Requesting a point count that is not
tabulated left the value-initialised containers untouched, so the caller
received a grid whose weights were all zero and which therefore
integrates every function to zero. No exception, assertion or warning was
raised, and the README claims the opposite ("will fail if the grid order
is incompatible").

Add a terminal else to each family that reports the requested size and
points at npts_by_algebraic_order().

Two data problems surfaced once the check was in place:

* lebedev_laikov_5294.hpp is shipped and included by the grids header but
  had no dispatch branch and no entry in the algebraic order tables, so
  the 5294-point Lebedev-Laikov grid was unreachable. Wire it up (order
  125, between 4802/119 and 5810/131).

* The tabulated 552-point Ahrens-Beylkin grid is not a valid order-39
  rule. Its weights sum to 0.9632 * 4pi, and after renormalisation it
  still integrates degree-8 monomials to only 5.7e-02 relative accuracy,
  while every neighbouring size is exact to machine precision. Withdraw
  it from the dispatch and the order tables until the table can be
  regenerated from the original Ahrens-Beylkin data; next_algebraic_order
  now rounds order 39 up to 41.

The new test checks both directions: unsupported sizes must throw, and
every size advertised by the order tables must be dispatched and must
integrate even monomials exactly up to its claimed algebraic order. That
second check is what found the Ahrens-Beylkin defect; a sum-of-weights
check alone would not have, since the grid is exact at low order once
renormalised.

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01FDTFYJMQ76iujDFNHzZyXF
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