File Coverage

blib/lib/Math/PlanePath/DragonMidpoint.pm
Criterion Covered Total %
statement 125 161 77.6
branch 25 34 73.5
condition 14 26 53.8
subroutine 15 22 68.1
pod 9 9 100.0
total 188 252 74.6


line stmt bran cond sub pod time code
1             # Copyright 2011, 2012, 2013, 2014, 2015, 2016, 2017, 2018, 2019 Kevin Ryde
2              
3             # This file is part of Math-PlanePath.
4             #
5             # Math-PlanePath is free software; you can redistribute it and/or modify
6             # it under the terms of the GNU General Public License as published by the
7             # Free Software Foundation; either version 3, or (at your option) any later
8             # version.
9             #
10             # Math-PlanePath is distributed in the hope that it will be useful, but
11             # WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
12             # or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
13             # for more details.
14             #
15             # You should have received a copy of the GNU General Public License along
16             # with Math-PlanePath. If not, see .
17              
18              
19             # math-image --path=DragonMidpoint --lines --scale=20
20             # math-image --path=DragonMidpoint --all --output=numbers_dash
21              
22             # A006466 contfrac 2*sum( 1/2^(2^n)), 1 and 2 only
23             # a(5n) recurrence ...
24             # 1,1,1,1, 2,
25             # 1,1,1,1,1,1,1, 2,
26             # 1,1,1,1, 2,
27             # 1,1,1,1, 2,
28             # 1, 2,
29             # 1,1,1,1, 2,
30             # 1,1,1,1,1,1,1, 2,
31             # 1,1,1,1, 2,
32             # 1, 2,
33             # 1,1,1,1,1,1,1, 2,
34             # 1,1,1,1, 2,
35             # 1, 2,
36             # 1,1,1,1, 2,
37             # 1,1,1,1, 2,
38             # 1,1,1,1,1,1,1, 2,
39             # 1,1,1,1, 2,
40             # 1, 2,
41             # 1,1,1,1,1,1,1, 2,
42             # 1,1,1,1, 2,
43             # 1,1,1,1, 2,
44             # 1, 2
45             # A076214 in decimal
46             #
47             # A073097 number of 4s - 6s - 2s - 1 is -1,0,1
48             # A081769 positions of 2s
49             # A073088 cumulative total multiples of 4 roughly, hence (4n-3-cum)/2
50             #
51             # A088435 (contfrac+1)/2 of sum(k>=1,1/3^(2^k)).
52             # A007404 in decimal
53             #
54              
55              
56             package Math::PlanePath::DragonMidpoint;
57 5     5   9466 use 5.004;
  5         24  
58 5     5   31 use strict;
  5         10  
  5         132  
59 5     5   25 use List::Util 'min'; # 'max'
  5         9  
  5         456  
60             *max = \&Math::PlanePath::_max;
61              
62 5     5   38 use vars '$VERSION', '@ISA';
  5         11  
  5         326  
63             $VERSION = 127;
64 5     5   730 use Math::PlanePath;
  5         9  
  5         138  
65 5     5   988 use Math::PlanePath::Base::NSEW;
  5         11  
  5         198  
66             @ISA = ('Math::PlanePath::Base::NSEW',
67             'Math::PlanePath');
68              
69             use Math::PlanePath::Base::Generic
70 5         279 'is_infinite',
71 5     5   28 'round_nearest';
  5         10  
72             use Math::PlanePath::Base::Digits
73 5         305 'bit_split_lowtohigh',
74             'digit_join_lowtohigh',
75 5     5   512 'round_up_pow';
  5         9  
76             *_divrem_mutate = \&Math::PlanePath::_divrem_mutate;
77              
78             # uncomment this to run the ### lines
79             # use Smart::Comments;
80              
81              
82             # whole plane when arms==4
83 5     5   1333 use Math::PlanePath::DragonCurve;
  5         13  
  5         162  
84              
85              
86 5     5   42 use constant n_start => 0;
  5         11  
  5         417  
87 5         7189 use constant parameter_info_array => [ { name => 'arms',
88             share_key => 'arms_4',
89             display => 'Arms',
90             type => 'integer',
91             minimum => 1,
92             maximum => 4,
93             default => 1,
94             width => 1,
95             description => 'Arms',
96 5     5   34 } ];
  5         10  
