# Perl Weekly Challenge 174.

My solutions (task 1, task 1a, task 2, task 2a, and task 2b ) to the The Weekly Challenge - 174.

# Task 1: Disarium Numbers

```
Submitted by: Mohammad S Anwar
Write a script to generate first 19 Disarium Numbers.
A disarium number is an integer where the sum of each digit
raised to the power of its position in the number, is equal to
the number.
For example,
518 is a disarium number as (5 ** 1) + (1 ** 2) + (8 ** 3) =>
5 + 1 + 512 => 518
```

I try a straightforward solution examining numbers one by one and testing them until I find enough Disarium numbers. I use PDL for conciseness.

```
perl -MPDL -E '
$N=shift; $n=$i=0; while($n<$N){$e=($p=long [split "", $i])->sequence+1;
++$n, say "$n: $i" if ($p**$e)->sumover==$i; ++$i }' 19
```

Results:

```
1: 0
2: 1
3: 2
4: 3
5: 4
6: 5
7: 6
8: 7
9: 8
10: 9
11: 89
12: 135
13: 175
14: 518
15: 598
16: 1306
17: 1676
18: 2427
19: 2646798
```

It seems to work fine and is fast for the first 18 desarium numbers, but the 19-th took about 2.5 minutes in my laptop. I’m not sure I would have the patience for waiting for the 20-th. The full code is

```
1 # Perl weekly challenge 174
2 # Task 1: Disarium numbers. Brute force
3 #
4 # See https://wlmb.github.io/2022/07/18/PWC174/#task-1-disarium-numbers
5 use v5.12;
6 use PDL;
7 say("Usage: $0 N\nto obtain the first N Disarium numbers"), exit unless @ARGV;
8 my $how_many=shift;
9 die "Expected a positive number" unless $how_many>0;
10 my ($counter, $test)=(0, 0);
11 while($counter<$how_many){
12 my $exponents=(my $digits=long [split "", $test])->sequence+1;
13 ++$counter, say "$counter: $test"
14 if ($digits**$exponents)->sumover==$test;
15 ++$test
16 }
```

Example:

```
./ch-1.pl 19
```

Results:

```
1: 0
2: 1
3: 2
4: 3
5: 4
6: 5
7: 6
8: 7
9: 8
10: 9
11: 89
12: 135
13: 175
14: 518
15: 598
16: 1306
17: 1676
18: 2427
19: 2646798
```

After a twitter exchange with Flavio Poletti `@polettix`

I learned
that one can do better than a brute force approach. I saw his Raku
solution PWC174 - Disarium Numbers - ETOOBUSY not quite understanding
it, but it inspired the following solution.

First, note that a number given by a sequence of N digits,
y_{1}y_{2}…y_{N} with a positive first digit is Disarium if
p(y_{1}y_{2}…y_{N})=d(y_{1}y_{2}…y_{N}), where
p(y_{1}y_{2}…y_{N})=∑_{n=1}^{N} y_{n}^{n} is the sum of powers of its digits and
d(y_{1}y_{2}…y_{N})=∑ 10^{N-n} y_{n} is its decimal value.

Notice that
p(x_{1}…x_{M}y_{M+1}…y_{N})=p(x_{1}…x_{M})+p(0_{1}…0_{M}y_{M+1}…y_{N})
and d(x_{1}…x_{M}y_{M+1}…y_{N})=10^{N-M}
d(x_{1}…x_{M})+d(y_{M+1}…y_{N}), where 0_{1}…0_{M} represents M
zeroes, one after another. Thus,
d(x_{1}…x_{M})=(p(0_{1}…0_{M}y_{M+1}…y_{N})-d(y_{M+1}…y_{N})+p(x_{1}…x_{M}))/10^{N-M}.

The smallest possible value for p(x_{1}…x_{M}) is 1, corresponding to
1_{1}0_{2}…0_{M}, while the largest possible value is ∑_{n=1}^{M}
9^{n}=(9^{M+1}-9)/8. Thus,
(p(0_{1}…0_{M}y_{M+1}…y_{N})-d(y_{M+1}…y_{N})+1)/10^{N-M}
≤ d(x_{1}…x_{M})
≤
(p(0_{1}…0_{M}y_{M+1}…y_{N})-d(y_{M+1}…y_{N})+(9^{M+1}-9)/8)/10^{N-M}.

