-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathexplicit_is_omega_prime.sf
67 lines (53 loc) · 1.45 KB
/
explicit_is_omega_prime.sf
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
#!/usr/bin/ruby
# a(n) is the smallest prime p such that omega(p^n - 1) is equal to n, where omega = A001221.
# https://oeis.org/A379450
# Known terms:
# 3, 5, 7, 11, 31, 11, 79, 47, 211, 113, 2473, 47, 45841, 389, 1123, 1061
# Jinyuan Wang: I get a(17) > 729607, a(18) = 373, a(19) > 1600000, a(20)-a(22)=2141, 83071, 43541. Computed with isomega function in https://oeis.org/history?seq=A219019
# PARI/GP program:
# a(n) = forprime(p=2, oo, if(omega(p^n-1) == n, return(p)));
# Upper-bounds:
# a(17) <= 2487619
# Lower-bounds:
# a(17) > 732181
# a(19) > 1600000
# Conjectured lower-bounds:
# a(17) > 2126167
# a(19) > 7955417
# a(23) > 198193
local Num!VERBOSE = true
local Num!USE_FACTORDB = true
local Num!USE_CONJECTURES = true
func my_is_omega(n, k) {
var f = [n]
var r = 0
while (k > 1) {
var t = f.last.iroot(k)+1
if (r < t) {
r = t
f = factor_upto(f.last, r).uniq
if (f.len == 1) {
return false
}
if (f.last.is_prime) {
return (k == f.len)
}
k = (k - f.len + 1)
}
else {
return false
}
}
return false
}
func a(n, from=2) {
for k in (from..Inf) {
k.is_prime || next
var t = (k**n - 1)
say "[#{n}] Checking k = #{k}: #{t}"
if (t.is_omega_prime(n)) {
#if (my_is_omega(t, n)) {
return k
}
}
}