File Coverage

blib/lib/Math/PlanePath/TerdragonCurve.pm

Criterion	Covered	Total	%
statement	163	288	56.6
branch	33	100	33.0
condition	26	44	59.0
subroutine	28	44	63.6
pod	13	13	100.0
total	263	489	53.7

# Copyright 2011, 2012, 2013, 2014, 2015, 2016, 2017, 2018, 2019, 2020 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 .



# points singles A052548 2^n + 2
# points doubles A000918 2^n - 2
# points triples A028243 3^(n-1) - 2*2^(n-1) + 1     cf A[k] = 2*3^(k-1) - 2*2^(k-1)

# T(3*N)   = (w+1)*T(N)                dir(N)=w^(2*count1digits)
# T(3*N+1) = (w+1)*T(N) + 1*dir(N)
# T(3*N+2) = (w+1)*T(N) + w*dir(N)

# T(0*3^k + N)  =             T(N)
# T(1*3^k + N)  = 2^k   + w^2*T(N)    # rotate and offset
# T(2*3^k + N)  = w*2^k +     T(N)    # offset only



package Math::PlanePath::TerdragonCurve;
use 5.004;
use strict;
use List::Util 'first';
use List::Util 'min'; # 'max'
*max = \&Math::PlanePath::_max;

use Math::PlanePath;
*_divrem_mutate = \&Math::PlanePath::_divrem_mutate;

use Math::PlanePath::Base::Generic
  'is_infinite',
  'round_nearest',
  'xy_is_even';
use Math::PlanePath::Base::Digits
  'digit_split_lowtohigh',
  'digit_join_lowtohigh',
  'round_up_pow';

use vars '$VERSION', '@ISA';
$VERSION = 128;
@ISA = ('Math::PlanePath');

use Math::PlanePath::TerdragonMidpoint;

# uncomment this to run the ### lines
# use Smart::Comments;


use constant n_start => 0;
use constant parameter_info_array =>
  [ { name      => 'arms',
      share_key => 'arms_6',
      display   => 'Arms',
      type      => 'integer',
      minimum   => 1,
      maximum   => 6,
      default   => 1,
      width     => 1,
      description => 'Arms',
    } ];

{
  my @x_negative_at_n = (undef, 13, 5, 5, 6, 7, 8);
  sub x_negative_at_n {
    my ($self) = @_;
    return $x_negative_at_n[$self->{'arms'}];
  }
}
{
  my @y_negative_at_n = (undef, 159, 75, 20, 11, 9, 10);
  sub y_negative_at_n {
    my ($self) = @_;
    return $y_negative_at_n[$self->{'arms'}];
  }
}
sub dx_minimum {
  my ($self) = @_;
  return ($self->{'arms'} == 1 ? -1 : -2);
}
use constant dx_maximum => 2;
use constant dy_minimum => -1;
use constant dy_maximum => 1;

sub _UNDOCUMENTED__dxdy_list {
  my ($self) = @_;
  return ($self->{'arms'} == 1
          ? Math::PlanePath::_UNDOCUMENTED__dxdy_list_three()
          : Math::PlanePath::_UNDOCUMENTED__dxdy_list_six());
}
{
  my @_UNDOCUMENTED__dxdy_list_at_n = (undef, 4, 9, 13, 7, 8, 5);
  sub _UNDOCUMENTED__dxdy_list_at_n {
    my ($self) = @_;
    return $_UNDOCUMENTED__dxdy_list_at_n[$self->{'arms'}];
  }
}
use constant absdx_minimum => 1;
use constant dsumxy_minimum => -2; # diagonals
use constant dsumxy_maximum => 2;
use constant ddiffxy_minimum => -2;
use constant ddiffxy_maximum => 2;

# arms=1 curve goes at 0,120,240 degrees
# arms=2 second +60 to 60,180,300 degrees
# so when arms==1 dir maximum is 240 degrees
sub dir_maximum_dxdy {
  my ($self) = @_;
  return ($self->{'arms'} == 1
          ? (-1,-1)    # 0,2,4 only           South-West
          : ( 1,-1));  # rotated to 1,3,5 too South-East
}

