File Coverage

blib/lib/Math/PlanePath/ComplexPlus.pm
Criterion Covered Total %
statement 35 132 26.5
branch 5 40 12.5
condition 6 15 40.0
subroutine 10 19 52.6
pod 10 10 100.0
total 66 216 30.5


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              
20             # math-image --path=ComplexPlus --all --scale=5
21             #
22             # math-image --path=ComplexPlus --expression='i<128?i:0' --output=numbers --size=132x40
23             #
24             # Realpart:
25             # math-image --path=ComplexPlus,realpart=2 --expression='i<50?i:0' --output=numbers --size=180
26             #
27             # Arms:
28             # math-image --path=ComplexPlus,arms=2 --expression='i<64?i:0' --output=numbers --size=79
29              
30              
31              
32             package Math::PlanePath::ComplexPlus;
33 1     1   9025 use 5.004;
  1         9  
34 1     1   6 use strict;
  1         2  
  1         49  
35             #use List::Util 'max';
36             *max = \&Math::PlanePath::_max;
37              
38 1     1   6 use vars '$VERSION', '@ISA';
  1         2  
  1         104  
39             $VERSION = 127;
40              
41 1     1   666 use Math::PlanePath;
  1         3  
  1         50  
42             @ISA = ('Math::PlanePath');
43             *_divrem_mutate = \&Math::PlanePath::_divrem_mutate;
44              
45             use Math::PlanePath::Base::Generic
46 1         46 'is_infinite',
47 1     1   6 'round_nearest';
  1         2  
48             use Math::PlanePath::Base::Digits
49 1         67 'round_up_pow',
50             'digit_split_lowtohigh',
51 1     1   461 'digit_join_lowtohigh';
  1         2  
52              
53             # uncomment this to run the ### lines
54             #use Smart::Comments;
55              
56              
57 1     1   8 use constant n_start => 0;
  1         2  
  1         94  
58 1         1243 use constant parameter_info_array =>
59             [ { name => 'realpart',
60             display => 'Real Part',
61             type => 'integer',
62             default => 1,
63             minimum => 1,
64             width => 2,
65             description => 'Real part r in the i+r complex base.',
66             },
67             { name => 'arms',
68             share_key => 'arms_2',
69             display => 'Arms',
70             type => 'integer',
71             minimum => 1,
72             maximum => 2,
73             default => 1,
74             width => 1,
75             description => 'Arms',
76             when_name => 'realpart',
77             when_value => '1',
78             },
79 1     1   7 ];
  1         2  
