Kronecker.periodic_discrepancy returns NaN under the default randomize='SHIFT'

[algo-hawk]

Summary

Kronecker.periodic_discrepancy returns NaN for the default randomize="SHIFT".
_square_periodic_discrepancies evaluates the kernel at x_{|a-b|}, which still
carries the random shift, rather than at the pairwise difference
x_a - x_b = (a-b)*alpha mod 1, where the shift cancels.

Reproducer

import numpy as np
from qmcpy import Kronecker
print(Kronecker(3, seed=7).periodic_discrepancy(32))  # contains NaN

Emits RuntimeWarning: invalid value encountered in sqrt from
qmcpy/discrete_distribution/kronecker.py:349.

Evidence

Negative squared discrepancies (n=32, default kernel, gamma=1):

d	randomize=False	randomize="SHIFT"
1	0 negatives	1 negative
2	0 negatives	0 negatives
3	0 negatives	25 negatives
5	0 negatives	1 negative

Two independent symptoms, consistent with one cause:

1. Shift-dependence. For a Kronecker sequence x_a - x_b = (a-b)*alpha mod 1,
   so the shift cancels and the discrepancy cannot depend on the seed. It currently
   does (seeds 1/7/42 give different values). Computing the same quantity from the
   pairwise definition gives an identical result for every seed, equal to the
   unshifted result.
2. Negative squared values. The second-Bernoulli-polynomial kernel is positive
   semi-definite, so the squared discrepancy should satisfy D^2 >= 0. Negative
   values make periodic_discrepancy take the square root of a negative number.

wssd_discrepancy consumes the same helper and is affected identically.

Impact

Silent wrong answers rather than a crash. Only masked today because
randomize=False happens to work (zero shift makes the two formulas coincide)
and the default shifted path had no unit test.

Status

Tests are in place on testing_env/brandon-unit-test-coverage, marked
@unittest.expectedFailure, in test/test_kronecker.py:

• test_squared_discrepancy_is_nonnegative_when_shifted
• test_squared_discrepancy_is_invariant_to_the_random_shift
• test_periodic_discrepancy_is_finite_when_shifted
• test_shifted_squared_discrepancy_matches_pairwise_definition

They pass as XFAIL today and will report unexpected success once the source is
corrected, at which point the markers should be removed. No source changes were made.

Note on the fix

The naive pairwise form is O(n^2). The O(n) structure can likely be kept by
evaluating the kernel at the unshifted (i*alpha) mod 1 instead of at x_i,
but that should be confirmed by a maintainer rather than assumed.

Reported by the testing team while raising unit-test coverage. Analysis was
AI-assisted; the numerical evidence above was reproduced locally.

[alegresor]

I think @fredhickernell and @AndersPride worked on this for unrandomized points

[AndersPride]

I think @fredhickernell and @AndersPride worked on this for unrandomized points

We have, this problem has been fixed in the lattice_kronecker branch