use constant turn_any_straight => 0; # never straight


#------------------------------------------------------------------------------

sub new {
  my $self = shift->SUPER::new(@_);
  $self->{'arms'} = max(1, min(6, $self->{'arms'} || 1));
  return $self;
}

my @dir6_to_si = (1,0,0, -1,0,0);
my @dir6_to_sj = (0,1,0, 0,-1,0);
my @dir6_to_sk = (0,0,1, 0,0,-1);

sub n_to_xy {
  my ($self, $n) = @_;
  ### TerdragonCurve n_to_xy(): $n

  if ($n < 0) { return; }
  if (is_infinite($n)) { return ($n, $n); }

  my $zero = ($n * 0);  # inherit bignum 0

  my $i = 0;
  my $j = 0;
  my $k = 0;
  my $si = $zero;
  my $sj = $zero;
  my $sk = $zero;

  # initial rotation from arm number
  {
    my $int = int($n);
    my $frac = $n - $int;  # inherit possible BigFloat
    $n = $int;             # BigFloat int() gives BigInt, use that

    my $rot = _divrem_mutate ($n, $self->{'arms'});

    my $s = $zero + 1;  # inherit bignum 1
    if ($rot >= 3) {
      $s = -$s;         # rotate 180
      $frac = -$frac;
      $rot -= 3;
    }
    if ($rot == 0)    { $i = $frac; $si = $s; } # rotate 0
    elsif ($rot == 1) { $j = $frac; $sj = $s; } # rotate +60
    else              { $k = $frac; $sk = $s; } # rotate +120
  }

  foreach my $digit (digit_split_lowtohigh($n,3)) {
    ### at: "$i,$j,$k   side $si,$sj,$sk"
    ### $digit

    if ($digit == 1) {
      ($i,$j,$k) = ($si-$j, $sj-$k, $sk+$i);  # rotate +120 and add
    } elsif ($digit == 2) {
      $i -= $sk;   # add rotated +60
      $j += $si;
      $k += $sj;
    }

    # add rotated +60
    ($si,$sj,$sk) = ($si - $sk,
                     $sj + $si,
                     $sk + $sj);
  }

  ### final: "$i,$j,$k   side $si,$sj,$sk"
  ### is: (2*$i + $j - $k).",".($j+$k)

  return (2*$i + $j - $k, $j+$k);
}


# all even points when arms==6
sub xy_is_visited {
  my ($self, $x, $y) = @_;
  if ($self->{'arms'} == 6) {
    return xy_is_even($self,$x,$y);
  } else {
    return defined($self->xy_to_n($x,$y));
  }
}

sub xy_to_n {
  return scalar((shift->xy_to_n_list(@_))[0]);
}
sub xy_to_n_list {
  my ($self, $x,$y) = @_;
  ### TerdragonCurve xy_to_n_list(): "$x, $y"

  $x = round_nearest($x);
  $y = round_nearest($y);
  {
    # nothing at an odd point, and trap overflows in $x+$y dividing out b
    my $sum = abs($x) + abs($y);
    if (is_infinite($sum)) { return $sum; }  # infinity
    if ($sum % 2) { return; }
  }

  if ($x==0 && $y==0) {
    return 0 .. $self->{'arms'}-1;
  }

  my $arms_count = $self->arms_count;
  my $zero = ($x * 0 * $y); # inherit bignum 0

  my @n_list;
  foreach my $d (0,1,2) {
    my ($ndigits,$arm) = _xy_d_to_ndigits_and_arm($x,$y,$d);
    next if $arm >= $arms_count;
    my $odd = ($arm & 1);
    if ($odd) {
      @$ndigits = (map {2-$_} @$ndigits);
      ### flip to: $ndigits
    }
    push @n_list,
      (digit_join_lowtohigh($ndigits, 3, $zero) + $odd) * $arms_count + $arm;
  }

  ### @n_list
  return sort {$a<=>$b} @n_list;
}

