View Single Post
Old 2020-09-10, 13:20   #1
Viliam Furik
 
"Viliam Furík"
Jul 2018
Martin, Slovakia

52×13 Posts
Default Python script for search for factors of M1277 using random k-intervals

Code:
import numpy as np
from random import randint

f = open("m1277.txt", "a")


def tf1277(k, r):
    expos = np.zeros((4, r), int)
    p = 1277
    y = 0
    z = 0
    e1 = 2555
    e2 = 3832
    e3 = 7663
    e4 = 8940
    pbinary = []
    q = p
    while q > 0:
        if q % 2 == 0:
            pbinary.append(0)
            q //= 2
        if q % 2 == 1:
            pbinary.append(1)
            q //= 2
    pbinary.reverse()
    for x in range(4):
        for prime in [3, 5, 7, 11, 13, 17, 19, 23, 29]:
            for a in range(prime + 5):
                if x == 0:
                    e = e1
                elif x == 1:
                    e = e2
                elif x == 2:
                    e = e3
                else:
                    e = e4
                if ((a + k) * 10216 + e) % prime == 0:
                    c = a
                    while c < r:
                        expos[x][c] = 1
                        c += prime
                    break
    for a in range(len(expos[0])):
        if expos[0][a] == 0:
            c = 1
            kp = (k + a) * 10216 + 2555
            for d in pbinary:
                c **= 2
                if d == 1:
                    c *= 2
                c %= kp
            if c == 1:
                print("-------------------------")
                print("factor = ", kp)
                print("k = ", (kp - 1) // 2 // 1277)
                f.write("\n")
                f.write("factor = " + str(kp))
                f.write("\n")
                f.write("k = " + str((kp - 1) // 2 // 1277))
                f.write("\n")
                quit()
            if y == 1000:
                y = 0
                z += 1
                print(z * 1000)
            y += 1
    for a in range(len(expos[1])):
        if expos[1][a] == 0:
            c = 1
            kp = (k + a) * 10216 + 3832
            for d in pbinary:
                c **= 2
                if d == 1:
                    c *= 2
                c %= kp
            if c == 1:
                print("-------------------------")
                print("factor = ", kp)
                print("k = ", (kp - 1) // 2 // 1277)
                f.write("\n")
                f.write("factor = " + str(kp))
                f.write("\n")
                f.write("k = " + str((kp - 1) // 2 // 1277))
                f.write("\n")
                quit()
            if y == 1000:
                y = 0
                z += 1
                print(z * 1000)
            y += 1
    for a in range(len(expos[2])):
        if expos[2][a] == 0:
            c = 1
            kp = (k + a) * 10216 + 7663
            for d in pbinary:
                c **= 2
                if d == 1:
                    c *= 2
                c %= kp
            if c == 1:
                print("-------------------------")
                print("factor = ", kp)
                print("k = ", (kp - 1) // 2 // 1277)
                f.write("\n")
                f.write("factor = " + str(kp))
                f.write("\n")
                f.write("k = " + str((kp - 1) // 2 // 1277))
                f.write("\n")
                quit()
            if y == 1000:
                y = 0
                z += 1
                print(z * 1000)
            y += 1
    for a in range(len(expos[3])):
        if expos[3][a] == 0:
            c = 1
            kp = (k + a) * 10216 + 8940
            for d in pbinary:
                c **= 2
                if d == 1:
                    c *= 2
                c %= kp
            if c == 1:
                print("-------------------------")
                print("factor = ", kp)
                print("k = ", (kp - 1) // 2 // 1277)
                f.write("\n")
                f.write("factor = " + str(kp))
                f.write("\n")
                f.write("k = " + str((kp - 1) // 2 // 1277))
                f.write("\n")
                quit()
            if y == 1000:
                y = 0
                z += 1
                print(z * 1000)
            y += 1


while True:
    tf1277(randint(58000000000000000, 2 ** (1277 / 2) / 1277 / 2), 1000000)
Hello, I have made the above Python script, that is called with random starting k value, and then tests about 1263000 sived-out k's above that value. It uses a grid of zeroes to sieve k's and combines the facts that factors must be 2*k*p+1 and also 1 or 7 (mod 8) - so it only tests factors that are 2555, 3832, 7663, 8940 (mod 10216).

If there is a mistake (either logical or mathematical) in my script, I'll be grateful for pointing out. I would also be grateful if somebody could help me with making the code runnable on a Nvidia GPU with CUDA (or in OpenCL, if it's easier), without necessarily using FFT or Barrett multiplication. (But if somebody is willing to explain, I'll be happy to listen, or read...).

And for the last, if you want, you can try to find the factor yourself with that code. It takes my CPU about 12 seconds on average to go through random range of size 10216000000 (that is difference between largest and smallest k tested). So if I was to randomly check the equivalent number of ranges (untested range from 2^67 to 2^639 divided by range size 10216000000), it would take me... about 2*10^172 years . But if more people engage in search, probability increases (unfortunately it's only slight change).

So if you want to help, you can run the script any time you like. If you want to help more, you can help with improving the speed of the script.
Viliam Furik is offline   Reply With Quote