123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081 |
- # Finding lucas pseudoprimes
- from math import sqrt, gcd
-
- def isqrt(n):
- """ Find the integer square root of n via newton's method, code via
- stackoverflow:
-
- (https://stackoverflow.com/questions/15390807/integer-square-root-in-python)
- """
- x = n
- y = (x + 1) // 2
- while y < x:
- x = y
- y = (x + n // x) // 2
- return x
-
- def hasIntSQRT(n):
- """Detect whether the square root of n is an integer,
- i.e. whether the isqrt(n) is the true square root.
- """
- isq = isqrt(n)
- return isq * isq == n
-
- def Dsequence():
- """Generate sequence 5, -7, 9, -11, 13, -15...
- """
- val = 5
- while True:
- if val % 4 == 1:
- yield val
- else:
- yield -val
- val = val + 2
-
- def Legendere(a, p):
- if (p % 2 == 0):
- raise ValueError("p must be odd, is {}".format(p))
- lv = pow(a, (p-1)//2, p)
- if lv == p - 1:
- lv = -1
- return lv
-
- def Jacobi(a, n):
- if (n % 2 == 0):
- raise ValueError("n must be odd, is {}".format(n))
- if a % n != a:
- return Jacobi(a%n, n)
- if a != 0 and a % 2 == 0:
- nm8 = n % 8
- if (nm8 == 3 or nm8 == 5):
- return -1 * Jacobi(a//2, n)
- else:
- return Jacobi(a//2, n)
- if a == 1:
- return 1
- if gcd(a, n) != 1:
- return 0
- if a == 0 and n == 1:
- return 1
- if n % 4 == 3 and a % 4 == 3:
- return -1 * Jacobi(n, a)
- return Jacobi(n, a)
-
- print(" ", end=" ")
- for k in range(1, 10):
- print("{}".format(k).rjust(2, " "), end=" ")
- print()
-
- for n in range(1, 21, 2):
- print("{}:".format(n).rjust(3, " "), end=" ")
- for k in range(1, 10):
- print("{}".format(Jacobi(k, n)).rjust(2, " "), end=" ")
- print()
-
-
-
-
-
-
-
|