I'm working on modernizing Rosetta Code's infrastructure. Starting with communications. Please accept this time-limited open invite to RC's Slack.. --Michael Mol (talk) 20:59, 30 May 2020 (UTC)

Brilliant numbers

From Rosetta Code
Task
Brilliant numbers
You are encouraged to solve this task according to the task description, using any language you may know.

Brilliant numbers are a subset of semiprime numbers. Specifically, they are numbers that are the product of exactly two prime numbers that both have the same number of digits when expressed in base 10.

Brilliant numbers are useful in cryptography and when testing prime factoring algorithms.


E.G.
  • 3 × 3 (9) is a brilliant number.
  • 2 × 7 (14) is a brilliant number.
  • 113 × 691 (78083) is a brilliant number.
  • 2 × 31 (62) is semiprime, but is not a brilliant number (different number of digits in the two factors).


Task
  • Find and display the first 100 brilliant numbers.
  • For the orders of magnitude 1 through 6, find and show the first brilliant number greater than or equal to the order of magnitude, and, its position in the series (or the count of brilliant numbers up to that point).


Stretch
  • Continue for larger orders of magnitude.


See also


ALGOL 68[edit]

BEGIN # Find Brilliant numbers - semi-primes whose two prime factors have    #
# the same number of digits #
PR read "primes.incl.a68" PR # include prime utilities #
INT max prime = 1 010; # maximum prime we will consider #
# should be enough to find the first brilliant number > 10^6 #
[]BOOL prime = PRIMESIEVE max prime; # sieve of primes to max prime #
# construct a table of brilliant numbers #
[ 1 : max prime * max prime ]BOOL brilliant;
FOR n FROM LWB brilliant TO UPB brilliant DO brilliant[ n ] := FALSE OD;
# brilliant numbers where one of the fators is 2 #
brilliant[ 4 ] := TRUE;
FOR p FROM 3 BY 2 TO 7 DO brilliant[ 2 * p ] := TRUE OD;
# brilliant numbers where both factors are odd #
INT p start := 1, p end := 9;
WHILE pstart < max prime DO
FOR p FROM p start BY 2 TO p end DO
IF prime[ p ] THEN
brilliant[ p * p ] := TRUE;
FOR q FROM p + 2 BY 2 TO p end DO
IF prime[ q ] THEN
brilliant[ p * q ] := TRUE
FI
OD
FI
OD;
p start := p end + 2;
p end := ( ( p start - 1 ) * 10 ) - 1;
IF p end > max prime THEN p end := max prime FI
OD;
# show the first 100 brilliant numbers #
INT b count := 0;
FOR n TO UPB brilliant WHILE b count < 100 DO
IF brilliant[ n ] THEN
print( ( whole( n, -6 ) ) );
IF ( b count +:= 1 ) MOD 10 = 0 THEN print( ( newline ) ) FI
FI
OD;
# first brilliant number >= 10^n, n = 1, 2, ..., 6 #
b count := 0;
INT power of 10 := 10;
FOR n TO UPB brilliant DO
IF brilliant[ n ] THEN
b count +:= 1;
IF n >= power of 10 THEN
print( ( "First brilliant number >= ", whole( power of 10, -8 )
, ": " , whole( n, -8 )
, " at position " , whole( b count, -6 )
, newline
)
);
power of 10 *:= 10
FI
FI
OD
END
Output:
     4     6     9    10    14    15    21    25    35    49
   121   143   169   187   209   221   247   253   289   299
   319   323   341   361   377   391   403   407   437   451
   473   481   493   517   527   529   533   551   559   583
   589   611   629   649   667   671   689   697   703   713
   731   737   767   779   781   793   799   803   817   841
   851   869   871   893   899   901   913   923   943   949
   961   979   989  1003  1007  1027  1037  1067  1073  1079
  1081  1121  1139  1147  1157  1159  1189  1207  1219  1241
  1247  1261  1271  1273  1333  1343  1349  1357  1363  1369
First brilliant number >=       10:       10 at position      4
First brilliant number >=      100:      121 at position     11
First brilliant number >=     1000:     1003 at position     74
First brilliant number >=    10000:    10201 at position    242
First brilliant number >=   100000:   100013 at position   2505
First brilliant number >=  1000000:  1018081 at position  10538

Arturo[edit]

brilliant?: function [x][
pf: factors.prime x
and? -> 2 = size pf
-> equal? size digits first pf
size digits last pf
]
 
brilliants: new []
i: 2
while [100 > size brilliants][
if brilliant? i -> 'brilliants ++ i
i: i + 1
]
 
print "First 100 brilliant numbers:"
loop split.every: 10 brilliants 'row [
print map to [:string] row 'item -> pad item 4
]
print ""
 
i: 4
nth: 0
order: 1
while [order =< 6] [
if brilliant? i [
nth: nth + 1
if i >= 10^order [
print ["First brilliant number >= 10 ^" order "is" i "at position" nth]
order: order + 1
]
]
 
i: i + 1
]
Output:
First 100 brilliant numbers:
   4    6    9   10   14   15   21   25   35   49 
 121  143  169  187  209  221  247  253  289  299 
 319  323  341  361  377  391  403  407  437  451 
 473  481  493  517  527  529  533  551  559  583 
 589  611  629  649  667  671  689  697  703  713 
 731  737  767  779  781  793  799  803  817  841 
 851  869  871  893  899  901  913  923  943  949 
 961  979  989 1003 1007 1027 1037 1067 1073 1079 
