In Raku, how does one calculate the sum of positive divisors from a prime factorization?

user_3134725 0
Viewed 133

In Raku, given a list of pairs (2 => 3, 3 => 2, 5 => 1, 7 => 4) ( representing the prime factorization of n = 2 3 · 3 2 · 5 1 · 7 4 ), how does construct a Raku expression for σ(n) = ( 2 0 + 2 1 + 2 2 + 2 3 ) · ( 3 0 + 3 1 + 3 2 ) · ( 5 0 + 5 1 ) · ( 7 0 + 7 1 + 7 2 + 7 3 + 7 4 ) ?

sub MAIN()
  {
  my $pairList = (2 => 3, 3 => 2, 5 => 1, 7 => 4) ;
  say '$pairList' ;
  say $pairList ;
  say $pairList.WHAT ;
  # Goal:
  #   from $pairList,
  #   the product (1 + 2 + 4 + 8) * (1 + 3 + 9) * (1 + 5) * (1 + 7 + 49 + 343 + 2401)
  #   = sigma ( 2^3 * 3^2 * 5^1 * 7^4 )
  } # end sub MAIN

Update 1

Based upon the answer of @raiph, the following program breaks the overall process into stages for the newcomer to Raku (such as me) …

sub MAIN()
  {
  my $pairList = (2 => 3, 3 => 2, 5 => 1, 7 => 4) ;
  say '$pairList' ;
  say $pairList ;
  say $pairList.WHAT ;

  # Goal:
  #   from $pairList,
  #   the product (1 + 2 + 4 + 8) * (1 + 3 + 9) * (1 + 5) * (1 + 7 + 49 + 343 + 2401)
  #   the product (15) * (13) * (6) * (2801)
  #   sigma ( 2^3 * 3^2 * 5^1 * 7^4 )
  #   3277170

  # Stage 1 : ((1 2 4 8) (1 3 9) (1 5) (1 7 49 343 2401))
  my $stage1 = $pairList.map: { (.key ** (my $++)) xx (.value + 1) } ;
  say '$stage1 : lists of powers' ;
  say $stage1 ;
  say $stage1.WHAT ;

  # Stage 2 : ((1 + 2 + 4 + 8) (1 + 3 + 9) (1 + 5) (1 + 7 + 49 + 343 + 2401))
  my $stage2 = $stage1.map: { sum $_ } ;
  say '$stage2 : sum each list' ;
  say $stage2 ;
  say $stage2.WHAT ;

  # Stage 3 : (1 + 2 + 4 + 8) * (1 + 3 + 9) * (1 + 5) * (1 + 7 + 49 + 343 + 2401)
  my $stage3 = $stage2.reduce( &infix:<*> ) ;
  say '$stage3 : product of list elements' ;
  say $stage3 ;
  say $stage3.WHAT ;
  } # end sub MAIN

A related post appears on Mathematics Stack Exchange.

Update 2

My original motivation had been to calculate aliquot sum s(n) = σ(n) - n. I found that prime factorization of each n is not necessary and seems inefficient. Raku and C++ programs calculating s(n) for n = 0 … 10 6 follow …

Raku
sub MAIN()
  {
  constant $limit = 1_000_000 ;

  my @s of Int = ( 1 xx ($limit + 1) ) ;
  @s[0] = 0 ;
  @s[1] = 0 ;

  loop ( my $column = 2; $column <= ($limit + 1) div 2; $column++ )
    {
    loop ( my $row = (2 * $column); $row <= $limit; $row += $column )
      {
      @s[$row] += $column ;
      } # end loop $row
    } # end loop $column

  say "s(", $limit, ") = ", @s[$limit] ; # s(1000000) = 1480437
  } # end sub MAIN
C++
(Observed to execute significantly faster than Raku)
#include <iostream>
#include <vector>
using namespace std ;

int main ( void )
  {
  const int LIMIT = 1000000 ;
  vector<int> s ( (LIMIT + 1), 1 ) ;
  s[0] = 0 ;
  s[1] = 0 ;
  for ( int col = 2 ; col <= (LIMIT + 1) / 2 ; col++ )
    for ( int row = (2 * col) ; row <= LIMIT ; row += col )
      s[row] += col ;
  cout << "s(" << LIMIT << ") = " << s[LIMIT] << endl ; // s(1000000) = 1480437
  } // end function main
key-value raku
2 Answers

There'll be bazillions of ways. I've ignored algorithmic efficiency. The first thing I wrote:

say [*] (2 => 3, 3 => 2, 5 => 1, 7 => 4) .map: { sum .key ** my $++ xx .value + 1 }

displays:

3277170

Explanation

 1 say
 2  [*]                       # `[op]` is a reduction. `[*] 6, 8, 9` is `432`.
 3    (2 => 3, 3 => 2, 5 => 1, 7 => 4)
 4      .map:
 5        {
 6          sum
 7            .key            # `.key` of `2 => 3` is `2`.
 8              **
 9                my          # `my` resets `$` for each call of enclosing `{...}`
10                  $++       # `$++` integer increments from `0` per thunk evaluation.
11                    xx      # `L xx R` forms list from `L` thunk evaluated `R` times
12                      .value + 1
13        }
user_10776720

It is unlikely that Raku is ever going to be faster than C++ for that kind of operation. It is still early in its life and there are lots of optimizations to be gained, but raw processing speed is not where it shines.

If you are trying to find the aliquot sum for all of the numbers in a continuous range then prime factorization is certainly less efficient than the method you arrived at. Sort of an inverse Sieve of Eratosthenes. There are a few things you could change to make it faster, though still probably much slower than C++

About twice as fast on my system:

constant $limit = 1_000_000;

my @s = 0,0;
@s.append: 1 xx $limit;

(2 .. $limit/2).race.map: -> $column {
    loop ( my $row = (2 * $column); $row <= $limit; $row += $column ) {
       @s[$row] += $column ;
    }
}

say "s(", $limit, ") = ", @s[$limit] ; # s(1000000) = 1480437

Where the prime factorization method really shines is for finding arbitrary aliquot sums.

This produces an answer in fractions of a second when the inverse sieve would likely take hours. Using the Prime::Factor module: from the Raku ecosystem

use Prime::Factor;
say "s(", 2**97-1, ") = ", (2**97-1).&proper-divisors.sum;
# s(158456325028528675187087900671) = 13842607235828485645777841
user_194800740
Related