NAME
Math::PlanePath::SierpinskiCurve  Sierpinski curve
SYNOPSIS
use Math::PlanePath::SierpinskiCurve;
my $path = Math::PlanePath::SierpinskiCurve>new (arms => 2);
my ($x, $y) = $path>n_to_xy (123);
DESCRIPTION
This is an integer version of the selfsimilar curve by Waclaw Sierpinski traversing the plane by right triangles. The default is a single arm of the curve in an eighth of the plane.
10  3132
 / \
9  30 33
  
8  29 34
 \ /
7  2526 28 35 3738
 / \ / \ / \
6  24 27 36 39
  
5  23 20 43 40
 \ / \ / \ /
4  78 2221 19 44 4241 55...
 / \ / \ /
3  6 9 18 45 54
     
2  5 10 17 46 53
 \ / \ / \
1  12 4 11 1314 16 47 4950 52
 / \ / \ / \ / \ / \ /
Y=0  . 0 3 12 15 48 51

+
^
X=0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16
The tiling it represents is
/
/\
/  \
/  \
/ 7 8 \
/ \  / \
/ \  / \
/ 6 \/ 9 \
/\
/\ 5 /\ 10 /\
/  \ /  \ /  \
/  \ /  \ /  \
/ 1 2 X 4 11 X 1314 \
/ \  / \  / \  / \ ...
/ \  / \  / \  / \
/ 0 \/ 3 \/ 12 \/ 15 \

The points are on a square grid with integer X,Y. 4 points are used in each 3x3 block. In general a point is used if
X%3==1 or Y%3==1 but not both
which means
((X%3)+(Y%3)) % 2 == 1
The X axis N=0,3,12,15,48,etc are all the integers which use only digits 0 and 3 in base 4. For example N=51 is 303 base4. Or equivalently the values all have doubled bits in binary, for example N=48 is 110000 binary. (Compare the CornerReplicate
which also has these values along the X axis.)
Level Ranges
Counting the N=0 point as level=0, and with each level being 4 copies of the previous, the levels end at
Nlevel = 4^level  1 = 0, 3, 15, ...
Xlevel = 3*2^level  2 = 1, 4, 10, ...
Ylevel = 0
For example level=2 is Nlevel = 2^(2*2)1 = 15 at X=3*2^22 = 10.
Doubling a level is the middle of the next level and is the top of the triangle in that next level.
Ntop = 2*4^level  1 = 1, 7, 31, ...
Xtop = 3*2^level  1 = 2, 5, 11, ...
Ytop = 3*2^level  2 = Xlevel = 1, 4, 10, ...
For example doubling level=2 is Ntop = 2*4^21 = 31 at X=3*2^21 = 11 and Y=3*2^22 = 10.
The factor of 3 arises from the three steps which make up the N=0,1,2,3 section. The Xlevel width grows as
Xlevel(1) = 3
Xlevel(level) = 2*Xwidth(level1) + 3
which dividing out the factor of 3 is 2*w+1, giving 2^k1 (in binary a left shift and bring in a new 1 bit).
Notice too the Nlevel points as a fraction of the triangular area Xlevel*(Xlevel1)/2 gives the 4 out of 9 points filled,
FillFrac = Nlevel / (Xlevel*(Xlevel1)/2)
> 4/9
Arms
The optional arms
parameter can draw multiple curves, each advancing successively. For example 2 arms,
arms => 2 ...

11  33 39 57 63
 / \ / \ / \ /
10  31 3537 41 55 5961 62...
 \ / \ /
9  29 43 53 60
    
8  27 45 51 58
 / \ / \
7  25 2119 4749 5052 56
 \ / \ / \ /
6  23 17 48 54
  
5  9 15 46 40
 / \ / \ / \
4  7 1113 1416 4442 38
 \ / \ /
3  5 12 18 36
    
2  3 10 20 34
 / \ / \
1  1 24 8 22 2628 32
 / \ / \ / \ /
Y=0  0 6 24 30

+
^
X=0 1 2 3 4 5 6 7 8 9 10 11
The N=0 point is at X=1,Y=0 (in all arms forms) so that the second arm is within the first quadrant.
1 to 8 arms can be done this way. For example 8 arms are
arms => 8
... ... 6
 
58 34 33 57 5
\ / \ / \ /
...59 5042 26 25 4149 56... 4
\ / \ /
51 18 17 48 3
   
43 10 9 40 2
/ \ / \
35 1911 2 1 816 32 1
\ / \ / \ /
27 3 . 0 24 < Y=0
28 4 7 31 1
/ \ / \ / \
36 2012 5 6 1523 39 2
\ / \ /
44 13 14 47 3
   
