optimization opportunities for sidef-scripts/Math /jacobi_symbol.sf

I've noticed that in the jacobi file

sidef-scripts/Math

/jacobi_symbol.sf

    var t = 1
    while (a %= n) {
        if (a.is_even) {
            var v = a.valuation(2)
            t *= (-1)**v if (n%8 ~~ [3,5])
            a >>= v
        }
        (a,n) = (n,a)
        t = -t if ([a%4, n%4] == [3,3])
    }

    n==1 ? t : 0

All the heavy duty is performed whenever it's even (which usually ends up performing a LOT of extra times). But due to Legendre / Jacobi / Kronecker symbols being completely multiplicative, all pairs of the same prime have no effect on the symbol as long as co-primality against the underlying prime could be asserted.

Of course it's a huge hassle to detect it all, but since this isn't a "square-free test", so let's ignore the hard stuff for that. Whenever the code reaches -

if ( a.is_even ) {   
 ....
}

You can save work by checking # of trailing zeros in "a" at that moment, and only go through the heavy duty stuff when the CTZ count is odd (and still get to discard all trailing zeros found that loop cycle).

Another optimization I use in my own jacobi function is that since jacobi symbols' biggest use case, by far, is primality testing, before I engage with the slowest loop to manually inverted the symbols, I have bespoke fast path code optimized for the "a" of "jacobi( a | n )" being a small prime under 100 for which the final symbol can directly be calculated without "inverting" -

          jacobi(29, 971)  

          => (971 % 29)^((29 - 1) / 2) % 29 --> 28 --> [-1]

The code also includes the logic to also deal with the extra inversion due to both being 3 (mod 4), but isn't applicable for this particular example since 29 is 1 (mod 4)

And this power-mod isn't the bit-by-bit variant either. Because I've intentionally limited the scope to primes under 100, and that the largest possible exponent is 48 (for 97), the double64 data type can either directly raise the entire exponent in one shot (for primes under 25), or mod twice, but raising the 1st mod's result to its 8th-power (the integer precision within double64 is guaranteed by modulus under 99)

jacobi(5, anything positive odd > 5)

is easiest to calculate. Actually nothing to calculate at all - just return the "signed boolean" condition of whether last decimal digit is 3 / 7 for non-residues or 1 / 9 for residues. Heck one can even write simple regex to calculate jacobi symbol for ( a := 5 )

Yes this made the function much larger, but since the super vast majority of inputs are being captured by it, the fast path is effectively the function itself.

Start with Math/jacobi_symbol.sf and read the existing Jacobi symbol loop and its handling of even values. Compare the proposed trailing-zero and small-prime fast paths against the current behavior and benchmark representative inputs; done means a measured optimization that preserves correct Legendre, Jacobi, and Kronecker results.