Wolstenholme numbers
A Wolstenholme number is a number that is the fully reduced numerator of a second order harmonic number:
Wolstenholme numbers 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.
Prime Wolstenholme numbers are Wolstenholme numbers that are prime.
(Oddly enough, there is also a class of numbers called Wolstenholme primes that are not Wolstenholme numbers that are prime. This is not them.)
- Task
- Find and display the first 20 Wolstenholme numbers. (Or as many as reasonably supported by your language if it is fewer.)
- Find and display the first 4 prime Wolstenholme numbers. (Or as many as reasonably supported by your language if it is fewer.)
- Stretch
- Find and display the digit count of the 500th, 1000th, 2500th, 5000th, 10000th Wolstenholme numbers.
- Find and display the digit count of the first 15 prime Wolstenholme numbers.
- See also
C
#include <stdio.h>
#include <string.h>
#include <gmp.h>
#include <locale.h>
const char *ord(unsigned long c) {
int m = c % 100;
if (m >= 4 && m <= 20) return "th";
m %= 10;
return (m == 1) ? "st" :
(m == 2) ? "nd" :
(m == 3) ? "rd" : "th";
}
void abbreviate(char a[], const char *s) {
size_t len = strlen(s);
if (len < 40) {
strcpy(a, s);
return;
}
strncpy(a, s, 20);
strcpy(a + 20, "...");
strncpy(a + 23, s + len - 20, 21);
}
int main() {
int i, pc = 0, si = 0;
unsigned long k, l[5] = {500, 1000, 2500, 5000, 10000};
char *s, a[44];
mpq_t w, h;
mpz_t n, primes[15];
mpq_init(w);
mpq_init(h);
mpz_init(n);
for (i = 0; i < 15; ++i) mpz_init(primes[i]);
printf("Wolstenholme numbers:\n");
setlocale(LC_NUMERIC, "");
for (k = 1; k <= 10000; ++k) {
mpq_set_ui(h, 1, k * k);
mpq_add(w, w, h);
mpq_get_num(n, w);
if (pc < 15 && mpz_probab_prime_p(n, 15) > 0) mpz_set(primes[pc++], n);
if (k <= 20) {
s = mpz_get_str(NULL, 10, n);
printf("%6ld%s: %s\n", k, ord(k), s);
} else if (k == l[si]) {
s = mpz_get_str(NULL, 10, n);
abbreviate(a, s);
printf("%'6ld%s: %s (digits: %ld)\n", k, ord(k), a, strlen(s));
++si;
}
}
printf("\nPrime Wolstenholme numbers:\n");
for (i = 0; i < 15; ++i) {
s = mpz_get_str(NULL, 10, primes[i]);
if (i < 4) {
printf("%6d%s: %s\n", i+1, ord(i+1), s);
} else {
abbreviate(a, s);
printf("%'6d%s: %s (digits: %ld)\n", i+1, ord(i+1), a, strlen(s));
}
}
mpq_clear(w);
mpq_clear(h);
mpz_clear(n);
for (i = 0; i < 15; ++i) mpz_clear(primes[i]);
return 0;
}
- Output:
Same as Wren example.
FreeBASIC
#include "isprime.bas"
Function GCD(n As Uinteger, d As Uinteger) As Uinteger
Return Iif(d = 0, n, GCD(d, n Mod d))
End Function
Function numArmonico(n As Uinteger, m As Uinteger) As Uinteger
Dim As Integer i, v = 0, f = 1
For i = 1 To n
f *= i ^ m
Next i
For i = 1 To n
v += f / i ^ m
Next i
v /= GCD(v, f)
Return v
End Function
Dim As Integer i, wols(12)
Print "First 12 Wolstenholme numbers:"
For i = 1 To 12
wols(i) = numArmonico (i,2)
Print wols(i)
Next i
Print !"\nTwo Wolstenholme primes:"
For i = 1 To 12
If isPrime(wols(i)) Then Print wols(i)
Next i
Sleep- Output:
First 12 Wolstenholme numbers: 1 5 49 205 5269 5369 266681 1077749 9778141 1968329 239437889 240505109 Two Wolstenholme primes: 5 266681
J
NB. the first y members of the Wholstenhome sequence
WolstenholmeSeq=: {{ 2 {."1@x:+/\%*:1x+i.y }}
WolstenholmeSeq 20
1 5 49 205 5269 5369 266681 1077749 9778141 1968329 239437889 240505109 40799043101 40931552621 205234915681 822968714749 238357395880861 238820721143261 86364397717734821 17299975731542641
(#~ 1 p: ])WolstenholmeSeq 20
5 266681 40799043101 86364397717734821
Phix
include mpfr.e atom t0 = time(), t1 = t0+1 sequence primes = {}, l5 = {500, 1000, 2500, 5000, 10000} mpq {h,w} = mpq_inits(2) mpz n = mpz_init() printf(1,"Wolstenholme numbers:\n"); for k=1 to iff(platform()=JS?1000:10000) do mpq_set_si(h, 1, k*k) mpq_add(w, w, h) mpq_get_num(n, w) -- if length(primes)<15 and mpz_prime(n) then if length(primes)<iff(platform()=JS?8:10) and mpz_prime(n) then primes = append(primes,mpz_get_short_str(n)) end if if k<=20 or k=l5[1] then printf(1,"%6d%s: %s\n", {k, ord(k), mpz_get_short_str(n)}) if k>20 then l5 = l5[2..$] end if elsif time()>t1 then progress("checking k=%d, length(primes)=%d...\r",{k,length(primes)}) t1 = time()+1 end if end for progress("") printf(1,"\nPrime Wolstenholme numbers:\n"); for i,pw in primes do printf(1,"%6d%s: %s\n", {i, ord(i), pw}) end for ?elapsed(time()-t0)
- Output:
Same as C, but capped at 1000/8 primes when run in a browser, to keep it under 10s.
Capping it to the first 10 primes makes it 25* faster, and also keeps it under 10s on desktop/Phix.
Python
""" rosettacode.orgwiki/Magnanimous_numbers """
from fractions import Fraction
from sympy import isprime
def wolstenholme(k):
""" Wolstenholme numbers """
return sum(Fraction(1, i*i) for i in range(1, k+1)).numerator
def process_wolstenholmes(nmax):
""" Run the tasks at rosettacode.org/wiki/Wolstenholme_numbers """
wols = [wolstenholme(n) for n in range(1, nmax+1)]
print('First 20 Wolstenholme numbers:')
for i in wols[:20]:
print(i)
print('\nFour Wolstenholme primes:')
for j in filter(isprime, wols):
print(j)
process_wolstenholmes(20)
- Output:
First 20 Wolstenholme numbers: 1 5 49 205 5269 5369 266681 1077749 9778141 1968329 239437889 240505109 40799043101 40931552621 205234915681 822968714749 238357395880861 238820721143261 86364397717734821 17299975731542641 Four Wolstenholme primes: 5 266681 40799043101 86364397717734821
Raku
use Lingua::EN::Numbers;
sub abbr ($_) { .chars < 41 ?? $_ !! .substr(0,20) ~ '..' ~ .substr(*-20) ~ " (digits: {.chars})" }
my @wolstenholme = lazy ([\+] (1..∞).hyper.map: {FatRat.new: 1, .²}).map: *.numerator;
say 'Wolstenholme numbers:';
printf "%8s: %s\n", .&ordinal-digit(:c), abbr @wolstenholme[$_-1] for flat 1..20, 5e2, 1e3, 2.5e3, 5e3, 1e4;
say "\nPrime Wolstenholme numbers:";
printf "%8s: %s\n", .&ordinal-digit(:c), @wolstenholme.grep(&is-prime)[$_-1]».&abbr for 1..15;
- Output:
Wolstenholme numbers:
1st: 1
2nd: 5
3rd: 49
4th: 205
5th: 5269
6th: 5369
7th: 266681
8th: 1077749
9th: 9778141
10th: 1968329
11th: 239437889
12th: 240505109
13th: 40799043101
14th: 40931552621
15th: 205234915681
16th: 822968714749
17th: 238357395880861
18th: 238820721143261
19th: 86364397717734821
20th: 17299975731542641
500th: 40989667509417020364..48084984597965892703 (digits: 434)
1,000th: 83545938483149689478..99094240207812766449 (digits: 866)
2,500th: 64537911900230612090..12785535902976933153 (digits: 2164)
5,000th: 34472086597488537716..22525144829082590451 (digits: 4340)
10,000th: 54714423173933343999..49175649667700005717 (digits: 8693)
Prime Wolstenholme numbers:
1st: 5
2nd: 266681
3rd: 40799043101
4th: 86364397717734821
5th: 36190908596780862323..79995976006474252029 (digits: 104)
6th: 33427988094524601237..48446489305085140033 (digits: 156)
7th: 22812704758392002353..84405125167217413149 (digits: 218)
8th: 28347687473208792918..45794572911130248059 (digits: 318)
9th: 78440559440644426017..30422337523878698419 (digits: 520)
10th: 22706893975121925531..02173859396183964989 (digits: 649)
11th: 27310394808585898968..86311385662644388271 (digits: 935)
12th: 13001072736642048751..08635832246554146071 (digits: 984)
13th: 15086863305391456002..05367804007944918693 (digits: 1202)
14th: 23541935187269979100..02324742766220468879 (digits: 1518)
15th: 40306783143871607599..58901192511859288941 (digits: 1539)
Wren
import "./gmp" for Mpq
import "./fmt" for Fmt
var w = Mpq.new()
var h = Mpq.new()
var primes = []
var l = [500, 1000, 2500, 5000, 10000]
System.print("Wolstenholme numbers:")
for (k in 1..10000) {
h.set(1, k * k)
w.add(h)
var n = w.num
if (primes.count < 15 && n.probPrime(15) > 0) primes.add(n)
if (k <= 20) {
Fmt.print("$8r: $i", k, n)
} else if (l.contains(k)) {
var ns = n.toString
Fmt.print("$,8r: $20a (digits: $d)", k, ns, ns.count)
}
}
System.print("\nPrime Wolstenholme numbers:")
for (i in 0...primes.count) {
if (i < 4) {
Fmt.print("$8r: $i", i+1, primes[i])
} else {
var ps = primes[i].toString
Fmt.print("$8r: $20a (digits: $d)", i+1, ps, ps.count)
}
}- Output:
Wolstenholme numbers:
1st: 1
2nd: 5
3rd: 49
4th: 205
5th: 5269
6th: 5369
7th: 266681
8th: 1077749
9th: 9778141
10th: 1968329
11th: 239437889
12th: 240505109
13th: 40799043101
14th: 40931552621
15th: 205234915681
16th: 822968714749
17th: 238357395880861
18th: 238820721143261
19th: 86364397717734821
20th: 17299975731542641
500th: 40989667509417020364...48084984597965892703 (digits: 434)
1,000th: 83545938483149689478...99094240207812766449 (digits: 866)
2,500th: 64537911900230612090...12785535902976933153 (digits: 2164)
5,000th: 34472086597488537716...22525144829082590451 (digits: 4340)
10,000th: 54714423173933343999...49175649667700005717 (digits: 8693)
Prime Wolstenholme numbers:
1st: 5
2nd: 266681
3rd: 40799043101
4th: 86364397717734821
5th: 36190908596780862323...79995976006474252029 (digits: 104)
6th: 33427988094524601237...48446489305085140033 (digits: 156)
7th: 22812704758392002353...84405125167217413149 (digits: 218)
8th: 28347687473208792918...45794572911130248059 (digits: 318)
9th: 78440559440644426017...30422337523878698419 (digits: 520)
10th: 22706893975121925531...02173859396183964989 (digits: 649)
11th: 27310394808585898968...86311385662644388271 (digits: 935)
12th: 13001072736642048751...08635832246554146071 (digits: 984)
13th: 15086863305391456002...05367804007944918693 (digits: 1202)
14th: 23541935187269979100...02324742766220468879 (digits: 1518)
15th: 40306783143871607599...58901192511859288941 (digits: 1539)