TWC 146: 10K Prime and CW Trees

In which we leap tall primes in a single bound, mis-take a tree, test percussion, and find the limits of a Curious Module.

Task 1: 10001st Prime Number - One-liners (expanded) in Raku and Perl.

Task 2: Curious Fraction Tree - Solutions in Raku and Perl (with 200+ tests), and another Perl solution using a CPAN module.

TWC Task #1, 10001st Prime Number

Perl

My Perl program is a 3-line version of this one-liner:

perl -Mntheory=nth_prime -wE 'say nth_prime(10_001)'

The ntheory (Number Theory) module has many routines that would solve the problem. nth_prime is the most direct.

Raku

constant @primes = (2, 3, 5, 7 ... *).grep(&is-prime);
say @primes.skip(10_000).head;

That first line creates a "constant" lazy infinte array of primes. It is constant in the sense that the program code cannot alter a value, but it is lazy, so the elements of @primes[0..N] are not generated until something tries to use @primes[N]. All the prior elements are cached and delivered on request without recalculation. If a higher element is requested later in the program's execution, then the generation of @primes elements will resume where it left off.

TWC Task #2, Curious Fraction Tree

I studied the patterns of parent->children and child->parent, slung some code to solve the task, wrote more to extend the tree, and searched OEIS and the Web.