my @x_to_digit = (0, 2, 1);  # digit = -X mod 3
my @digit_to_x = ([0,2,1],  [0,-1,-2],  [0,-1, 1]);
my @digit_to_y = ([0,0,1],  [0, 1, 0],  [0,-1,-1]);

# $d = 0,1,2 for segment leaving $x,$y at direction $d*120 degrees.
# For odd arms the digits are 0<->2 reversals.
sub _xy_d_to_ndigits_and_arm {
  my ($x,$y, $d) = @_;
  my @ndigits;
  my $arm;
  for (;;) {
    ### at: "$x,$y d=$d"
    if ($x==0 && $y==0) { $arm = 2*$d; last; }
    if ($d==0 && $x==-2 && $y==0) { $arm = 3; last; }
    if ($d==2 && $x==1  && $y==1) { $arm = 1; last; }
    if ($d==1 && $x==1  && $y==-1) { $arm = 5; last; }

    my $digit = $x_to_digit[$x%3];
    push @ndigits, $digit;

    if ($digit == 1) { $d = ($d-1) % 3; }
    $x -= $digit_to_x[$d]->[$digit];
    $y -= $digit_to_y[$d]->[$digit];

    ### $digit
    ### new d: $d
    ### subtract: "$digit_to_x[$d]->[$digit],$digit_to_y[$d]->[$digit] to $x,$y"

    # ### assert: ($x+$y) % 2 == 0
    # ### assert: $x % 3 == 0
    # ### assert: ($y-$x/3) % 2 == 0
    ($x,$y) = (($x+$y)/2,    # divide b = w6+1
               ($y-$x/3)/2);
  }
  ### $arm
  ### @ndigits
  return (\@ndigits, $arm);
}
# x+y*w3
# (x-y)+y*w3
# x/2 + y*sqrt3i/2
# sqrt3i/2 = w3+1/2
# x/2 + y*(w3+1/2) == 1/2*(x+y) + y*w3
# a = x+y = (x+3*y)/2
# GP-Test  my(x=0,y=0); (-x)%3 == 0
# GP-Test  my(x=2,y=0); (-x)%3 == 1
# GP-Test  my(x=1,y=1); (-x)%3 == 2

# minimum  -- no, not quite right
#
#                *----------*
#                 \
#                  \   *
#               *   \
#                    \
#          *----------*
#
# width = side/2
# minimum = side*sqrt(3)/2 - width
#         = side*(sqrt(3)/2 - 1)
#
# minimum 4/9 * 2.9^level roughly
# h = 4/9 * 2.9^level
# 2.9^level = h*9/4
# level = log(h*9/4)/log(2.9)
# 3^level = 3^(log(h*9/4)/log(2.9))
#         = h*9/4, but big bigger for log
#
# not exact
sub rect_to_n_range {
  my ($self, $x1,$y1, $x2,$y2) = @_;
  ### TerdragonCurve rect_to_n_range(): "$x1,$y1  $x2,$y2"
  my $xmax = int(max(abs($x1),abs($x2)));
  my $ymax = int(max(abs($y1),abs($y2)));
  return (0,
          ($xmax*$xmax + 3*$ymax*$ymax + 1)
          * 2
          * $self->{'arms'});
}

# direction
#
my @dir6_to_dx   = (2, 1,-1,-2, -1, 1);
my @dir6_to_dy   = (0, 1, 1, 0, -1,-1);
my @digit_to_nextturn = (2,-2);
sub n_to_dxdy {
  my ($self, $n) = @_;
  ### n_to_dxdy(): $n

  if ($n < 0) {
    return;  # first direction at N=0
  }
  if (is_infinite($n)) {
    return ($n,$n);
  }

  my $int = int($n);  # integer part
  $n -= $int;         # fraction part

  # initial direction from arm
  my $dir6 = _divrem_mutate ($int, $self->{'arms'});

  my @ndigits = digit_split_lowtohigh($int,3);
  $dir6 += 2 * scalar(grep {$_==1} @ndigits);  # count 1s for total turn
  $dir6 %= 6;
  my $dx = $dir6_to_dx[$dir6];
  my $dy = $dir6_to_dy[$dir6];

  if ($n) {
    # fraction part

    # find lowest non-2 digit, or zero if all 2s or no digits at all
    $dir6 += $digit_to_nextturn[ first {$_!=2} @ndigits, 0];
    $dir6 %= 6;
    $dx += $n*($dir6_to_dx[$dir6] - $dx);
    $dy += $n*($dir6_to_dy[$dir6] - $dy);
  }
  return ($dx, $dy);
}


