Math::PlanePath::PythagoreanTree -- primitive Pythagorean triples by tree
use Math::PlanePath::PythagoreanTree; my $path = Math::PlanePath::PythagoreanTree->new (tree_type => 'UAD', coordinates => 'AB'); my ($x, $y) = $path->n_to_xy (123);
This path enumerates primitive Pythagorean triples by a breadth-first traversal of a ternary tree, either a "UAD" or "FB" tree. Each point is an integer X,Y=A,B with integer hypotenuse, and is primitive in the sense that A and B have no common factor.
A^2 + B^2 = C^2 gcd(A,B)=1 ie. no common factor X=A, Y=B
A primitive triple always has one of A,B odd and the other even. The trees here give them ordered as A odd and B even.
The breadth-first traversal goes out to rather large A,B values while yet to complete smaller ones. The UAD tree goes out further than the FB.
The UAD tree by Berggren (1934) and later independently by Barning (1963), Hall (1970), and several others, uses three matrices U, A and D which can be multiplied onto an existing primitive triple to form three new primitive triples. See "UAD Matrices" below for details of the descent.
tree_type => "UAD" (the default Y=40 | 14 | | | | 7 Y=24 | 5 | Y=20 | 3 | Y=12 | 2 13 | | 4 Y=4 | 1 | +-------------------------------------------------- X=3 X=15 X=20 X=35 X=45
The starting point is N=1 at X=3,Y=4 which is the well-known 3^2 + 4^2 = 5^2. From there further N=2, N=3, N=4 are derived, then three more from each of those, etc,
depth=0 depth=1 depth=2 depth=3 N=1 N=2..4 N=5..13 N=14... +-> 7,24 +-> 5,12 --+-> 55,48 | +-> 45,28 | | +-> 39,80 3,4 --+-> 21,20 --+-> 119,120 | +-> 77,36 | | +-> 33,56 +-> 15,8 --+-> 65,72 +-> 35,12
Counting N=1 as depth=0, each level has 3^depth many points and the first N of a level, which is tree_depth_to_n(), is at
tree_depth_to_n()
Nstart = 1 + (1 + 3 + 3^2 + ... + 3^(depth-1)) = (3^depth + 1) / 2
These N levels are like a mixed-radix representation of N where the high digit is binary and the digits below are ternary.
+--------+---------+---------+-- --+---------+ N = | binary | ternary | ternary | ... | ternary | +--------+---------+---------+-- --+---------+ 1 0,1,2 0,1,2 0,1,2
The high digit is non-zero and so is always 1. The number of ternary digits is the "depth" and their value (dropping the high binary 1) is the position within that level.
Taking the middle "A" matrix repeatedly gives 3,4 -> 21,20 -> 119,120 -> 697,696 etc which are the triples with legs A,B differing by 1 and therefore just above or below the X=Y leading diagonal. The N values are 1,3,9,27,etc = 3^depth.
Taking the lower "D" matrix repeatedly gives 3,4 -> 15,8 -> 35,12 -> 63,16, etc which is the primitives among a sequence of triples known to the ancients,
A = k^2-1 B = 2*k C = k^2+1 so C=A+2
When k is even these are primitive. (If k is odd then A and B are both even, ie. a common factor of 2, so not primitive.) These points are the last of each level, so at N=(3^(depth+1)-1)/2 which is tree_depth_to_n_end().
tree_depth_to_n_end()
Taking the upper "U" matrix repeatedly gives 3.4 -> 5,12 -> 7,24 -> 9,40 etc with C=B+1. These are the first of each level so at Nstart described above. The resulting triples are a sequence known to Pythagoras
x^2 + ((x^2-1)/2)^2 = ((x^2+1)/2)^2 so A=x, B=(x^2-1)/2, C=(x^2+1)/2
and described by Fibonacci in his Liber Quadratorum (Book of Squares) in terms of sums of odd numbers
g = any odd square = A^2 B^2 = 1 + 3 + 5 + ... + g-2 = ((g-1)/2)^2 C^2 = 1 + 3 + 5 + ... + g-2 + g = ((g+1)/2)^2 so C^2 = A^2 + B^2 eg. g=25=A^2 B^2=((25-1)/2)^2=144 so A=5,B=12
The geometric interpretation is that an existing square of side B is extended by a "gnomon" around two sides (cf Math::PlanePath::Corner) making a new larger square of side C=B+1. If the length of the gnomon is a square then the total area is the sum of two squares.
*****gnomon******** gnomon length an odd square = A^2 +---------------+ * | | * so new bigger square area | square | * C^2 = A^2 + B^2 | side B | * | | * +---------------+ *
Option tree_type => "FB" selects the Fibonacci boxes tree by H. Lee Price
tree_type => "FB"
"The Pythagorean Tree: A New Species", 2008, http://arxiv.org/abs/0809.4324
This is based on expressing triples in certain "Fibonacci boxes" with a box of four values q',q,p,p' having p=q+q' and p'=p+q so each is the sum of the preceding two in a fashion similar to the Fibonacci sequence. A box where p and q have no common factor corresponds to a primitive triple. See "PQ Coordinates" and "FB Transformations" below.
tree_type => "FB" Y=40 | 5 | | | | 17 Y=24 | 4 | | 8 | Y=12 | 2 6 | | 3 Y=4 | 1 | +---------------------------------------------- X=3 X=15 x=21 X=35
For a given box three transformations can be applied to go to new boxes corresponding to new primitive triples. This visits all and only primitive triples, but in a different order to the UAD above.
The first point N=1 is again at X=3,Y=4, from which three further points N=2,3,4 are derived, then three more from each of those, etc.
N=1 N=2..4 N=5..13 N=14... +-> 9,40 +-> 5,12 --+-> 35,12 | +-> 11,60 | | +-> 21,20 3,4 --+-> 15,8 --+-> 55,48 | +-> 39,80 | | +-> 13,84 +-> 7,24 --+-> 63,16 +-> 15,112
Option coordinates => 'AC' gives the A and C legs of each triple as X=A,Y=C.
coordinates => 'AC'
coordinates => "AC" 85 | 122 10 | | 73 | 6 | 65 | 11 40 61 | 41 | | 7 | | 41 | 14 | 13 35 | | 3 25 | 5 | 17 | 4 13 | 2 | Y=5 | 1 | +------------------------------------------- X=3 7 9 21 35 45 55 63 77
Since A<C the coordinates are X<Y so all above the X=Y diagonal. The repeated "D" triples described with C=A+2 are just above the diagonal.
For the FB tree the set of points visited is the same, but with a different N numbering.
tree_type => "FB", coordinates => "AC" 85 | 11 35 | | 73 | 9 | 65 | 23 12 61 | 7 | | 17 | | 41 | 5 | 6 35 | | 8 25 | 4 | 17 | 3 13 | 2 | Y=5 | 1 | +------------------------------------------- X=3 7 9 21 35 45 55 63 77
Option coordinates => 'BC' gives the B and C legs of each triple as X=B,Y=C. This is the B=even and C=long legs of all primitive triples. This combination makes 45-degree straight lines.
coordinates => 'BC'
coordinates => "BC" 101 | 121 97 | 12 | 89 | 8 85 | 10 122 | | 73 | 6 | 65 | 40 11 61 | 41 | | 7 | | 41 | 14 | 13 35 | | 3 25 | 5 | 17 | 4 13 | 2 | Y=5 | 1 | +-------------------------------------------------- X=4 12 24 40 60 84
Since B<C the coordinates are X<Y and therefore above the X=Y leading diagonal. N=1,2,5,14,41,etc along the X=Y-1 diagonal are the repeated "U" matrix triples C=B+1 which are the start of each depth level described above.
tree_type => "FB", coordinates => "BC" 101 | 15 97 | 50 | 89 | 10 85 | 35 11 | | 73 | 9 | 65 | 12 23 61 | 7 | | 17 | | 41 | 5 | 6 35 | | 8 25 | 4 | 17 | 3 13 | 2 | Y=5 | 1 | +---------------------------------------------- X=4 12 24 40 60 84
These B,C points fall on 45-degree straight lines going up from X=Y-1. This occurs because a primitive triple A,B,C with A odd and B even can be written
A^2 = C^2 - B^2 A^2 = (C+B)*(C-B) gcd(C+B,C-B)=1 because gcd(A,B)=1 and hence gcd(B,C)=1 so C+B = s^2 C = (s^2 + t^2)/2 C-B = t^2 B = (s^2 - t^2)/2 s = odd integer >= 3 t = odd integer, and s > t >= 1 gcd(s,t)=1 so that gcd(C+B,C-B)=1
If t=1 then any odd integer s gives a primitive triple. These are the C=B+1 repeated "U" matrix described above. Further values of t give primitive triples so long as they have no common factor with s. As t increases the B,C coordinate combination makes a line at a 45-degree angle upwards,
C = B - t^2 so X+Y = t^2 opposite diagonal
All primitive triples start from a C=B+1 for C=(s^2+1)/2, half an odd square, and go up from there for gcd(s,t)=1.
Primitive Pythagorean triples can be parameterized as follows for A odd and B even. (As per Diophantus, and anonymous Arabic manuscript for constraining to primitive triples.)
A = P^2 - Q^2 B = 2*P*Q C = P^2 + Q^2 with P > Q >= 1, one odd, one even, and no common factor P = sqrt((C+A)/2) Q = sqrt((C-A)/2)
The first P=2,Q=1 is the triple A=3,B=4,C=5.
Option coordinates => 'PQ' gives these X=P,Y=Q as the returned X,Y coordinates (for either tree_type). Because P>Q>=1 the values fall in the eighth of the plane below the X=Y diagonal,
coordinates => 'PQ'
tree_type
coordinates => "PQ" 11 | * 10 | * 9 | * 8 | * * 7 | * * * 6 | * * 5 | * * * 4 | * * * * 3 | * * * 2 | * * * * * 1 | * * * * * * Y=0 | +------------------------ X=0 1 2 3 4 5 6 7 8 9 ...
The correspondence between P,Q and A,B means the trees visit all P,Q pairs with no common factor and one of them even. There's other ways to iterate through such coprime P,Q, and those methods would generate triples too, in a different order from the trees here.
The letters P and Q here are a little bit arbitrary. They're often written m,n or u,v, but don't want "n" to be confused that with the N point numbering or "u" to be confused with the U matrix in UAD.
See "FUNCTIONS" in Math::PlanePath for behaviour common to all path classes.
$path = Math::PlanePath::PythagoreanTree->new ()
$path = Math::PlanePath::PythagoreanTree->new (tree_type => $str, coordinates => $str)
Create and return a new path object. The tree_type option can be
"UAD" (the default) "FB"
The coordinates option can be
coordinates
"AB" (the default) "AC" "BC" "PQ"
$n = $path->n_start()
Return 1, the first N in the path.
($x,$y) = $path->n_to_xy ($n)
Return the X,Y coordinates of point number $n on the path. Points begin at 1 and if $n<1 then the return is an empty list.
$n
$n<1
$n = $path->xy_to_n ($x,$y)
Return the point number for coordinates $x,$y. If there's nothing at $x,$y then return undef.
$x,$y
undef
The return is undef if $x,$y is not a primitive Pythagorean triple, or for the PQ option if $x,$y doesn't satisfy the PQ constraints described above ("PQ Coordinates").
($n_lo, $n_hi) = $path->rect_to_n_range ($x1,$y1, $x2,$y2)
Return a range of N values which occur in a rectangle with corners at $x1,$y1 and $x2,$y2. The range is inclusive.
$x1
$y1
$x2
$y2
Both trees go off into large X,Y coordinates while yet to finish values close to the origin which means the N range for a rectangle can be quite large. For UAD $n_hi is roughly 3**max(x/2), or for FB smaller at roughly 3**log2(x).
$n_hi
3**max(x/2)
3**log2(x)
@n_children = $path->tree_n_children($n)
Return the three children of $n, or an empty list if $n < 1 (ie. before the start of the path).
$n < 1
This is simply 3*$n-1, 3*$n, 3*$n+1. This is like appending an extra ternary digit to go to the next level, but onto N+1 rather than N and then adjusting back.
3*$n-1, 3*$n, 3*$n+1
$num = $path->tree_n_num_children($n)
Return 3, since every node has three children, or return undef if $n<1 (ie. before the start of the path).
$n_parent = $path->tree_n_parent($n)
Return the parent node of $n, or undef if $n <= 1 (the top of the tree).
$n <= 1
This is simply floor(($n+1)/3), reversing the tree_n_children() calculation.
floor(($n+1)/3)
tree_n_children()
$depth = $path->tree_n_to_depth($n)
Return the depth of node $n, or undef if there's no point $n. The top of the tree at N=1 is depth=0, then its children depth=1, etc.
The structure of the tree with 3 nodes per point means the depth is floor(log3(2N-1)), so for example N=5 through N=13 all have depth=2.
$n = $path->tree_depth_to_n($depth)
$n = $path->tree_depth_to_n_end($depth)
Return the first or last N at tree level $depth in the path, or undef if nothing at that depth or not a tree. The top of the tree is depth=0.
$depth
The UAD matrices are
/ 1 2 2 \ U = | -2 -1 -2 | \ 2 2 3 / / 1 2 2 \ A = | 2 1 2 | \ 2 3 3 / / -1 -2 -2 \ D = | 2 1 2 | \ 2 2 3 /
They're multiplied on the right of an (A,B,C) vector, for example
(3, 4, 5) * U = (5, 12, 13)
Internally the code uses P,Q and calculates A,B at the end as necessary. The UAD transformations in P,Q coordinates are
U P -> 2P-Q Q -> P A P -> 2P+Q Q -> P D P -> P+2Q Q -> Q unchanged
The advantage of P,Q for the calculation is that it's 2 values instead of 3. The transformations could be written as 2x2 matrix multiplications if desired, but explicit steps are enough for the code.
Repeatedly applying "U" gives
U 2P-Q, P U^2 3P-2Q, 2P-Q U^3 4P-3Q, 3P-2Q ... U^k (k+1)P-kQ, kP-(k-1)Q = P+k(P-Q), Q+k*(P-Q)
If there's a run of k many high zeros in the Nrem = N-Nstart position in the level then they can be applied to the initial P=2,Q=1 as
U^k P=k+2, Q=k+1 start for k high zeros
The FB tree is calculated in P,Q and converted to A,B at the end as necessary. Its three transformations are
K1 P -> P+Q Q -> 2Q K2 P -> 2P Q -> P-Q K3 P -> 2P Q -> P+Q
Price's paper shows rearrangements of a set of four values q',q,p,p', but just the p and q are enough for the calculation.
xy_to_n() works in P,Q coordinates. An A,B input is converted per the formulas in "PQ Coordinates" above. A P,Q point can be reversed up the UAD tree to its parent point
xy_to_n()
if P > 3Q reverse "D" P -> P-2Q Q -> unchanged if P > 2Q reverse "A" P -> Q Q -> P-2Q otherwise reverse "U" P -> Q Q -> 2Q-P
This gives a ternary digit 2, 1, 0 respectively from low to high. Those plus a high "1" bit make N. The number of steps is the "depth" level.
If at any stage P,Q doesn't satisfy P>Q, one odd, the other even, then it means the original point, whether it was an A,B or a P,Q, was not a primitive triple. For a primitive triple the endpoint is always P=2,Q=1.
After converting to P,Q as necessary, a P,Q point can be reversed up the FB tree to its parent
if P odd reverse K1 P -> P-Q (so Q even) Q -> Q/2 if Q < P/2 reverse K2 P -> P/2 Q -> P/2 - Q otherwise reverse K3 P -> P/2 Q -> Q - P/2
This is a little like the binary greatest common divisor algorithm, but designed for one value odd and the other even. Like the UAD ascent above if at any stage P,Q doesn't satisfy P>Q, one odd, the other even, then the initial point wasn't a primitive triple.
For the UAD tree, the smallest A,B within each level is found at the topmost "U" steps for the smallest A or the bottom-most "D" steps for the smallest B. For example in the table above of level=2 N=5..13 the smallest A is the top A=7,B=24, and the smallest B is in the bottom A=35,B=12. In general
Amin = 2*level + 1 Bmin = 4*level
In P,Q coordinates the same topmost line is the smallest P and bottom-most the smallest Q. The values are
Pmin = level+1 Qmin = 1
The fixed Q=1 arises from the way the "D" transformation sends Q->Q unchanged, so every level includes a Q=1. This means if you ask what range of N is needed to cover all Q < someQ then there isn't one, only a P < someP has an N to go up to.
For the FB tree, the smallest A,B within each level is found in the topmost two final positions. For example in the table above of level=2 N=5..13 the smallest A is in the top A=9,B=40, and the smallest B is in the next row A=35,B=12. In general,
Amin = 2^level + 1 Bmin = 2^level + 4
In P,Q coordinates a Q=1 is found in that second row which is the minimum B, and the smallest P is found by taking K1 steps half-way then a K2 step, then K1 steps for the balance. This is a slightly complicated
Pmin = / 3*2^(k-1) + 1 if even level = 2*k \ 2^(k+1) + 1 if odd level = 2*k+1 Q = 1
The fixed Q=1 arises from the K1 steps giving
P = 2 + 1+2+4+8+...+2^(level-2) = 2 + 2^(level-1) - 1 = 2^(level-1) + 1 Q = 2^(level-1) followed by K2 step Q -> P-Q = 1
As for the UAD above this means small Q's always remain no matter how big N gets, only a P range determines an N range.
Entries in Sloane's Online Encyclopedia of Integer Sequences related to this path include,
http://oeis.org/A007051 (etc) A007051 N start of depth=n, (3^n+1)/2, ie. tree_depth_to_n() A003462 N end of depth=n-1, (3^n-1)/2, ie. tree_depth_to_n_end()
Math::PlanePath, Math::PlanePath::Hypot, Math::PlanePath::RationalsTree, Math::PlanePath::CoprimeColumns
http://user42.tuxfamily.org/math-planepath/index.html
Copyright 2011, 2012 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/>.
To install Math::PlanePath, copy and paste the appropriate command in to your terminal.
cpanm
cpanm Math::PlanePath
CPAN shell
perl -MCPAN -e shell install Math::PlanePath
For more information on module installation, please visit the detailed CPAN module installation guide.