Observations:

  • The "Curious Fraction Tree" is properly named the "Calkin-Wilf Tree", but is easily confused with the 150+year-earlier "Stern-Brocot tree".
    (by "easily confused", I mean that I confused them. Also confused C-W with Kepler's tree of harmonic fractions.)

  • Coolest illustration of the C-W tree: https://www.jasondavies.com/calkin-wilf-tree/

  • C-W tree covers every rational number.
    (Infinite in two dimensions? No big deal, I can do that easily via diagonals; we just get lots of duplicates.)

    • C-W tree does not need diagonals, and never has duplicates.
      • ...I can't do that.
        • Apparently, I can't even understand the proof!
  • Raku has rational numbers as a built-in type: Rat. Rat is auto-normalizing (3/6 becomes 1/2), has .numerator and .denominator methods, and a .nude method to get a list of both (Nu/De) at once.

  • If you follow a fraction's C-W lineage all the way back to (but not including) the 1/1 root, reverse that lineage to be root-to-descendant, and map each fraction along-the-way as proper-fraction==L, improper-fraction==R (i.e. more than 1 or less than 1), then you get a navigation tree of which branch to follow, all the way to the original number. e.g. 8/19 has the navigation of RLRRLL.

    • To a percussionist, those navigations look just like snare drum rudiments and exercises: 'RLRRLL' is a single paradiddle-diddle
      • Those navigations can translate back-and-forth (if we are careful) to binary numbers, and so to plain decimal integers. (More on this when I discover the Curious Module below.)
  • Exploring the tree via regular patterns of L(s) and R(s) finds some interesting properties:

    • 'LLLLLL'... and 'RRRRRR'... are the left and right edges of the tree, with top or bottom == 1.
      1/101 # 'L' x 100
      101/1 # 'R' x 100
    • R,'LLLLLL'... and L,'RRRRRR'... are the two center-most numbers ("framing" the half-way point) on any level after the first. Top or bottom will be '2'.
    • LR,'LLLLLL'... and LL,'RRRRRR'... frame the quarter-way point;
      RR,'LLLLLL'... and RL,'RRRRRR'... frame the three-quarter-way point on any level after the second. Top or bottom will be '3'.
    • (and so on, for row divisions by 8, 16, 32, ...)
    • Single-bit set (or unset) produce the top (or bottom) sequence 2,3,4...(size+1) as we slide the position of the bit:
      2/31 # RLLLLLLLLLLLLLLL
      3/44 # LRLLLLLLLLLLLLLL
      4/55 # LLRLLLLLLLLLLLLL
      ...
      31/2 # LRRRRRRRRRRRRRRR
      44/3 # RLRRRRRRRRRRRRRR
      55/4 # RRLRRRRRRRRRRRRR
    • Alternating L and R reaches the 1/3-point and 2/3-point in a row, producing approximations of the Golden Ratio φ ("The most irrational number"; 1.618033988749...), and its inverse|conjugate.
      These numbers are all ratios of successive Fibonacci numbers.
      e.g. LRLRLRLRLRLRLRLRLRL => 6765/10946 == F(20)/F(21) == 0.618034

After satisfying my urge for percussive exploration, I formalized the findings into a file of tests that the Perl and Raku solutions could both access, then wrote utility routines to help anyone who wants to explore more deeply.

I later found mention in a paper (URL?) that, unlike the Stern-Brocot tree, the C-W tree can be navigated directly across a row, also moving from the end of one row to the start of the next row. I added the calculation to the utility routines as `next-Calkin-Wilf-neighbor.

After all that, I found a Perl module on CPAN that solves the problem: Math::PlanePath::RationalsTree. This lead to a proof-of-concept one-liner, that calculates the parent but not the grandparent:

perl -MMath::PlanePath::RationalsTree -wE 'my $CW = Math::PlanePath::RationalsTree->new( tree_type => "CW" ); say join "/", $CW->n_to_xy( $CW->tree_n_parent( $CW->xy_to_n(@ARGV) ) );' 4817 5410

I cloned my original Perl solution and replaced my algorithm with calls to the module, as ch-2viamodule.pl, so I could be sure that the module could play in my drumming circle.

It worked perfectly, except where it didn't. Of 202 tests, only 166 passed.

Looking deeper into the 36 test failures, the problem became obvious. The module handles all types of rational trees with the efficient method of translating the LR navigation into zeros and ones to make (when decoded from binary) an integer. That integer represents the linear position of the rational.

$ perl -MMath::PlanePath::RationalsTree -wE 'my $CW = Math::PlanePath::RationalsTree->new( tree_type => "CW" ); say "N=$_ : tree_element=", join "/", $CW->n_to_xy( $_ ) for 1..5;'
N=1 : tree_element=1/1
N=2 : tree_element=1/2
N=3 : tree_element=2/1
N=4 : tree_element=1/3
N=5 : tree_element=3/2

The tree gets twice as big for every row (level?) you go down. Exponential growth! Around row 64, N will be around 2^64, and overflow the integer size of 64-bit Perl. All the tests below N=2^64 passed, as did some beyond that point (vagaries of rounding?), but you cannot rely on accuracy past that point. E.g.: finding the parent of "1/65" gives the wrong answer using the "convert to N" method, while 1/10001 and beyond all work with the direct calculations in Perl and Raku. (To be fair, Raku starts having problems past level 90 when drilling down along the Golden Ratio, due to Rat auto-collapsing into a more efficient Num. This could be prevented by using FatRat.)

Raku

Of note in my Raku solution :

  • Clean separation of types: Calkin-Wilf-tree-parent() takes and returns Rat, and task2() takes and returns Str.
  • Symmetry: task2() splits on slash at the start, and joins on slash at the end.
  • Sequence: @lineage is a lazy Seq going all the way back to the root, generated by repeated calls to get the parent of the last result. This "keep feeding the result back in as the next argument" (iterated function) is a common-enough pattern that Raku directly supports it via the ... or "sequence" operator. Since the Seq is lazy, and we only ask for the first two ancestors (via .head), only two are calculated.
  • Contain to Name: Placing the difference in its own variable allows not only DRY in the return statement, but also the symmetry between the two possible results. I particularly like how diff and -diff read, compared to n-d and d-n. (If I had not already committed, I might change diff/d to +diff/d to call attention the sign change).
  • multi sub MAIN clearly shows the two ways to run the program.

Perl

Of note in my Perl solution :

  • Clean separation of types: CalkinWilftree_parent() takes and returns a two-element [numerator,denominator] arrayref, and task2() takes and returns a "numerator/denominator" string.
  • map==grep++ : map is also serving as grep via the empty-list trick; if no non-whitespace characters exist in the line, then map produces nothing for that line, not even undef; the line is skipped.
  • /r : Using the r modifier (part of s{...}{...}msxr) on the substitution lets me return a changed copy, instead of modifying the original.

Kudos to manwar, and thanks to the whole team at The Weekly Challenge!

Elwood: What kind of music do you usually have here?
Claire: Oh, we got both kinds. We got Country and Western.

Rawhide!

Leave a comment

About Bruce Gray

user-pic "Util" on IRC and PerlMonks. Frequent speaker on Perl and Raku, but infrequent blogger.