2021-03-09, 22:44 | #1 |
Mar 2021
2^{2}×11 Posts |
Alternative Math for Lucas Lehmer Testing and Perfect Numbers
Alternative Math for Lucas Lehmer Testing and Perfect Numbers
I'm posting (for the first time) some alternative math i've discovered from studying primes that i wanted to share with others in case it's of interest to anyone. The first is an alternative to the Lucas Lehmer Test with comparable time results. It is a pow (with 3) in this form: pow(3, ((2**x-1)-1)//(x), 2**x-1) where the pow result is a powers of two number only for Mersenne Numbers M2 and above, and the result is 0 for M1. To use with the code simply use checkp2(mersenne_lr_k_ary(3, ((2**x-1)-1)//(x), 2**x-1, k=8) where k can be any number up to the limit you can calculate it too. The larger the k value, less multiplications are required. You can use debug and withstats ( without checkp2 ) to see the results. Here is example of it's use: Code:
# In [717]: x=mersenne[26] # In [718]: start = time.time() # ...: x=gmpy2.mpz(x) # ...: xx=checkp2(mersenne_lr_k_ary(3, ((2**x-1)-1)//(x), 2**x-1,8)) # ...: end = time.time() # ...: print(end-start) # 17.363231420516968 # This is python, so a c version of this code that uses a graphics cpu # # should be much faster, but this equivilent to lucas lehmer. # In [719]: xx # Out[719]: True # In [720]: x=89 # In [721]: mersenne_lr_k_ary(3, ((2**x-1)-1)//(x), 2**x-1,8) # Out[721]: mpz(68719476736) # In [722]: y = mersenne_lr_k_ary(3, ((2**x-1)-1)//(x), 2**x-1,8) # In [723]: checkp2(y) # Out[723]: True for this pow would make it faster than lucas lehmer, which is interesting, so i wanted to show the math here. The second is a corollary to to Perfect Numbers. Instead of testing the additive of the factors, we test for the bitlength of the oddfactor of the number ( gmpy2.mpz(temp//(two**(ffs(temp)))) ) where temp is the number you are testing. An example of this test for the perfect number for 127 which is 8128 ( altpn(127) ). If the prime is 2 return 6, this test works only for powers of 2 results. 2 returns 0 which isn't a powers of two. The first step is to check for the bit_length of 127: Code:
# temp=8128 # first we extract the odd prime: # In [712]: gmpy2.mpz(temp//(two**(ffs(temp)))) # Out[712]: mpz(127) # x = gmpy2.mpz(temp//(two**(ffs(temp)))) # Then we divide by 127 from temp. # # In [710]: gmpy2.mpz(temp//(gmpy2.mpz(temp//(two**(ffs(temp)))) )) # Out[710]: mpz(64) # y = gmpy2.mpz(temp//(gmpy2.mpz(temp//(two**(ffs(temp)))) )) # The actual test is that y=64 and x=127 have the same bit_length and that # y and x+1 are powers of two numbers. If both of these is true then we only need to # test with a primetest for the odd factor to test it's primality. # Here below, is an example of it's use. This alternative math will only test positive for # perfect numbers, making it sort of a corollary to the actual proof. # In [713]: levelx(8128) # Out[713]: True Here is the code for the above for anyone interested in it: Code:
# Here is code for the math provided: import gmpy2 def mer_rem(x, bits, mask): # Below is same as: return x % mask while True: r = gmpy2.mpz(0) while x != 0: r += x & mask x >>= bits if r == mask: r = gmpy2.mpz(0) if r < mask: return gmpy2.mpz(r) x = r def mersenne_test(x): # This pow test is equivilant to a Lucas Lehmer Test and its # speed would be increased with any modular exponentiation # techiniques or discoveries that reduce the number of multiplications # to get to the result return checkp2(mersenne_lr_k_ary(3,((2**x-1)-1)//x, 2**x-1, 5)) def mersenne_modp2(n, b): return gmpy2.mpz(n&(b-1)) def checkp2(n): return (n & (n-1) == 0) and n not in [0,1] def altpn(N): N = gmpy2.mpz(N) return pow(N,3) + -N * pow(N,2) + ((N+1)//2) * N + 0 def ffs(x): x=gmpy2.mpz(x) return gmpy2.bit_length(x&-x)-1 def mersenne_modp2(n, b): return gmpy2.mpz(n&(b-1)) def mersenne_lr_k_ary(a, b, n, k=5, withstats=False, debug=False): """ Compute a ** b (mod n) K-ary LR method, with a customizable 'k'. k changes work for every mersenne ( 5 seems optimal ) """ base = 2 << (k - 1) n = gmpy2.mpz(n) nbl = gmpy2.bit_length(n) if debug==True: count=gmpy2.mpz(0) # Precompute the table of exponents table = [1] * base for i in range(1, base): table[i] = mer_rem(table[i - 1] * a, nbl, n) #table[i - 1] * a % nbl #mer_rem(table[i - 1] * a, nbl, n) # Just like the binary LR method, just with a # different base # r = gmpy2.mpz(1) for digit in reversed(mersenne_digits_of_n(b, base)): for i in range(k): r = mer_rem(r * r, nbl, n) if debug == True: count+=1 if withstats == True: print(f"In k: {i}, {r}, {checkp2(r)}") if digit: r = mer_rem(r * table[digit], nbl, n) if debug == True: count+=1 if withstats == True: print(f"In digit: {table[digit]}, {r}, {checkp2(r)}") if debug == True: return r, count return r def mersenne_digits_of_n(n, b): """ Return the list of the digits in the base 'b' representation of n, from LSB to MSB """ digits = [] while n: digits.append(mersenne_modp2(n,b)) #n % b) n //= b return digits # Alternate Odd Factor ( not odd number ) Perfect Number Test (levelx) # Since 2 has no oddfactor return True for 6 # Step 1, check for powers of two for non odd factor # Step 2, check that bit_length offset of odd factor # equals that of the non odd factor # Step 3, check that odd factor + 1 is a p2 number # Step 4, check primality of odd factor # Result: Perfect Number is found or not found def levelx(N, withstats=False): if N <= 1: return False N = gmpy2.mpz(N) zero, one, two = gmpy2.mpz(0), gmpy2.mpz(1), gmpy2.mpz(2) temp = gmpy2.mpz(N) newlevel = gmpy2.bit_length(temp) primetest = gmpy2.mpz(temp//(two**(ffs(temp)))) offset = gmpy2.mpz(temp//primetest) nextlevel = newlevel // two check = temp // (two**nextlevel) prevtemp = two**nextlevel if N == 6: return True if withstats == True: print (newlevel, newlevel, temp) print (newlevel, nextlevel+one, prevtemp, primetest) if (prevtemp & (prevtemp-one) == zero) and prevtemp != zero: if gmpy2.bit_length(offset) == gmpy2.bit_length(primetest): if ((primetest+one) & ((primetest+one)-one) == zero): s = checkp2(mersenne_lr_k_ary(3,primetest//gmpy2.bit_length(primetest), primetest, 5)) if withstats == True: print(s) if s == True: return True return False else: return False else: return False else: return False |
Thread Tools | |
Similar Threads | ||||
Thread | Thread Starter | Forum | Replies | Last Post |
Lucas-Lehmer Test for Fermat Numbers. | amenezes | Math | 8 | 2020-09-26 15:46 |
Lucas-Lehmer Primes | henryzz | And now for something completely different | 42 | 2019-06-03 14:09 |
lucas-lehmer theorem | Robot2357 | Math | 6 | 2013-06-15 03:10 |
lucas lehmer outstretch | science_man_88 | Miscellaneous Math | 7 | 2010-07-14 12:35 |
Lucas-Lehmer | Dougal | Information & Answers | 9 | 2009-02-06 10:25 |