Math::PlanePath::KochCurve -- horizontal Koch curve
use Math::PlanePath::KochCurve; my $path = Math::PlanePath::KochCurve->new; my ($x, $y) = $path->n_to_xy (123);
This is an integer version of the self-similar curve by Helge von Koch going along the X axis and making triangular excursions upwards.
8 3 / \ 6---- 7 9----10 18-... 2 \ / \ 2 5 11 14 17 1 / \ / \ / \ / 0----1 3---- 4 12----13 15----16 <- Y=0 ^ X=0 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19
The replicating shape is the initial N=0 to N=4,
* / \ *---* *---*
which is rotated and repeated 3 times in the same pattern to give sections N=4 to N=8, N=8 to N=12, and N=12 to N=16. Then that N=0 to N=16 is itself replicated three times at the angles of the base pattern, and so on infinitely.
The X,Y coordinates are arranged on a square grid using every second point, per "Triangular Lattice" in Math::PlanePath. The result is flattened triangular segments with diagonals at a 45 degree angle.
Each replication adds 3 copies of the existing points and is thus 4 times bigger, so if N=0 to N=4 is reckoned as level 1 then a given replication level goes from
Nstart = 0 Nlevel = 4^level (inclusive)
Each replication is 3 times the width. The initial N=0 to N=4 figure is 6 wide and in general a level runs from
Xstart = 0 Xlevel = 2*3^level at N=Nlevel
The highest Y is 3 times greater at each level similarly. The peak is at the midpoint of each level,
Npeak = (4^level)/2 Ypeak = 3^level Xpeak = 3^level
It can be seen that the N=6 point backtracks horizontally to the same X as the start of its section N=4 to N=8. This happens in the further replications too and is the maximum extent of the backtracking.
The Nlevel is multiplied by 4 to get the end of the next higher level. The same 4*N can be applied to all points N=0 to N=Nlevel to get the same shape but a factor of 3 bigger X,Y coordinates. The in-between points 4*N+1, 4*N+2 and 4*N+3 are then new finer structure in the higher level.
Koch conceived the curve as having a fixed length and infinitely fine structure, making it continuous everywhere but differentiable nowhere. The code here can be pressed into use for that sort of construction for a given level of granularity by scaling
X/3^level Y/3^level
which makes it a fixed 2 wide by 1 high. Or for unit-side equilateral triangles then apply further factors 1/2 and sqrt(3)/2, as noted in "Triangular Lattice" in Math::PlanePath.
(X/2) / 3^level (Y*sqrt(3)/2) / 3^level
See "FUNCTIONS" in Math::PlanePath for behaviour common to all path classes.
$path = Math::PlanePath::KochCurve->new ()
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.
$n
$n < 0
Fractional positions give an X,Y position along a straight line between the integer positions.
($n_lo, $n_hi) = $path->rect_to_n_range ($x1,$y1, $x2,$y2)
The returned range is exact, meaning $n_lo and $n_hi are the smallest and biggest in the rectangle.
$n_lo
$n_hi
$n = $path->n_start()
Return 0, the first N in the path.
The sequence of turns made by the curve is straightforward. The curve always turns either +60 degrees or -120 degrees, it never goes straight ahead. In the base 4 representation of N the lowest non-zero digit gives the turn
low digit turn --------- ------------ 1 +60 degrees (left) 2 -120 degrees (right) 3 +60 degrees (left)
For example N=8 is 20 base 4, so lowest nonzero "2" means turn -120 degrees for the next segment.
When the least significant digit is non-zero it determines the turn, making the base N=0 to N=4 shape. When the low digit is zero then the next level up is in control, eg. N=0,4,8,12,16, making a turn where the base shape repeats.
The cumulative turn at a given N can be found by counting digits 1 and 2 in base 4.
direction = 60 * ((count of 1 digits in base 4) - (count of 2 digits in base 4)) degrees
For example N=11 is 23 in base 4, so 60*(0-1) = -60 degrees.
In this formula the count of 1s and 2s can go past 360 degrees, representing a spiralling around which occurs at progressively higher replication levels. The direction can be taken mod 360 degrees, or the count mod 6, for a direction 0 to 5 or as desired.
An easy over-estimate of the N values in a rectangle can be had from the Xlevel formula above. If Xlevel>rectangleX then Nlevel is past the rectangle extent.
X = 2*3^level floorlevel = floor log_base_3(X/2) Nhi = 4^(floorlevel+1) - 1
For example a rectangle extending to X=13 has floorlevel = floor(log3(13/2))=1 and Nhi=4^(1+1)-1=15.
The rounding-down of the log3 ensures a point such as X=18 which is the first in the next Nlevel will give that next level. So floorlevel=log3(18/2)=2 and Nhi=4^(2+1)-1=63.
The worst case for this over-estimate is when rectangleX==Xlevel, ie. just into the next level. In that case Nhi is almost a factor of 4 bigger than it needs to be.
The exact Nlo and Nhi in a rectangle can be found by searching along the curve. Nlo searches forward from the origin N=0. Nhi searches backward from the Nlevel over-estimate described above.
At a given digit position in the prospective N the sub-part of the curve comprising the lower digits has a certain triangular extent. If it's outside the target rectangle then step to the next digit value, and to the next of the digit above when past digit=3 (or below digit=0 when searching backwards).
There's six possible rotations for the curve sub-part. In the following "o" is the start and the surrounding lines show the triangular extent. There's just four curve parts shown in each, but these triangles bound a sub-curve of any level.
rot=0 -+- +-----------------+ -- -- - .-+-* *-+-o - -- * -- -- \ / -- -- / \ -- -- * -- - o-+-* *-+-. - -- -- +-----------------+ rot=3 -+- rot=1 +---------+ rot=4 /+ | . / / | | / / / o| |*-+-* / / / | | \ / / * | | * / / \ | | / / / *-+-*| |o / / / | | / / . | +/ +---------+ +\ rot=2 +---------+ | \ \ o | |. \ \ \ | | \ \ \ *-+-*| | * \ \ / | | / \ \ * | |*-+-* \ \ \ | | \ \ \ .| | o \ rot=5 \ | +---------+ \+
The "." is the start of the next sub-curve. It belongs to the next digit value and so can be excluded if desired. For rot=0 and rot=3 this means simply shortening the X range permitted, or for rot=1 and rot=4 similarly the Y range. For rot=2 and rot=5 it would require a separate test and doesn't matter very much.
Tight sub-part extent checking reduces the sub-parts which are examined, but it works perfectly well with a looser check, such as a square box for the sub-curve extents. Doing that might be easier if the target region was not a rectangle but some trickier shape.
The Koch curve is in Sloane's Online Encyclopedia of Integer Sequences in various forms,
http://oeis.org/A035263 (etc) A035263 -- turn 1=left,0=right, by morphism A096268 -- turn 0=left,1=right, by morphism A029883 -- turn +/-1=left,0=right, Thue-Morse first differences A089045 -- turn +/-1=left,0=right, by +/- something A003159 -- N positions of left turns, ending even number 0 bits A036554 -- N positions of right turns, ending odd number 0 bits
Math::PlanePath, Math::PlanePath::PeanoCurve, Math::PlanePath::HilbertCurve, Math::PlanePath::KochPeaks, Math::PlanePath::KochSnowflakes, Math::PlanePath::CCurve
Math::Fractal::Curve
http://user42.tuxfamily.org/math-planepath/index.html
Copyright 2011, 2012 Kevin Ryde
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.