laurent_r

These are some answers to the Week 169 of the Perl Weekly Challenge organized by Mohammad S. Anwar.

Spoiler Alert: This weekly challenge deadline is due in a couple of days from now (on June 19, 2022 at 23:59). This blog post offers some solutions to this challenge, please don’t read on if you intend to complete the challenge on your own.

Task 1: Brilliant Numbers

Write a script to generate first 20 Brilliant Numbers.

Brilliant numbers are numbers with two prime factors of the same length.

The number should have exactly two prime factors, i.e. it’s the product of two primes of the same length.

For example:

24287 = 149 x 163
24289 = 107 x 227

Therefore 24287 and 24289 are 2-brilliant numbers.

These two brilliant numbers happen to be consecutive as there are no even brilliant numbers greater than 14.

Output:

4, 6, 9, 10, 14, 15, 21, 25, 35, 49, 121, 143, 169, 187, 209, 221, 247, 253, 289, 299

There are essentially two ways to solve this type of problems: loop through all integers and filter out those that don’t satisfy the condition, and try to construct brilliant numbers from a list of primes. I prefer this second solution, because it requires much less computing resources, but, since brilliant numbers are not created in ascending order and the process can generate duplicates, we have to generate a few more numbers than what we strictly need to make sure that we obtain really the first 20 brilliant numbers.

Brilliant Numbers in Raku

In Raku, we’ll use the cross-product X operator to generate all pairs of primes. More precisely, we will use it as a metaoperator together with the * multiplication operator, X*, to generate directly the prime products. Since the primes have to have the same length, we’ll do it in two steps, one for primes of length 1 and one for primes of length 2. Since we need only a short list of small primes, I’ll simply hard code the list of the first ten primes. Task 2 of this challenge will show code to generate a longer list of primes when needed.

my @small-primes = 2, 3, 5, 7, 11, 13, 17, 19, 23, 29;
my @result = (@small-primes[0..3] X* @small-primes[0..3]).sort.squish;
append @result, (@small-primes[4..9] X* @small-primes[4..9]).sort.squish;
say @result[0..19];

This script displays the following output:

$ raku ./brilliant.raku
(4 6 9 10 14 15 21 25 35 49 121 143 169 187 209 221 247 253 289 299)

Brilliant Numbers in Perl

This is essentially a port to Perl of the Raku program above, except that we need to implement our own combine subroutine to replace the cross-product operator. We store the result in a hash to remove any duplicate.

use strict;
use warnings;
use feature "say";

my @small_primes = (2, 3, 5, 7, 11, 13, 17, 19, 23, 29);

sub combine {
    my @primes = @_;
    my %part_result;
    for my $i (0..$#primes) {
        for my $j ($i..$#primes) {
            $part_result{$primes[$i] * $primes[$j]} = 1;
        }
    }
    return sort {$a <=> $b} keys %part_result;
}
my @result = combine @small_primes[0..3];
push @result, combine @small_primes[4..9];
say "@result[0..19]";

This script displays the following output:

$ perl ./brilliant.pl
4 6 9 10 14 15 21 25 35 49 121 143 169 187 209 221 247 253 289 299

Task 2: Achilles Numbers

Write a script to generate first 20 Achilles Numbers. Please checkout wikipedia for more information.

An Achilles number is a number that is powerful but imperfect (not a perfect power). Named after Achilles, a hero of the Trojan war, who was also powerful but imperfect.

A positive integer n is a powerful number if, for every prime factor p of n, p^2 is also a divisor.

A number is a perfect power if it has any integer roots (square root, cube root, etc.).

For example 36 factors to (2, 2, 3, 3) - every prime factor (2, 3) also has its square as a divisor (4, 9). But 36 has an integer square root, 6, so the number is a perfect power.

But 72 factors to (2, 2, 2, 3, 3); it similarly has 4 and 9 as divisors, but it has no integer roots. This is an Achilles number.

Output:

72, 108,  200,  288,  392,  432,  500,  648,  675,  800,  864, 968, 972, 1125, 1152, 1323, 1352, 1372, 1568, 1800

Our first step will be to factorize the input candidates. An Achilles number is a powerful number, which means that every prime appearing in its prime factorization must have a power of at least 2. For a powerful number to be imperfect, the greatest common divisor (GCD) of its powers must be 1.

Achilles Numbers in Raku

The prime-factorssubroutine returns a BagHash containing a list of prime numbers together with their frequency in the prime factorization of the input integer.

