[seqfan] Re: Testing potential squares

Giovanni Resta giovanni.resta at iit.cnr.it
Thu Aug 23 10:28:38 CEST 2018


The drawback of a big table (you say 1GB) is that you risk to
disrupt data locality and the CPU cannot store all your table in the 
faster cache memory inside CPU and has to access frequently external RAM.

This clearly depends on your CPU and the overall program, because
it may be the case you have anyway to make a lot of memory accesses.

Since usually I do not have large tables of potential squares stored in 
memory, but more often I generate the numbers on-the-fly and then I test 
them one-by-one (or in parallel), I prefer a solution that use little 
memory.

In particular, for numbers of 64 and 128 bits (natively available in 
GCC), I use 3 tables with modules (in this order):

9360 = 2^4 * 3^2 * 5 * 13  (1 number every 27.86 passes)

9163 = 7^2 * 11 * 17 (1 number every 7.713)

8303 = 19^2 * 23 (1 number every 4.02)

The total is 27.86 * 7.7 * 4.02 = 864 which means I have to take
a square root quite seldom.

Extracting bits or even bytes on can be costly on modern architecture 
that are oriented to larger words, so my tables are stored as arrays of 
short unsigned integers (2 byte = 16 bit). I do not remember since
I've done this some years ago. Maybe on my CPU and typical program
this was faster than using 1 byte, or 4 byte integers for my tables.
Probably faster than extracting single bits.

I've not made extensive tests, and surely I'm not a expert.

You just to have to make some experiments to find the best
solution for your own case.

Giovanni


