Numbers which are the cube roots of the product of their proper divisors
Numbers which are the cube roots of the product of their proper divisors is a draft programming task. It is not yet considered ready to be promoted as a complete task, for reasons that should be found in its talk page.
- Example
Consider the number 24. Its proper divisors are: 1, 2, 3, 4, 6, 8 and 12. Their product is 13,824 and the cube root of this is 24. So 24 satisfies the definition in the task title.
- Task
Compute and show here the first 50 positive integers which are the cube roots of the product of their proper divisors.
Also show the 500th and 5,000th such numbers.
- Stretch
Compute and show the 50,000th such number.
- Reference
- Note
OEIS considers 1 to be the first number in this sequence even though, strictly speaking, it has no proper divisors. Please therefore do likewise.
ALGOL 68
As with the second Wren sample, uses the observation on the OEIS page to reduce the problem to finding numbers that are 1 or have 8 divisors (or 7 proper divisors).
BEGIN # find some numbers which are the cube roots of the product of their #
# proper divisors #
# the Online Encyclopedia of Integer Sequences states that these #
# numbers are 1 and those with eight divisors #
# NB: numbers with 8 divisors have 7 proper divisors #
INT max number = 500 000; # maximum number we will consider #
# form a table of proper divisor counts - assume the pdc of 1 is 7 #
[ 1 : max number ]INT pdc; FOR i TO UPB pdc DO pdc[ i ] := 1 OD;
pdc[ 1 ] := 7;
FOR i FROM 2 TO UPB pdc DO
FOR j FROM i + i BY i TO UPB pdc DO pdc[ j ] +:= 1 OD
OD;
# show the numbers which are the cube root of their proper divisor #
# product - equivalent to finding the numbers with a proper divisor #
# count of 7 ( we have "cheated" and set the pdc of 1 to 7 ) #
INT next show := 500;
INT c count := 0;
FOR n TO UPB pdc DO
IF pdc[ n ] = 7 THEN
# found a suitable number #
IF ( c count +:= 1 ) <= 50 THEN
print( ( " ", whole( n, -3 ) ) );
IF c count MOD 10 = 0 THEN print( ( newline ) ) FI
ELIF c count = next show THEN
print( ( whole( c count, -9 ), "th: ", whole( n, 0 ), newline ) );
next show *:= 10
FI
FI
OD
END- Output:
1 24 30 40 42 54 56 66 70 78
88 102 104 105 110 114 128 130 135 136
138 152 154 165 170 174 182 184 186 189
190 195 222 230 231 232 238 246 248 250
255 258 266 273 282 285 286 290 296 297
500th: 2526
5000th: 23118
50000th: 223735
BASIC
BASIC256
arraybase 1
limite = 500000
dim pdc(limite) fill 1
pdc[1] = 7
for i = 2 to pdc[?]
for j = i + i to pdc[?] step i
pdc[j] += 1
next j
next i
n5 = 500
cont = 0
print "First 50 numbers which are the cube roots"
print "of the products of their proper divisors:"
for i = 1 to pdc[?]
if pdc[i] = 7 then
cont += 1
if cont <= 50 then
print rjust(string(i),5);
if cont mod 10 = 0 then print
else
if cont = n5 then
print
print rjust(string(cont),9); "th: "; i;
n5 *= 10
end if
end if
end if
next i- Output:
Same as FreeBASIC entry.
True BASIC
LET limite = 500000
DIM pdc(1 to 500000)
FOR i = 1 to ubound(pdc)
LET pdc(i) = 1
NEXT i
LET pdc(1) = 7
FOR i = 2 to ubound(pdc)
FOR j = i+i to ubound(pdc) step i
LET pdc(j) = pdc(j)+1
NEXT j
NEXT i
LET n5 = 500
LET count = 0
PRINT "First 50 numbers which are the cube roots"
PRINT "of the products of their proper divisors:"
FOR i = 1 to ubound(pdc)
IF pdc(i) = 7 then
LET count = count + 1
IF count <= 50 THEN
PRINT using "####": i;
IF remainder(count, 10) = 0 THEN PRINT
ELSE
IF count = n5 THEN
PRINT
PRINT USING "#########th:": count;
PRINT i;
LET n5 = n5*10
END IF
END IF
END IF
NEXT i
END
- Output:
Same as FreeBASIC entry.
XBasic
PROGRAM "progname"
VERSION "0.0000"
DECLARE FUNCTION Entry ()
FUNCTION Entry ()
limite = 500000
DIM pdc[limite] '(1 TO limite)
FOR i = 1 TO UBOUND(pdc[])
pdc[i] = 1
NEXT i
pdc[1] = 7
FOR i = 2 TO UBOUND(pdc[])
FOR j = i + i TO UBOUND(pdc[]) STEP i
INC pdc[j]
NEXT j
NEXT i
n5 = 500
cont = 0
PRINT "First 50 numbers which are the cube roots"
PRINT "of the products of their proper divisors:"
FOR i = 1 TO UBOUND(pdc[])
IF pdc[i] = 7 THEN
INC cont
IF cont <= 50 THEN
PRINT RJUST$ (STRING$(i), 4);
IF cont MOD 10 = 0 THEN PRINT
ELSE
IF cont = n5 THEN
PRINT "\n"; FORMAT$("#########", cont); "th: "; i;
n5 = n5 * 10
END IF
END IF
END IF
NEXT i
END FUNCTION
END PROGRAM- Output:
Same as FreeBASIC entry.
Yabasic
limite = 500000
dim pdc(limite)
for i = 1 to arraysize(pdc(), 1)
pdc(i) = 1
next i
pdc(1) = 7
for i = 2 to arraysize(pdc(), 1)
for j = i + i to arraysize(pdc(), 1) step i
pdc(j) = pdc(j) + 1
next j
next i
n5 = 500
cont = 0
print "First 50 numbers which are the cube roots"
print "of the products of their proper divisors:"
for i = 1 to arraysize(pdc(), 1)
if pdc(i) = 7 then
cont = cont + 1
if cont <= 50 then
print i using("###");
if mod(cont, 10) = 0 print
else
if cont = n5 then
print "\n", cont using("#########"), "th: ", i;
n5 = n5 * 10
end if
end if
end if
next i- Output:
Same as FreeBASIC entry.
C++
#include <iomanip>
#include <iostream>
unsigned int divisor_count(unsigned int n) {
unsigned int total = 1;
for (; (n & 1) == 0; n >>= 1)
++total;
for (unsigned int p = 3; p * p <= n; p += 2) {
unsigned int count = 1;
for (; n % p == 0; n /= p)
++count;
total *= count;
}
if (n > 1)
total *= 2;
return total;
}
int main() {
std::cout.imbue(std::locale(""));
std::cout << "First 50 numbers which are the cube roots of the products of "
"their proper divisors:\n";
for (unsigned int n = 1, count = 0; count < 50000; ++n) {
if (n == 1 || divisor_count(n) == 8) {
++count;
if (count <= 50)
std::cout << std::setw(3) << n
<< (count % 10 == 0 ? '\n' : ' ');
else if (count == 500 || count == 5000 || count == 50000)
std::cout << std::setw(6) << count << "th: " << n << '\n';
}
}
}
- Output:
First 50 numbers which are the cube roots of the products of their proper divisors: 1 24 30 40 42 54 56 66 70 78 88 102 104 105 110 114 128 130 135 136 138 152 154 165 170 174 182 184 186 189 190 195 222 230 231 232 238 246 248 250 255 258 266 273 282 285 286 290 296 297 500th: 2,526 5,000th: 23,118 50,000th: 223,735
FreeBASIC
Dim As Single limite = 500000
Dim As Integer pdc(1 To limite)
Dim As Integer i, j
For i = 1 To Ubound(pdc)
pdc(i) = 1
Next i
pdc(1) = 7
For i = 2 To Ubound(pdc)
For j = i + i To Ubound(pdc) Step i
pdc(j) += 1
Next j
Next i
Dim As Integer n5 = 500, cont = 0
Print "First 50 numbers which are the cube roots"
Print "of the products of their proper divisors:"
For i = 1 To Ubound(pdc)
If pdc(i) = 7 Then
cont += 1
If cont <= 50 Then
Print Using "####"; i;
If cont Mod 10 = 0 Then Print
Elseif cont = n5 Then
Print Using !"\n#########th: &"; cont; i;
n5 *= 10
End If
End If
Next i
Sleep- Output:
First 50 numbers which are the cube roots
of the products of their proper divisors:
1 24 30 40 42 54 56 66 70 78
88 102 104 105 110 114 128 130 135 136
138 152 154 165 170 174 182 184 186 189
190 195 222 230 231 232 238 246 248 250
255 258 266 273 282 285 286 290 296 297
500th: 2526
5000th: 23118
50000th: 223735
J
Note that the cube root of the product of the proper divisors is the fourth root of the product of all divisors of a positive integer. That said, we do not need to find roots here -- we only need to inspect the powers of the prime factors of the number:
F=: 1 8 e.~_ */@:>:@q:"0 ]
Task examples:
N=: 1+I.F 1+i.2^18
5 10$N
1 24 30 40 42 54 56 66 70 78
88 102 104 105 110 114 128 130 135 136
138 152 154 165 170 174 182 184 186 189
190 195 222 230 231 232 238 246 248 250
255 258 266 273 282 285 286 290 296 297
499{N
2526
4999{N
23118
49999{N
223735
Phix
with javascript_semantics sequence n50 = {} integer count = 0, n = 1, n5 = 500 atom t0 = time() printf(1,"First 50 numbers which are the cube roots\n") printf(1," of the products of their proper divisors:\n") while count<500000 do -- if product(factors(n))=n*n*n then if n=1 or length(factors(n))=6 then -- safer/smidge faster count += 1 if count<=50 then n50 &= n if count=50 then printf(1,"%s\n",join_by(n50,1,10,"",fmt:="%4d")) end if elsif count=n5 then printf(1,"%,8dth: %,d (%s)\n",{n5,n,elapsed(time()-t0)}) n5 *= 10 end if end if n += 1 end while
By default factors() does not include 1 and n, or I could use length(factors(n,1))=8, both 25% faster than using product(), which exceeds precision limits on 32-bit for n=180, and on 64bit for n=240, though since you'll get exactly the same precision error on the n*n*n it kinda "worked by chance".
- Output:
First 50 numbers which are the cube roots
of the products of their proper divisors:
1 24 30 40 42 54 56 66 70 78
88 102 104 105 110 114 128 130 135 136
138 152 154 165 170 174 182 184 186 189
190 195 222 230 231 232 238 246 248 250
255 258 266 273 282 285 286 290 296 297
500th: 2,526 (0.0s)
5,000th: 23,118 (0.0s)
50,000th: 223,735 (0.6s)
500,000th: 2,229,229 (14.1s)
Python
''' Rosetta code rosettacode.org/wiki/Numbers_which_are_the_cube_roots_of_the_product_of_their_proper_divisors '''
from functools import reduce
from sympy import divisors
FOUND = 0
for num in range(1, 1_000_000):
divprod = reduce(lambda x, y: x * y, divisors(num)[:-1])if num > 1 else 1
if num * num * num == divprod:
FOUND += 1
if FOUND <= 50:
print(f'{num:5}', end='\n' if FOUND % 10 == 0 else '')
if FOUND == 500:
print(f'\nFive hundreth: {num:,}')
if FOUND == 5000:
print(f'\nFive thousandth: {num:,}')
if FOUND == 50000:
print(f'\nFifty thousandth: {num:,}')
break
- Output:
1 24 30 40 42 54 56 66 70 78 88 102 104 105 110 114 128 130 135 136 138 152 154 165 170 174 182 184 186 189 190 195 222 230 231 232 238 246 248 250 255 258 266 273 282 285 286 290 296 297 Five hundreth: 2,526 Five thousandth: 23,118 Fifty thousandth: 223,735
OEIS algorithm (see talk pages)
from sympy import divisors
numfound = 0
for num in range(1, 1_000_000):
if num == 1 or len(divisors(num)) == 8:
numfound += 1
if numfound <= 50:
print(f'{num:5}', end='\n' if numfound % 10 == 0 else '')
if numfound == 500:
print(f'\nFive hundreth: {num:,}')
if numfound == 5000:
print(f'\nFive thousandth: {num:,}')
if numfound == 50000:
print(f'\nFifty thousandth: {num:,}')
break
Output same as first algorithm.
Raku
use Prime::Factor;
use Lingua::EN::Numbers;
my @cube-div = lazy 1, |(2..∞).hyper.grep: { .³ == [×] .&proper-divisors }
put "First 50 numbers which are the cube roots of the products of their proper divisors:\n" ~
@cube-div[^50]».fmt("%3d").batch(10).join: "\n";
printf "\n%16s: %s\n", .Int.&ordinal.tc, comma @cube-div[$_ - 1] for 5e2, 5e3, 5e4;
- Output:
First 50 numbers which are the cube roots of the products of their proper divisors: 1 24 30 40 42 54 56 66 70 78 88 102 104 105 110 114 128 130 135 136 138 152 154 165 170 174 182 184 186 189 190 195 222 230 231 232 238 246 248 250 255 258 266 273 282 285 286 290 296 297 Five hundredth: 2,526 Five thousandth: 23,118 Fifty thousandth: 223,735
Wren
import "./math" for Int, Nums
import "./long" for ULong, ULongs
import "./fmt" for Fmt
var numbers50 = []
var count = 0
var n = 1
var ln
var maxSafe = Num.maxSafeInteger.cbrt.floor
System.print("First 50 numbers which are the cube roots of the products of their proper divisors:")
while (true) {
var pd = Int.properDivisors(n)
if ((n <= maxSafe && Nums.prod(pd) == n * n * n) ||
(ULongs.prod(pd.map { |f| ULong.new(f) }) == (ln = ULong.new(n)) * ln * ln )) {
count = count + 1
if (count <= 50) {
numbers50.add(n)
if (count == 50) Fmt.tprint("$3d", numbers50, 10)
} else if (count == 500) {
Fmt.print("\n500th : $,d", n)
} else if (count == 5000) {
Fmt.print("5,000th : $,d", n)
} else if (count == 50000) {
Fmt.print("50,000th: $,d", n)
break
}
}
n = n + 1
}- Output:
First 50 numbers which are the cube roots of the products of their proper divisors: 1 24 30 40 42 54 56 66 70 78 88 102 104 105 110 114 128 130 135 136 138 152 154 165 170 174 182 184 186 189 190 195 222 230 231 232 238 246 248 250 255 258 266 273 282 285 286 290 296 297 500th : 2,526 5,000th : 23,118 50,000th: 223,735
Alternatively and a bit quicker, inspired by the C++ entry and the OEIS comment that (apart from 1) n must have exactly 8 divisors:
import "./fmt" for Fmt
var divisorCount = Fn.new { |n|
var i = 1
var k = (n%2 == 0) ? 1 : 2
var count = 0
while (i <= n.sqrt) {
if (n%i == 0) {
count = count + 1
var j = (n/i).floor
if (j != i) count = count + 1
}
i = i + k
}
return count
}
var numbers50 = []
var count = 0
var n = 1
System.print("First 50 numbers which are the cube roots of the products of their proper divisors:")
while (true) {
var dc = divisorCount.call(n)
if (n == 1 || dc == 8) {
count = count + 1
if (count <= 50) {
numbers50.add(n)
if (count == 50) Fmt.tprint("$3d", numbers50, 10)
} else if (count == 500) {
Fmt.print("\n500th : $,d", n)
} else if (count == 5000) {
Fmt.print("5,000th : $,d", n)
} else if (count == 50000) {
Fmt.print("50,000th: $,d", n)
break
}
}
n = n + 1
}- Output:
Same as first version.