From 09f7f9432bb8bb5dfe4f47e8d3ba00a2fe461f53 Mon Sep 17 00:00:00 2001 From: James Smith Date: Sun, 3 Apr 2022 22:31:58 +0100 Subject: Update README.md --- challenge-158/james-smith/README.md | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/challenge-158/james-smith/README.md b/challenge-158/james-smith/README.md index 342202ccbf..a81bb61a28 100644 --- a/challenge-158/james-smith/README.md +++ b/challenge-158/james-smith/README.md @@ -94,7 +94,7 @@ sub additive_primes_split_sum0 { * `split` - 75% * `sum0 split` - 45% -# Challenge 2 - First series buban primes +# Challenge 2 - First series cuban primes ***Write a script to compute first series cuban primes <= 1000. (First series cuban primes have the form `((x+1)^3-x^3)/(x+1-x)` = `3x^2+3x+1`)*** -- cgit From 9bd3a76992e2c016e85625c4d94969ec05ca7a43 Mon Sep 17 00:00:00 2001 From: drbaggy Date: Mon, 4 Apr 2022 07:32:24 +0100 Subject: 159 --- challenge-159/james-smith/README.md | 130 +++++++++------------------------ challenge-159/james-smith/blog.txt | 1 + challenge-159/james-smith/perl/ch-1.pl | 20 +++++ challenge-159/james-smith/perl/ch-2.pl | 39 ++++++++++ 4 files changed, 94 insertions(+), 96 deletions(-) create mode 100644 challenge-159/james-smith/blog.txt create mode 100644 challenge-159/james-smith/perl/ch-1.pl create mode 100644 challenge-159/james-smith/perl/ch-2.pl diff --git a/challenge-159/james-smith/README.md b/challenge-159/james-smith/README.md index 342202ccbf..cb08f40942 100644 --- a/challenge-159/james-smith/README.md +++ b/challenge-159/james-smith/README.md @@ -1,6 +1,6 @@ -[< Previous 157](https://github.com/drbaggy/perlweeklychallenge-club/tree/master/challenge-157/james-smith) | -[Next 159 >](https://github.com/drbaggy/perlweeklychallenge-club/tree/master/challenge-159/james-smith) -# The Weekly Challenge 158 +[< Previous 156](https://github.com/drbaggy/perlweeklychallenge-club/tree/master/challenge-156/james-smith) | +[Next 160 >](https://github.com/drbaggy/perlweeklychallenge-club/tree/master/challenge-160/james-smith) +# The Weekly Challenge 159 You can find more information about this weeks, and previous weeks challenges at: @@ -12,120 +12,58 @@ submit solutions in whichever language you feel comfortable with. You can find the solutions here on github at: -https://github.com/drbaggy/perlweeklychallenge-club/tree/master/challenge-158/james-smith +https://github.com/drbaggy/perlweeklychallenge-club/tree/master/challenge-159/james-smith -# Challenge 1 - Additive primes +# Challenge 1 - Farey Sequence -***Additive primes are prime numbers for which the sum of their decimal digits are also primes.*** +*** You are given a positive number, `$n`. Write a script to compute Farey Sequence of the order `$n`. All rational numbers less than with with denominators less than `$n` ## The solution -We loop through each prime p, work out the digit sum (by repeated modulo 10/divide by 10) and check that it is prime. -We craft this as two loops - an outer `for` loop and an inner `do {} while`. - -```perl -use Math::Prime::Util qw(next_prime is_prime); - -sub additive_primes { - my @res; - for( - my $p = 2 ; - my $s = 0, ( my $q = $p ) <= $N; - (is_prime $s) && (push @res, $p), $p = next_prime $p - ) { - do { $s += $q % 10 } while $q = int $q / 10; - } - @res; -} -``` -## Output -``` -2, 3, 5, 7, 11, 23, 29, 41, 43, 47, 61, 67, 83, 89 -``` -## Notes - * We use a C-style `for` loop, with it's three parts *initialization*, *condition* and *increment*. - * The *condition* and *increment* statements are complicated, each with two parts separated by `,`s. - * The *condition* block is called at the start of each loop, and so we use it to initialise the variables as the loop, as well as checking the condition. - * The *increment* block is called at the end of the loop and so stores the value of `$p` if it is an additive prime, then it increments the loop with the next prime. - * Rather than doing a split and sum we use repeated dividing and summing, as it is more efficient around 20-30% more efficient. The increased performance is probably due to avoiding the "duality" of perl variables storing numbers as numbers/strings. - -## Extra code - -Rewritten with single line `for` ... (the original version of the code) - -```perl -sub additive_primes_div { - my @res; - for(my$p=2; my$s=0,(my $q=$p)<=$N;(is_prime$s)&&(push@res,$p),$p=next_prime$p) { - do { $s += $q % 10 } while $q = int $q / 10; - } - @res; -} -``` - -Alternative form with `for split`: - ```perl -sub additive_primes_split { - my @res; - for( my $p = 2; my $s = 0, $p <= $N ; $p = next_prime $p ) { - $s+=$_ for split //, $p; - push @res, $p if is_prime $s; +sub farley { + my($n,%v) = shift; + for my $d (1..$n) { + $v{$_/$d}||="$_/$d" for 1..$d; } - @res; + map{$v{$_}}sort{$a<=>$b}keys%v; } ``` -And an alternate using `sum0`... - -```perl -sub additive_primes_split_sum0 { - my @res; - for( my $p = 2; $p <= $N ; $p = next_prime $p ) { - push @res, $p if is_prime sum0 split //, $p; - } - @res; -} -``` +# Challenge 2 - Moebius Number -### Relative performance +***You are given a positive number `$n`. Write a script to generate the Moebius Number for the given number (see defn below) - * `div`/`mod` method - 100% - * `split` - 75% - * `sum0 split` - 45% +## Definition -# Challenge 2 - First series buban primes +The value of the Moebius number is: -***Write a script to compute first series cuban primes <= 1000. (First series cuban primes have the form `((x+1)^3-x^3)/(x+1-x)` = `3x^2+3x+1`)*** + * `0` if `$n` has a prime factor + * `1` if the number has an even number of prime factors + * `-1` otherwise ## The solution -The solution is rather short. We try each of the value of the sequence `3*$x**2+3*$x+1` in turn up to the value of 1000 (`$N`). -We output each value which in turn is prime. - ```perl -(is_prime $_) && say while $N >= ($_ = 3*$x*++$x+1); -``` -## Output -``` -7, 19, 37, 61, 127, 271, 331, 397, 547, 631, 919 +sub moebius { + my ($n,$p,$r) = (shift,1,1); + return -1 if is_prime($n); + $n%($p**2) ? ($n%$p && ($r=-$r)) : return 0 while ($p = next_prime $p) < $n; + $r; +} ``` -## Notes - * As we use `$_` as our temporary variable we can use `say` by itself to output it. - * We use our common trick of `(condition) && (command)` to work as an `command if(condition)` which can be embedded in a postfix loop. - * There is a "trick" as we increment `$x` half way through the calculation of `$_`. - * `$_ = 3*$x**2 + 3*$x + 1` => `$_ = 3 * $x * ($x+1) + 1` => `$_ = 3 * $x * ++$x + 1` - * We can replace the `$x+1` by the pre increment operator `++$x` so this becomes `3 * $x * ++$x + 1`. -## Aside +Expanding this out gives the more readable (but longer)... -Second series cuban primes have the form `((x+2)^3-x^3)/(x+2-x)` = `3x^2+3.2x+4`. We can tweak the code to give: - -```perl -(is_prime $_) && say while $N >= ($_ = 3 * $x * (2 + $x++) + 4); ``` +sub moebius_exp { + my ($n,$p,$r) = (shift,1,1); + return -1 if is_prime($n); + while( ($p = next_prime $p ) < $n ) { + return 0 unless $n%($p**2); + $r=-$r unless $n%$p; + } + $r; +} -which outputs: -``` -13, 109, 193, 433, 769 ``` diff --git a/challenge-159/james-smith/blog.txt b/challenge-159/james-smith/blog.txt new file mode 100644 index 0000000000..e631de6584 --- /dev/null +++ b/challenge-159/james-smith/blog.txt @@ -0,0 +1 @@ +https://github.com/drbaggy/perlweeklychallenge-club/tree/master/challenge-159/james-smith diff --git a/challenge-159/james-smith/perl/ch-1.pl b/challenge-159/james-smith/perl/ch-1.pl new file mode 100644 index 0000000000..5f6402ec41 --- /dev/null +++ b/challenge-159/james-smith/perl/ch-1.pl @@ -0,0 +1,20 @@ +#!/usr/local/bin/perl + +use strict; + +use warnings; +use feature qw(say); +use Test::More; +use Benchmark qw(cmpthese timethis); +use Data::Dumper qw(Dumper); + +say join ', ', farley( $_ ) for 1..10; + +sub farley { + my($n,%v) = shift; + for my $d (1..$n) { + $v{$_/$d}||="$_/$d" for 1..$d; + } + map{$v{$_}}sort{$a<=>$b}keys%v; +} + diff --git a/challenge-159/james-smith/perl/ch-2.pl b/challenge-159/james-smith/perl/ch-2.pl new file mode 100644 index 0000000000..5c86128a4d --- /dev/null +++ b/challenge-159/james-smith/perl/ch-2.pl @@ -0,0 +1,39 @@ +#!/usr/local/bin/perl + +use strict; + +use warnings; +use feature qw(say); +use Test::More; +use Benchmark qw(cmpthese timethis); +use Data::Dumper qw(Dumper); +use Math::Prime::Util qw(next_prime is_prime); + +my @TESTS = ( + [ 5, -1 ], + [ 10, 1 ], + [ 20, 0 ], +); + +is( moebius( $_->[0]), $_->[1] ) foreach @TESTS; +is( moebius_exp($_->[0]), $_->[1] ) foreach @TESTS; + +done_testing(); + +sub moebius { + my ($n,$p,$r) = (shift,1,1); + return -1 if is_prime($n); + $n%($p**2) ? ($n%$p && ($r=-$r)) : return 0 while ($p = next_prime $p) < $n; + $r; +} + +sub moebius_exp { + my ($n,$p,$r) = (shift,1,1); + return -1 if is_prime($n); + while( ($p = next_prime $p ) < $n ) { + return 0 unless $n%($p**2); + $r=-$r unless $n%$p; + } + $r; +} + -- cgit From d0335bccd0294b93ca8539b22390c0064baa3193 Mon Sep 17 00:00:00 2001 From: drbaggy Date: Mon, 4 Apr 2022 19:41:37 +0100 Subject: faster code --- challenge-159/james-smith/perl/ch-2.pl | 53 +++++++++++++++++++++++++++++----- 1 file changed, 46 insertions(+), 7 deletions(-) diff --git a/challenge-159/james-smith/perl/ch-2.pl b/challenge-159/james-smith/perl/ch-2.pl index 5c86128a4d..7d22c8df6d 100644 --- a/challenge-159/james-smith/perl/ch-2.pl +++ b/challenge-159/james-smith/perl/ch-2.pl @@ -13,26 +13,65 @@ my @TESTS = ( [ 5, -1 ], [ 10, 1 ], [ 20, 0 ], +## The following are good examples where the second algorithm works well! + [ 2, -1 ], + [ 6, 1 ], + [ 30, -1 ], + [ 210, 1 ], + [ 2310, -1 ], + [ 30030, 1 ], + [ 510510, -1 ], + [ 9699690, 1 ], + [ 6895531, -1 ], + [ 206865930, 1 ], + [ 223092870, -1 ], + [ 6469693230, 1 ], ); -is( moebius( $_->[0]), $_->[1] ) foreach @TESTS; -is( moebius_exp($_->[0]), $_->[1] ) foreach @TESTS; +warn "exp"; +is( moebius_exp( $_->[0] ), $_->[1] ) for @TESTS; +warn "opt"; +is( moebius_div_opt( $_->[0] ), $_->[1] ) for @TESTS; +timethis( 500_000, sub { moebius_div_opt( $_->[0] ) for @TESTS } ); +warn "div"; +is( moebius_div( $_->[0] ), $_->[1] ) for @TESTS; +timethis( 100, sub { moebius_div( $_->[0] ) for @TESTS } ); +warn "dumb"; +is( moebius( $_->[0] ), $_->[1] ) for @TESTS; +timethis( 5, sub { moebius( $_->[0] ) for @TESTS } ); done_testing(); sub moebius { my ($n,$p,$r) = (shift,1,1); return -1 if is_prime($n); - $n%($p**2) ? ($n%$p && ($r=-$r)) : return 0 while ($p = next_prime $p) < $n; + $n%($p**2) ? ($n%$p || ($r=-$r)) : return 0 while ($p = next_prime $p) <= $n; + $r; +} + +sub moebius_div { + my ($n,$p,$r) = (shift,1,1); + return -1 if is_prime $n; + $n%($p**2) ? ( $n%$p || ($r=-$r,$n/=$p)) : return 0 while ($p = next_prime $p) && $n>1; + $r; +} + +sub moebius_div_opt { + my ($n,$p,$r) = (shift,1,1); + return -1 if is_prime $n; + $n%$p || ($r=-$r,$n/=$p) && ( $n%$p ? ( (is_prime $n) && (return -$r) ) : return 0 ) while ($p = next_prime $p)**2 <= $n; $r; } sub moebius_exp { my ($n,$p,$r) = (shift,1,1); - return -1 if is_prime($n); - while( ($p = next_prime $p ) < $n ) { - return 0 unless $n%($p**2); - $r=-$r unless $n%$p; + return -1 if is_prime $n; + while( ($p = next_prime $p )**2 <= $n ) { + next if $n%$p; + $r=-$r; + $n/=$p; + return 0 unless $n%$p; + return -$r if is_prime $n; } $r; } -- cgit From e0ccb368c09d324df30fb8581b95bfda6bd5011f Mon Sep 17 00:00:00 2001 From: drbaggy Date: Mon, 4 Apr 2022 19:59:01 +0100 Subject: another tweak --- challenge-159/james-smith/perl/ch-2.pl | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/challenge-159/james-smith/perl/ch-2.pl b/challenge-159/james-smith/perl/ch-2.pl index 7d22c8df6d..fb9ee8a893 100644 --- a/challenge-159/james-smith/perl/ch-2.pl +++ b/challenge-159/james-smith/perl/ch-2.pl @@ -32,7 +32,8 @@ warn "exp"; is( moebius_exp( $_->[0] ), $_->[1] ) for @TESTS; warn "opt"; is( moebius_div_opt( $_->[0] ), $_->[1] ) for @TESTS; -timethis( 500_000, sub { moebius_div_opt( $_->[0] ) for @TESTS } ); +#timethis( 500_000, sub { moebius_div_opt( $_->[0] ) for @TESTS } ); +exit; warn "div"; is( moebius_div( $_->[0] ), $_->[1] ) for @TESTS; timethis( 100, sub { moebius_div( $_->[0] ) for @TESTS } ); @@ -59,7 +60,7 @@ sub moebius_div { sub moebius_div_opt { my ($n,$p,$r) = (shift,1,1); return -1 if is_prime $n; - $n%$p || ($r=-$r,$n/=$p) && ( $n%$p ? ( (is_prime $n) && (return -$r) ) : return 0 ) while ($p = next_prime $p)**2 <= $n; + $n%$p || ($n/=$p) && $n%$p ? ( is_prime $n ? (return $r) : ($r=-$r) ) : return 0 while ($p = next_prime $p)**2 <= $n; $r; } -- cgit From 103cda8b443a2cfe130adcd3839d652359f8b31f Mon Sep 17 00:00:00 2001 From: James Smith Date: Mon, 4 Apr 2022 19:59:12 +0100 Subject: Update README.md --- challenge-159/james-smith/README.md | 66 +++++++++++++++++++++++++++++++------ 1 file changed, 56 insertions(+), 10 deletions(-) diff --git a/challenge-159/james-smith/README.md b/challenge-159/james-smith/README.md index cb08f40942..31d01d0ee1 100644 --- a/challenge-159/james-smith/README.md +++ b/challenge-159/james-smith/README.md @@ -16,10 +16,14 @@ https://github.com/drbaggy/perlweeklychallenge-club/tree/master/challenge-159/ja # Challenge 1 - Farey Sequence -*** You are given a positive number, `$n`. Write a script to compute Farey Sequence of the order `$n`. All rational numbers less than with with denominators less than `$n` +***You are given a positive number, `$n`. Write a script to compute Farey Sequence of the order `$n`. All rational numbers less than with with denominators less than `$n`*** ## The solution +This is a relatively simple piece of code - we loop through the denominators `1`..`$n` and compute all the fractions greater than 0 and less than or equal to 1. +We use the value as the key to the hash and "`$_/$d`" as the values. +We start with the lower denominators so if we find the same fraction twice (e.g. `1/2` & `2/4`) the value we always store has the lowest denominator and so the fraction is in its simplest form. A simple (key based) sort at the end and returning the list of strings of fractions... + ```perl sub farley { my($n,%v) = shift; @@ -32,7 +36,7 @@ sub farley { # Challenge 2 - Moebius Number -***You are given a positive number `$n`. Write a script to generate the Moebius Number for the given number (see defn below) +***You are given a positive number `$n`. Write a script to generate the Moebius Number for the given number (see defn below)*** ## Definition @@ -42,7 +46,7 @@ The value of the Moebius number is: * `1` if the number has an even number of prime factors * `-1` otherwise -## The solution +## First pass... *naive* approach ```perl sub moebius { @@ -53,17 +57,59 @@ sub moebius { } ``` -Expanding this out gives the more readable (but longer)... +Here we just look for prime factors less than `$n` - for large `$n` - the number of primes checked is `$n/log $n`. + +## Second pass... *first-pass* optimization + +We can reduce the number of primes checked by reducing `$n` each time a prime is found, by dividing `$n` by the prime factor. +```perl +sub moebius_div { + my ($n,$p,$r) = (shift,1,1); + return -1 if is_prime $n; + $n%($p**2) ? ( $n%$p || ($r=-$r,$n/=$p)) : return 0 while ($p = next_prime $p) && $n>1; + $r; +} ``` -sub moebius_exp { + +## Third pass... *improved* performance + +Finally we ideally want to restrict the number of primes checked further. We know that if `$n/$p` is prime we can stop there when looking for factors. It then also reduces the number primes to check to the square root of `$n`, so the maximum number of primes needed to check is `sqrt($n)/log $n/2`. + +```perl +sub moebius_div_opt { my ($n,$p,$r) = (shift,1,1); - return -1 if is_prime($n); - while( ($p = next_prime $p ) < $n ) { - return 0 unless $n%($p**2); - $r=-$r unless $n%$p; - } + return -1 if is_prime $n; + $n%$p || ($n/=$p) && $n%$p ? ( is_prime $n ? (return $r) : ($r=-$r) ) : return 0 while ($p = next_prime $p)**2 <= $n; $r; } +``` + +## Relative performance + +On our extended test set the relative performance of the three methods is: + +| Method | Calls/second | Relative performance | +| :-----: | -----------: | -------------------: | +| simple | 0.0025/s | 1 x | +| divide | 3.0/s | 1,200 x | +| div opt | 38,000/s | 16,000,000 x | + +## Expanded version of optimal solution.. + +The format of the optimized code is compact - here is an expanded version of the code to show the logic. ``` +sub moebius { + my ($n,$p,$r) = (shift,1,1); + return -1 if is_prime $n; + while( ($p = next_prime $p )**2 <= $n ) { + next if $n%$p; + $n/=$p; + return 0 unless $n%$p; + return $r if is_prime $n; + $r=-$r; + } + $r; +} +``` -- cgit