1081 1121 1139 1147 1157 1159 1189 1207 1219 1241 
1247 1261 1271 1273 1333 1343 1349 1357 1363 1369 

First brilliant number >= 10 ^ 1 is 10 at position 4 
First brilliant number >= 10 ^ 2 is 121 at position 11 
First brilliant number >= 10 ^ 3 is 1003 at position 74 
First brilliant number >= 10 ^ 4 is 10201 at position 242 
First brilliant number >= 10 ^ 5 is 100013 at position 2505 
First brilliant number >= 10 ^ 6 is 1018081 at position 10538

C++[edit]

Library: Primesieve
#include <algorithm>
#include <chrono>
#include <iomanip>
#include <iostream>
#include <locale>
#include <vector>
 
#include <primesieve.hpp>
 
auto get_primes_by_digits(uint64_t limit) {
primesieve::iterator pi;
std::vector<std::vector<uint64_t>> primes_by_digits;
std::vector<uint64_t> primes;
for (uint64_t p = 10; p < limit;) {
uint64_t prime = pi.next_prime();
if (prime > p) {
primes_by_digits.push_back(std::move(primes));
p *= 10;
}
primes.push_back(prime);
}
return primes_by_digits;
}
 
int main() {
std::cout.imbue(std::locale(""));
 
auto start = std::chrono::high_resolution_clock::now();
 
auto primes_by_digits = get_primes_by_digits(1000000000);
 
std::cout << "First 100 brilliant numbers:\n";
std::vector<uint64_t> brilliant_numbers;
for (const auto& primes : primes_by_digits) {
for (auto i = primes.begin(); i != primes.end(); ++i)
for (auto j = i; j != primes.end(); ++j)
brilliant_numbers.push_back(*i * *j);
if (brilliant_numbers.size() >= 100)
break;
}
std::sort(brilliant_numbers.begin(), brilliant_numbers.end());
for (size_t i = 0; i < 100; ++i) {
std::cout << std::setw(5) << brilliant_numbers[i]
<< ((i + 1) % 10 == 0 ? '\n' : ' ');
}
 
std::cout << '\n';
uint64_t power = 10;
size_t count = 0;
for (size_t p = 1; p < 2 * primes_by_digits.size(); ++p) {
const auto& primes = primes_by_digits[p / 2];
size_t position = count + 1;
uint64_t min_product = 0;
for (auto i = primes.begin(); i != primes.end(); ++i) {
uint64_t p1 = *i;
auto j = std::lower_bound(i, primes.end(), (power + p1 - 1) / p1);
if (j != primes.end()) {
uint64_t p2 = *j;
uint64_t product = p1 * p2;
if (min_product == 0 || product < min_product)
min_product = product;
position += std::distance(i, j);
if (p1 >= p2)
break;
}
}
std::cout << "First brilliant number >= 10^" << p << " is "
<< min_product << " at position " << position << '\n';
power *= 10;
if (p % 2 == 1) {
size_t size = primes.size();
count += size * (size + 1) / 2;
}
}
 
auto end = std::chrono::high_resolution_clock::now();
std::chrono::duration<double> duration(end - start);
std::cout << "\nElapsed time: " << duration.count() << " seconds\n";
}
Output:
First 100 brilliant numbers:
    4     6     9    10    14    15    21    25    35    49
  121   143   169   187   209   221   247   253   289   299
  319   323   341   361   377   391   403   407   437   451
  473   481   493   517   527   529   533   551   559   583
  589   611   629   649   667   671   689   697   703   713
  731   737   767   779   781   793   799   803   817   841
  851   869   871   893   899   901   913   923   943   949
  961   979   989 1,003 1,007 1,027 1,037 1,067 1,073 1,079
1,081 1,121 1,139 1,147 1,157 1,159 1,189 1,207 1,219 1,241
1,247 1,261 1,271 1,273 1,333 1,343 1,349 1,357 1,363 1,369

First brilliant number >= 10^1 is 10 at position 4
First brilliant number >= 10^2 is 121 at position 11
First brilliant number >= 10^3 is 1,003 at position 74
First brilliant number >= 10^4 is 10,201 at position 242
First brilliant number >= 10^5 is 100,013 at position 2,505
First brilliant number >= 10^6 is 1,018,081 at position 10,538
First brilliant number >= 10^7 is 10,000,043 at position 124,364
First brilliant number >= 10^8 is 100,140,049 at position 573,929
First brilliant number >= 10^9 is 1,000,000,081 at position 7,407,841
First brilliant number >= 10^10 is 10,000,600,009 at position 35,547,995
First brilliant number >= 10^11 is 100,000,000,147 at position 491,316,167
First brilliant number >= 10^12 is 1,000,006,000,009 at position 2,409,600,866
First brilliant number >= 10^13 is 10,000,000,000,073 at position 34,896,253,010
First brilliant number >= 10^14 is 100,000,380,000,361 at position 174,155,363,187
First brilliant number >= 10^15 is 1,000,000,000,000,003 at position 2,601,913,448,897