The strategy is thus, choose a value for y_{N}, and find the bounds
for the rest of the digits x_{1}…x_{N-1}. If the bounds are
sufficiently tight, then explore all intermediate numbers looking for
Disarium numbers. If the bounds are too far apart, then choose also a value
for y_{N-1} and find bounds for x_{1}…x_{N-2}. Recourse until the
bounds are tight enough and then explore exhaustively the given
range. For each M explored, starting downwards from N-1 and stopping
when tight enough bounds are found, Choose all possible digits 0…9 for
Y_{M+1}.

These ideas lead to the following recursive code.

```
1 # Perl weekly challenge 174
2 # Task 1: Disarium numbers. Not so brute force
3 #
4 # See https://wlmb.github.io/2022/07/18/PWC174/#task-1-disarium-numbers
5 use v5.36;
6 use List::Util qw(sum min);
7 use Memoize;
8 memoize "power";
9 die "Usage: $0 N\nto obtain the first N Disarium numbers\n" unless @ARGV;
10 my $want=shift;
11 die "There are not that many Disarium numbers\n" if $want > 20;
12 die "We are not that patient today\n" if $want==20;
13 my $count=1;
14 say $count++, ": $_" for(0..min 9, $want);
15 for my $N(2..20){ # for each possible length
16 search($N, $N, "", 0, 0);
17 }
18 sub search ($N, $M, $rightmost, $p_old, $d_old){
19 --$M;
20 my $minmin=power(10,$M-1);
21 for my $right(0..9){ # for each possible rightmost digit
22 last if $count > $want;
23 my $p_new=power($right, $M+1)+$p_old;
24 my $d_new=$right*power(10, $N-$M-1)+$d_old;
25 my $min=($p_new-$d_new+1)/power(10, $N-$M);
26 my $max=($p_new-$d_new+(power(9,$M+1)-9)/8)/power(10, $N-$M);
27 $min=$minmin if $min<$minmin;
28 next if $max<=$min;
29 my $newright="$right$rightmost";
30 if($max-$min < 10){
31 for($min..$max){
32 my $candidate="$_$newright";
33 say $count++, ": $candidate" if disarium($candidate);
34 }
35 } else {
36 search($N, $M, $newright, $p_new, $d_new);
37 }
38 }
39 }
40 sub power($x, $y){ $x**$y }
41 sub disarium($candidate){
42 my @digits=split "", $candidate;
43 my $powersum=sum map {$digits[$_-1]**$_} (1..@digits);
44 return $powersum==$candidate;
45 }
```

I `memoize`

the `sub power`

to avoid recalculating the same powers
again and again. Line 10 deals with the 1 digit case. The recursive
subroutine `search`

takes as parameters the number of digits N, the
numbers of unknown digits M, a string with the rightmost N-M candidate
digits y_{M+1}…Y_{N}, and the corresponding values of p(y_{M+1}…Y_{N})
and d(y_{M+1}…Y_{N}). Then it decreases M and starts exploring values
of the new y_{M+1}, calling itself recursively if necessary.

Example:

```
./ch-1a.pl 19
```

Results:

```
1: 0
2: 1
3: 2
4: 3
5: 4
6: 5
7: 6
8: 7
9: 8
10: 9
11: 89
12: 135
13: 175
14: 518
15: 598
16: 1306
17: 1676
18: 2427
19: 2646798
```

So it works. The above run took 0.026s instead of 2.5m for the previous code, so it is about six thousand times faster.

To get the
last (20-th) Disarium number I would require `bigint`

and lots of
patience. Using `bigint`

and `Time::HiRes`

I timed the time required
for a few lengths (up to 14) and extrapolated. My estimate is that it
would take about 142 hours to get the 20-th Disarium number.

Curiously, the brute force solution but without using PDL is only about ten times slower for the first 19 Disarium numbers than the last, somewhat involved program, so it seems that in this case PDL is not so helpful.

# Task 2: Permutation Ranking

```
Submitted by: Mohammad S Anwar
You are given a list of integers with no duplicates, e.g. [0,
1, 2].
Write two functions, permutation2rank() which will take the
list and determine its rank (starting at 0) in the set of
possible permutations arranged in lexicographic order, and
rank2permutation() which will take the list and a rank number
and produce just that permutation.
Please checkout this post for more informations and
algorithm.
Given the list [0, 1, 2] the ordered permutations are:
0: [0, 1, 2]
1: [0, 2, 1]
2: [1, 0, 2]
3: [1, 2, 0]
4: [2, 0, 1]
5: [2, 1, 0]
and therefore:
permutation2rank([1, 0, 2]) = 2
rank2permutation([0, 1, 2], 1) = [0, 2, 1]
```

A straightforward solution may be easily constructed by producing all permutations in order, which yields the following 3-liner:

