Numbers whose count of divisors is prime
- Task
Find positive integers n which count of divisors is prime, but not equal to 2, where n < 1,000.
Stretch goal: (as above), but where n < 100,000.
11l
F is_prime(a)
I a == 2
R 1B
I a < 2 | a % 2 == 0
R 0B
L(i) (3 .. Int(sqrt(a))).step(2)
I a % i == 0
R 0B
R 1B
print(‘Numbers which count of divisors is prime are:’)
V row = 0
L(n) 1..99999
V num = 0
L(m) 1 .. n
I n % m == 0
num++
I is_prime(num) & num != 2
print(‘#6’.format(n), end' ‘ ’)
row++
I row % 5 == 0
print()
print("\n\nFound "row‘ numbers’)
- Output:
Numbers which count of divisors is prime are: 4 9 16 25 49 64 81 121 169 289 361 529 625 729 841 961 1024 1369 1681 1849 2209 2401 2809 3481 3721 4096 4489 5041 5329 6241 6889 7921 9409 10201 10609 11449 11881 12769 14641 15625 16129 17161 18769 19321 22201 22801 24649 26569 27889 28561 29929 32041 32761 36481 37249 38809 39601 44521 49729 51529 52441 54289 57121 58081 59049 63001 65536 66049 69169 72361 73441 76729 78961 80089 83521 85849 94249 96721 97969 Found 79 numbers
Action!
INCLUDE "H6:SIEVE.ACT"
INT FUNC CountDivisors(INT x)
INT i,max,count
count=2 i=2 max=x
WHILE i<max
DO
IF x MOD i=0 THEN
max=x/i
IF i<max THEN
count==+2
ELSEIF i=max THEN
count==+1
FI
FI
i==+1
OD
RETURN (count)
PROC Main()
DEFINE MAXNUM="999"
BYTE ARRAY primes(MAXNUM+1)
INT i,count
Put(125) PutE() ;clear the screen
Sieve(primes,MAXNUM+1)
FOR i=2 TO MAXNUM
DO
IF primes(i)=0 THEN
count=CountDivisors(i)
IF count>2 AND primes(count)=1 THEN
PrintF("%I has %I divisors%E",i,count)
FI
FI
OD
RETURN
- Output:
Screenshot from Atari 8-bit computer
4 has 3 divisors 9 has 3 divisors 16 has 5 divisors 25 has 3 divisors 49 has 3 divisors 64 has 7 divisors 81 has 5 divisors 121 has 3 divisors 169 has 3 divisors 289 has 3 divisors 361 has 3 divisors 529 has 3 divisors 625 has 5 divisors 729 has 7 divisors 841 has 3 divisors 961 has 3 divisors
ALGOL 60
begin
comment - return n mod m;
integer procedure mod(n, m);
value n, m; integer n, m;
begin
mod := n - m * entier(n / m);
end;
comment - return count of divisors of n;
integer procedure divcount(n);
value n; integer n;
begin
integer count, f1, f2;
comment - every integer has at least 2 divisors (1 and n);
count := 2;
f1 := 2;
for f1 := f1 while (f1 * f1) <= n do
begin
if mod(n, f1) = 0 then
begin
f2 := n / f1;
if f2 > f1 then
count := count + 2
else
count := count + 1;
end;
f1 := f1 + 1;
end;
divcount := count;
end;
comment - return true if n is prime;
boolean procedure isprime(n);
value n; integer n;
begin
if n < 2 then
isprime := false
else if mod(n, 2) = 0 then
isprime := (n = 2)
else
begin
comment - check odd divisors up to sqrt(n);
integer i, limit;
boolean divisible;
i := 3;
limit := entier(sqrt(n));
divisible := false;
for i := i while i <= limit and not divisible do
begin
if mod(n, i) = 0 then
divisible := true;
i := i + 2
end;
isprime := not divisible;
end;
end;
integer i, sq, limit, found;
limit := 100000;
outstring(1,"Integers up to");
outinteger(1, limit);
outstring(1, "with a count of divisors that is prime:\n");
found := 0;
comment - we only need to check perfect squares;
i := 2;
for sq := i * i while sq < limit do
begin
if isprime(divcount(sq)) then
begin
outinteger(1, sq);
found := found + 1;
if mod(found, 10) = 0 then outstring(1,"\n");
end;
i := i + 1;
end;
outstring(1,"\n");
outinteger(1,found);
outstring(1,"were found\n");
end
- Output:
Integers up to 100000 with a count of divisors that is prime: 4 9 16 25 49 64 81 121 169 289 361 529 625 729 841 961 1024 1369 1681 1849 2209 2401 2809 3481 3721 4096 4489 5041 5329 6241 6889 7921 9409 10201 10609 11449 11881 12769 14641 15625 16129 17161 18769 19321 22201 22801 24649 26569 27889 28561 29929 32041 32761 36481 37249 38809 39601 44521 49729 51529 52441 54289 57121 58081 59049 63001 65536 66049 69169 72361 73441 76729 78961 80089 83521 85849 94249 96721 97969 79 were found
ALGOL 68
Counts the divisors without using division.
BEGIN # find numbers with prime divisor counts #
INT max number := 1 000;
TO 2 DO
# construct a table of the divisor counts #
[ 1 : max number ]INT ndc; FOR i FROM 1 TO UPB ndc DO ndc[ i ] := 1 OD;
FOR i FROM 2 TO UPB ndc DO
FOR j FROM i BY i TO UPB ndc DO ndc[ j ] +:= 1 OD
OD;
# show the numbers with prime divisor counts #
print( ( "Numbers up to ", whole( max number, 0 ), " with odd prime divisor counts:", newline ) );
INT p count := 0;
FOR i TO UPB ndc DO
INT divisor count = ndc[ i ];
IF ODD divisor count AND ndc[ divisor count ] = 2 THEN
print( ( whole( i, -8 ) ) );
IF ( p count +:= 1 ) MOD 10 = 0 THEN print( ( newline ) ) FI
FI
OD;
print( ( newline ) );
max number := 100 000
OD
END
- Output:
Numbers up to 1000 with odd prime divisor counts: 4 9 16 25 49 64 81 121 169 289 361 529 625 729 841 961 Numbers up to 100000 with odd prime divisor counts: 4 9 16 25 49 64 81 121 169 289 361 529 625 729 841 961 1024 1369 1681 1849 2209 2401 2809 3481 3721 4096 4489 5041 5329 6241 6889 7921 9409 10201 10609 11449 11881 12769 14641 15625 16129 17161 18769 19321 22201 22801 24649 26569 27889 28561 29929 32041 32761 36481 37249 38809 39601 44521 49729 51529 52441 54289 57121 58081 59049 63001 65536 66049 69169 72361 73441 76729 78961 80089 83521 85849 94249 96721 97969
AppleScript
Only squares checked, as per the Discussion page.
on countFactors(n) -- Positive ns only.
if (n < 4) then return 2 - ((n = 1) as integer)
set factorCount to 2
set sqrt to n ^ 0.5
set limit to sqrt div 1
if (limit = sqrt) then
set factorCount to 3
set limit to limit - 1
end if
repeat with i from 2 to limit
if (n mod i = 0) then set factorCount to factorCount + 2
end repeat
return factorCount
end countFactors
on join(lst, delim)
set astid to AppleScript's text item delimiters
set AppleScript's text item delimiters to delim
set txt to lst as text
set AppleScript's text item delimiters to astid
return txt
end join
on task()
set limit to 100000 - 1
set nWidth to (count (limit as text))
set inset to " "'s text 1 thru (nWidth - 1)
set output to {"Positive integers < " & (limit + 1) & " whose divisor count is an odd prime:"}
set row to {}
set counter to 0
repeat with sqrt from 2 to (limit ^ 0.5 div 1)
set n to sqrt * sqrt
if (countFactors(countFactors(n)) = 2) then
set counter to counter + 1
set row's end to (inset & n)'s text -nWidth thru -1
if ((count row) = 10) then
set output's end to join(row, " ")
set row to {}
end if
end if
end repeat
if (row ≠ {}) then set output's end to join(row, " ")
set output's end to linefeed & counter & " such integers"
return join(output, linefeed)
end task
task()
- Output:
"Positive integers < 100000 whose divisor count is an odd prime:
4 9 16 25 49 64 81 121 169 289
361 529 625 729 841 961 1024 1369 1681 1849
2209 2401 2809 3481 3721 4096 4489 5041 5329 6241
6889 7921 9409 10201 10609 11449 11881 12769 14641 15625
16129 17161 18769 19321 22201 22801 24649 26569 27889 28561
29929 32041 32761 36481 37249 38809 39601 44521 49729 51529
52441 54289 57121 58081 59049 63001 65536 66049 69169 72361
73441 76729 78961 80089 83521 85849 94249 96721 97969
79 such integers"
Arturo
numsWithPrimeNofDivisors: select 1..100000 'x [
nofDivisors: size factors x
and? [prime? nofDivisors]
[nofDivisors <> 2]
]
loop split.every: 5 numsWithPrimeNofDivisors 'x ->
print map x 's -> pad to :string s 6
- Output:
4 9 16 25 49 64 81 121 169 289 361 529 625 729 841 961 1024 1369 1681 1849 2209 2401 2809 3481 3721 4096 4489 5041 5329 6241 6889 7921 9409 10201 10609 11449 11881 12769 14641 15625 16129 17161 18769 19321 22201 22801 24649 26569 27889 28561 29929 32041 32761 36481 37249 38809 39601 44521 49729 51529 52441 54289 57121 58081 59049 63001 65536 66049 69169 72361 73441 76729 78961 80089 83521 85849 94249 96721 97969
AWK
# syntax: GAWK -f NUMBERS_WHOSE_COUNT_OF_DIVISORS_IS_PRIME.AWK
BEGIN {
start = 2
stop = 99999
stop2 = 999
for (i=start; i*i<=stop; i++) {
n = count_divisors(i*i)
if (n>2 && is_prime(n)) {
printf("%6d%1s",i*i,++count%10?"":"\n")
if (i*i <= stop2) {
count2++
}
}
}
printf("\nNumbers with odd prime divisor counts %d-%d: %d\n",start,stop2,count2)
printf("Numbers with odd prime divisor counts %d-%d: %d\n",start,stop,count)
exit(0)
}
function count_divisors(n, count,i) {
for (i=1; i*i<=n; i++) {
if (n % i == 0) {
count += (i == n / i) ? 1 : 2
}
}
return(count)
}
function is_prime(x, i) {
if (x <= 1) {
return(0)
}
for (i=2; i<=int(sqrt(x)); i++) {
if (x % i == 0) {
return(0)
}
}
return(1)
}
- Output:
4 9 16 25 49 64 81 121 169 289 361 529 625 729 841 961 1024 1369 1681 1849 2209 2401 2809 3481 3721 4096 4489 5041 5329 6241 6889 7921 9409 10201 10609 11449 11881 12769 14641 15625 16129 17161 18769 19321 22201 22801 24649 26569 27889 28561 29929 32041 32761 36481 37249 38809 39601 44521 49729 51529 52441 54289 57121 58081 59049 63001 65536 66049 69169 72361 73441 76729 78961 80089 83521 85849 94249 96721 97969 Numbers with odd prime divisor counts 2-999: 16 Numbers with odd prime divisor counts 2-99999: 79
BASIC
BASIC256
function isPrime(v)
if v < 2 then return False
if (v mod 2) = 0 then return v = 2
if (v mod 3) = 0 then return v = 3
d = 5
while d * d <= v
if (v mod d) = 0 then return False else d = d + 2
end while
return True
end function
row = 0
print "Numbers which count of divisors is prime are:"
for n = 1 to 1000
num = 0
for m = 1 to n
if (n mod m) = 0 then num = num + 1
next m
if isPrime(num) and num <> 2 then
print ""; n; " ";
row = row + 1
end if
next n
print
print "Found "; row; " numbers"
end
FreeBASIC
Function isPrime(Byval ValorEval As Integer) As Boolean
If ValorEval < 2 Then Return False
If ValorEval Mod 2 = 0 Then Return ValorEval = 2
If ValorEval Mod 3 = 0 Then Return ValorEval = 3
Dim d As Integer = 5
While d * d <= ValorEval
If ValorEval Mod d = 0 Then Return False Else d += 2
Wend
Return True
End Function
Dim As Integer row = 0
Dim As Uinteger n, num, m
Print "Numbers which count of divisors is prime are:"
For n = 1 To 100000
num = 0
For m = 1 To n
If n Mod m = 0 Then num += 1 : End If
Next m
If isPrime(num) And num <> 2 Then
Print Using " ##### "; n;
row += 1
If row Mod 5 = 0 Then Print : End If
End If
Next n
Print !"\n\nFound"; row; " numbers"
Sleep
- Output:
Numbers which count of divisors is prime are: 4 9 16 25 49 64 81 121 169 289 361 529 625 729 841 961 1024 1369 1681 1849 2209 2401 2809 3481 3721 4096 4489 5041 5329 6241 6889 7921 9409 10201 10609 11449 11881 12769 14641 15625 16129 17161 18769 19321 22201 22801 24649 26569 27889 28561 29929 32041 32761 36481 37249 38809 39601 44521 49729 51529 52441 54289 57121 58081 59049 63001 65536 66049 69169 72361 73441 76729 78961 80089 83521 85849 94249 96721 97969 Found 79 numbers
QBasic
row = 0
PRINT "Numbers which count of divisors is prime are:"
FOR n = 1 TO 1000
num = 0
FOR m = 1 TO n
IF n MOD m = 0 THEN num = num + 1
NEXT m
IF isPrime(num) AND num <> 2 THEN
PRINT USING " ##### "; n;
row = row + 1
IF row MOD 5 = 0 THEN PRINT
END IF
NEXT n
PRINT : PRINT "Found"; row; " numbers"
Sleep
END
FUNCTION isPrime (ValorEval)
IF ValorEval < 2 THEN isPrime = False
IF ValorEval MOD 2 = 0 THEN isPrime = 2
IF ValorEval MOD 3 = 0 THEN isPrime = 3
d = 5
WHILE d * d <= ValorEval
IF ValorEval MOD d = 0 THEN isPrime = False ELSE d = d + 2
WEND
isPrime = True
END FUNCTION
PureBasic
Procedure isPrime(v.i)
If v <= 1 : ProcedureReturn #False
ElseIf v < 4 : ProcedureReturn #True
ElseIf v % 2 = 0 : ProcedureReturn #False
ElseIf v < 9 : ProcedureReturn #True
ElseIf v % 3 = 0 : ProcedureReturn #False
Else
Protected r = Round(Sqr(v), #PB_Round_Down)
Protected f = 5
While f <= r
If v % f = 0 Or v % (f + 2) = 0
ProcedureReturn #False
EndIf
f + 6
Wend
EndIf
ProcedureReturn #True
EndProcedure
OpenConsole()
fila.i = 0
PrintN("Numbers which count of divisors is prime are:")
For n.i = 1 To 100000
num.i = 0
For m.i = 1 To n
If Mod(n, m) = 0
num + 1
EndIf
Next
If isPrime(num) And num <> 2
Print(" " + Str(n) + " ")
fila + 1
If Mod(fila, 5) = 0
PrintN("")
EndIf
EndIf
Next
PrintN(#CRLF$ + "Found " + Str(fila) + " numbers")
Input()
CloseConsole()
End
Yabasic
row = 0
print "Numbers which count of divisors is prime are:"
for n = 1 to 1000
num = 0
for m = 1 to n
if mod(n, m) = 0 then num = num + 1 : fi
next m
if isPrime(num) and num <> 2 then
print n using "#####",
row = row + 1
if mod(row, 5) = 0 then print : fi
end if
next n
print "\n\nFound ", row, " numbers"
end
sub isPrime(v)
if v < 2 then return False : fi
if mod(v, 2) = 0 then return v = 2 : fi
if mod(v, 3) = 0 then return v = 3 : fi
d = 5
while d * d <= v
if mod(v, d) = 0 then return False else d = d + 2 : fi
end while
return True
end sub
- Output:
Igual que la entrada de FreeBASIC.
C++
#include <cmath>
#include <cstdlib>
#include <iomanip>
#include <iostream>
int divisor_count(int n) {
int total = 1;
for (; (n & 1) == 0; n >>= 1)
++total;
for (int p = 3; p * p <= n; p += 2) {
int count = 1;
for (; n % p == 0; n /= p)
++count;
total *= count;
}
if (n > 1)
total *= 2;
return total;
}
bool is_prime(int n) {
if (n < 2)
return false;
if (n % 2 == 0)
return n == 2;
if (n % 3 == 0)
return n == 3;
for (int p = 5; p * p <= n; p += 4) {
if (n % p == 0)
return false;
p += 2;
if (n % p == 0)
return false;
}
return true;
}
int main(int argc, char** argv) {
int limit = 1000;
switch (argc) {
case 1:
break;
case 2:
limit = std::strtol(argv[1], nullptr, 10);
if (limit <= 0) {
std::cerr << "Invalid limit\n";
return EXIT_FAILURE;
}
break;
default:
std::cerr << "usage: " << argv[0] << " [limit]\n";
return EXIT_FAILURE;
}
int width = static_cast<int>(std::ceil(std::log10(limit)));
int count = 0;
for (int i = 1;; ++i) {
int n = i * i;
if (n >= limit)
break;
int divisors = divisor_count(n);
if (divisors != 2 && is_prime(divisors))
std::cout << std::setw(width) << n << (++count % 10 == 0 ? '\n' : ' ');
}
std::cout << "\nCount: " << count << '\n';
return EXIT_SUCCESS;
}
- Output:
Default input:
4 9 16 25 49 64 81 121 169 289 361 529 625 729 841 961 Count: 16
Stretch goal:
4 9 16 25 49 64 81 121 169 289 361 529 625 729 841 961 1024 1369 1681 1849 2209 2401 2809 3481 3721 4096 4489 5041 5329 6241 6889 7921 9409 10201 10609 11449 11881 12769 14641 15625 16129 17161 18769 19321 22201 22801 24649 26569 27889 28561 29929 32041 32761 36481 37249 38809 39601 44521 49729 51529 52441 54289 57121 58081 59049 63001 65536 66049 69169 72361 73441 76729 78961 80089 83521 85849 94249 96721 97969 Count: 79
CLU
% Find the amount of divisors for 1..N
div_counts = proc (n: int) returns (sequence[int])
divs: array[int] := array[int]$fill(1,n,1)
for d: int in int$from_to(2, n) do
for m: int in int$from_to_by(d, n, d) do
divs[m] := divs[m] + 1
end
end
return(sequence[int]$a2s(divs))
end div_counts
% Find maximum element of sequence
max = proc (seq: sequence[int]) returns (int)
n: int := 0
for i: int in sequence[int]$elements(seq) do
if i>n then n:=i end
end
return(n)
end max
% Sieve primes up to N
sieve = proc (n: int) returns (sequence[bool])
prime: array[bool] := array[bool]$fill(1,n,true)
prime[1] := false
for p: int in int$from_to(2, n/2) do
for c: int in int$from_to_by(p*p, n, p) do
prime[c] := false
end
end
return(sequence[bool]$a2s(prime))
end sieve
start_up = proc ()
MAX_N = 100000
po: stream := stream$primary_output()
div_count: sequence[int] := div_counts(MAX_N)
prime: sequence[bool] := sieve(max(div_count))
count: int := 0
for i: int in int$from_to(1, MAX_N) do
dc: int := div_count[i]
if dc ~= 2 cand prime[dc] then
stream$putright(po, int$unparse(i), 8)
count := count + 1
if count//10 = 0 then stream$putl(po, "") end
end
end
stream$putl(po, "\nFound " || int$unparse(count) || " numbers.")
end start_up
- Output:
4 9 16 25 49 64 81 121 169 289 361 529 625 729 841 961 1024 1369 1681 1849 2209 2401 2809 3481 3721 4096 4489 5041 5329 6241 6889 7921 9409 10201 10609 11449 11881 12769 14641 15625 16129 17161 18769 19321 22201 22801 24649 26569 27889 28561 29929 32041 32761 36481 37249 38809 39601 44521 49729 51529 52441 54289 57121 58081 59049 63001 65536 66049 69169 72361 73441 76729 78961 80089 83521 85849 94249 96721 97969 Found 79 numbers.
Delphi
procedure GetUniqueFactors(N: integer; var IA: TIntegerDynArray);
{Get unique factors of a number}
var I: integer;
begin
SetLength(IA,1);
for I:=2 to N do
if (N mod I)=0 then
begin
SetLength(IA,Length(IA)+1);
IA[High(IA)]:=I;
end;
end;
procedure ShowCountPrimes(Memo: TMemo);
{Count the number of Unique factors that are prime}
var I,C,Cnt: integer;
var IA: TIntegerDynArray;
var S: string;
begin
Cnt:=0; S:='';
for I:=1 to 100000-1 do
begin
GetUniqueFactors(I,IA);
C:=Length(IA);
if (C=2) or (not IsPrime(C)) then continue;
Inc(Cnt);
S:=S+Format('%8D',[I]);
If (Cnt mod 5)=0 then S:=S+CRLF;
end;
Memo.Lines.Add('Count='+IntToStr(Cnt));
Memo.Lines.Add(S);
end;
- Output:
Count=79 4 9 16 25 49 64 81 121 169 289 361 529 625 729 841 961 1024 1369 1681 1849 2209 2401 2809 3481 3721 4096 4489 5041 5329 6241 6889 7921 9409 10201 10609 11449 11881 12769 14641 15625 16129 17161 18769 19321 22201 22801 24649 26569 27889 28561 29929 32041 32761 36481 37249 38809 39601 44521 49729 51529 52441 54289 57121 58081 59049 63001 65536 66049 69169 72361 73441 76729 78961 80089 83521 85849 94249 96721 97969 Elapsed Time: 17.921 Sec.
EasyLang
fastfunc isprim num .
if num mod 2 = 0
return 0
.
i = 3
while i <= sqrt num
if num mod i = 0
return 0
.
i += 2
.
return 1
.
for n to 999
cnt = 0
for m to n
cnt += if n mod m = 0
.
if cnt > 2 and isprim cnt = 1
write n & " "
.
.
- Output:
4 9 16 25 49 64 81 121 169 289 361 529 625 729 841 961
F#
This task uses Extensible Prime Generator (F#)
// Numbers whose divisor count is prime. Nigel Galloway: July 13th., 2021
primes64()|>Seq.takeWhile(fun n->n*n<100000L)|>Seq.collect(fun n->primes32()|>Seq.skip 1|>Seq.map(fun g->pown n (g-1))|>Seq.takeWhile((>)100000L))|>Seq.sort|>Seq.iter(printf "%d "); printfn ""
- Output:
4 9 16 25 49 64 81 121 169 289 361 529 625 729 841 961 1024 1369 1681 1849 2209 2401 2809 3481 3721 4096 4489 5041 5329 6241 6889 7921 9409 10201 10609 11449 11881 12769 14641 15625 16129 17161 18769 19321 22201 22801 24649 26569 27889 28561 29929 32041 32761 36481 37249 38809 39601 44521 49729 51529 52441 54289 57121 58081 59049 63001 65536 66049 69169 72361 73441 76729 78961 80089 83521 85849 94249 96721 97969
Factor
USING: formatting grouping io kernel math math.primes
math.primes.factors math.ranges sequences sequences.extras ;
FROM: math.extras => integer-sqrt ;
: odd-prime? ( n -- ? ) dup 2 = [ drop f ] [ prime? ] if ;
: pdc-upto ( n -- seq )
integer-sqrt [1,b]
[ sq ] [ divisors length odd-prime? ] map-filter ;
100,000 pdc-upto 10 group [ [ "%-8d" printf ] each nl ] each
- Output:
4 9 16 25 49 64 81 121 169 289 361 529 625 729 841 961 1024 1369 1681 1849 2209 2401 2809 3481 3721 4096 4489 5041 5329 6241 6889 7921 9409 10201 10609 11449 11881 12769 14641 15625 16129 17161 18769 19321 22201 22801 24649 26569 27889 28561 29929 32041 32761 36481 37249 38809 39601 44521 49729 51529 52441 54289 57121 58081 59049 63001 65536 66049 69169 72361 73441 76729 78961 80089 83521 85849 94249 96721 97969
Go
package main
import (
"fmt"
"rcu"
)
func countDivisors(n int) int {
count := 0
i := 1
k := 1
if n%2 == 1 {
k = 2
}
for ; i*i <= n; i += k {
if n%i == 0 {
count++
j := n / i
if j != i {
count++
}
}
}
return count
}
func main() {
const limit = 1e5
var results []int
for i := 2; i * i < limit; i++ {
n := countDivisors(i * i)
if n > 2 && rcu.IsPrime(n) {
results = append(results, i * i)
}
}
climit := rcu.Commatize(limit)
fmt.Printf("Positive integers under %7s whose number of divisors is an odd prime:\n", climit)
under1000 := 0
for i, n := range results {
fmt.Printf("%7s", rcu.Commatize(n))
if (i+1)%10 == 0 {
fmt.Println()
}
if n < 1000 {
under1000++
}
}
fmt.Printf("\n\nFound %d such integers (%d under 1,000).\n", len(results), under1000)
}
- Output:
Positive integers under 100,000 whose number of divisors is an odd prime: 4 9 16 25 49 64 81 121 169 289 361 529 625 729 841 961 1,024 1,369 1,681 1,849 2,209 2,401 2,809 3,481 3,721 4,096 4,489 5,041 5,329 6,241 6,889 7,921 9,409 10,201 10,609 11,449 11,881 12,769 14,641 15,625 16,129 17,161 18,769 19,321 22,201 22,801 24,649 26,569 27,889 28,561 29,929 32,041 32,761 36,481 37,249 38,809 39,601 44,521 49,729 51,529 52,441 54,289 57,121 58,081 59,049 63,001 65,536 66,049 69,169 72,361 73,441 76,729 78,961 80,089 83,521 85,849 94,249 96,721 97,969 Found 79 such integers (16 under 1,000).
J
To find the divisors of a number we can use any of the implementations defined at in the Factors of an integer task: foi
, factrs
, factrst
, factorsOfNumber
or factors
. Here, we will use factors
as it's the most efficient:
(#~ (2&< * 1&p:)@#@factors"0) }. i. 100000
4 9 16 25 49 64 81 121 169 289 361 529 625 729 841 961 1024 1369 1681 1849 2209 2401 2809 3481 3721 4096 4489 5041 5329 6241 6889 7921 9409 10201 10609 11449 11881 12769 14641 15625 16129 17161 18769 19321 22201 22801 24649 26569 27889 28561 29929 32041 32761 36481 37249 38809 39601 44521 49729 51529 52441 54289 57121 58081 59049 63001 65536 66049 69169 72361 73441 76729 78961 80089 83521 85849 94249 96721 97969
(Note that we first conditioned J's environment to show longer lines, using 9!:37(10*9!:36)
.)
Java
import java.util.ArrayList;
import java.util.List;
import java.util.TreeSet;
import java.util.stream.Collectors;
public final class NumbersWhoseCountOfDivisorsIsPrime {
public static void main(String[] args) {
List<Integer> primes = listPrimeNumbers((int) Math.sqrt(100_000));
primes.stream().flatMap( p -> primes.stream().skip(1).map( i -> (int) Math.pow(p, i - 1) ))
.filter( i -> i < 100_000)
.collect(Collectors.toCollection(TreeSet::new))
.forEach( i -> System.out.print(i + " ") );
}
private static List<Integer> listPrimeNumbers(int maximum) {
final int halfMaximum = ( maximum + 1 ) / 2;
boolean[] composite = new boolean[halfMaximum];
for ( int i = 1, p = 3; i < halfMaximum; p += 2, i++ ) {
if ( ! composite[i] ) {
for ( int j = i + p; j < halfMaximum; j += p ) {
composite[j] = true;
}
}
}
List<Integer> primes = new ArrayList<Integer>(List.of( 2 ));
for ( int i = 1, p = 3; i < halfMaximum; p += 2, i++ ) {
if ( ! composite[i] ) {
primes.add(p);
}
}
return primes;
}
}
- Output:
4 9 16 25 49 64 81 121 169 289 361 529 625 729 841 961 1024 1369 1681 1849 2209 2401 2809 3481 3721 4096 4489 5041 5329 6241 6889 7921 9409 10201 10609 11449 11881 12769 14641 15625 16129 17161 18769 19321 22201 22801 24649 26569 27889 28561 29929 32041 32761 36481 37249 38809 39601 44521 49729 51529 52441 54289 57121 58081 59049 63001 65536 66049 69169 72361 73441 76729 78961 80089 83521 85849 94249 96721 97969
jq
Works with gojq, the Go implementation of jq
For a suitable definition of `is_prime`, see Erdős-primes#jq.
def add(s): reduce s as $x (null; .+$x);
def count_divisors:
add(
if . == 1 then 1
else . as $n
| label $out
| range(1; $n) as $i
| ($i * $i) as $i2
| if $i2 > $n then break $out
else if $i2 == $n
then 1
elif ($n % $i) == 0
then 2
else empty
end
end
end);
1000, 100000
| "\nn with odd prime divisor counts, 1 < n < \(.):",
(range(1;.) | select(count_divisors | (. > 2 and is_prime)))
- Output:
n with odd prime divisor counts, 1 < n < 1000: 4 9 16 25 49 64 81 121 169 289 361 529 625 729 841 961 n with odd prime divisor counts, 1 < n < 100000: 4 9 .... 85849 94249 96721 97969
Julia
using Primes
ispdc(n) = (ndivs = prod(collect(values(factor(n))).+ 1); ndivs > 2 && isprime(ndivs))
foreach(p -> print(rpad(p[2], 8), p[1] % 10 == 0 ? "\n" : ""), enumerate(filter(ispdc, 1:100000)))
- Output:
4 9 16 25 49 64 81 121 169 289 361 529 625 729 841 961 1024 1369 1681 1849 2209 2401 2809 3481 3721 4096 4489 5041 5329 6241 6889 7921 9409 10201 10609 11449 11881 12769 14641 15625 16129 17161 18769 19321 22201 22801 24649 26569 27889 28561 29929 32041 32761 36481 37249 38809 39601 44521 49729 51529 52441 54289 57121 58081 59049 63001 65536 66049 69169 72361 73441 76729 78961 80089 83521 85849 94249 96721 97969
Mathematica / Wolfram Language
max = 100000;
maxPrime = NextPrime[Sqrt@max, -1];
maxPower = NextPrime[Log[2, max], -1];
base = NestWhileList[NextPrime, 2, # < maxPrime &];
g = NestWhileList[NextPrime, 3, # < maxPower &] - 1;
ans = Sort@Select[Flatten@Table[base^n, {n, g}], # < max &];
Labeled[Partition[Select[ans, # < 1000 &], UpTo[8]] //
TableForm, "Numbers up to 1000 with prime divisor counts:", Top]
Labeled[Partition[ans, UpTo[8]] //
TableForm, "Numbers up to 100,000 with prime divisor counts:", Top]
- Output:
Numbers up to 1000 with prime divisor counts: 4 9 16 25 49 64 81 121 169 289 361 529 625 729 841 961
Numbers up to 100,000 with prime divisor counts: 4 9 16 25 49 64 81 121 169 289 361 529 625 729 841 961 1024 1369 1681 1849 2209 2401 2809 3481 3721 4096 4489 5041 5329 6241 6889 7921 9409 10201 10609 11449 11881 12769 14641 15625 16129 17161 18769 19321 22201 22801 24649 26569 27889 28561 29929 32041 32761 36481 37249 38809 39601 44521 49729 51529 52441 54289 57121 58081 59049 63001 66049 69169 72361 73441 76729 78961 80089 83521 85849 94249 96721 97969
Nim
Checking only divisors of squares (see discussion).
import math, sequtils, strformat, strutils
func divCount(n: Positive): int =
var n = n
for d in 1..n:
if d * d > n: break
if n mod d == 0:
inc result
if n div d != d:
inc result
func isOddPrime(n: Positive): bool =
if n < 3 or n mod 2 == 0: return false
if n mod 3 == 0: return n == 3
var d = 5
while d <= sqrt(n.toFloat).int:
if n mod d == 0: return false
inc d, 2
if n mod d == 0: return false
inc d, 4
result = true
iterator numWithOddPrimeDivisorCount(lim: Positive): int =
for k in 1..sqrt(lim.toFloat).int:
let n = k * k
if n.divCount().isOddPrime():
yield n
var list = toSeq(numWithOddPrimeDivisorCount(1000))
echo &"Found {list.len} numbers between 1 and 999 whose number of divisors is an odd prime:"
echo list.join(" ")
echo()
list = toSeq(numWithOddPrimeDivisorCount(100_000))
echo &"Found {list.len} numbers between 1 and 99_999 whose number of divisors is an odd prime:"
for i, n in list:
stdout.write &"{n:5}", if (i + 1) mod 10 == 0: '\n' else: ' '
echo()
- Output:
Found 16 numbers between 1 and 999 whose number of divisors is an odd prime: 4 9 16 25 49 64 81 121 169 289 361 529 625 729 841 961 Found 79 numbers between 1 and 99_999 whose number of divisors is an odd prime: 4 9 16 25 49 64 81 121 169 289 361 529 625 729 841 961 1024 1369 1681 1849 2209 2401 2809 3481 3721 4096 4489 5041 5329 6241 6889 7921 9409 10201 10609 11449 11881 12769 14641 15625 16129 17161 18769 19321 22201 22801 24649 26569 27889 28561 29929 32041 32761 36481 37249 38809 39601 44521 49729 51529 52441 54289 57121 58081 59049 63001 65536 66049 69169 72361 73441 76729 78961 80089 83521 85849 94249 96721 97969
PARI/GP
ispdc(n) = {
local(factors, ndivs);
factors = factor(n);
ndivs = numdiv(n);
ndivs > 2 && isprime(ndivs)
};
{
cnt=0;
for (i=1, 100000,
if (ispdc(i),
print1(i, " ");
cnt=cnt+1;
if(cnt%10==0,print(""));
);
);
}
- Output:
4 9 16 25 49 64 81 121 169 289 361 529 625 729 841 961 1024 1369 1681 1849 2209 2401 2809 3481 3721 4096 4489 5041 5329 6241 6889 7921 9409 10201 10609 11449 11881 12769 14641 15625 16129 17161 18769 19321 22201 22801 24649 26569 27889 28561 29929 32041 32761 36481 37249 38809 39601 44521 49729 51529 52441 54289 57121 58081 59049 63001 65536 66049 69169 72361 73441 76729 78961 80089 83521 85849 94249 96721 97969
Pascal
program FacOfInteger;
{$IFDEF FPC}
// {$R+,O+} //debuging purpose
{$MODE DELPHI}
{$Optimization ON,ALL}
{$ELSE}
{$APPTYPE CONSOLE}
{$ENDIF}
uses
sysutils;
//#############################################################################
//Prime decomposition
type
tPot = record
potSoD : Uint64;
potPrim,
potMax :Uint32;
end;
tprimeFac = record
pfPrims : array[0..13] of tPot;
pfSumOfDivs : Uint64;
pfCnt,
pfNum,
pfDivCnt: Uint32;
end;
tSmallPrimes = array[0..6541] of Word;
tItem = NativeUint;
tDivisors = array of tItem;
tpDivisor = pNativeUint;
var
SmallPrimes: tSmallPrimes;
procedure InsertSort(pDiv:tpDivisor; Left, Right : NativeInt );
var
I, J: NativeInt;
Pivot : tItem;
begin
for i:= 1 + Left to Right do
begin
Pivot:= pDiv[i];
j:= i - 1;
while (j >= Left) and (pDiv[j] > Pivot) do
begin
pDiv[j+1]:=pDiv[j];
Dec(j);
end;
pDiv[j+1]:= pivot;
end;
end;
procedure InitSmallPrimes;
var
pr,testPr,j,maxprimidx,delta: Uint32;
isPrime : boolean;
Begin
SmallPrimes[0] := 2;
SmallPrimes[1] := 3;
delta := 2;
maxprimidx := 1;
pr := 5;
repeat
isprime := true;
j := 0;
repeat
testPr := SmallPrimes[j];
IF testPr*testPr > pr then
break;
If pr mod testPr = 0 then
Begin
isprime := false;
break;
end;
inc(j);
until false;
if isprime then
Begin
inc(maxprimidx);
SmallPrimes[maxprimidx]:= pr;
end;
inc(pr,delta);
delta := 2+4-delta;
until pr > 1 shl 16 -1;
end;
function isPrime(n:Uint32):boolean;
var
pr,idx: NativeInt;
begin
result := n in [2,3];
if NOT(result) AND (n>4) AND (n AND 1 <> 0 ) then
begin
idx := 1;
repeat
pr := SmallPrimes[idx];
result := (n mod pr) <>0;
inc(idx);
until NOT(result) or (sqr(pr)>n) or (idx > High(SmallPrimes));
end;
end;
procedure PrimeFacOut(const primeDecomp:tprimeFac;proper:Boolean=true);
var
i,k : Int32;
begin
with primeDecomp do
Begin
write(pfNum,' = ');
k := pfCnt-1;
For i := 0 to k-1 do
with pfPrims[i] do
If potMax = 1 then
write(potPrim,'*')
else
write(potPrim,'^',potMax,'*');
with pfPrims[k] do
If potMax = 1 then
write(potPrim)
else
write(potPrim,'^',potMax);
if proper then
writeln(' got ',pfDivCnt-1,' proper divisors with sum : ',pfSumOfDivs-pfNum)
else
writeln(' got ',pfDivCnt,' divisors with sum : ',pfSumOfDivs);
end;
end;
procedure PrimeDecomposition(var res:tprimeFac;n:Uint32);
var
DivSum,fac:Uint64;
i,pr,cnt,DivCnt,quot{to minimize divisions} : NativeUint;
Begin
if SmallPrimes[0] <> 2 then
InitSmallPrimes;
res.pfNum := n;
cnt := 0;
DivCnt := 1;
DivSum := 1;
i := 0;
if n <= 1 then
Begin
with res.pfPrims[0] do
Begin
potPrim := n;
potMax := 1;
end;
cnt := 1;
end
else
repeat
pr := SmallPrimes[i];
IF pr*pr>n then
Break;
quot := n div pr;
IF pr*quot = n then
with res do
Begin
with pfPrims[Cnt] do
Begin
potPrim := pr;
potMax := 0;
fac := pr;
repeat
n := quot;
quot := quot div pr;
inc(potMax);
fac *= pr;
until pr*quot <> n;
DivCnt *= (potMax+1);
DivSum *= (fac-1)DIV (pr-1);
end;
inc(Cnt);
end;
inc(i);
until false;
//a big prime left over?
IF n > 1 then
with res do
Begin
with pfPrims[Cnt] do
Begin
potPrim := n;
potMax := 1;
end;
inc(Cnt);
DivCnt *= 2;
DivSum *= n+1;
end;
with res do
Begin
pfCnt:= cnt;
pfDivCnt := DivCnt;
pfSumOfDivs := DivSum;
end;
end;
function isAbundant(const pD:tprimeFac):boolean;inline;
begin
with pd do
result := pfSumOfDivs-pfNum > pfNum;
end;
function DivCount(const pD:tprimeFac):NativeUInt;inline;
begin
result := pD.pfDivCnt;
end;
function SumOfDiv(const primeDecomp:tprimeFac):NativeUInt;inline;
begin
result := primeDecomp.pfSumOfDivs;
end;
procedure GetDivs(var pD:tprimeFac;var Divs:tDivisors);
var
pDivs : tpDivisor;
i,len,j,l,p,pPot,k: NativeInt;
Begin
i := DivCount(pD);
IF i > Length(Divs) then
setlength(Divs,i);
pDivs := @Divs[0];
pDivs[0] := 1;
len := 1;
l := len;
For i := 0 to pD.pfCnt-1 do
with pD.pfPrims[i] do
Begin
//Multiply every divisor before with the new primefactors
//and append them to the list
k := potMax-1;
p := potPrim;
pPot :=1;
repeat
pPot *= p;
For j := 0 to len-1 do
Begin
pDivs[l]:= pPot*pDivs[j];
inc(l);
end;
dec(k);
until k<0;
len := l;
end;
//Sort. Insertsort much faster than QuickSort in this special case
InsertSort(pDivs,0,len-1);
end;
Function GetDivisors(var pD:tprimeFac;n:Uint32;var Divs:tDivisors):Int32;
var
i:Int32;
Begin
if pD.pfNum <> n then
PrimeDecomposition(pD,n);
i := DivCount(pD);
IF i > Length(Divs) then
setlength(Divs,i+1);
GetDivs(pD,Divs);
result := DivCount(pD);
end;
procedure AllFacsOut(var pD:tprimeFac;n: Uint32;Divs:tDivisors;proper:boolean=true);
var
k,j: Int32;
Begin
k := GetDivisors(pD,n,Divs)-1;// zero based
PrimeFacOut(pD,proper);
IF proper then
dec(k);
IF k > 0 then
Begin
For j := 0 to k-1 do
write(Divs[j],',');
writeln(Divs[k]);
end;
end;
//Prime decomposition
//#############################################################################
procedure SpeedTest(var pD: tprimeFac;Limit:Uint32);
var
Ticks : Int64;
number,numSqr,Cnt: UInt32;
Begin
Ticks := GetTickCount64;
Cnt := 0;
number := 1;
numSqr:=1;
repeat
number += 1;
numSqr := sqr(number);
PrimeDecomposition(pD,numSqr);
IF DivCount(pD)>2 then
if isPrime(DivCount(pD)) then
inc(cnt);//writeln(number:5,numSqr:10,DivCount(pD):5);
until numSqr>= Limit;
writeln('SpeedTest ',(GetTickCount64-Ticks)/1000:0:3,' secs for 1..',Limit,' found ',Cnt);
writeln;
end;
var
pD: tprimeFac;
Divisors : tDivisors;
numroot,num,cnt : Uint32;
BEGIN
InitSmallPrimes;
setlength(Divisors,1);
write('':4);
for cnt := 1 to 10 do
write(cnt:7);
writeln;
cnt := 0;
write(cnt:3,':');
For numroot := 2 to 1000 do
begin
num := sqr(numroot);
PrimeDecomposition(pD,num);
IF DivCount(pD)>2 then
if isPrime(DivCount(pD)) then
begin
write(num:7);
inc(cnt);
if cnt MOD 10 =0 then
Begin
writeln;write(cnt:3,':');
end;
end;
end;
if cnt MOD 8 <>0 then
writeln;
writeln;
SpeedTest(pD,4000*1000*1000);
END.
- Output:
1 2 3 4 5 6 7 8 9 10 0: 4 9 16 25 49 64 81 121 169 289 10: 361 529 625 729 841 961 1024 1369 1681 1849 20: 2209 2401 2809 3481 3721 4096 4489 5041 5329 6241 30: 6889 7921 9409 10201 10609 11449 11881 12769 14641 15625 40: 16129 17161 18769 19321 22201 22801 24649 26569 27889 28561 50: 29929 32041 32761 36481 37249 38809 39601 44521 49729 51529 60: 52441 54289 57121 58081 59049 63001 65536 66049 69169 72361 70: 73441 76729 78961 80089 83521 85849 94249 96721 97969 100489 80: 109561 113569 117649 120409 121801 124609 128881 130321 134689 139129 90: 143641 146689 151321 157609 160801 167281 175561 177241 185761 187489 100: 192721 196249 201601 208849 212521 214369 218089 229441 237169 241081 110: 249001 253009 259081 262144 271441 273529 279841 292681 299209 310249 120: 316969 323761 326041 332929 344569 351649 358801 361201 368449 375769 130: 380689 383161 398161 410881 413449 418609 426409 434281 436921 452929 140: 458329 466489 477481 491401 502681 516961 528529 531441 537289 546121 150: 552049 564001 573049 579121 591361 597529 619369 635209 654481 657721 160: 674041 677329 683929 687241 703921 707281 727609 734449 737881 744769 170: 769129 776161 779689 786769 822649 829921 844561 863041 877969 885481 180: 896809 908209 923521 935089 942841 954529 966289 982081 994009 SpeedTest 0.230 secs for 1..4000000000 found 6417
Perl
use strict;
use warnings;
use ntheory <is_prime divisors>;
push @matches, $_**2 for grep { is_prime divisors $_**2 } 1..int sqrt 1e5;
print @matches . " matching:\n" . (sprintf "@{['%6d' x @matches]}", @matches) =~ s/(.{72})/$1\n/gr;
- Output:
79 matching: 4 9 16 25 49 64 81 121 169 289 361 529 625 729 841 961 1024 1369 1681 1849 2209 2401 2809 3481 3721 4096 4489 5041 5329 6241 6889 7921 9409 10201 10609 11449 11881 12769 14641 15625 16129 17161 18769 19321 22201 22801 24649 26569 27889 28561 29929 32041 32761 36481 37249 38809 39601 44521 49729 51529 52441 54289 57121 58081 59049 63001 65536 66049 69169 72361 73441 76729 78961 80089 83521 85849 94249 96721 97969
Phix
atom t0 = time()
function pd(integer n) n = length(factors(n,1)) return n!=2 and is_prime(n) end function
for n in {1e3,1e5} do
sequence res = filter(tagset(n),pd)
printf(1,"%d < %,d found: %v\n",{length(res),n,shorten(res,"",5)})
end for
?elapsed(time()-t0)
- Output:
16 < 1,000 found: {4,9,16,25,49,"...",529,625,729,841,961} 79 < 100,000 found: {4,9,16,25,49,"...",83521,85849,94249,96721,97969} "0.4s"
smarter, faster
As per Nigel's analysis on the talk page, valid numbers must be of the form a^(b-1) where a is prime and b is an odd prime, with no further checking required.
atom t1 = time()
for n in {1e3,1e5,4e9} do
sequence r = {}
for a in get_primes_le(floor(sqrt(n))) do
integer b = 2
while true do
atom c = power(a,get_prime(b)-1)
if c>n then exit end if
r &= c
b += 1
end while
end for
r = sort(deep_copy(r))
printf(1,"%d < %,d found: %v\n",{length(r),n,shorten(r,"",5)})
end for
?elapsed(time()-t1)
- Output:
16 < 1,000 found: {4,9,16,25,49,"...",529,625,729,841,961} 79 < 100,000 found: {4,9,16,25,49,"...",83521,85849,94249,96721,97969} 6417 < 4,000,000,000 found: {4,9,16,25,49,"...",3991586041.0,3993860809.0,3994113601.0,3995630521.0,3999424081.0} "0.0s"
Python
from sympy import isprime
def main():
num_count = 1
for num in range(4, 100_000):
count = 1
for divisor in range(1, (num // 2) + 1):
if num % divisor == 0:
count += 1
if count > 2 and isprime(count):
print("%6d" % num, end=' ')
if num_count % 5 == 0:
print()
num_count += 1
if __name__ == '__main__':
main()
- Output:
4 9 16 25 49 64 81 121 169 289 361 529 625 729 841 961 1024 1369 1681 1849 2209 2401 2809 3481 3721 4096 4489 5041 5329 6241 6889 7921 9409 10201 10609 11449 11881 12769 14641 15625 16129 17161 18769 19321 22201 22801 24649 26569 27889 28561 29929 32041 32761 36481 37249 38809 39601 44521 49729 51529 52441 54289 57121 58081 59049 63001 65536 66049 69169 72361 73441 76729 78961 80089 83521 85849 94249 96721 97969
Quackery
Noting that only perfect squares can fit the criteria (see discussion).
factors
and isqrt
are defined at Factors of an integer#Quackery.
[ dip number$
over size -
space swap of
swap join echo$ ] is recho ( n n --> )
[ []
100000 isqrt drop times
[ i^ dup * factors size
dup 2 = iff drop done
factors size 2 = while
i^ dup * join ] ]
[ dup [] != while
10 split swap
witheach [ 6 recho ]
cr again ]
drop
- Output:
4 9 16 25 49 64 81 121 169 289 361 529 625 729 841 961 1024 1369 1681 1849 2209 2401 2809 3481 3721 4096 4489 5041 5329 6241 6889 7921 9409 10201 10609 11449 11881 12769 14641 15625 16129 17161 18769 19321 22201 22801 24649 26569 27889 28561 29929 32041 32761 36481 37249 38809 39601 44521 49729 51529 52441 54289 57121 58081 59049 63001 65536 66049 69169 72361 73441 76729 78961 80089 83521 85849 94249 96721 97969
Raku
use Prime::Factor;
my $ceiling = ceiling sqrt 1e5;
say display :10cols, :fmt('%6d'), (^$ceiling)».² .grep: { .&divisors.is-prime };
sub display ($list, :$cols = 10, :$fmt = '%6d', :$title = "{+$list} matching:\n" ) {
cache $list;
$title ~ $list.batch($cols)».fmt($fmt).join: "\n"
}
- Output:
79 matching: 4 9 16 25 49 64 81 121 169 289 361 529 625 729 841 961 1024 1369 1681 1849 2209 2401 2809 3481 3721 4096 4489 5041 5329 6241 6889 7921 9409 10201 10609 11449 11881 12769 14641 15625 16129 17161 18769 19321 22201 22801 24649 26569 27889 28561 29929 32041 32761 36481 37249 38809 39601 44521 49729 51529 52441 54289 57121 58081 59049 63001 65536 66049 69169 72361 73441 76729 78961 80089 83521 85849 94249 96721 97969
REXX
/*REXX pgm finds positive integers N whose # of divisors is prime (& ¬=2), where N<1000.*/
parse arg hi cols . /*obtain optional arguments from the CL*/
if hi=='' | hi=="," then hi= 1000 /*Not specified? Then use the defaults*/
if cols=='' | cols=="," then cols= 10 /* " " " " " " */
call genP /*build array of semaphores for primes.*/
w= 10 /*W: the maximum width of any column. */
title= ' positive integers N whose number of divisors is prime (and not equal to 2), ' ,
"where N < " commas(hi)
say ' index │'center(title, 1 + cols*(w+1) )
say '───────┼'center("" , 1 + cols*(w+1), '─')
finds= 0; idx= 1; $=
do j=2; jj= j*j; if jj>=hi then leave /*process positive square ints in range*/
n= nDivs(jj); if n==2 then iterate /*get number of divisors of composite J*/
if \!.n then iterate /*Number divisors prime? No, then skip*/
finds= finds + 1 /*bump the number of found numbers. */
$= $ right( commas(j), w) /*add a positive integer ──► $ list. */
if finds//cols\==0 then iterate /*have we populated a line of output? */
say center(idx, 7)'│' substr($, 2); $= /*display what we have so far (cols). */
idx= idx + cols /*bump the index count for the output*/
end /*j*/ /* [↑] process a range of integers. */
if $\=='' then say center(idx, 7)"│" substr($, 2) /*possible display residual output.*/
say '───────┴'center("" , 1 + cols*(w+1), '─')
say
say 'Found ' commas(finds) title
exit 0 /*stick a fork in it, we're all done. */
/*──────────────────────────────────────────────────────────────────────────────────────*/
commas: parse arg ?; do jc=length(?)-3 to 1 by -3; ?=insert(',', ?, jc); end; return ?
/*──────────────────────────────────────────────────────────────────────────────────────*/
nDivs: procedure; parse arg x; d= 2; if x==1 then return 1 /*handle special case of 1*/
odd= x // 2 /* [↓] process EVEN or ODD ints. ___*/
do j=2+odd by 1+odd while j*j<x /*divide by all the integers up to √ x */
if x//j==0 then d= d + 2 /*÷? Add two divisors to the total. */
end /*j*/ /* [↑] % ≡ integer division. */
if j*j==x then return d + 1 /*Was X a square? Then add 1 to total.*/
return d /*return the total. */
/*──────────────────────────────────────────────────────────────────────────────────────*/
genP: @.1=2; @.2=3; @.3=5; @.4=7; @.5=11 /*define some low primes. */
!.=0; !.2=1; !.3=1; !.5=1; !.7=1; !.11=1 /* " " " " semaphores. */
#=5; s.#= @.# **2 /*number of primes so far; prime². */
do j=@.#+2 by 2 to hi-1 /*find odd primes from here on. */
parse var j '' -1 _; if _==5 then iterate /*J divisible by 5? (right dig)*/
if j// 3==0 then iterate /*" " " 3? */
if j// 7==0 then iterate /*" " " 7? */
do k=5 while s.k<=j /* [↓] divide by the known odd primes.*/
if j // @.k == 0 then iterate j /*Is J ÷ X? Then not prime. ___ */
end /*k*/ /* [↑] only process numbers ≤ √ J */
#= #+1; @.#= j; s.#= j*j; !.j= 1 /*bump # of Ps; assign next P; P²; P# */
end /*j*/; return
- output when using the default inputs:
index │ positive integers N whose number of divisors is prime (and not equal to 2), where N < 1,000 ───────┼─────────────────────────────────────────────────────────────────────────────────────────────────────────────── 1 │ 4 9 16 25 49 64 81 121 169 289 11 │ 361 529 625 729 841 961 ───────┴─────────────────────────────────────────────────────────────────────────────────────────────────────────────── Found 16 positive integers N whose number of divisors is prime (and not equal to 2), where N < 1,000
- output when using the input of: 100000
index │ positive integers N whose number of divisors is prime (and not equal to 2), where N < 100,000 ───────┼─────────────────────────────────────────────────────────────────────────────────────────────────────────────── 1 │ 4 9 16 25 49 64 81 121 169 289 11 │ 361 529 625 729 841 961 1,024 1,369 1,681 1,849 21 │ 2,209 2,401 2,809 3,481 3,721 4,096 4,489 5,041 5,329 6,241 31 │ 6,889 7,921 9,409 10,201 10,609 11,449 11,881 12,769 14,641 15,625 41 │ 16,129 17,161 18,769 19,321 22,201 22,801 24,649 26,569 27,889 28,561 51 │ 29,929 32,041 32,761 36,481 37,249 38,809 39,601 44,521 49,729 51,529 61 │ 52,441 54,289 57,121 58,081 59,049 63,001 65,536 66,049 69,169 72,361 71 │ 73,441 76,729 78,961 80,089 83,521 85,849 94,249 96,721 97,969 ───────┴─────────────────────────────────────────────────────────────────────────────────────────────────────────────── Found 79 positive integers N whose number of divisors is prime (and not equal to 2), where N < 100,000
Ring
load "stdlib.ring"
row = 0
see "working..." + nl
see "Numbers which count of divisors is prime are:" + nl
for n = 1 to 1000
num = 0
for m = 1 to n
if n%m = 0
num++
ok
next
if isprime(num) and num != 2
see "" + n + " "
row++
if row%5 = 0
see nl
ok
ok
next
see nl + "Found " + row + " numbers" + nl
see "done..." + nl
- Output:
working... Numbers which count of divisors is prime are: 4 9 16 25 49 64 81 121 169 289 361 529 625 729 841 961 Found 16 numbers done...
RPL
≪ { }
1 999 FOR j
j DIVIS SIZE
IF DUP 2 == THEN DROP ELSE
IF ISPRIME? THEN j + END
END
NEXT
≫ 'TASK' STO
- Output:
1; { 4 9 16 25 49 64 81 121 169 289 361 529 625 729 841 961 }
Runs in 2 minutes 16 seconds on a HP-50g.
Rust
fn count_divisors(number : u64 ) -> u64 {
let mut divisors : u64 = 0 ;
for n in 1..=number {
if number % n == 0 {
divisors += 1 ;
}
}
divisors
}
fn main() {
let mut the_numbers : Vec<u64> = Vec::new( ) ;
let mut count : i32 = 0 ;
for n in 1..100000 {
let divis = count_divisors( n as u64 ) ;
if divis > 2 && primal::is_prime( divis ) {
the_numbers.push( n ) ;
}
}
for n in the_numbers {
count += 1 ;
print!(" {:6}" , n ) ;
if count % 8 == 0 {
println!( ) ;
}
}
println!( ) ;
println!("count is {}" , count) ;
}
- Output:
4 9 16 25 49 64 81 121 169 289 361 529 625 729 841 961 1024 1369 1681 1849 2209 2401 2809 3481 3721 4096 4489 5041 5329 6241 6889 7921 9409 10201 10609 11449 11881 12769 14641 15625 16129 17161 18769 19321 22201 22801 24649 26569 27889 28561 29929 32041 32761 36481 37249 38809 39601 44521 49729 51529 52441 54289 57121 58081 59049 63001 65536 66049 69169 72361 73441 76729 78961 80089 83521 85849 94249 96721 97969 count is 79
Ruby
Testing squares only, according to observation on the discussion page.
require 'prime'
def tau(n) = n.prime_division.inject(1){|res, (d, exp)| res *= exp+1}
res = (1..Integer.sqrt(100_000)).filter_map{|n| sqr = n*n; sqr if tau(sqr).prime? }
res.each_slice(10){|slice| puts "%10d"*slice.size % slice}
- Output:
4 9 16 25 49 64 81 121 169 289 361 529 625 729 841 961 1024 1369 1681 1849 2209 2401 2809 3481 3721 4096 4489 5041 5329 6241 6889 7921 9409 10201 10609 11449 11881 12769 14641 15625 16129 17161 18769 19321 22201 22801 24649 26569 27889 28561 29929 32041 32761 36481 37249 38809 39601 44521 49729 51529 52441 54289 57121 58081 59049 63001 65536 66049 69169 72361 73441 76729 78961 80089 83521 85849 94249 96721 97969
Sidef
var limit = 100_000
say "Positive integers under #{limit.commify} whose number of divisors is an odd prime:"
1..limit -> grep { !.is_prime && .sigma0.is_prime }.each_slice(10, {|*a|
say a.map{'%6s' % _}.join(' ')
})
- Output:
Positive integers under 100,000 whose number of divisors is an odd prime: 4 9 16 25 49 64 81 121 169 289 361 529 625 729 841 961 1024 1369 1681 1849 2209 2401 2809 3481 3721 4096 4489 5041 5329 6241 6889 7921 9409 10201 10609 11449 11881 12769 14641 15625 16129 17161 18769 19321 22201 22801 24649 26569 27889 28561 29929 32041 32761 36481 37249 38809 39601 44521 49729 51529 52441 54289 57121 58081 59049 63001 65536 66049 69169 72361 73441 76729 78961 80089 83521 85849 94249 96721 97969
Wren
import "./math" for Int
import "./fmt" for Fmt
var limit = 1e5
var results = []
var i = 2
while (i * i < limit) {
var n = Int.divisors(i * i).count
if (n > 2 && Int.isPrime(n)) results.add(i * i)
i = i + 1
}
Fmt.print("Positive integers under $,7d whose number of divisors is an odd prime:", limit)
Fmt.tprint("$,7d", results, 10)
var under1000 = results.count { |r| r < 1000 }
System.print("\nFound %(results.count) such integers (%(under1000) under 1,000).")
- Output:
Positive integers under 100,000 whose number of divisors is an odd prime: 4 9 16 25 49 64 81 121 169 289 361 529 625 729 841 961 1,024 1,369 1,681 1,849 2,209 2,401 2,809 3,481 3,721 4,096 4,489 5,041 5,329 6,241 6,889 7,921 9,409 10,201 10,609 11,449 11,881 12,769 14,641 15,625 16,129 17,161 18,769 19,321 22,201 22,801 24,649 26,569 27,889 28,561 29,929 32,041 32,761 36,481 37,249 38,809 39,601 44,521 49,729 51,529 52,441 54,289 57,121 58,081 59,049 63,001 65,536 66,049 69,169 72,361 73,441 76,729 78,961 80,089 83,521 85,849 94,249 96,721 97,969 Found 79 such integers (16 under 1,000).
XPL0
func IsPrime(N); \Return 'true' if N is a prime number
int N, I;
[if N <= 1 then return false;
for I:= 2 to sqrt(N) do
if rem(N/I) = 0 then return false;
return true;
];
func Divisors(N); \Return number of unique divisors of N
int N, SN, Count, D;
[SN:= sqrt(N); \N must be a perfect square to get an odd (prime>2) count
if SN*SN # N then return 0;
Count:= 3; \SN, 1 and N are unique divisors of N >= 4
for D:= 2 to SN-1 do
if rem(N/D) = 0 then Count:= Count+2;
return Count;
];
int N, Count;
[Count:= 0;
for N:= 4 to 100_000-1 do
if IsPrime(Divisors(N)) then
[Count:= Count+1;
IntOut(0, N);
if rem(Count/10) = 0 then CrLf(0) else ChOut(0, 9\tab\);
];
]
- Output:
4 9 16 25 49 64 81 121 169 289 361 529 625 729 841 961 1024 1369 1681 1849 2209 2401 2809 3481 3721 4096 4489 5041 5329 6241 6889 7921 9409 10201 10609 11449 11881 12769 14641 15625 16129 17161 18769 19321 22201 22801 24649 26569 27889 28561 29929 32041 32761 36481 37249 38809 39601 44521 49729 51529 52441 54289 57121 58081 59049 63001 65536 66049 69169 72361 73441 76729 78961 80089 83521 85849 94249 96721 97969
- Draft Programming Tasks
- Prime Numbers
- 11l
- Action!
- Action! Sieve of Eratosthenes
- ALGOL 60
- ALGOL 68
- AppleScript
- Arturo
- AWK
- BASIC
- BASIC256
- FreeBASIC
- QBasic
- PureBasic
- Yabasic
- C++
- CLU
- Delphi
- SysUtils,StdCtrls
- EasyLang
- F Sharp
- Factor
- Go
- Go-rcu
- J
- Java
- Jq
- Julia
- Mathematica
- Wolfram Language
- Nim
- PARI/GP
- Pascal
- Perl
- Ntheory
- Phix
- Python
- SymPy
- Quackery
- Raku
- REXX
- Ring
- RPL
- Rust
- Ruby
- Sidef
- Wren
- Wren-math
- Wren-fmt
- XPL0