Elapsed time: 0.167788 seconds

Go[edit]

Translation of: Wren
Library: Go-rcu
package main
 
import (
"fmt"
"math"
"rcu"
"sort"
)
 
var primes = rcu.Primes(1e8 - 1)
 
type res struct {
bc interface{}
next int
}
 
func getBrilliant(digits, limit int, countOnly bool) res {
var brilliant []int
count := 0
pow := 1
next := math.MaxInt
for k := 1; k <= digits; k++ {
var s []int
for _, p := range primes {
if p >= pow*10 {
break
}
if p > pow {
s = append(s, p)
}
}
for i := 0; i < len(s); i++ {
for j := i; j < len(s); j++ {
prod := s[i] * s[j]
if prod < limit {
if countOnly {
count++
} else {
brilliant = append(brilliant, prod)
}
} else {
if next > prod {
next = prod
}
break
}
}
}
pow *= 10
}
if countOnly {
return res{count, next}
}
return res{brilliant, next}
}
 
func main() {
fmt.Println("First 100 brilliant numbers:")
brilliant := getBrilliant(2, 10000, false).bc.([]int)
sort.Ints(brilliant)
brilliant = brilliant[0:100]
for i := 0; i < len(brilliant); i++ {
fmt.Printf("%4d ", brilliant[i])
if (i+1)%10 == 0 {
fmt.Println()
}
}
fmt.Println()
for k := 1; k <= 13; k++ {
limit := int(math.Pow(10, float64(k)))
r := getBrilliant(k, limit, true)
total := r.bc.(int)
next := r.next
climit := rcu.Commatize(limit)
ctotal := rcu.Commatize(total + 1)
cnext := rcu.Commatize(next)
fmt.Printf("First >= %18s is %14s in the series: %18s\n", climit, ctotal, cnext)
}
}
Output:
First 100 brilliant numbers:
   4    6    9   10   14   15   21   25   35   49 
 121  143  169  187  209  221  247  253  289  299 
 319  323  341  361  377  391  403  407  437  451 
 473  481  493  517  527  529  533  551  559  583 
 589  611  629  649  667  671  689  697  703  713 
 731  737  767  779  781  793  799  803  817  841 
 851  869  871  893  899  901  913  923  943  949 
 961  979  989 1003 1007 1027 1037 1067 1073 1079 
1081 1121 1139 1147 1157 1159 1189 1207 1219 1241 
1247 1261 1271 1273 1333 1343 1349 1357 1363 1369 

First >=                 10 is              4 in the series:                 10
First >=                100 is             11 in the series:                121
First >=              1,000 is             74 in the series:              1,003
First >=             10,000 is            242 in the series:             10,201
First >=            100,000 is          2,505 in the series:            100,013
First >=          1,000,000 is         10,538 in the series:          1,018,081
First >=         10,000,000 is        124,364 in the series:         10,000,043
First >=        100,000,000 is        573,929 in the series:        100,140,049
First >=      1,000,000,000 is      7,407,841 in the series:      1,000,000,081
First >=     10,000,000,000 is     35,547,995 in the series:     10,000,600,009
First >=    100,000,000,000 is    491,316,167 in the series:    100,000,000,147
First >=  1,000,000,000,000 is  2,409,600,866 in the series:  1,000,006,000,009
First >= 10,000,000,000,000 is 34,896,253,010 in the series: 10,000,000,000,073

J[edit]

oprimes=: {{ NB. all primes of order y
p:(+i.)/-/\ p:inv +/\1 9*10^y
}}
 
obrill=: {{ NB. all brilliant numbers of order y primes
~.,*/~oprimes y
}}
 
brillseq=: {{ NB. sequences of brilliant numbers up through order y-1 primes
/:~;obrill each i.y
}}

Task examples:

   10 10 $brillseq 2
