45 lines
960 B
Text
45 lines
960 B
Text
get "libhdr"
|
|
|
|
let calcsigmas(sig, n) be
|
|
$( sig!0 := 0
|
|
for i = 0 to n do sig!i := 0
|
|
for i = 1 to n/2 do
|
|
$( let j = i
|
|
while 0 < j <= n do
|
|
$( sig!j := sig!j + i
|
|
j := j + i
|
|
$)
|
|
$)
|
|
$)
|
|
|
|
let gcd(m, n) = n=0 -> m, gcd(n, m rem n)
|
|
|
|
let duff(sig, n) = sig!n > n+1 & gcd(n, sig!n) = 1
|
|
let triple(sig, n) = duff(sig, n) & duff(sig, n+1) & duff(sig, n+2)
|
|
|
|
let first(sig, f, max, cb) be
|
|
$( let n = 0
|
|
for i = 1 to max
|
|
$( n := n+1 repeatuntil f(sig, n)
|
|
cb(i, n)
|
|
$)
|
|
$)
|
|
|
|
let start() be
|
|
$( let showsingle(i, n) be
|
|
$( writef("%I4", n)
|
|
if i rem 10=0 then wrch('*N')
|
|
$)
|
|
|
|
let showtriple(i, n) be writef("%I2: %I6 %I6 %I6*N", i, n, n+1, n+2)
|
|
|
|
let sig = getvec(20000)
|
|
calcsigmas(sig, 20000)
|
|
|
|
writes("First 50 Duffinian numbers:*N")
|
|
first(sig, duff, 50, showsingle)
|
|
|
|
writes("*NFirst 15 Duffinian triples:*N")
|
|
first(sig, triple, 15, showtriple)
|
|
freevec(sig)
|
|
$)
|