Math::...tionals(3pm)
💖 Math::PlanePath::GcdRationals(3pm) 💖
okay bestie, gather round, because this one is genuinely iconic — it's a Perl module page from the User Contributed Perl Documentation, and it's literally called Math::PlanePath::GcdRationals. it lives at Math::PlanePath::GcdRationals(3pm) and it's about enumerating rationals with a triangular GCD trick. no cap, this is the good stuff. 🦋
✨ NAME ✨
Math::PlanePath::GcdRationals — rationals by triangular GCD 💅
🌸 SYNOPSIS 🌸
use Math::PlanePath::GcdRationals;
my $path = Math::PlanePath::GcdRationals->new;
my ($x, $y) = $path->n_to_xy (123);👑 DESCRIPTION 👑
so this path enumerates X/Y rationals using a method by Lance Fortnow that pulls a greatest common divisor out of a triangular position. obsessed, honestly.
<http://blog.computationalcomplexity.org/2004/03/counting-rationals-quickly.html>
the whole vibe of this approach is that it's both efficient to calculate and it visits blocks of X/Y rationals using a modest range of N values — roughly a square, N=2*max(num,den)^2 in the default rows style. lowkey perfect. 🎀
13 | 79 80 81 82 83 84 85 86 87 88 89 90
12 | 67 71 73 77 278
11 | 56 57 58 59 60 61 62 63 64 65 233 235
10 | 46 48 52 54 192 196
9 | 37 38 40 41 43 44 155 157 161
8 | 29 31 33 35 122 126 130
7 | 22 23 24 25 26 27 93 95 97 99 101 103
6 | 16 20 68 76 156
5 | 11 12 13 14 47 49 51 53 108 111 114
4 | 7 9 30 34 69 75 124
3 | 4 5 17 19 39 42 70 74 110
2 | 2 8 18 32 50 72 98
1 | 1 3 6 10 15 21 28 36 45 55 66 78 91
Y=0 |
--------------------------------------------------------
X=0 1 2 3 4 5 6 7 8 9 10 11 12 13the mapping from N to rational is
N = i + j*(j-1)/2 for upper triangle 1 <= i <= j
gcd = GCD(i,j)
rational = i/j + gcd-1which means X=numerator, Y=denominator are
X = (i + j*(gcd-1))/gcd = j + (i-j)/gcd
Y = j/gcdthe i,j position is a numbering of points above the X=Y diagonal by rows, in the style of Math::PlanePath::PyramidRows with step=1, but starting from i=1,j=1. 🧚
j=4 | 7 8 9 10
j=3 | 4 5 6
j=2 | 2 3
j=1 | 1
+-------------
i=1 2 3 4if GCD(i,j)=1 then X/Y is simply X=i,Y=j unchanged. this means fractions X/Y < 1 are numbered by rows with increasing numerator, but skipping positions where i,j have a common factor. so picky, so valid. 💅
the skipped positions where i,j have a common factor become rationals X/Y>1, ie. below the X=Y diagonal. the integer part is GCD(i,j)-1 so rational = gcd-1 + i/j. for example
N=51 is at i=6,j=10 by rows
common factor gcd(6,10)=2
so rational R = 2-1 + 6/10 = 1+3/5 = 8/5
ie. X=8,Y=5if j is prime then gcd(i,j)=1 and so X=i,Y=j. this means that in rows with prime Y are numbered by consecutive N across to the X=Y diagonal. for example in row Y=7 above N=22 to N=27. 🔥
🎀 Triangular Numbers 🎀
N=1,3,6,10,etc along the bottom Y=1 row is the triangular numbers N=k*(k-1)/2. such an N is at i=k,j=k and has gcd(i,j)=k which divides out to Y=1.
N=k*(k-1)/2 i=k,j=k
Y = j/gcd
= 1 on the bottom row
X = (i + j*(gcd-1)) / gcd
= (k + k*(k-1)) / k
= k-1 successive points on that bottom rowN=1,2,4,7,11,etc in the column at X=1 immediately follows each of those bottom row triangulars, ie. N+1. the way that works out is so satisfying fr fr.
N in X=1 column = Y*(Y-1)/2 + 1💜 Primes 💜
if N is prime then it's above the sloping line X=2*Y. if N is composite then it might be above or below, but the primes are always above. here's the table with dots ... marking the X=2*Y line.
primes and composites above
|
6 | 16 20 68
| .... X=2*Y
5 | 11 12 13 14 47 49 51 53 ....
| ....
4 | 7 9 30 34 .... 69
| ....
3 | 4 5 17 19 .... 39 42 70 only
| .... composite
2 | 2 8 .... 18 32 50 below
| ....
1 | 1 ..3. 6 10 15 21 28 36 45 55
| ....
Y=0 | ....
---------------------------------------------
X=0 1 2 3 4 5 6 7 8 9 10values below X=2*Y such as 39 and 42 are always composite. values above such as 19 and 30 are either prime or composite. only X=2,Y=1 is exactly on the line, which is prime N=3 as it happens. the rest of the line X=2*k,Y=k is not visited since common factor k would mean X/Y is not a rational in least terms. periodt. 👑
this pattern of primes and composites occurs because N is a multiple of gcd(i,j) when that gcd is odd, or a multiple of gcd/2 when that gcd is even.
N = i + j*(j-1)/2
gcd = gcd(i,j)
N = gcd * (i/gcd + j/gcd * (j-1)/2) when gcd odd
gcd/2 * (2i/gcd + j/gcd * (j-1)) when gcd evenif gcd odd then either j/gcd or j-1 is even, to take the /2 divisor. if gcd even then only gcd/2 can come out as a factor, since taking out the full gcd might leave both j/gcd and j-1 odd and so the /2 not an integer. that happens for example to N=70 😭
N = 70
i = 4, j = 12 for 4 + 12*11/2 = 70 = N
gcd(i,j) = 4
but N is not a multiple of 4, only of 4/2=2of course knowing gcd or gcd/2 is a factor of N is only useful when that factor is 2 or more, so
odd gcd >= 2 means gcd >= 3
even gcd with gcd/2 >= 2 means gcd >= 4
so N composite when gcd(i,j) >= 3if gcd<3 then the "factor" coming out is only 1 and says nothing about whether N is prime or composite. there are both prime and composite N with gcd<3, as can be seen among the values above the X=2*Y line in the table above. ✨
🦋 Rows Reverse 🦋
option "pairs_order => "rows_reverse"" reverses the order of points within the rows of i,j pairs,
j=4 | 10 9 8 7
j=3 | 6 5 4
j=2 | 3 2
j=1 | 1
+------------
i=1 2 3 4the X,Y numbering becomes
pairs_order => "rows_reverse"
11 | 66 65 64 63 62 61 60 59 58 57
10 | 55 53 49 47 209
9 | 45 44 42 41 39 38 170 168
8 | 36 34 32 30 135 131
7 | 28 27 26 25 24 23 104 102 100 98
6 | 21 17 77 69
5 | 15 14 13 12 54 52 50 48 118
4 | 10 8 35 31 76 70
3 | 6 5 20 18 43 40 75 71
2 | 3 9 19 33 51 73
1 | 1 2 4 7 11 16 22 29 37 46 56
Y=0 |
------------------------------------------------
X=0 1 2 3 4 5 6 7 8 9 10 11the triangular numbers, per "Triangular Numbers" above, are now in the X=1 column, ie. at the left rather than at the Y=1 bottom row. that bottom row is now the next after each triangular, ie. T(X)+1. 🌸
🌈 Diagonals 🌈
option "pairs_order => "diagonals_down"" takes the i,j pairs by diagonals down from the Y axis. "pairs_order => "diagonals_up"" likewise but upwards from the X=Y centre up to the Y axis. (these numberings are in the style of Math::PlanePath::DiagonalsOctant.) iykyk. 💖
diagonals_down diagonals_up
j=7 | 13 j=7 | 16
j=6 | 10 14 j=6 | 12 15
j=5 | 7 11 15 j=5 | 9 11 14
j=4 | 5 8 12 16 j=4 | 6 8 10 13
j=3 | 3 6 9 j=3 | 4 5 7
j=2 | 2 4 j=2 | 2 3
j=1 | 1 j=1 | 1
+------------ +------------
i=1 2 3 4 i=1 2 3 4the resulting path becomes
pairs_order => "diagonals_down"
9 | 21 27 40 47 63 72
8 | 17 28 41 56 74
7 | 13 18 23 29 35 42 58 76
6 | 10 30 44
5 | 7 11 15 20 32 46 62 80
4 | 5 12 22 48 52
3 | 3 6 14 24 33 55
2 | 2 8 19 34 54
1 | 1 4 9 16 25 36 49 64 81
Y=0 |
--------------------------------
X=0 1 2 3 4 5 6 7 8 9the Y=1 bottom row is the perfect squares which are at i=j in the "DiagonalsOctant" and have gcd(i,j)=i thus becoming X=i,Y=1.
pairs_order => "diagonals_up"
9 | 25 29 39 45 58 65
8 | 20 28 38 50 80
7 | 16 19 23 27 32 37 63 78
6 | 12 26 48
5 | 9 11 14 17 35 46 59 74
4 | 6 10 24 44 54
3 | 4 5 15 22 34 51
2 | 2 8 18 33 52
1 | 1 3 7 13 21 31 43 57 73
Y=0 |
--------------------------------
X=0 1 2 3 4 5 6 7 8 9N=1,2,4,6,9 etc in the X=1 column is the perfect squares k*k and the pronics k*(k+1) interleaved, also called the quarter-squares. N=2,5,10,17,etc on Y=X+1 above the leading diagonal are the squares+1, and N=3,8,15,24,etc below on Y=X-1 below the diagonal are the squares-1. slay. 🔥
the GCD division moves points downwards and shears them across horizontally. the effect on diagonal lines of i,j points is as follows
| 1
| 1 gcd=1 slope=-1
| 1
| 1
| 1
| 1
| 1
| 1
| 1
| . gcd=2 slope=0
| . 2
| . 2 3 gcd=3 slope=1
| . 2 3 gcd=4 slope=2
| . 2 3 4
| . 3 4 5 gcd=5 slope=3
| . 4 5
| . 4 5
| . 5
+-------------------------------the line of "1"s is the diagonal with gcd=1 and thus X,Y=i,j unchanged.
the line of "2"s is when gcd=2 so X=(i+j)/2,Y=j/2. since i+j=d is constant within the diagonal this makes X=d fixed, ie. vertical.
then gcd=3 becomes X=(i+2j)/3 which slopes across by +1 for each i, or gcd=4 has X=(i+3j)/4 slope +2, etc. 💅
...and that's only part 1, bestie — there's more coming, don't go anywhere 🫶
Of course only some of the points on an i,j diagonal have a given gcd, but the ones that do get transformed this way, and that's the whole vibe. The effect? For N up to a given diagonal row, all the * points below get traversed — plus some extras hanging out in wedge-shaped arms off to the side. it's giving very "main character plus entourage" 💅
| *
| * * up to a given diagonal points "*"
| * * * all visited, plus some wedges out
| * * * * to the right
| * * * * *
| * * * * * /
| * * * * * / --
| * * * * * --
| * * * * *--
+--------------And in terms of the rationals X/Y, bestie, the tea is: up to N=d^2 with diagonal d=2j, the fractions you enumerate are
N=d^2
enumerates all num/den where num <= d and num+den <= 2*d✨ FUNCTIONS ✨
Go peep "FUNCTIONS" in Math::PlanePath for the behaviour that's common to all the path classes, okay? 🫶
$path = Math::PlanePath::GcdRationals->new ()$path = Math::PlanePath::GcdRationals->new (pairs_order => $str)— makes you a shiny new path object, no cap. Thepairs_orderoption can be any of these cuties:
"rows" (default)
"rows_reverse"
"diagonals_down"
"diagonals_up"🧮 FORMULAS 🧮
X,Y to N — Rows
The defining formula up there for X,Y can be inverted to hand you back i,j and N. Important tho: this calculation does not notice if X,Y share a common factor, so if you need a coprime(X,Y) check you gotta do it separately (for "xy_to_n()" it does get done).
X/Y = g-1 + (i/g)/(j/g)The g-1 integer part gets recovered with a division of X by Y:
X = quot*Y + rem division by Y rounded towards 0
where 0 <= rem < Y
unless Y=1 in which case use quot=X-1, rem=1
g-1 = quot
g = quot+1The Y=1 special case can instead be left as the regular kind of division quot=X,rem=0, so 0<=rem<Y. That gives i=0, which is outside the intended 1<=i<=j range, but j is 1 bigger and the combo still lands you the correct N. It's like the i=g,j=g point at the end of a row got moved to i=0,j=g+1 right before the start of the next row. So if you only care about N and not i,j, you can just leave rem=0. slay ✨
Equating the denominators in that X/Y formula gives you j by:
Y = j/g the definition above
j = g*Y
= (quot+1)*Y
j = X+Y-rem per the division X=quot*Y+remAnd equating the numerators gives you i by:
X = (g-1)*Y + i/g the definition above
i = X*g - (g-1)*Y*g
= X*g - quot*Y*g
= X*g - (X-rem)*g per the division X=quot*Y+rem
i = rem*g
i = rem*(quot+1)Then N from i,j, using the definition above:
N = i + j*(j-1)/2For example, X=11, Y=4 divides X/Y as 11=4*2+3 for quot=2, rem=3, so i=3*(2+1)=9, j=11+4-3=12, and so N=9+12*11/2=75 (which is exactly what's in that first table up above, bestie 💖).
You can get away with using only the quotient p and not the remainder rem, by taking j=(quot+1)*Y instead of j=X+Y-rem. But usually a division op hands you the remainder at no extra cost, or a cost small enough that it's worth swapping a multiply for an add or two. efficiency queen 👑
The gcd g can also be recovered by rounding up in the division, instead of rounding down then incrementing with g=quot+1.
g = ceil(X/Y)
= cquot for division X=cquot*Y - cremBut heads up: division in most programming languages rounds towards 0 or towards -infinity, not upwards towards +infinity. so like, watch out 🫠
X,Y to N — Rows Reverse
For pairs_order="rows_reverse", the horizontal i gets reversed to j-i+1. You can work that right into the triangular part of the N formula like so:
Nrrev = (j-i+1) + j*(j-1)/2 for 1<=i<=j
= j*(j+1)/2 - i + 1The Y=1 case described above cannot be left to go through with rem=0 giving i=0 and j+1, because the reversal j-i+1 is then not correct. So either use rem=1 as described, or if you don't, compensate at the end:
if r=0 then j -= 2 adjust
Nrrev = j*(j+1)/2 - i + 1 same Nrrev as aboveFor example, X=5, Y=1 is quot=5, rem=0, which gives i=0*(5+1)=0, j=5+1-0=6. Without the adjustment it'd be Nrrev=6*7/2-0+1=22, which is wrong 😭. But adjusting j-=2 so that j=6-2=4 gives the desired Nrrev=4*5/2-0+1=11 (per the table in "Rows Reverse" above). we love a correction arc 💅
🔍 OEIS 🔍
This enumeration of rationals lives in Sloane's Online Encyclopedia of Integer Sequences in the following forms:
<http://oeis.org/A054531> (etc)
pairs_order="rows" (the default)
A226314 X coordinate
A054531 Y coordinate, being N/GCD(i,j)
A000124 N in X=1 column, triangular+1
A050873 ceil(X/Y), gcd by rows
A050873-A023532 floor(X/Y)
gcd by rows and subtract 1 unless i=j
pairs_order="diagonals_down"
A033638 N in X=1 column, quartersquares+1 and pronic+1
A000290 N in Y=1 row, perfect squares
pairs_order="diagonals_up"
A002620 N in X=1 column, squares and pronics
A002061 N in Y=1 row, central polygonals (extra initial 1)
A002522 N at Y=X+1 above leading diagonal, squares+1👀 SEE ALSO 👀
Go read these next, bestie, they're the whole friend group:
Math::PlanePathMath::PlanePath::DiagonalRationalsMath::PlanePath::RationalsTreeMath::PlanePath::CoprimeColumnsMath::PlanePath::DiagonalsOctant
🏠 HOME PAGE 🏠
<http://user42.tuxfamily.org/math-planepath/index.html>📜 LICENSE 📜
Copyright 2011, 2012, 2013, 2014, 2015, 2016, 2017, 2018, 2019, 2020 Kevin Ryde
This file is part of Math-PlanePath.
Math-PlanePath is free software; you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation; either version 3, or (at your option) any later version.
Math-PlanePath is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.
You should have received a copy of the GNU General Public License along with Math-PlanePath. If not, see <http://www.gnu.org/licenses/>.
and that's the whole page, bestie 🎀 this was perl v5.32.0, dated 2021-01-23, and yeah the title got clipped to Math::...tionals(3pm) but you know exactly which one it is. mwah 💋
📦 where this page comes from
- debian bookworm: libmath-planepath-perl 129-1 · GPL-3.0-or-later
- debian trixie: libmath-planepath-perl 129-1 · GPL-3.0-or-later
- ubuntu jammy: libmath-planepath-perl 129-1 · GPL-3.0-or-later
- ubuntu noble: libmath-planepath-perl 129-1 · GPL-3.0-or-later
- ubuntu resolute: libmath-planepath-perl 129-1 · GPL-3.0-or-later
- ubuntu stonking: libmath-planepath-perl 129-1 · GPL-3.0-or-later
- 💖 girlified by deepseek-ai/DeepSeek-V4.1-Flash on : a language model's rewrite of the original, under its terms