52 21 22 55 4
/ \ / \
...60 5345 29 30 4654 63... 5
/ \ / \ / \
61 37 38 62 6
 
... ... 7
^
7 6 5 4 3 2 1 X=0 1 2 3 4 5 6
The middle "." is the origin X=0,Y=0. It would be more symmetrical to make the origin the middle of the eight arms, at X=0.5,Y=0.5 in the above, but that would give fractional X,Y values. Apply an offset X+0.5,Y+0.5 to centre it if desired.
Spacing
The optional diagonal_spacing
and straight_spacing
can increase the space between points diagonally or vertically+horizontally. The default for each is 1.
straight_spacing => 2
diagonal_spacing => 1
7  8
/ \
6 9
 
 
 
5 10 ...
\ / \
1  2 4 11 13  14 16
/ \ / \ / \ /
0 3 12 15
X=0 1 2 3 4 5 6 7 8 9 10 11 12 13 ...
The effect is only to spread the points. The straight lines are both horizontal and vertical so when they're stretched the curve remains on a 45 degree angle in an eighth of the plane.
In the level formulas above the "3" factor becomes 2*d+s, effectively being the N=0 to N=3 section sized as d+s+d.
d = diagonal_spacing
s = straight_spacing
Xlevel = (2d+s)*(2^level  1) + 1
Xtop = (2d+s)*2^(level1)  d  s + 1
Ytop = (2d+s)*2^(level1)  d  s
Closed Curve
Sierpinski's original conception was a closed curve filling a unit square by ever greater selfsimilar detail,
/\_/\ /\_/\ /\_/\ /\_/\
\ / \ / \ / \ /
       
/ _ \_/ _ \ / _ \_/ _ \
\/ \ / \/ \/ \ / \/
   
/\_/ _ \_/\ /\_/ _ \_/\
\ / \ / \ / \ /
       
/ _ \ / _ \_/ _ \ / _ \
\/ \/ \/ \ / \/ \/ \/
 
/\_/\ /\_/ _ \_/\ /\_/\
\ / \ / \ / \ /
       
/ _ \_/ _ \ / _ \_/ _ \
\/ \ / \/ \/ \ / \/
   
/\_/ _ \_/\ /\_/ _ \_/\
\ / \ / \ / \ /
       
/ _ \ / _ \ / _ \ / _ \
\/ \/ \/ \/ \/ \/ \/ \/
The code here might be pressed into use for this by drawing a mirror image of the curve N=0 through Nlevel. Or using the arms=>2
form N=0 to N=4^level  1, inclusive, and joining up the ends.
The curve is also usually conceived as scaling down by quarters. This can be had with straight_spacing => 2
and then an offset to X+1,Y+1 to centre in a 4*2^level square
Koch Curve Midpoints
The replicating structure is the same as the Koch curve (Math::PlanePath::KochCurve) in that the curve repeats four times to make the next level.
The Sierpinski curve points are midpoints of a Koch curve of 90 degree angles with a unit gap between verticals.
Koch Curve Koch Curve
90 degree angles, unit gap
/\  
/ \  
/ \  
   
Sierpinski curve points "*" as midpoints
 
7 8
 
6 9
5 10
     
1 2 4 11 13 14
     
0 3 12 15
Koch Curve Rounded
The Sierpinski curve in mirror image across the X=Y diagonal and rotated 45 degrees is pairs of points on the lines of the Koch curve 90 degree angles unit gap from above.
Sierpinski curve mirror image and turn 45 degrees
two points on each Koch line segment
15 16
 
14 17
1213 . . 1819
1110 . . 2120
3 4 9 22 27 28
     
2 5 8 23 26 29
01 . . 67 . . 2425 . . 3031
This is a kind of "rounded" form of the 90degree Koch, similar what DragonRounded
does for the DragonCurve
. Each 90 turn of the Koch curve is done by two turns of 45 degrees in the Sierpinski curve here, and each 180 degree turn in the Koch is two 90 degree turns here. So the Sierpinski turn sequence is pairs of the Koch turn sequence, as follows. The mirroring means a swap left<>right between the two.
N=1 2 3 4 5 6 7 8
Koch L R L L L R L R ...
N=1,2 3,4 5,6 7,8 9,10 11,12 13,14 15,16
Sierp R R L L R R R R R R L L R R L L ...
FUNCTIONS
See "FUNCTIONS" in Math::PlanePath for behaviour common to all path classes.
$path = Math::PlanePath::SierpinskiCurve>new ()
$path = Math::PlanePath::SierpinskiCurve>new (arms => $integer, diagonal_spacing => $integer, straight_spacing => $integer)