#-----------------------------------------------------------------------------
# eg. arms=5 0 .. 5*3^k    step by 5s
#            1 .. 5*3^k+1  step by 5s
#            4 .. 5*3^k+4  step by 5s
#
sub level_to_n_range {
  my ($self, $level) = @_;
  return (0,  (3**$level + 1) * $self->{'arms'} - 1);
}
sub n_to_level {
  my ($self, $n) = @_;
  if ($n < 0) { return undef; }
  if (is_infinite($n)) { return $n; }
  $n = round_nearest($n);
  _divrem_mutate ($n, $self->{'arms'});
  my ($pow, $exp) = round_up_pow ($n, 3);
  return $exp;
}

#-----------------------------------------------------------------------------
# right boundary N

# mixed radix binary, ternary
# no 11, 12, 20
# 11 -> 21, including low digit
# run of 11111 becomes 22221
# low to high 1 or 0 <- 0   cannot 20 can 10 00
#             2 or 0 <- 1   cannot 11 can 21 01
#             2 or 0 <- 2   cannot 12 can 02 22
sub _UNDOCUMENTED__right_boundary_i_to_n {
  my ($self, $i) = @_;
  my @digits = _digit_split_mix23_lowtohigh($i);
  for ($i = $#digits; $i >= 1; $i--) {   # high to low
    if ($digits[$i] == 1 && $digits[$i-1] != 0) {
      $digits[$i] = 2;
    }
  }
  return digit_join_lowtohigh(\@digits, 3, $i*0);

  # {
  #   for (my $i = 0; $i < $#digits; $i++) {   # low to high
  #     if ($digits[$i+1] == 1 && ($digits[$i] == 1 || $digits[$i] == 2)) {
  #       $digits[$i+1] = 2;
  #     }
  #   }
  #   return digit_join_lowtohigh(\@digits,3);
  # }
}

# Return a list of digits, low to high, which is a mixed radix
# representation low digit ternary and the rest binary.
sub _digit_split_mix23_lowtohigh {
  my ($n) = @_;
  if ($n == 0) {
    return ();
  }
  my $low = _divrem_mutate($n,3);
  return ($low, digit_split_lowtohigh($n,2));
}

{
  # disallowed digit pairs $disallowed[high][low]
  my @disallowed;
  $disallowed[1][1] = 1;
  $disallowed[1][2] = 1;
  $disallowed[2][0] = 1;

  sub _UNDOCUMENTED__n_segment_is_right_boundary {
    my ($self, $n) = @_;
    if (is_infinite($n)) { return 0; }
    unless ($n >= 0) { return 0; }
    $n = int($n);

    # no boundary when arms=6, right boundary is only in arm 0
    {
      my $arms = $self->{'arms'};
      if ($arms == 6) { return 0; }
      if (_divrem_mutate($n,$arms)) { return 0; }
    }

    my $prev = _divrem_mutate($n,3);
    while ($n) {
      my $digit = _divrem_mutate($n,3);
      if ($disallowed[$digit][$prev]) {
        return 0;
      }
      $prev = $digit;
    }
    return 1;
  }
}

#-----------------------------------------------------------------------------
# left boundary N