```
perl -MAlgorithm::Combinatorics=permutations -MPDL -MPDL::NiceSlice -E '
$p=pdl "[".shift."]"; @o=$p->qsort->list; $a=pdl[permutations(\@o)];
for(0..$a->dim(1)-1){say("rank=$_"), last if all($a(,$_)==$p);}
say "$_->",$a(:,($_)) foreach(@ARGV);
' "1 0 2" 0 2 4 1 3 5
```

Results:

```
rank=2
0->[0 1 2]
2->[1 0 2]
4->[2 0 1]
1->[0 2 1]
3->[1 2 0]
5->[2 1 0]
```

Similar example with a different set of numbers:

```
perl -MAlgorithm::Combinatorics=permutations -MPDL -MPDL::NiceSlice -E '
$p=pdl "[".shift."]"; @o=$p->qsort->list; $a=pdl[permutations(\@o)];
for(0..$a->dim(1)-1){say("rank=$_"), last if all($a(,$_)==$p);}
say "$_->",$a(:,($_)) foreach(@ARGV);
' "10 5 15" 0 2 4 1 3 5
```

Results:

```
rank=2
0->[5 10 15]
2->[10 5 15]
4->[15 5 10]
1->[5 15 10]
3->[10 15 5]
5->[15 10 5]
```

It works but for larger sets it becomes extremely inefficient, as the number of permutations grows as the factorial of the number of elements. Thus, for the full code I change to an algorithm based on the page The rank of a permutation, TryAlgo.

Consider a permutation P_{0}P_{1}…P_{n-1} of the ordered set
O_{0}O_{1}…O_{n-1}. There are
(n-1)! permutations that start with O_{0}, (n-2)! that start with O_{1}
and so on. Thus, the rank r(P_{0}P_{1}…P_{n-1}) is simply
r_{0}×(n-1)!+r(P_{1}…P_{n-1}), where r_{0} is the sequence number of
P_{0} in O_{0}O_{1}…O_{n-1}, counting up from 0, and r(P_{1}…P_{n-1})
is the rank of the permutation P_{1}…P_{n-1} of the ordered
set O_{0}O_{1}…O_{n-1} from which we have removed P_{0}. We
can proceed recursively to obtain finally
r(P_{0}P_{1}…P_{n-1})=r_{0}×(n-1)!+r_{1}×(n-2)!+…r_{n-1}×0!,
where r_{i} is the position of P_{i} in the set
O_{0}O_{1}…O_{n-1} from which we have removed P_{0},
P_{1}…P_{i-1}. Notice that r_{n-1}=0 as there is no choice
left for the last element.

Conversely, given a rank r, we can obtain r_{0} by dividing it
by (n-1)!, we can obtain r_{1} by dividing the residue by
(n-2)!, and so on. Then, we choose P_{0} as the r_{0}-th
element of the set O_{0}O_{1}…O_{n-1}, we remove P_{0} from
the set and choose P_{1} as the r_{1}-th element of the
remaining set, and so on. This algorithm is embodied in the
following code:

```
1 # Perl weekly challenge 174
2 # Task 2: Permutation ranking
3 #
4 # See https://wlmb.github.io/2022/07/18/PWC174/#task-2-permutation-ranking
5 use v5.36;
6 use PDL;
7 use PDL::NiceSlice;
8 my $usage=<<"END";
9 Usage: $0 "P1 P2...Pn" R1 R2...
10 to obtain rank of permutation P1 P2...Pn
11 and the rank-Ri permutations. Note quotations.
12 END
13 say($usage), exit unless @ARGV>0;
14 my $permutation=long "[".shift."]";
15 my $size=$permutation->nelem;
16 die "Elements should be distinct" unless $permutation->uniq->nelem==$size;
17 my @ordered=$permutation->qsort->list;
18 my %element2index = map { ($ordered[$_], $_) } (0..$size-1);
19 my $permuted_indices=long [@element2index{$permutation->list}];
20 my $factorials=factorials($size);
21 say "permutation2rank($permutation)=", permutation2rank($permuted_indices);
22 say "rank2permutation($permutation, $_)=", rank2permutation($_) for(@ARGV);
23
24 sub permutation2rank($permutation){ # ranks a permutation of 0...N-1
25 return $factorials->inner(ranks($permutation));
26 }
27 sub ranks($permutation){
28 my $cmp=$permutation(*1)>$permutation; #P_i>P_j
29 $cmp->inner($cmp->xvals>=$cmp->yvals); #r_i=sum_{j>=i}(P_i>P_j)
30 }
31 sub factorials($size){
32 my $f=1;
33 (long [1, map {$f=$_*$f} (1..$size-1)])->(-1:0);
34 }
35 sub rank2permutation($rank){
36 my @indices=map {my $index=floor($rank/$_); $rank%=$_; $index} $factorials->list;
37 my @copy=@ordered;
38 long [ map {splice @copy, $_, 1} @indices];
39 }
40
```

