Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Improve isprime performance for 32-bit integers #16349

Closed
wants to merge 1 commit into from

Conversation

pabloferz
Copy link
Contributor

As discussed here #11594 (mainly by @VicDrastik and @danaj). isprime has room for improvement. In particular this PR takes care only of the case of integers < 2^32. A version for 64-bit integers would require another strategy and should be done in another PR.

This version is much faster and uses less memory than the current implementation.

There is another possibility that was mentioned also in #11594, that is to use the algorithm due to Forišek and Jančina (which can be found here http://ceur-ws.org/Vol-1326/020-Forisek.pdf). @VicDrastik didn't seem to be fond of that solution but it is certainly faster than the one proposed here (although it uses a little more of memory 512 bytes, but we were already pre computing al primes up to 2^16, so this doesn't seem like much more.

I leave the code for that option below in case someone wants to consider it.

const bases = UInt16[15591,2018,166,7429,8064,16045,10503,4399,1949,1295,2776,3620,560,3128,5212,
2657,2300,2021,4652,1471,9336,4018,2398,20462,10277,8028,2213,6219,620,3763,4852,5012,3185,
1333,6227,5298,1074,2391,5113,7061,803,1269,3875,422,751,580,4729,10239,746,2951,556,2206,
3778,481,1522,3476,481,2487,3266,5633,488,3373,6441,3344,17,15105,1490,4154,2036,1882,1813,
467,3307,14042,6371,658,1005,903,737,1887,7447,1888,2848,1784,7559,3400,951,13969,4304,177,41,
19875,3110,13221,8726,571,7043,6943,1199,352,6435,165,1169,3315,978,233,3003,2562,2994,10587,
10030,2377,1902,5354,4447,1555,263,27027,2283,305,669,1912,601,6186,429,1930,14873,1784,1661,
524,3577,236,2360,6146,2850,55637,1753,4178,8466,222,2579,2743,2031,2226,2276,374,2132,813,
23788,1610,4422,5159,1725,3597,3366,14336,579,165,1375,10018,12616,9816,1371,536,1867,10864,
857,2206,5788,434,8085,17618,727,3639,1595,4944,2129,2029,8195,8344,6232,9183,8126,1870,3296,
7455,8947,25017,541,19115,368,566,5674,411,522,1027,8215,2050,6544,10049,614,774,2333,3007,
35201,4706,1152,1785,1028,1540,3743,493,4474,2521,26845,8354,864,18915,5465,2447,42,4511,
1660,166,1249,6259,2553,304,272,7286,73,6554,899,2816,5197,13330,7054,2818,3199,811,922,350,
7514,4452,3449,2663,4708,418,1621,1171,3471,88,11345,412,1559,194]

function witnesses(n::Union{UInt8,Int8,UInt16,Int16,UInt32,Int32})
    i = UInt64(n)
    i = ((i >> 16) $ i) * 0x45d9f3b
    i = ((i >> 16) $ i) * 0x45d9f3b
    i = ((i >> 16) $ i) & 255 + 1
    bases[i]
end

@kmsquire
Copy link
Member

Is primes used in base? Should it be moved to a package?

@pabloferz
Copy link
Contributor Author

Just for reference, the time it takes to check that all 32-bit the primes are primes is

635.517114 seconds (6.64 G allocations: 150.434 GB, 2.44% gc time) (master)
410.274658 seconds (2.73 G allocations: 89.169 GB, 2.15% gc time) (this PR)
309.419580 seconds (2.74 G allocations: 89.359 GB, 2.66% gc time) (Forišek and Jančina)

Also, checking whether or not each element of an array of 10^7 random UInt32s is prime, takes

4.734542 seconds (54.87 M allocations: 913.643 MB, 1.23% gc time) (master)
1.464444 seconds (1.85 M allocations: 28.153 MB, 0.07% gc time) (this PR)
1.188251 seconds (2.46 M allocations: 37.486 MB, 0.15% gc time) (Forišek and Jančina)

@pabloferz
Copy link
Contributor Author

As far as i can see, none of the functions in primes.jl is used anywhere in base and I think all of them could be moved into a package. I don't know if there is a particular interest to keep them here.

for m in (3,5,7,11,13,17,19,23)
n % m == 0 && return false
end
n < 841 && return n > 1
s = trailing_zeros(n-1)
d = (n-1) >>> s
for a in witnesses(n)
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think you can give a type hint

 for a in witnesses(n)::Tuple{Vararg{Int}}

making the witnesses stack allocated.

@JeffBezanson
Copy link
Member

This should now be moved to https://github.com/JuliaMath/Primes.jl

@pabloferz Great work, thanks for this.

@pabloferz
Copy link
Contributor Author

No problem. Moved over JuliaMath/Primes.jl#9

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.

4 participants