4 6 9 10 14 15 21 25 35 49
121 143 169 187 209 221 247 253 289 299
319 323 341 361 377 391 403 407 437 451
473 481 493 517 527 529 533 551 559 583
589 611 629 649 667 671 689 697 703 713
731 737 767 779 781 793 799 803 817 841
851 869 871 893 899 901 913 923 943 949
961 979 989 1003 1007 1027 1037 1067 1073 1079
1081 1121 1139 1147 1157 1159 1189 1207 1219 1241
1247 1261 1271 1273 1333 1343 1349 1357 1363 1369
NB. order, index, value
(brillseq 4) (],.(I. 10^]) ([,.{) [) 1 2 3 4 5 6
1 3 10
2 10 121
3 73 1003
4 241 10201
5 2504 100013
6 10537 1018081

Stretch goal (results are order, index, value):

   (brillseq 4) (],.(I. 10^]) ([,.{) [) ,7
7 124363 10000043
(brillseq 5) (],.(I. 10^]) ([,.{) [) 8 9
8 573928 100140049
9 7407840 1000000081

Julia[edit]

 
using Primes
 
function isbrilliant(n)
p = factor(n).pe
return (length(p) == 1 && p[1][2] == 2) ||
length(p) == 2 && ndigits(p[1][1]) == ndigits(p[2][1]) && p[1][2] == p[2][2] == 1
end
 
function testbrilliants()
println("First 100 brilliant numbers:")
foreach(p -> print(lpad(p[2], 5), p[1] % 20 == 0 ? "\n" : ""),
enumerate(filter(isbrilliant, 1:1370)))
bcount, results, positions = 0, zeros(Int, 9), zeros(Int, 9)
for n in 1:10^10
if isbrilliant(n)
bcount += 1
for i in 1:9
if n >= 10^i && results[i] == 0
results[i] = n
positions[i] = bcount
println("First >=", lpad(10^i, 12), " is", lpad(bcount, 8),
" in the series: $n")
end
end
end
end
return results, positions
end
 
testbrilliants()
 
Output:
First 100 brilliant numbers:
    4    6    9   10   14   15   21   25   35   49  121  143  169  187  209  221  247  253  289  299
  319  323  341  361  377  391  403  407  437  451  473  481  493  517  527  529  533  551  559  583
  589  611  629  649  667  671  689  697  703  713  731  737  767  779  781  793  799  803  817  841
  851  869  871  893  899  901  913  923  943  949  961  979  989 1003 1007 1027 1037 1067 1073 1079
 1081 1121 1139 1147 1157 1159 1189 1207 1219 1241 1247 1261 1271 1273 1333 1343 1349 1357 1363 1369
First >=          10 is       4 in the series: 10
First >=         100 is      11 in the series: 121
First >=        1000 is      74 in the series: 1003
First >=       10000 is     242 in the series: 10201
First >=      100000 is    2505 in the series: 100013
First >=     1000000 is   10538 in the series: 1018081
First >=    10000000 is  124364 in the series: 10000043
First >=   100000000 is  573929 in the series: 100140049
First >=  1000000000 is 7407841 in the series: 1000000081

Perl[edit]

Library: ntheory
use strict;
use warnings;
use feature 'say';
use List::AllUtils <max head firstidx uniqint>;
use ntheory <primes is_semiprime forsetproduct>;
 
sub table { my $t = shift() * (my $c = 1 + length max @_); ( sprintf( ('%'.$c.'d')x@_, @_) ) =~ s/.{1,$t}\K/\n/gr }
sub comma { reverse ((reverse shift) =~ s/(.{3})/$1,/gr) =~ s/^,//r }
 
my(@B,@Br);
for my $oom (1..5) {
my @P = grep { $oom == length } @{primes(10**$oom)};
forsetproduct { is_semiprime($_[0] * $_[1]) and push @B, $_[0] * $_[1] } \@P, \@P;
@Br = uniqint sort { $a <=> $b } @Br, @B;
}
 
say "First 100 brilliant numbers:\n" . table 10, head 100, @Br;
 
for my $oom (1..9) {
my $key = firstidx { $_ > 10**$oom } @Br;
printf "First >= %13s is position %9s in the series: %13s\n", comma(10**$oom), comma($key), comma $Br[$key];
}
Output:
First 100 brilliant numbers:
    4    6    9   10   14   15   21   25   35   49
  121  143  169  187  209  221  247  253  289  299
  319  323  341  361  377  391  403  407  437  451
  473  481  493  517  527  529  533  551  559  583
  589  611  629  649  667  671  689  697  703  713
  731  737  767  779  781  793  799  803  817  841
  851  869  871  893  899  901  913  923  943  949
  961  979  989 1003 1007 1027 1037 1067 1073 1079
 1081 1121 1139 1147 1157 1159 1189 1207 1219 1241
 1247 1261 1271 1273 1333 1343 1349 1357 1363 1369

First >=            10 is position         4 in the series:            14
First >=           100 is position        10 in the series:           121
First >=         1,000 is position        73 in the series:         1,003
First >=        10,000 is position       241 in the series:        10,201
First >=       100,000 is position     2,504 in the series:       100,013
First >=     1,000,000 is position    10,537 in the series:     1,018,081
First >=    10,000,000 is position    124364 in the series:    10,000,043
First >=   100,000,000 is position    573929 in the series:   100,140,049
First >= 1,000,000,000 is position   7407841 in the series: 1,000,000,081

Faster approach[edit]

Translation of: Sidef
use 5.020;
use strict;
use warnings;
 
use ntheory qw(:all);
use experimental qw(signatures);
 
sub is_briliant_number ($n) {
is_semiprime($n) || return;
my @f = factor($n);
length($f[0]) == length($f[1]);
}
 
sub next_brilliant_number ($n) {
++$n while not is_briliant_number($n);
$n;
}
 
sub brilliant_numbers_count ($n) {
 
use integer;
 
my $count = 0;
my $len = length(sqrtint($n));
 
foreach my $k (1 .. $len - 1) {
my $pi = prime_count(10**($k - 1), 10**$k - 1);
$count += binomial($pi, 2) + $pi;
}
 
my $min = 10**($len - 1);
my $max = 10**$len - 1;
 
my $pi_min = prime_count($min);
my $pi_max = prime_count($max);
 
my $j = -1;
 
forprimes {
if ($_*$_ <= $n) {
$count += (($max <= $n/$_) ? $pi_max : prime_count($n/$_)) - $pi_min - ++$j;
}
else {
lastfor;
}
} $min, $max;
 
return $count;
}
 
say "First 100 brilliant numbers:";
 
my @nums;
 
for (my $k = 1 ; scalar(@nums) < 100 ; ++$k) {
push(@nums, $k) if is_briliant_number($k);
}
 
while (@nums) {
my @slice = splice(@nums, 0, 10);
say join ' ', map { sprintf("%4s", $_) } @slice;
}
 
say '';
 
foreach my $n (1 .. 13) {
my $v = next_brilliant_number(vecprod((10) x $n));
printf("First brilliant number >= 10^%d is %s", $n, $v);
printf(" at position %s\n", brilliant_numbers_count($v));
}
Output:
First 100 brilliant numbers:
   4    6    9   10   14   15   21   25   35   49
 121  143  169  187  209  221  247  253  289  299
 319  323  341  361  377  391  403  407  437  451
 473  481  493  517  527  529  533  551  559  583
 589  611  629  649  667  671  689  697  703  713
 731  737  767  779  781  793  799  803  817  841
 851  869  871  893  899  901  913  923  943  949
 961  979  989 1003 1007 1027 1037 1067 1073 1079
1081 1121 1139 1147 1157 1159 1189 1207 1219 1241
1247 1261 1271 1273 1333 1343 1349 1357 1363 1369

First brilliant number >= 10^1 is 10 at position 4
First brilliant number >= 10^2 is 121 at position 11
First brilliant number >= 10^3 is 1003 at position 74
First brilliant number >= 10^4 is 10201 at position 242
First brilliant number >= 10^5 is 100013 at position 2505
First brilliant number >= 10^6 is 1018081 at position 10538
First brilliant number >= 10^7 is 10000043 at position 124364
First brilliant number >= 10^8 is 100140049 at position 573929
First brilliant number >= 10^9 is 1000000081 at position 7407841
First brilliant number >= 10^10 is 10000600009 at position 35547995
First brilliant number >= 10^11 is 100000000147 at position 491316167
First brilliant number >= 10^12 is 1000006000009 at position 2409600866
First brilliant number >= 10^13 is 10000000000073 at position 34896253010

Phix[edit]

Translation of: C++
Library: Phix/online

Replaced with C++ translation; much faster and now goes comfortably to 1e15 even on 32 bit. You can run this online here.

--
-- demo\rosetta\BrilliantNumbers.exw
-- =================================
--
with javascript_semantics
requires("1.0.2") -- (for in)
atom t0 = time()

function get_primes_by_digits(integer limit)
    sequence primes = get_primes_le(power(10,limit)),
             primes_by_digits = {}
    integer p = 10
    while length(primes) do
        integer pi = abs(binary_search(p,primes))-1
        primes_by_digits &= {primes[1..pi]}
        primes = primes[pi+1..$]
        p*= 10
    end while
    return primes_by_digits
end function
sequence primes_by_digits = get_primes_by_digits(8)
 
procedure first100()
    sequence brilliant_numbers = {}
    for primes in primes_by_digits do
        for i=1 to length(primes) do
--see talk page
--          for j=i to length(primes) do 
            for j=1 to i do 
                brilliant_numbers &= primes[i]*primes[j]
            end for
        end for
        if length(brilliant_numbers)>=100 then exit end if
    end for
    brilliant_numbers = sort(brilliant_numbers)[1..100]
    sequence j100 = join_by(brilliant_numbers,1,10," ","\n","%,5d")
    printf(1,"First 100 brilliant numbers:\n%s\n\n",{j100})
end procedure
first100()

atom pwr = 10, count = 0
for p=1 to 2*length(primes_by_digits)-1 do
    sequence primes = primes_by_digits[floor(p/2)+1]
    atom pos = count+1,
         min_product = 0
    for i=1 to length(primes) do
        integer p1 = primes[i],
                j = abs(binary_search(floor((pwr+p1-1)/p1),primes,i))
        if j<=length(primes) then -- (always is, I think)
            integer p2 = primes[j]
            atom prod = p1*p2
            if min_product=0 or prod<min_product then
                min_product = prod
            end if          
            pos += j-i
            if p1>=p2 then exit end if
        end if
    end for
    printf(1,"First brilliant number >= 10^%d is %,d at position %,d\n", {p, min_product, pos})
    pwr *= 10;
    if odd(p) then
        integer size = length(primes)
        count += size * (size + 1) / 2;
    end if
end for
?elapsed(time()-t0)
{} = wait_key()
Output:
First 100 brilliant numbers:
    4     6     9    10    14    15    21    25    35    49
  121   143   169   187   209   221   247   253   289   299
  319   323   341   361   377   391   403   407   437   451
  473   481   493   517   527   529   533   551   559   583
  589   611   629   649   667   671   689   697   703   713
  731   737   767   779   781   793   799   803   817   841
  851   869   871   893   899   901   913   923   943   949
  961   979   989 1,003 1,007 1,027 1,037 1,067 1,073 1,079
1,081 1,121 1,139 1,147 1,157 1,159 1,189 1,207 1,219 1,241
1,247 1,261 1,271 1,273 1,333 1,343 1,349 1,357 1,363 1,369


First brilliant number >= 10^1 is 10 at position 4
First brilliant number >= 10^2 is 121 at position 11
First brilliant number >= 10^3 is 1,003 at position 74
First brilliant number >= 10^4 is 10,201 at position 242
First brilliant number >= 10^5 is 100,013 at position 2,505
First brilliant number >= 10^6 is 1,018,081 at position 10,538
First brilliant number >= 10^7 is 10,000,043 at position 124,364
First brilliant number >= 10^8 is 100,140,049 at position 573,929
First brilliant number >= 10^9 is 1,000,000,081 at position 7,407,841
First brilliant number >= 10^10 is 10,000,600,009 at position 35,547,995
First brilliant number >= 10^11 is 100,000,000,147 at position 491,316,167
First brilliant number >= 10^12 is 1,000,006,000,009 at position 2,409,600,866
First brilliant number >= 10^13 is 10,000,000,000,073 at position 34,896,253,010
First brilliant number >= 10^14 is 100,000,380,000,361 at position 174,155,363,187
First brilliant number >= 10^15 is 1,000,000,000,000,003 at position 2,601,913,448,897
"3.3s"

Raku[edit]

1 through 7 are fast. 8 and 9 take a bit longer.

use Lingua::EN::Numbers;
 
# Find an abundance of primes to use to generate brilliants
my %primes = (2..100000).grep( &is-prime ).categorize: { .chars };
 
# Generate brilliant numbers
my @brilliant = lazy flat (1..*).map: -> $digits {
sort flat (^%primes{$digits}).race.map: { %primes{$digits}[$_] X× (flat %primes{$digits}[$_ .. *]) }
};
 
# The task
put "First 100 brilliant numbers:\n" ~ @brilliant[^100].batch(10)».fmt("%4d").join("\n") ~ "\n" ;
 
for 1 .. 7 -> $oom {
my $threshold = exp $oom, 10;
my $key = @brilliant.first: :k, * >= $threshold;
printf "First >= %13s is %9s in the series: %13s\n", comma($threshold), ordinal-digit(1 + $key, :u), comma @brilliant[$key];
}
Output:
First 100 brilliant numbers:
   4    6    9   10   14   15   21   25   35   49
 121  143  169  187  209  221  247  253  289  299
 319  323  341  361  377  391  403  407  437  451
 473  481  493  517  527  529  533  551  559  583
 589  611  629  649  667  671  689  697  703  713
 731  737  767  779  781  793  799  803  817  841
 851  869  871  893  899  901  913  923  943  949
 961  979  989 1003 1007 1027 1037 1067 1073 1079
1081 1121 1139 1147 1157 1159 1189 1207 1219 1241
1247 1261 1271 1273 1333 1343 1349 1357 1363 1369

First >=            10 is       4ᵗʰ in the series:            10
First >=           100 is      11ᵗʰ in the series:           121
First >=         1,000 is      74ᵗʰ in the series:         1,003
First >=        10,000 is     242ⁿᵈ in the series:        10,201
First >=       100,000 is    2505ᵗʰ in the series:       100,013
First >=     1,000,000 is   10538ᵗʰ in the series:     1,018,081
First >=    10,000,000 is  124364ᵗʰ in the series:    10,000,043
First >=   100,000,000 is  573929ᵗʰ in the series:   100,140,049
First >= 1,000,000,000 is 7407841ˢᵗ in the series: 1,000,000,081

Rust[edit]

Translation of: C++
// [dependencies]
// primal = "0.3"
// indexing = "0.4.1"
 
fn get_primes_by_digits(limit: usize) -> Vec<Vec<usize>> {
let mut primes_by_digits = Vec::new();
let mut power = 10;
let mut primes = Vec::new();
for prime in primal::Primes::all().take_while(|p| *p < limit) {
if prime > power {
primes_by_digits.push(primes);
primes = Vec::new();
power *= 10;
}
primes.push(prime);
}
primes_by_digits
}
 
fn main() {
use indexing::algorithms::lower_bound;
use std::time::Instant;
 
let start = Instant::now();
 
let primes_by_digits = get_primes_by_digits(1000000000);
 
println!("First 100 brilliant numbers:");
let mut brilliant_numbers = Vec::new();
for primes in &primes_by_digits {
for i in 0..primes.len() {
let p1 = primes[i];
for j in i..primes.len() {
let p2 = primes[j];
brilliant_numbers.push(p1 * p2);
}
}
if brilliant_numbers.len() >= 100 {
break;
}
}
brilliant_numbers.sort();
for i in 0..100 {
let n = brilliant_numbers[i];
print!("{:4}{}", n, if (i + 1) % 10 == 0 { '\n' } else { ' ' });
}
 
println!();
let mut power = 10;
let mut count = 0;
for p in 1..2 * primes_by_digits.len() {
let primes = &primes_by_digits[p / 2];
let mut position = count + 1;
let mut min_product = 0;
for i in 0..primes.len() {
let p1 = primes[i];
let n = (power + p1 - 1) / p1;
let j = lower_bound(&primes[i..], &n);
let p2 = primes[i + j];
let product = p1 * p2;
if min_product == 0 || product < min_product {
min_product = product;
}
position += j;
if p1 >= p2 {
break;
}
}
println!("First brilliant number >= 10^{p} is {min_product} at position {position}");
power *= 10;
if p % 2 == 1 {
let size = primes.len();
count += size * (size + 1) / 2;
}
}
 
let time = start.elapsed();
println!("\nElapsed time: {} milliseconds", time.as_millis());
}
Output:
First 100 brilliant numbers:
   4    6    9   10   14   15   21   25   35   49
 121  143  169  187  209  221  247  253  289  299
 319  323  341  361  377  391  403  407  437  451
 473  481  493  517  527  529  533  551  559  583
 589  611  629  649  667  671  689  697  703  713
 731  737  767  779  781  793  799  803  817  841
 851  869  871  893  899  901  913  923  943  949
 961  979  989 1003 1007 1027 1037 1067 1073 1079
1081 1121 1139 1147 1157 1159 1189 1207 1219 1241
1247 1261 1271 1273 1333 1343 1349 1357 1363 1369

First brilliant number >= 10^1 is 10 at position 4
First brilliant number >= 10^2 is 121 at position 11
First brilliant number >= 10^3 is 1003 at position 74
First brilliant number >= 10^4 is 10201 at position 242
First brilliant number >= 10^5 is 100013 at position 2505
First brilliant number >= 10^6 is 1018081 at position 10538
First brilliant number >= 10^7 is 10000043 at position 124364
First brilliant number >= 10^8 is 100140049 at position 573929
First brilliant number >= 10^9 is 1000000081 at position 7407841
First brilliant number >= 10^10 is 10000600009 at position 35547995
First brilliant number >= 10^11 is 100000000147 at position 491316167
First brilliant number >= 10^12 is 1000006000009 at position 2409600866
First brilliant number >= 10^13 is 10000000000073 at position 34896253010
First brilliant number >= 10^14 is 100000380000361 at position 174155363187
First brilliant number >= 10^15 is 1000000000000003 at position 2601913448897

Elapsed time: 676 milliseconds

Sidef[edit]

func is_briliant_number(n) {
n.is_semiprime && (n.factor.map{.len}.uniq.len == 1)
}
 
func brilliant_numbers_count(n) {
 
var count = 0
var len = n.isqrt.len
 
for k in (1 .. len-1) {
var pi = prime_count(10**(k-1), 10**k - 1)
count += binomial(pi, 2)+pi
}
 
var min = (10**(len - 1))
var max = (10**len - 1)
 
each_prime(min, max, {|p|
count += prime_count(p, max `min` idiv(n, p))
})
 
return count
}
 
say "First 100 brilliant numbers:"
 
100.by(is_briliant_number).each_slice(10, {|*a|
say a.map { '%4s' % _}.join(' ')
})
 
say ''
 
for n in (1 .. 12) {
var v = (10**n .. Inf -> first_by(is_briliant_number))
printf("First brilliant number >= 10^%d is %s", n, v)
printf(" at position %s\n", brilliant_numbers_count(v))
}
Output:
First 100 brilliant numbers:
   4    6    9   10   14   15   21   25   35   49
 121  143  169  187  209  221  247  253  289  299
 319  323  341  361  377  391  403  407  437  451
 473  481  493  517  527  529  533  551  559  583
 589  611  629  649  667  671  689  697  703  713
 731  737  767  779  781  793  799  803  817  841
 851  869  871  893  899  901  913  923  943  949
 961  979  989 1003 1007 1027 1037 1067 1073 1079
1081 1121 1139 1147 1157 1159 1189 1207 1219 1241
1247 1261 1271 1273 1333 1343 1349 1357 1363 1369

First brilliant number >= 10^1 is 10 at position 4
First brilliant number >= 10^2 is 121 at position 11
First brilliant number >= 10^3 is 1003 at position 74
First brilliant number >= 10^4 is 10201 at position 242
First brilliant number >= 10^5 is 100013 at position 2505
First brilliant number >= 10^6 is 1018081 at position 10538
First brilliant number >= 10^7 is 10000043 at position 124364
First brilliant number >= 10^8 is 100140049 at position 573929
First brilliant number >= 10^9 is 1000000081 at position 7407841
First brilliant number >= 10^10 is 10000600009 at position 35547995
First brilliant number >= 10^11 is 100000000147 at position 491316167
First brilliant number >= 10^12 is 1000006000009 at position 2409600866

Wren[edit]

Library: Wren-math
Library: Wren-seq
Library: Wren-fmt
import "./math" for Int
import "./seq" for Lst
import "./fmt" for Fmt
 
var primes = Int.primeSieve(1e7-1)
 
var getBrilliant = Fn.new { |digits, limit, countOnly|
var brilliant = []
var count = 0
var pow = 1
var next = Num.maxSafeInteger
for (k in 1..digits) {
var s = primes.where { |p| p > pow && p < pow * 10 }.toList
for (i in 0...s.count) {
for (j in i...s.count) {
var prod = s[i] * s[j]
if (prod < limit) {
if (countOnly) {
count = count + 1
} else {
brilliant.add(prod)
}
} else {
next = next.min(prod)
break
}
}
}
pow = pow * 10
}
return countOnly ? [count, next] : [brilliant, next]
}
 
System.print("First 100 brilliant numbers:")
var brilliant = getBrilliant.call(2, 10000, false)[0]
brilliant.sort()
brilliant = brilliant[0..99]
for (chunk in Lst.chunks(brilliant, 10)) Fmt.print("$4d", chunk)
System.print()
for (k in 1..12) {
var limit = 10.pow(k)
var res = getBrilliant.call(k, limit, true)
var total = res[0]
var next = res[1]
Fmt.print("First >= $,17d is $,15r in the series: $,17d", limit, total + 1, next)
}
Output:
First 100 brilliant numbers:
   4    6    9   10   14   15   21   25   35   49
 121  143  169  187  209  221  247  253  289  299
 319  323  341  361  377  391  403  407  437  451
 473  481  493  517  527  529  533  551  559  583
 589  611  629  649  667  671  689  697  703  713
 731  737  767  779  781  793  799  803  817  841
 851  869  871  893  899  901  913  923  943  949
 961  979  989 1003 1007 1027 1037 1067 1073 1079
1081 1121 1139 1147 1157 1159 1189 1207 1219 1241
1247 1261 1271 1273 1333 1343 1349 1357 1363 1369

First >=                10 is             4th in the series:                10
First >=               100 is            11th in the series:               121
First >=             1,000 is            74th in the series:             1,003
First >=            10,000 is           242nd in the series:            10,201
First >=           100,000 is         2,505th in the series:           100,013
First >=         1,000,000 is        10,538th in the series:         1,018,081
First >=        10,000,000 is       124,364th in the series:        10,000,043
First >=       100,000,000 is       573,929th in the series:       100,140,049
First >=     1,000,000,000 is     7,407,841st in the series:     1,000,000,081
First >=    10,000,000,000 is    35,547,995th in the series:    10,000,600,009
First >=   100,000,000,000 is   491,316,167th in the series:   100,000,000,147
First >= 1,000,000,000,000 is 2,409,600,866th in the series: 1,000,006,000,009

XPL0[edit]

 
func NumDigits(N); \Return number of digits in N
int N, Cnt;
[Cnt:= 0;
repeat N:= N/10;
Cnt:= Cnt+1;
until N = 0;
return Cnt;
];
 
func Brilliant(N); \Return 'true' if N is a brilliant number
int N, Limit, Cnt, F;
int A(3);
[Limit:= sqrt(N);
Cnt:= 0; F:= 2;
loop [if rem(N/F) = 0 then
[A(Cnt):= F;
Cnt:= Cnt+1;
if Cnt > 2 then quit;
N:= N/F;
]
else F:= F+1;
if F > N then quit;
if F > Limit then
[A(Cnt):= N;
Cnt:= Cnt+1;
quit;
];
];
if Cnt # 2 then return false;
return NumDigits(A(0)) = NumDigits(A(1));
];
 
int Cnt, N, Mag;
[Format(5, 0);
Cnt:= 0; N:= 4;
loop [if Brilliant(N) then
[RlOut(0, float(N));
Cnt:= Cnt+1;
if Cnt >= 100 then quit;
if rem(Cnt/10) = 0 then CrLf(0);
];
N:= N+1;
];
CrLf(0); CrLf(0);
Format(7, 0);
Cnt:= 0; N:= 4; Mag:= 10;
loop [if Brilliant(N) then
[Cnt:= Cnt+1;
if N >= Mag then
[Text(0, "First >= ");
RlOut(0, float(Mag));
Text(0, " is ");
RlOut(0, float(Cnt));
Text(0, " in series: ");
RlOut(0, float(N));
CrLf(0);
if Mag >= 1_000_000 then quit;
Mag:= Mag*10;
];
];
N:= N+1;
];
]
Output:
    4    6    9   10   14   15   21   25   35   49
  121  143  169  187  209  221  247  253  289  299
  319  323  341  361  377  391  403  407  437  451
  473  481  493  517  527  529  533  551  559  583
  589  611  629  649  667  671  689  697  703  713
  731  737  767  779  781  793  799  803  817  841
  851  869  871  893  899  901  913  923  943  949
  961  979  989 1003 1007 1027 1037 1067 1073 1079
 1081 1121 1139 1147 1157 1159 1189 1207 1219 1241
 1247 1261 1271 1273 1333 1343 1349 1357 1363 1369

First >=      10 is       4 in series:      10
First >=     100 is      11 in series:     121
First >=    1000 is      74 in series:    1003
First >=   10000 is     242 in series:   10201
First >=  100000 is    2505 in series:  100013
First >= 1000000 is   10538 in series: 1018081