80              
81             # b=i+r
82             # theta = atan(1/r)
83             sub x_negative_at_n {
84 0     0 1 0 my ($self) = @_;
85 0 0       0 if ($self->{'realpart'} == 1) { return 8; }
  0         0  
86 0         0 return $self->{'norm'} ** _ceil((2*atan2(1,1)) / atan2(1,$self->{'realpart'}));
87             }
88             sub y_negative_at_n {
89 0     0 1 0 my ($self) = @_;
90 0 0       0 if ($self->{'realpart'} == 1) { return 32; }
  0         0  
91 0         0 return $self->{'norm'} ** _ceil((4*atan2(1,1)) / atan2(1,$self->{'realpart'}));
92             }
93             sub _ceil {
94 0     0   0 my ($x) = @_;
95 0         0 my $int = int($x);
96 0 0       0 return ($x > $int ? $int+1 : $int);
97             }
98              
99             sub absdx_minimum {
100 0     0 1 0 my ($self) = @_;
101 0 0       0 return ($self->{'realpart'} == 1
102             ? 0 # i+1 N=1 dX=0,dY=1
103             : 1); # i+r otherwise always diff
104             }
105             # use constant dir_maximum_dxdy => (0,0); # supremum, almost full way
106              
107             sub turn_any_straight {
108 0     0 1 0 my ($self) = @_;
109 0         0 return ($self->{'realpart'} != 1); # realpart=1 never straight
110             }
111              
112              
113             #------------------------------------------------------------------------------
114             sub new {
115 5     5 1 1308 my $self = shift->SUPER::new(@_);
116              
117 5         11 my $realpart = $self->{'realpart'};
118 5 100 66     21 if (! defined $realpart || $realpart < 1) {
119 4         11 $self->{'realpart'} = $realpart = 1;
120             }
121 5         11 $self->{'norm'} = $realpart*$realpart + 1;
122              
123 5         8 my $arms = $self->{'arms'};
124 5 100 66     23 if (! defined $arms || $arms <= 0 || $realpart != 1) { $arms = 1; }
  4 50 66     6  
125 0         0 elsif ($arms > 2) { $arms = 2; }
126 5         9 $self->{'arms'} = $arms;
127              
128 5         14 return $self;
129             }
130              
131             sub n_to_xy {
132 0     0 1 0 my ($self, $n) = @_;
133             ### ComplexPlus n_to_xy(): $n
134              
135 0 0       0 if ($n < 0) { return; }
  0         0  
136 0 0       0 if (is_infinite($n)) { return ($n,$n); }
  0         0  
137              
138             {
139 0         0 my $int = int($n);
  0         0  
140             ### $int
141             ### $n
142 0 0       0 if ($n != $int) {
143 0         0 my ($x1,$y1) = $self->n_to_xy($int);
144 0         0 my ($x2,$y2) = $self->n_to_xy($int+$self->{'arms'});
145 0         0 my $frac = $n - $int; # inherit possible BigFloat
146 0         0 my $dx = $x2-$x1;
147 0         0 my $dy = $y2-$y1;
148 0         0 return ($frac*$dx + $x1, $frac*$dy + $y1);
149             }
150 0         0 $n = $int; # BigFloat int() gives BigInt, use that
151             }
152              
153 0         0 my $realpart = $self->{'realpart'};
154 0         0 my $norm = $self->{'norm'};
155             ### $norm
156             ### $realpart
157              
158             # for i+1, arm=0 start X=0,Y=0, arm=1 start X=0,Y=1
159 0         0 my $x = 0;
160 0         0 my $y = _divrem_mutate ($n, $self->{'arms'});
161              
162             # for i+1, arm=0 start dX=1,dY=0, arm=1 start dX=-1,dY=0
163 0         0 my $dy = ($n * 0); # 0, inheriting bignum from $n
164 0 0       0 my $dx = ($y ? -1 : 1) + $dy; # inheriting bignum from $n
165              
166 0         0 foreach my $digit (digit_split_lowtohigh($n,$norm)) {
167             ### at: "$x,$y digit=$digit dxdy=$dx,$dy"
168              
169 0         0 $x += $dx * $digit;
170 0         0 $y += $dy * $digit;
171              
172             # multiply i+r, ie. (dx,dy) = (dx + i*dy)*(i+$realpart)
173 0         0 ($dx,$dy) = ($realpart*$dx - $dy, $dx + $realpart*$dy);
174             }
175              
176             ### final: "$x,$y"
177 0         0 return ($x,$y);
178             }
179              
180             sub xy_to_n {
181 0     0 1 0 my ($self, $x, $y) = @_;
182             ### ComplexPlus xy_to_n(): "$x, $y"
183              
184 0         0 $x = round_nearest ($x);
185 0         0 $y = round_nearest ($y);
186              
187 0         0 my $realpart = $self->{'realpart'};
188             {
189 0         0 my $rx = $realpart*$x;
  0         0  
190 0         0 my $ry = $realpart*$y;
191 0         0 foreach my $overflow ($rx+$ry, $rx-$ry) {
192 0 0       0 if (is_infinite($overflow)) { return $overflow; }
  0         0  
193             }
194             }
195              
196 0         0 my $orig_x = $x;
197 0         0 my $orig_y = $y;
198              
199 0         0 my $norm = $self->{'norm'};
200 0         0 my $zero = ($x * 0 * $y); # inherit bignum 0
201 0         0 my @n; # digits low to high
202              
203 0         0 my $prev_x = 0;
204 0         0 my $prev_y = 0;
205 0   0     0 while ($x || $y) {
206 0         0 my $neg_y = $x - $y*$realpart;
207 0         0 my $digit = $neg_y % $norm;
208             ### at: "$x,$y n=$n digit $digit"
209              
210 0         0 push @n, $digit;
211 0         0 $x -= $digit;
212 0         0 $neg_y -= $digit;
213              
214             ### assert: ($neg_y % $norm) == 0
215             ### assert: (($x * $realpart + $y) % $norm) == 0
216              
217             # divide i+r = mul (i-r)/(i^2 - r^2)
218             # = mul (i-r)/-norm
219             # is (i*y + x) * (i-realpart)/-norm
220             # x = [ x*-realpart - y ] / -norm
221             # = [ x*realpart + y ] / norm
222             # y = [ y*-realpart + x ] / -norm
223             # = [ y*realpart - x ] / norm
224             #
225 0         0 ($x,$y) = (($x*$realpart+$y)/$norm, -$neg_y/$norm);
226              
227 0 0 0     0 if ($x == $prev_x && $y == $prev_y) {
228 0         0 last;
229             }
230 0         0 $prev_x = $x;
231 0         0 $prev_y = $y;
232             }
233              
234             ### final: "$x,$y n=$n cf arms $self->{'arms'}"
235              
236 0 0       0 if ($y) {
237 0 0       0 if ($self->{'arms'} > 1) {
238             ### not on first arm, re-run as: -$orig_x, 1-$orig_y
239 0         0 local $self->{'arms'} = 1;
240 0         0 my $n = $self->xy_to_n(-$orig_x,1-$orig_y);
241 0 0       0 if (defined $n) {
242 0         0 return 1 + 2*$n; # 180 degrees
243             }
244             }
245             ### X,Y not visited
246 0         0 return undef;
247             }
248              
249 0         0 my $n = digit_join_lowtohigh (\@n, $norm, $zero);
250 0 0       0 if ($self->{'arms'} > 1) {
251 0         0 $n *= 2;
252             }
253 0         0 return $n;
254             }
255              
256             # not exact
257             sub rect_to_n_range {
258 0     0 1 0 my ($self, $x1,$y1, $x2,$y2) = @_;
259             ### ComplexPlus rect_to_n_range(): "$x1,$y1 $x2,$y2"
260              
261 0         0 my $xm = max(abs($x1),abs($x2));
262 0         0 my $ym = max(abs($y1),abs($y2));
263 0         0 my $n_hi = ($xm*$xm + $ym*$ym) * $self->{'arms'};
264 0 0       0 if ($self->{'realpart'} == 1) {
265 0         0 $n_hi *= 16; # 2**4
266             }
267 0         0 return (0, int($n_hi));
268             }
269              
270             #------------------------------------------------------------------------------
271             # levels
272              
273             sub level_to_n_range {
274 9     9 1 660 my ($self, $level) = @_;
275 9         32 return (0, $self->{'norm'}**$level * $self->{'arms'} - 1);
276             }
277             sub n_to_level {
278 0     0 1   my ($self, $n) = @_;
279 0 0         if ($n < 0) { return undef; }
  0            
280 0 0         if (is_infinite($n)) { return $n; }
  0            
281 0           $n = round_nearest($n);
282 0           _divrem_mutate ($n, $self->{'arms'});
283 0           my ($pow, $exp) = round_up_pow ($n+1, $self->{'norm'});
284 0           return $exp;
285             }
286              
287             #------------------------------------------------------------------------------
288             1;
289             __END__