pythran compromises scipy precision on armhf

[drew-parsons]

Issue #2002 reports that pythran breaks scipy on 32-bit architectures. Patch scipy/scipy#16646 in scipy addresses that problem (with two more points needing to be fixed).

However, even after the scipy PR#16646 patch is applied, pythran still causes more than a dozen packages to fail their tests on armhf. In this case the packages run "fine" with scipy+pythran, the problem is that the calculated results do not agree with expected results within the expected tolerance. By contrast, tests on other architectures (including 32-bit architectures armel and i386) are passing. That suggests the problem must in the way pythran uses the hardware floating point unit of the armhf architecture.

I will therefore disable pythran support from scipy in the debian build for armhf. But before that is in place, I'm documenting the problem here with a screenshot listing the affected packages:

That is, affected packages are
astropy/5.1-3
astroscrappy/1.1.0-2
ccdproc/2.3.1-1
healpy/1.15.0-2
libflame/5.2.0-3
libgpuarray/0.7.6-9
pandas/1.3.5+dfsg-5
phonopy/2.15.1-1
pyfai/0.21.3+dfsg1-1
python-bayespy/0.5.22-2
scipy/1.8.1-9
skimage/0.19.3-1
statsmodels/0.13.2+dfsg-4
xrayutilities/1.7.3-1

Also failing on armhf, though possibly not because of pythran:
silx/1.0.0+dfsg-4

The test logs for these packages can be found via https://ci.debian.net/
The relevant tests are those run against scipy/1.8.1-9
e.g. for astropy: https://ci.debian.net/packages/a/astropy/testing/armhf/ with failing log https://ci.debian.net/data/autopkgtest/testing/armhf/a/astropy/25113914/log.gz

[drew-parsons]

Here is a sample error (from astropy):

___________ test_nterms_methods[standard-4-none-True-False-fastchi2] ___________

method = 'fastchi2', center_data = False, fit_mean = True, errors = 'none'
nterms = 4, normalization = 'standard'
data = (array([12.73923375,  5.39573428,  0.81947048,  0.33055271, 16.26540478,
       18.25511155, 12.13271552, 14.58993122,...7, 1.02602225, 1.02343473, 0.58893564, 1.48194269,
       1.0713956 , 0.50640888, 1.2726492 , 1.47826571, 1.08987003]))

    @pytest.mark.parametrize('method', NTERMS_METHODS)
    @pytest.mark.parametrize('center_data', [True, False])
    @pytest.mark.parametrize('fit_mean', [True, False])
    @pytest.mark.parametrize('errors', ['none', 'partial', 'full'])
    @pytest.mark.parametrize('nterms', [0, 2, 4])
    @pytest.mark.parametrize('normalization', NORMALIZATIONS)
    def test_nterms_methods(method, center_data, fit_mean, errors,
                            nterms, normalization, data):
        t, y, dy = data
        frequency = 0.8 + 0.01 * np.arange(40)
    
        if errors == 'none':
            dy = None
        elif errors == 'partial':
            dy = dy[0]
        elif errors == 'full':
            pass
        else:
            raise ValueError(f"Unrecognized error type: '{errors}'")
    
        ls = LombScargle(t, y, dy, center_data=center_data,
                         fit_mean=fit_mean, nterms=nterms,
                         normalization=normalization)
    
        if nterms == 0 and not fit_mean:
            with pytest.raises(ValueError) as err:
                ls.power(frequency, method=method)
            assert 'nterms' in str(err.value) and 'bias' in str(err.value)
        else:
            P_expected = ls.power(frequency)
    
            # don't use fast fft approximations here
            kwds = {}
            if 'fast' in method:
                kwds['method_kwds'] = dict(use_fft=False)
            P_method = ls.power(frequency, method=method, **kwds)
    
>           assert_allclose(P_expected, P_method, rtol=1E-7, atol=1E-25)
E           AssertionError: 
E           Not equal to tolerance rtol=1e-07, atol=1e-25
E           
E           Mismatched elements: 2 / 40 (5%)
E           Max absolute difference: 5.61386531e-05
E           Max relative difference: 0.00088298
E            x: array([0.064194, 0.045562, 0.075394, 0.097386, 0.095926, 0.117525,
E                  0.088898, 0.066776, 0.050198, 0.02832 , 0.076128, 0.113121,
E                  0.058681, 0.125867, 0.133244, 0.134014, 0.062316, 0.206491,...
E            y: array([0.064194, 0.045562, 0.075394, 0.097386, 0.095926, 0.117525,
E                  0.088898, 0.066776, 0.050198, 0.02832 , 0.076128, 0.113121,
E                  0.058681, 0.125867, 0.133244, 0.134014, 0.062316, 0.206491,...

/usr/lib/python3/dist-packages/astropy/timeseries/periodograms/lombscargle/tests/test_lombscargle.py:185: AssertionError

[serge-sans-paille]

Let's try to debug that one. I unfortunately don't have access to an armhf box, but this comment

That suggests the problem must in the way pythran uses the hardware floating point unit of the armhf

caught my attention. Pythran is not doing anything special with that respect... I assume scipy is not compiled with xsimd support, correct?

[serge-sans-paille]

Just to confirm: these issues also appear on armhf when PR#16646 is not applied?

[drew-parsons]

The armhf build log for scipy 1.8.1-9 is https://buildd.debian.org/status/fetch.php?pkg=scipy&arch=armhf&ver=1.8.1-9&stamp=1661219183&raw=0

It was build against xsimd 8.1.0-7 - pythran still needs it to satisfy scalar function headers. Is it possible it's making general ARM assumptions that don't apply on armhf?

scipy 1.8.1-5 had armhf configured with pythran without PR#16646. astropy 5.0.4-1 passed tests at that point. Likewise phonopy 2.15.1-1 passed tests then. xrayutilities 1.7.3-1 gave the same SIGSEGV fault.

Other things may have changed between the scipy 1.8.1-5 and 1.8.1-9 builds, e.g. switch from gcc 11 to gcc 12

[drew-parsons]

Update after applying scipy PR#16646 (the latest revised version) with the debian scipy build 1.8.1-13.

Many packages are still failing on armhf (screenshot below). It's not absolutely clear why they're failing. The tolerance problems I reported above aren't showing now, perhaps they were fixed by PR#16646. What's left is a mix of memory overflows and segfaults. It's not clear that the errors are in fact caused by pythran, except that they only seem to occur when scipy is configured with pythran support on armhf. Perhaps pythran increases memory requirements exceeding what's available on the armhf machines.

Next I'll try deactivating pythran on armhf to confirm which packages fail tests without pythran.

Current screenshot of armhf package failures with scipy 1.8.1-13 configured to use pythran 0.11.0+ds-7: