Sequence: nth number with exactly n divisors: Difference between revisions
Sequence: nth number with exactly n divisors (view source)
Revision as of 11:36, 4 February 2024
, 3 months ago→{{header|Wren}}: Minor tidy
m (→{{header|Wren}}: Minor tidy) |
|||
| (7 intermediate revisions by 6 users not shown) | |||
Line 14:
:*[[Sequence: smallest number greater than previous term with exactly n divisors]]
:*[[Sequence: smallest number with exactly n divisors]]
<br/><br/>
=={{header|ALGOL 68}}==
{{works with|ALGOL 68G|Any - tested with release 2.8.3.win32}}
Uses Algol 68G's LONG LONG INT (with default precision) to handle the large numbers.<br/>
Builds a table of proper divisor counts for testing non-prime n (As with other samples uses the fact that for prime n, the required number is the nth prime to the power n - 1.<br/>
Using a table is practical for smallish n, (a table of 500 000 elements is used here, big enough for up to n = 24. For n = 25, a table of over 52 000 000 would be required.<br/>
For non-prime n, the divisors are also shown.
<syntaxhighlight lang="algol68">
BEGIN
INT max number = 500 000; # maximum number we will count the divisors of #
# form a table of proper divisor counts #
[ 1 : max number ]INT pdc; FOR i TO UPB pdc DO pdc[ i ] := 1 OD;
pdc[ 1 ] := 0;
FOR i FROM 2 TO UPB pdc DO
FOR j FROM i + i BY i TO UPB pdc DO pdc[ j ] +:= 1 OD
OD;
# find the first few primes #
[ 1 : 30 ]INT prime;
INT p count := 0;
FOR i WHILE p count < UPB prime DO
IF pdc[ i ] = 1 THEN prime[ p count +:= 1 ] := i FI
OD;
# show the nth number with n divisors #
INT w = -43; # width to print the numbers (negative means no leading +) #
print( ( " 1: ", whole( 1, w ), " | 1", newline ) );
FOR i FROM 2 TO 23 DO
print( ( whole( i, -2 ), ": " ) );
IF pdc( i ) = 1 THEN
print( ( whole( LENG LENG prime[ i ] ^ ( i - 1 ), w ), newline ) )
ELSE
INT c := 0;
FOR j TO UPB pdc WHILE c < i DO
IF pdc[ j ] = i - 1 THEN
c +:= 1;
IF c = i THEN
print( ( whole( j, w ), " | 1" ) );
FOR d FROM 2 TO j OVER 2 DO
IF j MOD d = 0 THEN print( ( " ", whole( d, 0 ) ) ) FI
OD;
print( ( " ", whole( j, 0 ), newline ) )
FI
FI
OD
FI
OD
END
</syntaxhighlight>
{{out}}
<pre style="font-size: 12px">
1: 1 | 1
2: 3
3: 25
4: 14 | 1 2 7 14
5: 14641
6: 44 | 1 2 4 11 22 44
7: 24137569
8: 70 | 1 2 5 7 10 14 35 70
9: 1089 | 1 3 9 11 33 99 121 363 1089
10: 405 | 1 3 5 9 15 27 45 81 135 405
11: 819628286980801
12: 160 | 1 2 4 5 8 10 16 20 32 40 80 160
13: 22563490300366186081
14: 2752 | 1 2 4 8 16 32 43 64 86 172 344 688 1376 2752
15: 9801 | 1 3 9 11 27 33 81 99 121 297 363 891 1089 3267 9801
16: 462 | 1 2 3 6 7 11 14 21 22 33 42 66 77 154 231 462
17: 21559177407076402401757871041
18: 1044 | 1 2 3 4 6 9 12 18 29 36 58 87 116 174 261 348 522 1044
19: 740195513856780056217081017732809
20: 1520 | 1 2 4 5 8 10 16 19 20 38 40 76 80 95 152 190 304 380 760 1520
21: 141376 | 1 2 4 8 16 32 47 64 94 188 376 752 1504 2209 3008 4418 8836 17672 35344 70688 141376
22: 84992 | 1 2 4 8 16 32 64 83 128 166 256 332 512 664 1024 1328 2656 5312 10624 21248 42496 84992
23: 1658509762573818415340429240403156732495289</pre>
=={{header|C}}==
{{trans|C++}}
<
#include <stdbool.h>
#include <stdint.h>
Line 135 ⟶ 210:
return 0;
}</
{{out}}
<pre>A073916(1) = 1
Line 155 ⟶ 230:
=={{header|C++}}==
{{trans|Java}}
<
#include <vector>
Line 246 ⟶ 321:
return 0;
}</
{{out}}
<pre>A073916(1) = 1
Line 263 ⟶ 338:
A073916(14) = 2752
A073916(15) = 9801</pre>
=={{header|D}}==
{{trans|Java}}
<syntaxhighlight lang="d">import std.bigint;
import std.math;
import std.stdio;
bool isPrime(long test) {
if (test == 2) {
return true;
}
if (test % 2 == 0) {
return false;
}
for (long d = 3 ; d * d <= test; d += 2) {
if (test % d == 0) {
return false;
}
}
return true;
}
int[] calcSmallPrimes(int numPrimes) {
int[] smallPrimes;
smallPrimes ~= 2;
int count = 0;
int n = 3;
while (count < numPrimes) {
if (isPrime(n)) {
smallPrimes ~= n;
count++;
}
n += 2;
}
return smallPrimes;
}
immutable MAX = 45;
immutable smallPrimes = calcSmallPrimes(MAX);
int getDivisorCount(long n) {
int count = 1;
while (n % 2 == 0) {
n /= 2;
count += 1;
}
for (long d = 3; d * d <= n; d += 2) {
long q = n / d;
long r = n % d;
int dc = 0;
while (r == 0) {
dc += count;
n = q;
q = n / d;
r = n % d;
}
count += dc;
}
if (n != 1) {
count *= 2;
}
return count;
}
BigInt OEISA073916(int n) {
if (isPrime(n) ) {
return BigInt(smallPrimes[n-1]) ^^ (n - 1);
}
int count = 0;
int result = 0;
for (int i = 1; count < n; i++) {
if (n % 2 == 1) {
// The solution for an odd (non-prime) term is always a square number
int root = cast(int) sqrt(cast(real) i);
if (root * root != i) {
continue;
}
}
if (getDivisorCount(i) == n) {
count++;
result = i;
}
}
return BigInt(result);
}
void main() {
foreach (n; 1 .. MAX + 1) {
writeln("A073916(", n, ") = ", OEISA073916(n));
}
}</syntaxhighlight>
{{out}}
<pre>A073916(1) = 1
A073916(2) = 3
A073916(3) = 25
A073916(4) = 14
A073916(5) = 14641
A073916(6) = 44
A073916(7) = 24137569
A073916(8) = 70
A073916(9) = 1089
A073916(10) = 405
A073916(11) = 819628286980801
A073916(12) = 160
A073916(13) = 22563490300366186081
A073916(14) = 2752
A073916(15) = 9801
A073916(16) = 462
A073916(17) = 21559177407076402401757871041
A073916(18) = 1044
A073916(19) = 740195513856780056217081017732809
A073916(20) = 1520
A073916(21) = 141376
A073916(22) = 84992
A073916(23) = 1658509762573818415340429240403156732495289
A073916(24) = 1170
A073916(25) = 52200625
A073916(26) = 421888
A073916(27) = 52900
A073916(28) = 9152
A073916(29) = 1116713952456127112240969687448211536647543601817400964721
A073916(30) = 6768
A073916(31) = 1300503809464370725741704158412711229899345159119325157292552449
A073916(32) = 3990
A073916(33) = 12166144
A073916(34) = 9764864
A073916(35) = 446265625
A073916(36) = 5472
A073916(37) = 11282036144040442334289838466416927162302790252609308623697164994458730076798801
A073916(38) = 43778048
A073916(39) = 90935296
A073916(40) = 10416
A073916(41) = 1300532588674810624476094551095787816112173600565095470117230812218524514342511947837104801
A073916(42) = 46400
A073916(43) = 635918448514386699807643535977466343285944704172890141356181792680152445568879925105775366910081
A073916(44) = 240640
A073916(45) = 327184</pre>
=={{header|Factor}}==
This makes use of most of the optimizations discussed in the Go example.
<
lists.lazy.examples literals math math.functions math.primes
math.primes.factors math.ranges sequences ;
Line 297 ⟶ 511:
: main ( -- ) 45 [1,b] [ dup fn "%2d : %d\n" printf ] each ;
MAIN: main</
{{out}}
<pre>
Line 346 ⟶ 560:
45 : 327184
</pre>
=={{header|FreeBASIC}}==
But it takes an "eternity" to resolve.
<syntaxhighlight lang="freebasic">Dim As Integer n, num, pnum
Dim As Ulongint m, p
Dim As Ulongint limit = 18446744073709551615
Print "The first 15 terms of OEIS:A073916 are:"
For n = 1 To 15
num = 0
For m = 1 To limit
pnum = 0
For p = 1 To limit
If (m Mod p = 0) Then pnum += 1
Next p
If pnum = n Then num += 1
If num = n Then
Print Using "## : &"; n; m
Exit For
End If
Next m
Next n
Sleep</syntaxhighlight>
=={{header|Go}}==
Line 351 ⟶ 590:
The remaining terms (up to the 33rd) are not particularly large and so are calculated by brute force.
<
import (
Line 433 ⟶ 672:
}
}
}</
{{out}}
Line 481 ⟶ 720:
3. While searching for the nth number with exactly n divisors, where feasible a record is kept of any numbers found to have exactly k divisors (k > n) so that the search for these numbers can start from a higher base.
<
import (
Line 601 ⟶ 840:
}
}
}</
{{out}}
Line 653 ⟶ 892:
</pre>
=={{header|Haskell}}==
<
import Math.NumberTheory.ArithmeticFunctions (divisorCount)
import Math.NumberTheory.Primes (Prime, unPrime)
Line 680 ⟶ 919:
main :: IO ()
main = mapM_ print nths</
{{out}}
<pre>(1,1)
Line 717 ⟶ 956:
(34,9764864)
(35,446265625)</pre>
=={{header|J}}==
Implementation:
<syntaxhighlight lang="j">
A073916=: {{
if.1 p: y do. (p:^x:)y-1 return.
elseif.1|y do.
f= *:
else.
f=. ]
end. r=.i.0
off=. 1
while. y>#r do.
r=. r,f off+I.y=*/|:1+_ q:f off+i.y
off=. off+y
end.
(y-1){r
}}"0
</syntaxhighlight>
Task example:
<syntaxhighlight lang="j">
(,. A073916) 1+i.20
1 1
2 3
3 25
4 14
5 14641
6 44
7 24137569
8 70
9 1089
10 405
11 819628286980801
12 160
13 22563490300366186081
14 2752
15 9801
16 462
17 21559177407076402401757871041
18 1044
19 740195513856780056217081017732809
20 1520
</syntaxhighlight>
=={{header|Java}}==
Line 722 ⟶ 1,008:
Replace translation with Java native implementation.
<
import java.math.BigInteger;
import java.util.ArrayList;
Line 811 ⟶ 1,097:
}
</syntaxhighlight>
{{out}}
Line 860 ⟶ 1,146:
A073916(44) = 240640
A073916(45) = 327184
</pre>
=={{header|jq}}==
'''Adapted from [[#Wren|Wren]]'''<br>
'''Works with gojq, the Go implementation of jq'''
See [[Erd%C5%91s-primes#jq]] for a suitable implementation of
`is_prime`.
The precision of the integer arithmetic of the C implementation of jq
is only precise enough
for computing the n-th value up to and including [16,462].
Accordingly gojq was used to produce the output shown below.
'''Preliminaries'''
<syntaxhighlight lang="jq">def count(stream): reduce stream as $i (0; .+1);
# To maintain precision:
def power($b): . as $in | reduce range(0;$b) as $i (1; . * $in);
def primes: 2, (range(3; infinite; 2) | select(is_prime));
# divisors as an unsorted stream
def divisors:
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 $i
elif ($n % $i) == 0
then $i, ($n/$i)
else empty
end
end
end;
</syntaxhighlight>
'''The Task'''
<syntaxhighlight lang="jq"># Emit [n, nth_with_n_divisors] for n in range(1; .+1)
def nth_with_n_divisors:
| [limit( .; primes)] as $primes
| range( 1; 1 + .) as $i
| if $i | is_prime
then [$i, ($primes[$i-1]|power($i-1))]
else {count: 0, j: 1}
| until(.count == $i ;
.cont = false
| if ($i % 2) == 1 then (.j|sqrt|floor) as $sq
| if ($sq * $sq) != .j then .j += 1 | .cont = true else . end
else .
end
| if .cont == false
then if (.j | count(divisors)) == $i
then .count += 1
else .
end
| if .count != $i then .j += 1 else . end
else .
end )
| [ $i, .j]
end;
"The first 33 terms in the sequence are:",
(33 | nth_with_n_divisors)</syntaxhighlight>
{{out}}
<pre>
The first 33 terms in the sequence are:
[1,1]
[2,3]
[3,25]
[4,14]
[5,14641]
[6,44]
[7,24137569]
[8,70]
[9,1089]
[10,405]
[11,819628286980801]
[12,160]
[13,22563490300366186081]
[14,2752]
[15,9801]
[16,462]
[17,21559177407076402401757871041]
[18,1044]
[19,740195513856780056217081017732809]
[20,1520]
[21,141376]
[22,84992]
[23,1658509762573818415340429240403156732495289]
[24,1170]
[25,52200625]
[26,421888]
[27,52900]
[28,9152]
[29,1116713952456127112240969687448211536647543601817400964721]
[30,6768]
[31,1300503809464370725741704158412711229899345159119325157292552449]
[32,3990]
[33,12166144]
</pre>
=={{header|Julia}}==
<
function countdivisors(n)
Line 892 ⟶ 1,281:
nthwithndivisors(35)
</
<pre>
1 : 1
Line 933 ⟶ 1,322:
=={{header|Kotlin}}==
{{trans|Go}}
<
import java.math.BigInteger
Line 1,011 ⟶ 1,400:
}
}
}</
{{output}}
Line 1,052 ⟶ 1,441:
=={{header|Mathematica}} / {{header|Wolfram Language}}==
<
Length[Divisors[n]], {n,
200000}]; t = {}; n = 0; ok = True; While[ok, n++;
Line 1,058 ⟶ 1,447:
c = Flatten[Position[d, n, 1, n]];
If[Length[c] >= n, AppendTo[t, c[[n]]], ok = False]]];
t</
{{out}}
<pre>{1,3,25,14,14641,44,24137569,70,1089,405,819628286980801,160,22563490300366186081,2752,9801,462,21559177407076402401757871041,1044,740195513856780056217081017732809,1520,141376,84992,1658509762573818415340429240403156732495289,1170}</pre>
Line 1,066 ⟶ 1,455:
{{libheader|bignum}}
This is a translation of the fast Go version. It runs in about 23s on our laptop.
<
import bignum
Line 1,144 ⟶ 1,533:
k > records[cd].num and (d == 1 or d == 2 and not cd.isOdd):
records[cd].num = k
inc records[cd].count</
{{out}}
Line 1,197 ⟶ 1,586:
{{libheader|ntheory}}
{{trans|Raku}}
<
use warnings;
use bigint;
Line 1,217 ⟶ 1,606:
}
}
}</
{{out}}
<pre>First 20 terms of OEIS:A073916
Line 1,227 ⟶ 1,616:
Certainly not the fastest way to do it, hence the relatively small limit of 24, which takes less than 0.4s,<br>
whereas a limit of 25 would need to invoke factors() 52 million times which would no doubt take a fair while.
<!--<
<span style="color: #008080;">with</span> <span style="color: #008080;">javascript_semantics</span>
<span style="color: #008080;">constant</span> <span style="color: #000000;">LIMIT</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">24</span>
Line 1,251 ⟶ 1,640:
<span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
<!--</
{{out}}
<pre>
Line 1,283 ⟶ 1,672:
===cheating slightly===
No real patterns that I could see here, but you can still identify and single out the troublemakers (of which there are about 30).
<!--<
<span style="color: #008080;">with</span> <span style="color: #008080;">javascript_semantics</span>
<span style="color: #008080;">include</span> <span style="color: #004080;">mpfr</span><span style="color: #0000FF;">.</span><span style="color: #000000;">e</span>
Line 1,348 ⟶ 1,737:
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"completed in %s\n"</span><span style="color: #0000FF;">,{</span><span style="color: #7060A8;">elapsed</span><span style="color: #0000FF;">(</span><span style="color: #7060A8;">time</span><span style="color: #0000FF;">()-</span><span style="color: #000000;">t0</span><span style="color: #0000FF;">)})</span>
<!--</
{{out}}
<pre>
Line 1,457 ⟶ 1,846:
=={{header|Python}}==
This implementation exploits the fact that terms corresponding to a prime value for n are always the nth prime to the (n-1)th power.
<syntaxhighlight lang="python">
def divisors(n):
divs = [1]
Line 1,521 ⟶ 1,910:
for item in sequence(15):
print(item)
</syntaxhighlight>
<b>Output:</b>
<syntaxhighlight lang="python">
1
3
Line 1,539 ⟶ 1,928:
2752
9801
</syntaxhighlight>
=={{header|Raku}}==
Line 1,547 ⟶ 1,936:
[https://tio.run/##dVLbTsJAEH3nKw4IpgW6ATQY2QAao4kvmsijEFPoVpr05nZr2hh@yk/wx@p0CwUf3JfuzJwzc/ZMYyH9cVEk6RqO92ltojRUMJaZia8G6EihUhliBM9FxrzEiqUXCK5rPde3CTwEY1RLPqRirh9F0mSBHU9gzbB09m3Kk4a@SBJk6IDSCHIsc0wppsFwOCiYUmU@p1uzCSPvwzGx0/xdg74NorR9L/AU0UYDXmVutKKEUu9SxNS4VoldH0PGuryxx7xW7BXHGSqE2grYUto5w@Lu6YU6xqlC68GTiTrMUkIGCSIXz/ePi8nt4OriejhucY00qH8FM9k2j4U0JqO1rTbbowftcO@BQbcZLmGHDiVrlSa9WNdFFpcQC8M@asE6XplkiMY40YmhpR0evXvA/6ZIsK0iSRWidzq0PPK0VNo1tbH6Vunr/nwfyTWTueX7JyejyhOKTB2WSI1pWey8/mf4v9BerxRZavmLab/VYb1jXhS/ Try it online!]
<syntaxhighlight lang="raku"
return 2 if x.is-prime;
+flat (1 .. x.sqrt.floor).map: -> \d {
Line 1,572 ⟶ 1,961:
}
}
};</
<pre>First 20 terms of OEIS:A073916
Line 1,580 ⟶ 1,969:
Programming note: this REXX version has minor optimization, and all terms of the sequence are determined (found) in order.
===little optimization===
<
parse arg N . /*obtain optional argument from the CL.*/
if N=='' | N=="," then N= 15 /*Not specified? Then use the default.*/
Line 1,639 ⟶ 2,028:
do j=11 by 6 until j*j>#; if # // j==0 | # // (J+2)==0 then return 0
end /*j*/ /* ___ */
return 1 /*Exceeded √ # ? Then # is prime. */</
{{out|output|text= when using the input: <tt> 45 </tt>}}
Line 1,710 ⟶ 2,099:
This REXX version (unlike the 1<sup>st</sup> version), only goes through the numbers once, instead
of looking for numbers that have specific number of divisors.
<
parse arg N . /*obtain optional argument from the CL.*/
if N=='' | N=="," then N= 15 /*Not specified? Then use the default.*/
Line 1,786 ⟶ 2,175:
if #//(i+2)==0 then return 0
end /*i*/ /* ___ */
return 1 /*Exceeded √ # ? Then # is prime. */</
{{out|output|text= is identical to the 1<sup>st</sup> REXX version.}} <br><br>
=={{header|Ring}}==
<
load "stdlib.ring"
Line 1,819 ⟶ 2,208:
see nl + "done..." + nl
</syntaxhighlight>
{{out}}
<pre>
Line 1,844 ⟶ 2,233:
=={{header|Ruby}}==
{{trans|Java}}
<
return false if n < 2
return n == 2 if n % 2 == 0
Line 1,932 ⟶ 2,321:
print "A073916(", n, ") = ", OEISA073916(n), "\n"
n = n + 1
end</
{{out}}
<pre>A073916(1) = 1
Line 1,951 ⟶ 2,340:
=={{header|Sidef}}==
<
n.prime**(n-1)
}
Line 1,959 ⟶ 2,348:
}
say 20.of { f(_+1) }</
{{out}}
<pre>
Line 1,970 ⟶ 2,359:
{{libheader|Wren-big}}
{{libheader|Wren-fmt}}
<
import "./big" for BigInt
import "./fmt" for Fmt
var MAX = 33
Line 2,005 ⟶ 2,394:
}
}
}</
{{out}}
Line 2,051 ⟶ 2,440:
[[Extensible prime generator#zkl]] could be used instead.
<
p:=BI(1);
primes:=pmax.pump(List(0), p.nextPrime, "copy"); //-->(0,3,5,7,11,13,17,19,...)
Line 2,090 ⟶ 2,479:
}
}
}</
{{out}}
<pre>
| |||