Solve PWC170

3 files changed, 40 insertions, 0 deletions

diff --git a/challenge-170/wlmb/blog.txt b/challenge-170/wlmb/blog.txt
new file mode 100644
index 0000000000..7bdf7a144e
--- /dev/null
+++ b/challenge-170/wlmb/blog.txt
@@ -0,0 +1 @@
+https://wlmb.github.io/2022/06/20/PWC170/
diff --git a/challenge-170/wlmb/perl/ch-1.pl b/challenge-170/wlmb/perl/ch-1.pl
new file mode 100755
index 0000000000..d41c8a92f6
--- /dev/null
+++ b/challenge-170/wlmb/perl/ch-1.pl
@@ -0,0 +1,18 @@
+#!/usr/bin/env perl
+# Perl weekly challenge 170
+# Task 1: Primorial numbers from scratch
+#
+# See https://wlmb.github.io/2022/06/20/PWC170/#task-1-primorial-numbers
+use v5.12;
+use warnings;
+use PDL;
+use PDL::NiceSlice;
+die "Usage: ./ch-1.pl N\nto obtain the first N Primorial numbers" unless @ARGV;
+my $N=shift; # desired primorials
+my $M=$N>6?1+$N*(log($N)+log(log($N))):14; #upper bound on N-th prime
+my $sieve=ones($M); # large enough Eratosthenes sieve
+$sieve(0:1).=0; # 0 and 1 are not primes
+$sieve($_**2:-1:$_).=0 for(2..sqrt($M)); # all non-trivial multiples are not primes
+my $primes=sequence($M)->where($sieve); # primes correspond to non-zeroed positions in sieve
+say "P($_)=", $_==0?1:$primes(0:$_-1)->prodover # multiply first primes to obtain primorials
+    foreach(0..$N-1);
diff --git a/challenge-170/wlmb/perl/ch-2.pl b/challenge-170/wlmb/perl/ch-2.pl
new file mode 100755
index 0000000000..2b2e12da3a
--- /dev/null
+++ b/challenge-170/wlmb/perl/ch-2.pl
@@ -0,0 +1,21 @@
+#!/usr/bin/env perl
+# Perl weekly challenge 170
+# Task 2: Kronecker product
+#
+# See https://wlmb.github.io/2022/06/20/PWC170/#task-2-kronecker-product
+use v5.12;
+use warnings;
+use PDL;
+use PDL::NiceSlice;
+die "Usage: ./ch-2.pl A B\nto obtain the Kronecker product of A and B" unless @ARGV==2;
+# The input matrices should be written as an array of rows, each row
+# as an array of numbers and within quotes, as in  "[[1,2,3],[4,5,6]]"
+# for a 2x3 matrix
+my ($A, $B)=map {pdl $_} @ARGV;
+my $C=$A(*1,*1,:,:)*$B(:,:,*1,*1); #use dummy indices to build tensor
+				#product # C_{ijkl}=A_{ij}B_{kl}
+# get size of each dimension
+my ($I, $J, $K, $L)=($A->dim(1), $A->dim(0), $B->dim(1), $B->dim(0));
+my $Kronecker=$C->mv(1,2) # change indices to ikjl
+    ->reshape($J*$L, $I*$K); # clump indices i and k, and j and l
+say "The Kronecker product of $A and $B is $Kronecker";