Solve PWC146

3 files changed, 80 insertions, 0 deletions

diff --git a/challenge-146/wlmb/blog.txt b/challenge-146/wlmb/blog.txt
new file mode 100644
index 0000000000..7919e7679b
--- /dev/null
+++ b/challenge-146/wlmb/blog.txt
@@ -0,0 +1 @@
+https://wlmb.github.io/2022/01/03/PWC146/
diff --git a/challenge-146/wlmb/perl/ch-1.pl b/challenge-146/wlmb/perl/ch-1.pl
new file mode 100755
index 0000000000..483b024a67
--- /dev/null
+++ b/challenge-146/wlmb/perl/ch-1.pl
@@ -0,0 +1,31 @@
+#!/usr/bin/env perl
+# Perl weekly challenge 146
+# Task 1: 10001st prime
+#
+# See https://wlmb.github.io/2022/01/03/PWC146/#task-1-10001st-prime
+use v5.12;
+use warnings;
+use PDL;
+use PDL::NiceSlice;
+for my $N(@ARGV?@ARGV:10001){
+    die "Argument should be positive" unless $N>=1;
+    # Estimate size $M of required sieve by solving $M/log($M) approx $N
+    # unless $N is too small
+    my $M=$N<4?6:find_zero(sub {my $x=shift; $N-$x/log($x)},
+			sub {my $l=log($_[0]); 1/$l**2-1/$l}, $N);
+    my $sieve=ones($M); # fill sieve with ones
+    $sieve(0:1).=0; # 0 and 1 are not primes
+    $sieve($_*$_:-1:$_).=0 foreach(2..sqrt($M-1)); # multiples of 'it' are not prime
+    my $primes=$sieve->xvals->where($sieve); # first primes
+    die "Short sieve" unless $N<=$primes->nelem; # shouldn't happen
+    say "$N-th prime is ", $primes(($N-1));
+}
+
+no PDL::NiceSlice; # NiceSlice destroys indirect function calls!
+sub find_zero { # Find zero of function using Newton's iteration
+    my ($f, $d, $x)=@_; # function, derivative, initial guess
+    my $y;
+    my $max=10; # guard against non-convergence
+    do{($y, $x)= ($x, $x-$f->($x)/$d->($x))} until approx($y,$x) or --$max<=0;
+    return $x;
+}
diff --git a/challenge-146/wlmb/perl/ch-2.pl b/challenge-146/wlmb/perl/ch-2.pl
new file mode 100755
index 0000000000..721d90f75d
--- /dev/null
+++ b/challenge-146/wlmb/perl/ch-2.pl
@@ -0,0 +1,48 @@
+#!/usr/bin/env perl
+# Perl weekly challenge 146
+# Task 2: Curious Fraction Tree
+#
+# See https://wlmb.github.io/2022/01/03/PWC146/#task-2-curious-fraction-tree
+use v5.12;
+use warnings;
+use Try::Tiny;
+foreach(@ARGV){
+    try {
+	my($n, $d)=($1,$2) if m{^\s*(\d+)\s*/\s*(\d+)\s*$};
+	die "Wrong argument $_\n" unless defined $n and defined $d;
+	die "Numerator and denominator in $_  should be positive\n"
+	    unless $n>0 and $d>0;
+	my $gcd=gcd($n, $d);
+	say "Warning: $_ not reduced" unless $gcd==1;
+	($n, $d)=map {$_/$gcd} ($n, $d);
+	try {
+	    my @parent=parent($n, $d);
+	    try {
+		my @grand_parent=parent(@parent);
+		say "Input: $n/$d\nParent:$parent[0]/$parent[1]\n",
+		    "Grand parent: $grand_parent[0]/$grand_parent[1]\n";
+	    }
+	    catch {
+		die "No grandparent of $n/$d\n";
+	    }
+
+	}
+	catch {
+	    die $_;
+	}
+    }
+    catch {
+	say $_;
+    }
+}
+
+sub gcd {
+    my ($n, $m)=@_;
+    ($n, $m)=($m, $n%$m) while ($m);
+    return $n;
+}
+sub parent {
+    my ($n, $d)=@_;
+    die "No parent of $n/$d\n" if $n==$d;
+    return $d>$n?($n,$d-$n):($n-$d,$d);
+}