Create and return a new path object.
($x,$y) = $path>n_to_xy ($n)

Return the X,Y coordinates of point number
$n
on the path. Points begin at 0 and if$n < 0
then the return is an empty list.Fractional positions give an X,Y position along a straight line between the integer positions.
$n = $path>n_start()

Return 0, the first N in the path.
Level Methods
($n_lo, $n_hi) = $path>level_to_n_range($level)

Return
(0, 4**$level  1)
, or for multiple arms return(0, $arms * 4**$level  1)
.There are 4^level points in a level, or arms*4^level when multiple arms, numbered starting from 0.
FORMULAS
N to dX,dY
The curve direction at N even can be calculated from the base4 digits of N/2 in a fashion similar to the Koch curve ("N to Direction" in Math::PlanePath::KochCurve). Counting direction in eighths so 0=East, 1=NorthEast, 2=North, etc,
digit direction
 
0 0
1 2
2 2
3 0
direction = 1 + sum direction[base4 digits of N/2]
for N even
For example the direction at N=10 has N/2=5 which is "11" in base4, so direction = 1+(2)+(2) = 3 = southwest.
The 1 in 1+sum is direction northeast for N=0, then 2 or +2 for the digits follow the curve. For an odd arm the curve is mirrored and the sign of each digit direction is flipped, so a subtract instead of add,
direction
mirrored = 1  sum direction[base4 digits of N/2]
for N even
For odd N=2k+1 the direction at N=2k is calculated and then also the turn which is made from N=2k to N=2(k+1). This is similar to the Koch curve next turn ("N to Next Turn" in Math::PlanePath::KochCurve).
lowest non3 next turn
digit of N/2 (at N=2k+1,N=2k+2)
 
0 1 (right)
1 +2 (left)
2 1 (right)
Again the turn is in eighths, so 1 means 45 degrees (to the right). For example at N=14 has N/2=7 which is "13" in base4 so lowest non3 is "1" which is turn +2, so at N=15 and N=16 turn by 90 degrees left.
direction = 1 + sum direction[base4 digits of k]
+ if N odd then nextturn[lownon3 of k]
for N=2k or 2k+1
dX,dY = direction to 1,0 1,1 0,1 etc
For fractional N the same nextturn is applied to calculate the direction of the next segment, and combined with the integer dX,dY as per "N to dX,dY  Fractional" in Math::PlanePath.
N=2k or 2k+1 + frac
direction = 1 + sum direction[base4 digits of k]
if (frac != 0 or N odd)
turn = nextturn[lownon3 of k]
if N odd then direction += turn
dX,dY = direction to 1,0 1,1 0,1 etc
if frac!=0 then
direction += turn
next_dX,next_dY = direction to 1,0 1,1 0,1 etc
dX += frac*(next_dX  dX)
dY += frac*(next_dY  dY)
For the straight_spacing
and diagonal_spacing
options the dX,dY values are not units like dX=1,dY=0 but instead are the spacing amount, either straight or diagonal so
direction delta with spacing
 
0 dX=straight_spacing, dY=0
1 dX=diagonal_spacing, dY=diagonal_spacing
2 dX=0, dY=straight_spacing
3 dX=diagonal_spacing, dY=diagonal_spacing
etc
As an alternative, it's possible to take just base4 digits of N, without separate handling for the lowbit of N, but it requires an adjustment on the low base4 digit, and the next turn calculation for fractional N becomes hairier. A little state table could encode the cumulative and lowest whatever if desired, to take N by base4 digits high to low, or equivalently by bits high to low with an initial state based on high bit at an odd or even bit position.
OEIS
The Sierpinski curve is in Sloane's Online Encyclopedia of Integer Sequences as,
http://oeis.org/A039963 (etc)
A039963 turn 1=right,0=left, doubling the KochCurve turns
A081706 N1 of left turn positions
(first values 2,3 whereas N=3,4 here)
A127254 abs(dY), so 0=horizontal, 1=vertical or diagonal,
except extra initial 1
A081026 X at N=2^k, being successively 3*2^j1, 3*2^j
A039963 is numbered starting n=0 for the first turn, which is at the point N=1 in the path here.
SEE ALSO
Math::PlanePath, Math::PlanePath::SierpinskiCurveStair, Math::PlanePath::SierpinskiArrowhead, Math::PlanePath::KochCurve
HOME PAGE
http://user42.tuxfamily.org/mathplanepath/index.html
LICENSE
Copyright 2011, 2012, 2013, 2014, 2015, 2016, 2017, 2018, 2019, 2020 Kevin Ryde
MathPlanePath 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.
MathPlanePath 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 MathPlanePath. If not, see <http://www.gnu.org/licenses/>.