Then, we apply the following rules: * Any prime appearing in the factorization must have a power of at last 2; * There must be at least two distinct prime factors; * The greatest common divisor (GCD) of the powers must be 1.

Note that the Raku built-in gcd routine is an infix operator, which can therefore be used only with 2 operands. However, the [] metaoperator makes it possible to use it with a list of more than 2 operands.

my @primes = (2..1000).grep({.is-prime});

sub prime-factors (UInt $num-in) {
    my $factors = BagHash.new;
    my $num = $num-in;
    for @primes -> $div {
        while ($num %% $div) {
            $factors{$div}++;
            $num div= $div;
        }
        return $factors if $num == 1;
    }
    $factors{$num}++ unless $num == $num-in;
    return $factors;
}

my $count = 0;
for 1..Inf -> $n { 
    my @powers = (prime-factors $n).values;
    if @powers.none < 2 and @powers.elems > 1 and ([gcd] @powers) == 1 {
        say $n.fmt("%4d"), " -> ", @powers;
        $count++;
        last if $count >= 20
    }
}

This script displays the following output:

$ raku ./Achilles_nums.raku
  72 -> [3 2]
 108 -> [3 2]
 200 -> [2 3]
 288 -> [5 2]
 392 -> [2 3]
 432 -> [3 4]
 500 -> [2 3]
 648 -> [4 3]
 675 -> [3 2]
 800 -> [5 2]
 864 -> [5 3]
 968 -> [2 3]
 972 -> [2 5]
1125 -> [2 3]
1152 -> [2 7]
1323 -> [3 2]
1352 -> [2 3]
1372 -> [2 3]
1568 -> [2 5]
1800 -> [2 3 2]

Achilles Numbers in Perl

This is a port to Perl of the Raku program above, with the same rules as above. Since there is no built-in gcd function in Perl, we implement our own. When it receives more than two arguments, it calls itself recursively until it ends up with 2 arguments.

#!/usr/bin/perl
use strict;
use warnings;
use feature qw/say/;
use Data::Dumper;
use constant MAX => 500;

my @primes = (2, 3, 5, 7);

# Greatest common divisor of two or more integers
sub gcd {
    if (@_ > 2) {
        return gcd( gcd($_[0], $_[1]), @_[2..$#_]);
    } else {
        my ($i, $j) = sort { $a <=> $b } @_;
        while ($j) {
            ($i, $j) = ($j, $i % $j);
        }
        return $i;
    }
}

# Creating a hash of prime factors (as keys) with their powers (as values) 
sub prime_factors {
    my $num = shift;
    my $origin_num = $num;
    my %factors;
    for my $div (@primes) {
        while ($num % $div == 0) {
            $factors{$div}++;
            $num /= $div;
        }
        return %factors if $num == 1;
    }
    $factors{$num}++ unless $num == $origin_num;
    return %factors;
}


# Populating an array of primes up to MAX
my $current = 9;
while (1) {
    my $prime = 1;
    for my $i (@primes) {
        my $i_sq = $i * $i;
        last if $i_sq > $current;
        $prime = 0, last if $current % $i == 0;
    }
    push @primes, $current if $prime;;
    $current += 2;
    last if $current > MAX;
}

my $count = 0;
for my $n (1..MAX*MAX) { 
    my %factors =  prime_factors $n;
    my @powers = map $factors{$_}, keys %factors;
    next if grep { $_ < 2} @powers;
    next if scalar @powers < 2;
    if (gcd(@powers) == 1) {
        printf "%4d -> %s\n", $n,  join " ", @powers;
        $count++;
        last if $count >= 20
    }

This script displays the following output:

$ perl  ./Achilles_nums.pl
  72 -> 2 3
 108 -> 2 3
 200 -> 3 2
 288 -> 5 2
 392 -> 3 2
 432 -> 3 4
 500 -> 3 2
 648 -> 4 3
 675 -> 2 3
 800 -> 2 5
 864 -> 3 5
 968 -> 3 2
 972 -> 5 2
1125 -> 2 3
1152 -> 2 7
1323 -> 2 3
1352 -> 2 3
1372 -> 2 3
1568 -> 5 2
1800 -> 2 2 3

Wrapping up

The next week Perl Weekly Challenge will start soon. If you want to participate in this challenge, please check https://perlweeklychallenge.org/ and make sure you answer the challenge before 23:59 BST (British summer time) on June 26, 2022. And, please, also spread the word about the Perl Weekly Challenge if you can.