41 lines
1,022 B
Nim
41 lines
1,022 B
Nim
import std/[algorithm, math, strformat]
|
|
|
|
const MaxNumber = 500_000
|
|
|
|
# Construct a table of the divisor counts.
|
|
var ds: array[1..MaxNumber, int]
|
|
ds.fill 1
|
|
for i in 2..MaxNumber:
|
|
for j in countup(i, MaxNumber, i):
|
|
ds[j] += i
|
|
|
|
# Set the divisor counts of non-Duffinian numbers to 0.
|
|
ds[1] = 0 # 1 is not Duffinian.
|
|
for n in 2..MaxNumber:
|
|
let nds = ds[n]
|
|
if nds == n + 1 or gcd(n, nds) != 1:
|
|
# "n" is prime or is not relatively prime to its divisor sum.
|
|
ds[n] = 0
|
|
|
|
# Show the first 50 Duffinian numbers.
|
|
echo "First 50 Duffinian numbers:"
|
|
var dcount = 0
|
|
var n = 1
|
|
while dcount < 50:
|
|
if ds[n] != 0:
|
|
stdout.write &" {n:3}"
|
|
inc dcount
|
|
if dcount mod 25 == 0:
|
|
echo()
|
|
inc n
|
|
echo()
|
|
|
|
# Show the Duffinian triplets below MaxNumber.
|
|
echo &"The Duffinian triplets up to {MaxNumber}:"
|
|
dcount = 0
|
|
for n in 3..MaxNumber:
|
|
if ds[n - 2] != 0 and ds[n - 1] != 0 and ds[n] != 0:
|
|
inc dcount
|
|
stdout.write &" {(n - 2, n - 1, n): ^24}"
|
|
stdout.write if dcount mod 4 == 0: '\n' else: ' '
|
|
echo()
|