On 08/21/2018 02:09 PM, Keith F. Lynch wrote:
> I recently wanted to test a large number of large numbers to see
> if they were squares.  I wanted to do so as quickly as possible.
> The odds of a positive integer around N being square are about
> 1/(2*sqrt(N)).  So if N is around a trillion (10^12), the odds
> are about 1 in 2 million.
> 
> The obvious way to test is to take the square root.  But that's slow.
> 
> To quickly rule out most non-squares, I created, just once, at the
> beginning of the program, a table of all quadratic residues mod 2^16
> (65536).  Then I tested every potential square by ANDing it with
> 2^16-1 (65535), and looking up the result in that table, both of which
> are very fast.  That rules out nearly 5/6 of all non-squares.  The
> proportions of quadratic residues modulo higher and higher powers of 2
> quickly approach 1/6.  (2^16 already gets 1/5.999268.  See A023105.)
> 
> To rule out most of the remaining non-squares, I searched for numbers
> that had a lower proportion of quadratic residues than any smaller
> number (e.g. 1441440, the first number such that fewer than one
> percent of numbers are quadratic residues modulo it), intending to
> create a table of residues mod such a number and test every potential
> square by calculating it mod that number and seeing if it was in the
> table.  After calculating the first few numbers that have a record
> low proportion of quadratic residues, I looked them up in OEIS, but
> couldn't find them.  Later I found they're at A085635; I had missed
> it because it lists the first term as 2 when I think it should be 1.
> 
> Should I correct that, or is there a reason it's 2?  Thanks.  I also
> have more terms to add to it and to its companion sequence A084848.
> They currently show 44 terms, the last of which is around a hundred
> thousand and has a ratio of about 50.  I've computed 193 terms so
> far, of which the last is about two thirds of a trillion and has a
> ratio of about 1588.
> 
> Later I decided that, having already weeded most non-squares out mod
> 2^16, I think I should instead look for odd numbers that have a lower
> proportion of quadratic residues than any smaller odd number.  That
> list is not in OEIS.  I'll add it as soon as someone tells me whether
> I should include 1 as its first term, and if not, why not.  Thanks.
> 
> Those two tests (mod 65536 and mod 740025 (the first odd number to
> have fewer than 3% of numbers modulo it to be residues)), together
> leave fewer than 1/200 of non-squares.  So instead of one in two
> million, squares are one in ten thousand.  So I only need try taking
> the square root of ten thousand numbers for every number that turns
> out to actually be a square.  Still not very good.
> 
> If I have one gigabyte to devote to the task, I can instead test my
> prospective square mod 7661478825, a larger record-setting odd number,
> which has about 1/233 residues.  That requires one bit per residue
> to get it to fit in one gigabyte.  If N is my prospective square mod
> 7661478825, I shift N right 3 bits to get the byte to look up in the
> table.  There will be about a 1/29 chance of a hit, i.e. at least one
> residue in the block of 8 containing the correct bit.  Only if there
> is a hit do I check the bit, which I do by ANDing N with 7, and using
> that as the index of which bit in that byte to look at.  So I've ruled
> out all but about 1/1400 non-squares.
> 
> I can do better.  Suppose I look at the record low proportions of
> residues, not for all odd numbers, but separately for numbers all of
> whose prime factors are congruent to 1 mod 4 and for numbers all of
> whose prime factors are congruent to 3 mod 4.  Since numbers in one of
> those sets will all be relatively prime to all numbers in the other of
> the two sets, and since for any two relatively prime numbers A and B,
> the number of residues mod AB will always be the product of the number
> of residues mod A with the number of residues mod B, the result will
> be the same by testing mod A followed by testing mod B as by testing
> mod A*B.  (A table of size A*B would require a prohibitive amount
> of memory.)
> 
> If, once again, I have one gigabyte for the purpose, I can test mod
> 3877273323 and (if it passes that test) mod 3833254945.  The former
> is the largest record-holding number less than 4*10^9 all of whose
> prime factors are congruent to 3 mod 4, and the latter is the largest
> record-holding number less than 4*10^9 all of whose prime factors are
> congruent to 1 mod 4.  The factorizations are 3^2 7^2 11 19 23 31 and
> 5 13 17 29 37 53 61 respectively.  And the proportions of quadratic
> residues are about 1/127 and about 1/85 respectively, or about
> 1/10800.  Combined with the first test, mod 65536, I've ruled out all
> but about 1/65000 squares.  And testing two mod tables rather than one
> doesn't slow things down much, since only if it passes the first one
> do I try the second one.  So I only need try taking the square root
> of about 30 numbers for every number that turns out to actually be a
> square (assuming the prospective squares are around a trillion).
> 
> I could do even better by having four tables, one for numbers all of
> whose prime factors are mod 1 mod 8, another for numbers all of whose
> prime factors are mod 3 mod 8, and likewise for 5 and 7.  Or eight
> tables, one for each possible remainder mod 16.  (Again, these are for
> all of the prime factors of the large number with the record low (for
> its size) proportion of quadratic residues, not for that number itself.)
> 
> I looked into instead checking mod some power of each successive
> prime.  After checking mod 2^16 (with an AND) I could check mod 3^16,
> then 5^16, then 7^16, etc.  But that would require more taking of mods,
> since in the limit of high powers, powers of 2 have 1/6 residues, and
> powers of odd primes have P/(2(P+1)) residues, i.e. 3/8, 5/12, 7/16,
> etc.  That's 1/2 in the limit of large power of large primes.
> 
> That's not a ridiculously bad approach, especially if you're short
> on memory.  The average number of mods, assuming that the potential
> square passed the 65536 test, would be the sum of the partial products
> of each odd prime over two more than the prime, i.e.
> 1 + 3/8 + (3/8 * 5/12) + (3/8 * 5/12 * 7/16) + (3/8 * 5/12 * 7/16 * 11/24) + ...
> which comes to about 1.658653437258581054.  For comparison, my two odd
> mods approach would, with the numbers given, require 1.0079 mods on
> average.
> 
> (Of course in every algorithm you'd stop doing mods and take the square
> root once the probability of the number being square got high enough.
> So only the first few digits of my 1.658... number above are meaningful
> unless you're testing truly stupendous potential squares.  But anything
> worth doing is worth overdoing.  Especially since that presumably-
> transcendental prime-related number doesn't seem to have ever appeared
> anywhere before.)
> 
> I'm curious whether anyone can think of a faster way to test large
> numbers for squareness.
> 
> I'm also curious what the ultimate behavior of the A085635/A084848
> proportions is.  I know they grow without limit.  But at what rate
> do the record-setting ratios increase?  Thanks.
> 
> --
> Seqfan Mailing list - http://list.seqfan.eu/
> 




More information about the SeqFan mailing list