The input to the program is a list of distinct numbers between
quotation marks, that define the permutation followed by a
list of ranks to examine. I use PDL to simplify some
expressions. In line 17 I produce the ordered set O_{0}O_{1}…O{n-1}, which
in lines 18 and 19 I map to the simpler set 0,1…n-1, for which the
indices and correponding values coincide. Line 20
prepares a vector of decreasing factorials by calling the
subroutine `factorials`

. The subroutine `ranks`

calculates a
vector of partial ranks r_{i} by counting how many P_{j}’s
obey P_{j}>P_{i}, with i≥j. To that end, in line 28 we use a
dummy index to build a matrix with entries P_{i}>P_{j} at row
i column j. Then, we use an inner product in line 29 to select
and sum its upper right triangle, yielding
r_{i}=∑_{j>=i}(P_{i}>P_{j}). The subroutine `rank2permutation`

maps a rank r to a permutation by iteratively obtaining an array of
indices r_{i}=r/f_{i} updating the rank r→r%f_{i}. I use the
`splice`

routine in line 38 to chose and remove the r_{i}-th
term P_{i} of the set, and return a vector with
P_{0}P_{1}…P_{n-1}.

Example:

```
./ch-2.pl "1 0 2" 0 1 2 3 4 5
```

Results:

```
permutation2rank([1 0 2])=2
rank2permutation([1 0 2], 0)=[0 1 2]
rank2permutation([1 0 2], 1)=[0 2 1]
rank2permutation([1 0 2], 2)=[1 0 2]
rank2permutation([1 0 2], 3)=[1 2 0]
rank2permutation([1 0 2], 4)=[2 0 1]
rank2permutation([1 0 2], 5)=[2 1 0]
```

Another example, similar to above but with different numbers:

```
./ch-2.pl "10 9 20" 0 1 2 3 4 5
```

Results:

```
permutation2rank([10 9 20])=2
rank2permutation([10 9 20], 0)=[9 10 20]
rank2permutation([10 9 20], 1)=[9 20 10]
rank2permutation([10 9 20], 2)=[10 9 20]
rank2permutation([10 9 20], 3)=[10 20 9]
rank2permutation([10 9 20], 4)=[20 9 10]
rank2permutation([10 9 20], 5)=[20 10 9]
```

The run-time of this algorithm is proportional to the square of the number of elements, so it is much faster than the oneliner above.

Notice that by replacing line 17 by

```
my @ordered=sort {$a cmp $b} $permutation->list;
```

I would obtain a lexicographic order instead of a numerical order (see task 2a).

Example:

```
echo "Numerical order"
./ch-2.pl "22 111 3" 0 1 2 3 4 5
echo; echo "Lexicographical order (with modified line 17)"
./ch-2a.pl "22 111 3" 0 1 2 3 4 5
```

Results:

```
Numerical order
permutation2rank([22 111 3])=3
rank2permutation([22 111 3], 0)=[3 22 111]
rank2permutation([22 111 3], 1)=[3 111 22]
rank2permutation([22 111 3], 2)=[22 3 111]
rank2permutation([22 111 3], 3)=[22 111 3]
rank2permutation([22 111 3], 4)=[111 3 22]
rank2permutation([22 111 3], 5)=[111 22 3]
Lexicographical order (with modified line 17)
permutation2rank([22 111 3])=2
rank2permutation([22 111 3], 0)=[111 22 3]
rank2permutation([22 111 3], 1)=[111 3 22]
rank2permutation([22 111 3], 2)=[22 111 3]
rank2permutation([22 111 3], 3)=[22 3 111]
rank2permutation([22 111 3], 4)=[3 111 22]
rank2permutation([22 111 3], 5)=[3 22 111]
```

Furthermore, I can also generalize to non-numerical inputs by replacing some ndarrays by perl arrays (see task 2b).

```
./ch-2b.pl "BB AAA C" 0 1 2 3 4 5
```

Results:

```
permutation2rank([BB AAA C])=2
rank2permutation([BB AAA C], 0)=[AAA BB C]
rank2permutation([BB AAA C], 1)=[AAA C BB]
rank2permutation([BB AAA C], 2)=[BB AAA C]
rank2permutation([BB AAA C], 3)=[BB C AAA]
rank2permutation([BB AAA C], 4)=[C AAA BB]
rank2permutation([BB AAA C], 5)=[C BB AAA]
```