# mixed 0,1, 2, 10, 11, 12, 100, 101, 102, 110, 111, 112, 1000, 1001, 1002, 1010, 1011, 1012, 1100, 1101, 1102,
# vals  0,1,12,120,121,122,1200,1201,1212,1220,1221,1222,12000,12001,12012,12120,12121,12122,12200,12201,12212,
{
  my @_UNDOCUMENTED__left_boundary_i_to_n = ([0,2],  # 0
                                             [0,2],  # 1
                                             [1,2]); # 2
  sub _UNDOCUMENTED__left_boundary_i_to_n {
    my ($self, $i, $level) = @_;
    ### _UNDOCUMENTED__left_boundary_i_to_n(): $i
    ### $level

    if (defined $level && $level < 0) {
      if ($i <= 2) {
        return $i;
      }
      $i += 2;
    }

    my @digits = _digit_split_mix23_lowtohigh($i);
    ### @digits

    if (defined $level) {
      if ($level >= 0) {
        if (@digits > $level) {
          ### beyond given level ...
          return undef;
        }
        # pad for $level, total $level many digits
        push @digits, (0) x ($level - scalar(@digits));
      } else {
        ### union all levels ...
        pop @digits;
        if ($digits[-1]) {
          push @digits, 0;     # high 0,1  or 0,2 when i=3
        } else {
          $digits[-1] = 1;     # high   1
        }
      }
    } else {
      ### infinite curve, an extra high 0 ...
      push @digits, 0;
    }
    ### @digits

    my $prev = $digits[0];
    foreach my $i (1 .. $#digits) {
      $prev = $digits[$i] = $_UNDOCUMENTED__left_boundary_i_to_n[$prev][$digits[$i]];
    }
    ### ternary: @digits
    return digit_join_lowtohigh(\@digits, 3, $i*0);
  }
}

{
  # disallowed digit pairs $disallowed[high][low]
  my @disallowed;
  $disallowed[0][2] = 1;
  $disallowed[1][0] = 1;
  $disallowed[1][1] = 1;

  sub _UNDOCUMENTED__n_segment_is_left_boundary {
    my ($self, $n, $level) = @_;
    ### _UNDOCUMENTED__n_segment_is_left_boundary(): $n
    ### $level

    if (is_infinite($n)) { return 0; }
    unless ($n >= 0) { return 0; }
    $n = int($n);

    if (defined $level && $level == 0) {
      ### level 0 curve, N=0 is only segment: ($n == 0)
      return ($n == 0);
    }

    {
      my $arms = $self->{'arms'};
      if ($arms == 6) {
        return 0;
      }
      my $arm = _divrem_mutate($n,$arms);
      if ($arm != $arms-1) {
        return 0;
      }
    }

    my $prev = _divrem_mutate($n,3);
    if (defined $level) { $level -= 1; }

    for (;;) {
      if (defined $level && $level == 0) {
        ### end of level many digits, must be N < 3**$level
        return ($n == 0);
      }
      last unless $n;

      my $digit = _divrem_mutate($n,3);
      if ($disallowed[$digit][$prev]) {
        return 0;
      }
      if (defined $level) { $level -= 1; }
      $prev = $digit;
    }

    return ((defined $level && $level < 0)   # union all levels
            || ($prev != 2));                # not high 2 otherwise
  }

  sub _UNDOCUMENTED__n_segment_is_any_left_boundary {
    my ($self, $n) = @_;
    my $prev = _divrem_mutate($n,3);
    while ($n) {
      my $digit = _divrem_mutate($n,3);
      if ($disallowed[$digit][$prev]) {
        return 0;
      }
      $prev = $digit;
    }
    return 1;
  }

  # sub left_boundary_n_pred {
  #   my ($n) = @_;
  #   my $n3 = '0' . Math::BaseCnv::cnv($n,10,3);
  #   return ($n3 =~ /02|10|11/ ? 0 : 1);
  # }
}
sub _UNDOCUMENTED__n_segment_is_boundary {
  my ($self, $n, $level) = @_;
  return $self->_UNDOCUMENTED__n_segment_is_right_boundary($n)
    || $self->_UNDOCUMENTED__n_segment_is_left_boundary($n,$level);
}

1;
__END__