97              
98             {
99             my @x_negative_at_n = (undef, 6,5,2,2);
100             sub x_negative_at_n {
101 0     0 1 0 my ($self) = @_;
102 0         0 return $x_negative_at_n[$self->{'arms'}];
103             }
104             }
105             {
106             my @y_negative_at_n = (undef, 27,19,11,7);
107             sub y_negative_at_n {
108 0     0 1 0 my ($self) = @_;
109 0         0 return $y_negative_at_n[$self->{'arms'}];
110             }
111             }
112             {
113             my @_UNDOCUMENTED__dxdy_list_at_n = (undef, 9, 9, 5, 3);
114             sub _UNDOCUMENTED__dxdy_list_at_n {
115 0     0   0 my ($self) = @_;
116 0         0 return $_UNDOCUMENTED__dxdy_list_at_n[$self->{'arms'}];
117             }
118             }
119              
120              
121             #------------------------------------------------------------------------------
122              
123             sub new {
124 15     15 1 2621 my $self = shift->SUPER::new(@_);
125 15   100     100 $self->{'arms'} = max(1, min(4, $self->{'arms'} || 1));
126 15         36 return $self;
127             }
128              
129             # sub n_to_xy {
130             # my ($self, $n) = @_;
131             # ### DragonMidpoint n_to_xy(): $n
132             #
133             # if ($n < 0) { return; }
134             # if (is_infinite($n)) { return ($n, $n); }
135             #
136             # {
137             # my $int = int($n);
138             # if ($n != $int) {
139             # my ($x1,$y1) = $self->n_to_xy($int);
140             # my ($x2,$y2) = $self->n_to_xy($int+$self->{'arms'});
141             # my $frac = $n - $int; # inherit possible BigFloat
142             # my $dx = $x2-$x1;
143             # my $dy = $y2-$y1;
144             # return ($frac*$dx + $x1, $frac*$dy + $y1);
145             # }
146             # $n = $int; # BigFloat int() gives BigInt, use that
147             # }
148             #
149             # my ($x1,$y1) = Math::PlanePath::DragonCurve->n_to_xy($n);
150             # my ($x2,$y2) = Math::PlanePath::DragonCurve->n_to_xy($n+1);
151             #
152             # my $dx = $x2-$x1;
153             # my $dy = $y2-$y1;
154             # return ($x1+$y1 + ($dx+$dy-1)/2,
155             # $y1-$x1 + ($dy-$dx+1)/2);
156             # }
157              
158             sub n_to_xy {
159 264     264 1 19132 my ($self, $n) = @_;
160             ### DragonMidpoint n_to_xy(): $n
161              
162 264 50       647 if ($n < 0) { return; }
  0         0  
163 264 50       649 if (is_infinite($n)) { return ($n, $n); }
  0         0  
164              
165 264         468 my $frac;
166             {
167 264         348 my $int = int($n);
  264         408  
168 264         412 $frac = $n - $int; # inherit possible BigFloat
169 264         394 $n = $int; # BigFloat int() gives BigInt, use that
170             }
171 264         380 my $zero = ($n * 0); # inherit bignum 0
172              
173             # arm as initial rotation
174 264         656 my $rot = _divrem_mutate ($n, $self->{'arms'});
175              
176             ### $arms
177             ### rot from arm: $rot
178             ### $n
179              
180             # ENHANCE-ME: sx,sy just from len,len
181 264         679 my @digits = bit_split_lowtohigh($n);
182 264         491 my @sx;
183             my @sy;
184              
185             {
186 264         367 my $sx = $zero + 1;
  264         392  
187 264         372 my $sy = -$sx;
188 264         518 foreach (@digits) {
189 1898         2672 push @sx, $sx;
190 1898         2457 push @sy, $sy;
191              
192             # (sx,sy) + rot+90(sx,sy)
193 1898         3047 ($sx,$sy) = ($sx - $sy,
194             $sy + $sx);
195             }
196             }
197              
198             ### @digits
199 264         388 my $rev = 0;
200 264         368 my $x = $zero;
201 264         345 my $y = $zero;
202 264         348 my $above_low_zero = 0;
203              
204 264         599 for (my $i = $#digits; $i >= 0; $i--) { # high to low
205 1898         2672 my $digit = $digits[$i];
206 1898         2487 my $sx = $sx[$i];
207 1898         2450 my $sy = $sy[$i];
208             ### at: "$x,$y $digit side $sx,$sy"
209             ### $rot
210              
211 1898 100       3163 if ($rot & 2) {
212 879         1185 $sx = -$sx;
213 879         1171 $sy = -$sy;
214             }
215 1898 100       3115 if ($rot & 1) {
216 939         1519 ($sx,$sy) = (-$sy,$sx);
217             }
218             ### rotated side: "$sx,$sy"
219              
220 1898 100       2990 if ($rev) {
221 929 100       1392 if ($digit) {
222 492         649 $x -= $sy;
223 492         950 $y += $sx;
224             ### rev add to: "$x,$y next is still rev"
225             } else {
226 437         676 $above_low_zero = $digits[$i+1];
227 437         550 $rot ++;
228 437         859 $rev = 0;
229             ### rev rot, next is no rev ...
230             }
231             } else {
232 969 100       1472 if ($digit) {
233 555         706 $rot ++;
234 555         761 $x += $sx;
235 555         685 $y += $sy;
236 555         1081 $rev = 1;
237             ### plain add to: "$x,$y next is rev"
238             } else {
239 414         817 $above_low_zero = $digits[$i+1];
240             }
241             }
242             }
243              
244             # Digit above the low zero is the direction of the next turn, 0 for left,
245             # 1 for right.
246             #
247             ### final: "$x,$y rot=$rot above_low_zero=".($above_low_zero||0)
248              
249 264 100       524 if ($rot & 2) {
250 126         206 $frac = -$frac; # rotate 180
251 126         199 $x -= 1;
252             }
253 264 100       579 if (($rot+1) & 2) {
254             # rot 1 or 2
255 155         208 $y += 1;
256             }
257 264 100 100     672 if (!($rot & 1) && $above_low_zero) {
258 52         91 $frac = -$frac;
259             }
260 264         464 $above_low_zero ^= ($rot & 1);
261 264 100       520 if ($above_low_zero) {
262 145         219 $y = $frac + $y;
263             } else {
264 119         175 $x = $frac + $x;
265             }
266              
267             ### rotated return: "$x,$y"
268 264         983 return ($x,$y);
269             }
270              
271             # or tables arithmetically,
272             #
273             # my $ax = ((($x+1) ^ ($y+1)) >> 1) & 1;
274             # my $ay = (($x^$y) >> 1) & 1;
275             # ### assert: $ax == - $yx_adj_x[$y%4]->[$x%4]
276             # ### assert: $ay == - $yx_adj_y[$y%4]->[$x%4]
277             #
278             my @yx_adj_x = ([0,1,1,0],
279             [1,0,0,1],
280             [1,0,0,1],
281             [0,1,1,0]);
282              
283             my @yx_adj_y = ([0,0,1,1],
284             [0,0,1,1],
285             [1,1,0,0],
286             [1,1,0,0]);
287              
288             # arm $x $y 2 | 1 Y=1
289             # 0 0 0 3 | 0 Y=0
290             # 1 0 1 ----+----
291             # 2 -1 1 X=-1 X=0
292             # 3 -1 0
293             my @xy_to_arm = ([0, # x=0,y=0
294             1], # x=0,y=1
295             [3, # x=-1,y=0
296             2]); # x=-1,y=1
297              
298             sub xy_to_n {
299 512     512 1 11789 my ($self, $x, $y) = @_;
300             ### DragonMidpoint xy_to_n(): "$x, $y"
301              
302 512         1096 $x = round_nearest($x);
303 512         1048 $y = round_nearest($y);
304              
305 512         770 { my $overflow = abs($x)+abs($y)+2;
  512         855  
306 512 50       1395 if (is_infinite($overflow)) { return $overflow; }
  0         0  
307             }
308 512         1041 my $zero = ($x * 0 * $y);
309 512         740 my @nbits; # low to high
310              
311 512   100     1696 while ($x < -1 || $x > 0 || $y < 0 || $y > 1) {
      100        
      100        
312 7018         5165 my $y4 = $y % 4;
313 7018         4735 my $x4 = $x % 4;
314 7018         5007 my $ax = $yx_adj_x[$y4]->[$x4];
315 7018         4870 my $ay = $yx_adj_y[$y4]->[$x4];
316              
317             ### at: "$x,$y n=$n axy=$ax,$ay bit=".($ax^$ay)
318              
319 7018         5353 push @nbits, $ax^$ay;
320              
321 7018         4677 $x -= $ax;
322 7018         4667 $y -= $ay;
323             ### assert: ($x+$y)%2 == 0
324 7018         29540 ($x,$y) = (($x+$y)/2, # rotate -45 and divide sqrt(2)
325             ($y-$x)/2);
326             }
327              
328             ### final: "xy=$x,$y"
329              
330 512         863 my $arm = $xy_to_arm[$x]->[$y];
331             ### $arm
332 512         1304 my $arms_count = $self->arms_count;
333 512 100       1114 if ($arm >= $arms_count) {
334 112         253 return undef;
335             }
336              
337 400 100       838 if ($arm & 1) {
338             ### flip ...
339 109         274 @nbits = map {$_^1} @nbits;
  729         1271  
340             }
341              
342 400         1119 return digit_join_lowtohigh(\@nbits, 2, $zero) * $arms_count + $arm;
343             }
344              
345             #------------------------------------------------------------------------------
346             # xy_is_visited()
347              
348             sub xy_is_visited {
349 0     0 1 0 my ($self, $x, $y) = @_;
350             return ($self->{'arms'} >= 4
351 0   0     0 || _xy_to_arm($x,$y) < $self->{'arms'});
352             }
353              
354             # return arm number 0,1,2,3
355             sub _xy_to_arm {
356 0     0   0 my ($x, $y) = @_;
357             ### DragonMidpoint _xy_to_arm(): "$x, $y"
358              
359 0         0 $x = round_nearest($x);
360 0         0 $y = round_nearest($y);
361              
362 0         0 { my $overflow = abs($x)+abs($y)+2;
  0         0  
363 0 0       0 if (is_infinite($overflow)) { return $overflow; }
  0         0  
364             }
365              
366 0   0     0 while ($x < -1 || $x > 0 || $y < 0 || $y > 1) {
      0        
      0        
367 0         0 my $y4 = $y % 4;
368 0         0 my $x4 = $x % 4;
369 0         0 $x -= $yx_adj_x[$y4]->[$x4];
370 0         0 $y -= $yx_adj_y[$y4]->[$x4];
371              
372             ### assert: ($x+$y)%2 == 0
373 0         0 ($x,$y) = (($x+$y)/2, # rotate -45 and divide sqrt(2)
374             ($y-$x)/2);
375             }
376 0         0 return $xy_to_arm[$x]->[$y];
377             }
378              
379             #------------------------------------------------------------------------------
380              
381             # not exact
382             sub rect_to_n_range {
383 94     94 1 8341 my ($self, $x1,$y1, $x2,$y2) = @_;
384             ### DragonMidpoint rect_to_n_range(): "$x1,$y1 $x2,$y2 arms=$self->{'arms'}"
385 94         160 $x1 = abs($x1);
386 94         126 $x2 = abs($x2);
387 94         131 $y1 = abs($y1);
388 94         147 $y2 = abs($y2);
389 94         225 my $xmax = int(max($x1,$x2));
390 94         186 my $ymax = int(max($y1,$y2));
391             return (0,
392 94         311 ($xmax*$xmax + $ymax*$ymax + 1) * $self->{'arms'} * 5);
393             }
394              
395             # sub rect_to_n_range {
396             # my ($self, $x1,$y1, $x2,$y2) = @_;
397             # ### DragonMidpoint rect_to_n_range(): "$x1,$y1 $x2,$y2"
398             #
399             # return Math::PlanePath::DragonCurve->rect_to_n_range
400             # (sqrt(2)*$x1, sqrt(2)*$y1, sqrt(2)*$x2, sqrt(2)*$y2);
401             # }
402              
403             #------------------------------------------------------------------------------
404              
405             sub level_to_n_range {
406 0     0 1   my ($self, $level) = @_;
407 0           return (0, 2**$level * $self->{'arms'} - 1);
408             }
409             sub n_to_level {
410 0     0 1   my ($self, $n) = @_;
411 0 0         if ($n < 0) { return undef; }
  0            
412 0 0         if (is_infinite($n)) { return $n; }
  0            
413 0           $n = round_nearest($n);
414 0           _divrem_mutate ($n, $self->{'arms'});
415 0           my ($pow, $exp) = round_up_pow ($n+1, 2);
416 0           return $exp;
417             }
418              
419             #------------------------------------------------------------------------------
420             1;
421             __END__