ENH: Use AVX for float32 implementation of np.sin & np.cos - #13368

I realize the first 6 commits are exactly same as in this PR: #13134. The only commit that is new is the last one which implements sine and cosine (but it relies on the previous commits).

For reviewers - if the other PR is merged, you can switch the base of this PR between <anything> and back to master, and it will recompute the diff and throw out all the merged commits.

#13134 is merged. Could you fix the conflicts/rebase/switch the base?

Thanks! Its updated now and has just the one commit.

How do current benchmarks compare before/after this PR?

As mentioned in the commit message, the performance gains are leveraged only if the input arguments to sine/cosine are between [-71476.0625f, 71476.0625f] for cosine and [-117435.992f,
117435.992f] for sine. I measured 3 different benchmarks.

(1) For micro-benchmarks:
Array size = 10000, Command = "%timeit np.cos([myarr])"

(2) Package ai.cs provides an API to transform spherical coordinates to
cartesean system:
Array size = 10000, Command = "%timeit ai.cs.sp2cart(r,theta,phi)"

(3) Package photutils provides an API to find the best fit of first and
second harmonic functions to a set of (angle, intensity) pairs:
Array size = 1000, Command = "%timeit fit_first_and_second_harmonics(E, data)"

Sorry for not being clear. I meant that I would expect that a performance PR would come with changes in benchmarks in numpy.benchmarks, run by

python runtests.py --bench

As far as I could tell, those benchmarks were running the float64 version of sine and cosine. Is there a way I could run the float32 version of it? Or do they already do that?

I will take a look at this today and let you know if I am able to run the benchmarks for float32.

I ran the ufunc benchmarks python runtests.py --bench ufunc and saw this:

But the thing is, the benchmarks for ufunc seems to reporting just one metric (sum of the total time, I assume) for a whole bunch of dtypes (int16, int32, float16, float32, float128, complex64, complex256), whereas my patch only improves performance for float32. So I am not sure if these numbers are representative of performance of this patch.

I added my own benchmarks which measures computing sine and cosine for just float32:

class Trignometry(Benchmark):                                                   
    def setup(self):                                                            
        self.f = np.ones(10000, dtype=np.float32)                           
                                                                                
    def time_sine(self):                                                        
        np.sin(self.f)                                                          
                                                                                
    def time_cosine(self):                                                      
        np.cos(self.f)  

and this is what I measured:

This makes _multiarray_umath.so on Ubuntu 18.04 with gcc-8 grow from 17271072 to 17299296, which is 27.5 kB or 0.16%.

Made some edits to handle NAN and added some tests as well.

Added a commit which fixes the bug reported in #13493 for sin and cos.

Updated with following fixes:

• Added tests to validate sin/cos
• using AVX for strided arrays (fixes the same issue reported for exp/log)
• sin/cos will require an fma instead of just avx2. Otherwise the results are in-consistent across avx2 and avx512 implementations.
• Rounded up the max tolerable ulp error to an int (the assert_array_max_ulp returns an int)

Unused variable warnings.

yeah, thanks :) It was in cpuid.c, just fixed it.

seems like some network issues resulted in the one test failing and shouldn't be related to the patch itself.

just checking if there is any more feedback for this patch?

@mattip @juliantaylor ping ..

Hi @mattip, just wanted to know what you are thinking with respect to this patch. I have more performance patches I am working on but they rely on some of the commits in this patch.

Re-based with master. The AVX functions for sin and cos pass the tests on my local skylake machine (though it is currently disabled in the CI).

@r-devulap what would it take to enable the tests in CI? Do we need to find a skylake machine?

If you have a cpu > Haswell then the AVX2/FMA algorithms should already be tested. SkylakeX will be needed for AVX512F. Additionally, I also need a way to enable the test only on these platforms and ignore on the rest.

@mattip ping, how do you want to proceed with PR? Is there a way to enable the tests only for AVX machines?

My inclination would be to put this in as-is, since you made sure tests pass on the relevant systems. I am not sure I understand the problem. I see all the checks pass sucessfully. Which tests are enabled for AVX-only?

The test test_validate_transcendentals in the file numpy/core/tests/test_umath_accuracy.py is currently disabled completely in the CI. The reason is that it is different platforms/compilers seem to generate different outputs for the same functions and the accuracy varies significantly. I would like to know if there is way to enable this test in the CI only for platforms that NumPy utilizes the AVX version of sin/cos/exp/log and ignore the test on other platforms.

I guess you could add a function to numpy/core/src/umath/_umath_tests.c which would export a C-level detection of AVX, then call it as part of the setup for those tests.

Would you like to do that as part of this PR or in a separate, follow-on one?

that might take me some time to implement, perhaps a separate follow up PR is better?

OK, let's put this in as-is, in the hope we will improve #14048 's testing over time

thanks @r-devulap

OK, let's put this in as-is, in the hope we will improve #14048 's testing over time

sounds good, thanks!

I've done the SkylakeX CI emulation testing upstream in OpenBLAS for AVX-related testing purposes: OpenMathLib/OpenBLAS#2134

Not sure that NumPy would really want to go that far though.