Duffinian numbers: Difference between revisions
m →{{header|ALGOL 68}}: notes and comments |
m →{{header|ALGOL 68}}: tweak |
||
Line 60: | Line 60: | ||
FOR n FROM 2 TO UPB ds DO |
FOR n FROM 2 TO UPB ds DO |
||
IF INT nds = ds[ n ]; |
IF INT nds = ds[ n ]; |
||
nds = n + 1 |
IF nds = n + 1 THEN TRUE ELSE gcd( n, nds ) /= 1 FI |
||
OR gcd( n, nds ) /= 1 |
|||
THEN |
THEN |
||
# n is prime or is not |
# n is prime or is not relatively prime to its divisor sum # |
||
ds[ n ] := 0 |
ds[ n ] := 0 |
||
FI |
FI |
Revision as of 19:22, 25 February 2022
A Duffinian number is a composite number k that is relatively prime to its sigma sum σ.
The sigma sum of k is the sum of the divisors of k.
- E.G.
161 is a Duffinian number.
- It is composite. (7 × 23)
- The sigma sum 192 (1 + 7 + 23 + 161) is relatively prime to 161.
Duffinian numbers are very common.
It is not uncommon for two consecutive integers to be Duffinian (a Duffinian twin) (8, 9), (35, 36), (49, 50), etc.
Less common are Duffinian triplets; three consecutive Duffinian numbers. (63, 64, 65), (323, 324, 325), etc.
Much, much less common are Duffinian quadruplets and quintuplets. The first Duffinian quintuplet is (202605639573839041, 202605639573839042, 202605639573839043, 202605639573839044, 202605639573839045).
It is not possible to have six consecutive Duffinian numbers
- Task
- Find and show the first 50 Duffinian numbers.
- Find and show at least the first 15 Duffinian triplets.
- See also
ALGOL 68
Constructs a table of divisor counts without doing any divisions. <lang algol68>BEGIN # find Duffinian numbers: non-primes relatively prime to their divisor count #
INT max number := 500 000; # largest number we will consider # # iterative Greatest Common Divisor routine, returns the gcd of m and n # PROC gcd = ( INT m, n )INT: BEGIN INT a := ABS m, b := ABS n; WHILE b /= 0 DO INT new a = b; b := a MOD b; a := new a OD; a END # gcd # ; # construct a table of the divisor counts # [ 1 : max number ]INT ds; FOR i TO UPB ds DO ds[ i ] := 1 OD; FOR i FROM 2 TO UPB ds DO FOR j FROM i BY i TO UPB ds DO ds[ j ] +:= i OD OD; # set the divisor counts of non-Duffinian numbers to 0 # ds[ 1 ] := 0; # 1 is not Duffinian # FOR n FROM 2 TO UPB ds DO IF INT nds = ds[ n ]; IF nds = n + 1 THEN TRUE ELSE gcd( n, nds ) /= 1 FI THEN # n is prime or is not relatively prime to its divisor sum # ds[ n ] := 0 FI OD; # show the first 50 Duffinian numbers # print( ( "The first 50 Duffinian numbers:", newline ) ); INT dcount := 0; FOR n WHILE dcount < 50 DO IF ds[ n ] /= 0 THEN # found a Duffinian number # print( ( " ", whole( n, -3) ) ); IF ( dcount +:= 1 ) MOD 25 = 0 THEN print( ( newline ) ) FI FI OD; print( ( newline ) ); # show the duffinian triplets below UPB ds # print( ( "The Duffinian triplets up to ", whole( UPB ds, 0 ), ":", newline ) ); dcount := 0; FOR n FROM 3 TO UPB ds DO IF ds[ n - 2 ] /= 0 AND ds[ n - 1 ] /= 0 AND ds[ n ] /= 0 THEN # found a Duffinian triplet # print( ( " (", whole( n - 2, -7 ), " ", whole( n - 1, -7 ), " ", whole( n, -7 ), ")" ) ); IF ( dcount +:= 1 ) MOD 4 = 0 THEN print( ( newline ) ) FI FI OD
END</lang>
- Output:
The first 50 Duffinian numbers: 4 8 9 16 21 25 27 32 35 36 39 49 50 55 57 63 64 65 75 77 81 85 93 98 100 111 115 119 121 125 128 129 133 143 144 155 161 169 171 175 183 185 187 189 201 203 205 209 215 217 The Duffinian triplets up to 500000: ( 63 64 65) ( 323 324 325) ( 511 512 513) ( 721 722 723) ( 899 900 901) ( 1443 1444 1445) ( 2303 2304 2305) ( 2449 2450 2451) ( 3599 3600 3601) ( 3871 3872 3873) ( 5183 5184 5185) ( 5617 5618 5619) ( 6049 6050 6051) ( 6399 6400 6401) ( 8449 8450 8451) ( 10081 10082 10083) ( 10403 10404 10405) ( 11663 11664 11665) ( 12481 12482 12483) ( 13447 13448 13449) ( 13777 13778 13779) ( 15841 15842 15843) ( 17423 17424 17425) ( 19043 19044 19045) ( 26911 26912 26913) ( 30275 30276 30277) ( 36863 36864 36865) ( 42631 42632 42633) ( 46655 46656 46657) ( 47523 47524 47525) ( 53137 53138 53139) ( 58563 58564 58565) ( 72961 72962 72963) ( 76175 76176 76177) ( 79523 79524 79525) ( 84099 84100 84101) ( 86527 86528 86529) ( 94177 94178 94179) ( 108899 108900 108901) ( 121103 121104 121105) ( 125315 125316 125317) ( 128017 128018 128019) ( 129599 129600 129601) ( 137287 137288 137289) ( 144399 144400 144401) ( 144721 144722 144723) ( 154567 154568 154569) ( 158403 158404 158405) ( 166463 166464 166465) ( 167041 167042 167043) ( 175231 175232 175233) ( 177607 177608 177609) ( 181475 181476 181477) ( 186623 186624 186625) ( 188497 188498 188499) ( 197191 197192 197193) ( 199711 199712 199713) ( 202499 202500 202501) ( 211249 211250 211251) ( 230399 230400 230401) ( 231199 231200 231201) ( 232561 232562 232563) ( 236195 236196 236197) ( 242063 242064 242065) ( 243601 243602 243603) ( 248003 248004 248005) ( 260099 260100 260101) ( 260641 260642 260643) ( 272483 272484 272485) ( 274575 274576 274577) ( 285155 285156 285157) ( 291599 291600 291601) ( 293763 293764 293765) ( 300303 300304 300305) ( 301087 301088 301089) ( 318095 318096 318097) ( 344449 344450 344451) ( 354481 354482 354483) ( 359551 359552 359553) ( 359999 360000 360001) ( 367235 367236 367237) ( 374543 374544 374545) ( 403201 403202 403203) ( 406801 406802 406803) ( 417697 417698 417699) ( 419903 419904 419905) ( 423199 423200 423201) ( 435599 435600 435601) ( 468511 468512 468513) ( 470449 470450 470451) ( 488071 488072 488073)
Factor
<lang factor>USING: combinators.short-circuit.smart grouping io kernel lists lists.lazy math math.primes math.primes.factors math.statistics prettyprint sequences sequences.deep ;
- duffinian? ( n -- ? )
{ [ prime? not ] [ dup divisors sum simple-gcd 1 = ] } && ;
- duffinians ( -- list ) 3 lfrom [ duffinian? ] lfilter ;
- triples ( -- list )
duffinians dup cdr dup cdr lzip lzip [ flatten ] lmap-lazy [ differences { 1 1 } = ] lfilter ;
"First 50 Duffinian numbers:" print 50 duffinians ltake list>array 10 group simple-table. nl
"First 15 Duffinian triplets:" print 15 triples ltake list>array simple-table.</lang>
- Output:
First 50 Duffinian numbers: 4 8 9 16 21 25 27 32 35 36 39 49 50 55 57 63 64 65 75 77 81 85 93 98 100 111 115 119 121 125 128 129 133 143 144 155 161 169 171 175 183 185 187 189 201 203 205 209 215 217 First 15 Duffinian triplets: 63 64 65 323 324 325 511 512 513 721 722 723 899 900 901 1443 1444 1445 2303 2304 2305 2449 2450 2451 3599 3600 3601 3871 3872 3873 5183 5184 5185 5617 5618 5619 6049 6050 6051 6399 6400 6401 8449 8450 8451
Julia
<lang julia>using Primes
function σ(n)
f = [one(n)] for (p,e) in factor(n) f = reduce(vcat, [f*p^j for j in 1:e], init=f) end return sum(f)
end
isDuffinian(n) = !isprime(n) && gcd(n, σ(n)) == 1
function testDuffinians()
println("First 50 Duffinian numbers:") foreach(p -> print(rpad(p[2], 4), p[1] % 25 == 0 ? "\n" : ""), enumerate(filter(isDuffinian, 2:217))) n, found = 2, 0 println("\nFifteen Duffinian triplets:") while found < 15 if isDuffinian(n) && isDuffinian(n + 1) && isDuffinian(n + 2) println(lpad(n, 6), lpad(n +1, 6), lpad(n + 2, 6)) found += 1 end n += 1 end
end
testDuffinians()
</lang>
- Output:
First 50 Duffinian numbers: 4 8 9 16 21 25 27 32 35 36 39 49 50 55 57 63 64 65 75 77 81 85 93 98 100 111 115 119 121 125 128 129 133 143 144 155 161 169 171 175 183 185 187 189 201 203 205 209 215 217 Fifteen Duffinian triplets: 63 64 65 323 324 325 511 512 513 721 722 723 899 900 901 1443 1444 1445 2303 2304 2305 2449 2450 2451 3599 3600 3601 3871 3872 3873 5183 5184 5185 5617 5618 5619 6049 6050 6051 6399 6400 6401 8449 8450 8451
Phix
with javascript_semantics sequence duffinian = {false} integer n = 2, count = 0, triplet = 0, triple_count = 0 while triple_count<50 do bool bDuff = not is_prime(n) and gcd(n,sum(factors(n,1)))=1 duffinian &= bDuff if bDuff then count += 1 if count=50 then sequence s50 = apply(true,sprintf,{{"%3d"},find_all(true,duffinian)}) printf(1,"First 50 Duffinian numbers:\n%s\n",join_by(s50,1,25," ")) end if triplet += 1 triple_count += (triplet>=3) else triplet = 0 end if n += 1 end while sequence s = apply(true,sq_add,{match_all({true,true,true},duffinian),{{0,1,2}}}), p = apply(true,pad_tail,{apply(true,sprintf,{{"[%d,%d,%d]"},s}),24}) printf(1,"First 50 Duffinian triplets:\n%s\n",{join_by(p,1,4," ")})
- Output:
First 50 Duffinian numbers: 4 8 9 16 21 25 27 32 35 36 39 49 50 55 57 63 64 65 75 77 81 85 93 98 100 111 115 119 121 125 128 129 133 143 144 155 161 169 171 175 183 185 187 189 201 203 205 209 215 217 First 50 Duffinian triplets: [63,64,65] [323,324,325] [511,512,513] [721,722,723] [899,900,901] [1443,1444,1445] [2303,2304,2305] [2449,2450,2451] [3599,3600,3601] [3871,3872,3873] [5183,5184,5185] [5617,5618,5619] [6049,6050,6051] [6399,6400,6401] [8449,8450,8451] [10081,10082,10083] [10403,10404,10405] [11663,11664,11665] [12481,12482,12483] [13447,13448,13449] [13777,13778,13779] [15841,15842,15843] [17423,17424,17425] [19043,19044,19045] [26911,26912,26913] [30275,30276,30277] [36863,36864,36865] [42631,42632,42633] [46655,46656,46657] [47523,47524,47525] [53137,53138,53139] [58563,58564,58565] [72961,72962,72963] [76175,76176,76177] [79523,79524,79525] [84099,84100,84101] [86527,86528,86529] [94177,94178,94179] [108899,108900,108901] [121103,121104,121105] [125315,125316,125317] [128017,128018,128019] [129599,129600,129601] [137287,137288,137289] [144399,144400,144401] [144721,144722,144723] [154567,154568,154569] [158403,158404,158405] [166463,166464,166465] [167041,167042,167043]
Raku
<lang perl6>use Prime::Factor;
my @duffinians = lazy (3..*).hyper.grep: { !.is-prime && $_ gcd .&divisors.sum == 1 };
put "First 50 Duffinian numbers:\n" ~ @duffinians[^50].batch(10)».fmt("%3d").join: "\n";
put "\nFirst 40 Duffinian triplets:\n" ~
((^∞).grep: -> $n { (@duffinians[$n] + 1 == @duffinians[$n + 1]) && (@duffinians[$n] + 2 == @duffinians[$n + 2]) })[^40]\ .map( { "({@duffinians[$_ .. $_+2].join: ', '})" } ).batch(4)».fmt("%-24s").join: "\n";</lang>
- Output:
First 50 Duffinian numbers: 4 8 9 16 21 25 27 32 35 36 39 49 50 55 57 63 64 65 75 77 81 85 93 98 100 111 115 119 121 125 128 129 133 143 144 155 161 169 171 175 183 185 187 189 201 203 205 209 215 217 First 40 Duffinian triplets: (63, 64, 65) (323, 324, 325) (511, 512, 513) (721, 722, 723) (899, 900, 901) (1443, 1444, 1445) (2303, 2304, 2305) (2449, 2450, 2451) (3599, 3600, 3601) (3871, 3872, 3873) (5183, 5184, 5185) (5617, 5618, 5619) (6049, 6050, 6051) (6399, 6400, 6401) (8449, 8450, 8451) (10081, 10082, 10083) (10403, 10404, 10405) (11663, 11664, 11665) (12481, 12482, 12483) (13447, 13448, 13449) (13777, 13778, 13779) (15841, 15842, 15843) (17423, 17424, 17425) (19043, 19044, 19045) (26911, 26912, 26913) (30275, 30276, 30277) (36863, 36864, 36865) (42631, 42632, 42633) (46655, 46656, 46657) (47523, 47524, 47525) (53137, 53138, 53139) (58563, 58564, 58565) (72961, 72962, 72963) (76175, 76176, 76177) (79523, 79524, 79525) (84099, 84100, 84101) (86527, 86528, 86529) (94177, 94178, 94179) (108899, 108900, 108901) (121103, 121104, 121105)
Wren
<lang ecmascript>import "./math" for Int, Nums import "./seq" for Lst import "./fmt" for Fmt
var limit = 200000 // say var d = Int.primeSieve(limit-1, false) d[1] = false for (i in 2...limit) {
if (!d[i]) continue if (i % 2 == 0 && !Int.isSquare(i) && !Int.isSquare(i/2)) { d[i] = false continue } var sigmaSum = Nums.sum(Int.divisors(i)) if (Int.gcd(sigmaSum, i) != 1) d[i] = false
}
var duff = (1...d.count).where { |i| d[i] }.toList System.print("First 50 Duffinian numbers:") for (chunk in Lst.chunks(duff[0..49], 10)) Fmt.print("$3d", chunk)
var triplets = [] for (i in 2...limit) {
if (d[i] && d[i-1] && d[i-2]) triplets.add([i-2, i-1, i])
} System.print("\nFirst 50 Duffinian triplets:") for (chunk in Lst.chunks(triplets[0..49], 4)) Fmt.print("$-25s", chunk)</lang>
- Output:
First 50 Duffinian numbers: 4 8 9 16 21 25 27 32 35 36 39 49 50 55 57 63 64 65 75 77 81 85 93 98 100 111 115 119 121 125 128 129 133 143 144 155 161 169 171 175 183 185 187 189 201 203 205 209 215 217 First 50 Duffinian triplets: [63, 64, 65] [323, 324, 325] [511, 512, 513] [721, 722, 723] [899, 900, 901] [1443, 1444, 1445] [2303, 2304, 2305] [2449, 2450, 2451] [3599, 3600, 3601] [3871, 3872, 3873] [5183, 5184, 5185] [5617, 5618, 5619] [6049, 6050, 6051] [6399, 6400, 6401] [8449, 8450, 8451] [10081, 10082, 10083] [10403, 10404, 10405] [11663, 11664, 11665] [12481, 12482, 12483] [13447, 13448, 13449] [13777, 13778, 13779] [15841, 15842, 15843] [17423, 17424, 17425] [19043, 19044, 19045] [26911, 26912, 26913] [30275, 30276, 30277] [36863, 36864, 36865] [42631, 42632, 42633] [46655, 46656, 46657] [47523, 47524, 47525] [53137, 53138, 53139] [58563, 58564, 58565] [72961, 72962, 72963] [76175, 76176, 76177] [79523, 79524, 79525] [84099, 84100, 84101] [86527, 86528, 86529] [94177, 94178, 94179] [108899, 108900, 108901] [121103, 121104, 121105] [125315, 125316, 125317] [128017, 128018, 128019] [129599, 129600, 129601] [137287, 137288, 137289] [144399, 144400, 144401] [144721, 144722, 144723] [154567, 154568, 154569] [158403, 158404, 158405] [166463, 166464, 166465] [167041, 167042, 167043]
XPL0
<lang XPL0>func IsPrime(N); \Return 'true' if N is prime int N, I; [if N <= 2 then return N = 2; if (N&1) = 0 then \even >2\ return false; for I:= 3 to sqrt(N) do
[if rem(N/I) = 0 then return false; I:= I+1; ];
return true; ];
func SumDiv(Num); \Return sum of proper divisors of Num int Num, Div, Sum, Quot; [Div:= 2; Sum:= 0; loop [Quot:= Num/Div;
if Div > Quot then quit; if rem(0) = 0 then [Sum:= Sum + Div; if Div # Quot then Sum:= Sum + Quot; ]; Div:= Div+1; ];
return Sum+1; ];
func GCD(A, B); \Return greatest common divisor of A and B int A, B; [while A#B do
if A>B then A:= A-B else B:= B-A;
return A; ];
func Duff(N); \Return 'true' if N is a Duffinian number int N; [if IsPrime(N) then return false; return GCD(SumDiv(N), N) = 1; ];
int C, N; [Format(4, 0); C:= 0; N:= 4; loop [if Duff(N) then
[RlOut(0, float(N)); C:= C+1; if C >= 50 then quit; if rem(C/20) = 0 then CrLf(0); ]; N:= N+1; ];
CrLf(0); CrLf(0); Format(5, 0); C:= 0; N:= 4; loop [if Duff(N) & Duff(N+1) & Duff(N+2) then
[RlOut(0, float(N)); RlOut(0, float(N+1)); RlOut(0, float(N+2)); CrLf(0); C:= C+1; if C >= 15 then quit; ]; N:= N+1; ];
]</lang>
- Output:
4 8 9 16 21 25 27 32 35 36 39 49 50 55 57 63 64 65 75 77 81 85 93 98 100 111 115 119 121 125 128 129 133 143 144 155 161 169 171 175 183 185 187 189 201 203 205 209 215 217 63 64 65 323 324 325 511 512 513 721 722 723 899 900 901 1443 1444 1445 2303 2304 2305 2449 2450 2451 3599 3600 3601 3871 3872 3873 5183 5184 5185 5617 5618 5619 6049 6050 6051 6399 6400 6401 8449 8450 8451