Sequence: nth number with exactly n divisors: Difference between revisions

Content added Content deleted
m (→‎{{header|REXX}}: added an optimized REXX version.)
(→‎{{header|Go}}: Added faster version.)
Line 139: Line 139:
32 : 3990
32 : 3990
33 : 12166144
33 : 12166144
</pre>

The following much faster version (runs in less than 30 seconds on my 1.6GHz Celeron) uses three further optimizations:

1. Apart from the 2nd and 10th terms, all the even terms are themselves even.

2. A sieve is used to generate all prime divisors needed. This doesn't take up much time or memory but speeds up the counting of all divisors considerably.

3. While searching for the nth number with exactly n divisors, where feasible a record is kept of any numbers found to have exactly k divisors (k > n) so that the search for these numbers can start from a higher base.

<lang go>package main

import (
"fmt"
"math"
"math/big"
)

type record struct{ num, count int }

var (
bi = new(big.Int)
primes = []int{2}
)

func isPrime(n int) bool {
bi.SetUint64(uint64(n))
return bi.ProbablyPrime(0)
}

func sieve(limit int) {
c := make([]bool, limit+1) // composite = true
// no need to process even numbers
p := 3
for {
p2 := p * p
if p2 > limit {
break
}
for i := p2; i <= limit; i += 2 * p {
c[i] = true
}
for {
p += 2
if !c[p] {
break
}
}
}
for i := 3; i <= limit; i += 2 {
if !c[i] {
primes = append(primes, i)
}
}
}

func countDivisors(n int) int {
count := 1
for i, p := 0, primes[0]; p*p <= n; i, p = i+1, primes[i+1] {
if n%p != 0 {
continue
}
n /= p
count2 := 1
for n%p == 0 {
n /= p
count2++
}
count *= (count2 + 1)
if n == 1 {
return count
}
}
if n != 1 {
count *= 2
}
return count
}

func isOdd(x int) bool {
return x%2 == 1
}

func main() {
sieve(22000)
const max = 45
records := [max + 1]record{}
z := new(big.Int)
p := new(big.Int)
fmt.Println("The first", max, "terms in the sequence are:")
for i := 1; i <= max; i++ {
if isPrime(i) {
z.SetUint64(uint64(primes[i-1]))
p.SetUint64(uint64(i - 1))
z.Exp(z, p, nil)
fmt.Printf("%2d : %d\n", i, z)
} else {
count := records[i].count
if count == i {
fmt.Printf("%2d : %d\n", i, records[i].num)
continue
}
odd := isOdd(i)
k := records[i].num
l := 1
if !odd && i != 2 && i != 10 {
l = 2
}
for j := k + l; ; j += l {
if odd {
sq := int(math.Sqrt(float64(j)))
if sq*sq != j {
continue
}
}
cd := countDivisors(j)
if cd == i {
count++
if count == i {
fmt.Printf("%2d : %d\n", i, j)
break
}
} else if cd > i && cd <= max && records[cd].count < cd &&
j > records[cd].num && (l == 1 || (l == 2 && !isOdd(cd))) {
records[cd].num = j
records[cd].count++
}
}
}
}
}</lang>

{{out}}
<pre>
The first 45 terms in the sequence are:
1 : 1
2 : 3
3 : 25
4 : 14
5 : 14641
6 : 44
7 : 24137569
8 : 70
9 : 1089
10 : 405
11 : 819628286980801
12 : 160
13 : 22563490300366186081
14 : 2752
15 : 9801
16 : 462
17 : 21559177407076402401757871041
18 : 1044
19 : 740195513856780056217081017732809
20 : 1520
21 : 141376
22 : 84992
23 : 1658509762573818415340429240403156732495289
24 : 1170
25 : 52200625
26 : 421888
27 : 52900
28 : 9152
29 : 1116713952456127112240969687448211536647543601817400964721
30 : 6768
31 : 1300503809464370725741704158412711229899345159119325157292552449
32 : 3990
33 : 12166144
34 : 9764864
35 : 446265625
36 : 5472
37 : 11282036144040442334289838466416927162302790252609308623697164994458730076798801
38 : 43778048
39 : 90935296
40 : 10416
41 : 1300532588674810624476094551095787816112173600565095470117230812218524514342511947837104801
42 : 46400
43 : 635918448514386699807643535977466343285944704172890141356181792680152445568879925105775366910081
44 : 240640
45 : 327184
</pre>
</pre>


Retrieved from "https://rosettacode.org/wiki/Sequence:_nth_number_with_exactly_n_divisors"