Math-PlanePath-Toothpick
view release on metacpan or search on metacpan
devel/lib/Math/PlanePath/ToothpickTreeByCells-oeis.t view on Meta::CPAN
}
return \@got;
});
# sub full_from_wedge {
# my ($n) = @_;
# return 2*wedge(n) + 2*a(n+1) - 4n - 1 for n>0. - N. J. A.
#
# }
# use Memoize;
# BEGIN { Memoize::memoize('wedge_formula'); }
#------------------------------------------------------------------------------
# A170890 - unwedge_down_W total cells
MyOEIS::compare_values
(anum => 'A170890',
max_count => $max_count,
func => sub {
my ($count) = @_;
my $path = Math::PlanePath::ToothpickTreeByCells->new (parts => 'unwedge_down_W');
devel/one-of-eight.pl view on Meta::CPAN
my $n = $c->tree_depth_to_n($depth);
my $n2 = Math::BaseCnv::cnv($n,10,2);
my $pn = $p->tree_depth_to_n($depth);
my $calc = (16*4**$k + 24*$k - 7) / 9;
my $delta = $n - $prev_n;
my $d2 = Math::BaseCnv::cnv($delta,10,2);
printf "%5d path=%8d formula=%8d cells=%8d %20s\n",
$depth, $pn, $calc, $n, $n2;
# printf "%5d %8d %20s\n", $depth, $delta, $d2;
$prev_n = $n;
}
exit 0;
}
{
# rect_to_n_range() on 2^k
require Math::PlanePath::OneOfEight;
devel/toothpick.pl view on Meta::CPAN
# = 2^k - 1 - k
#
require Math::PlanePath::ToothpickTreeByCells;
require Math::BaseCnv;
my $cells = Math::PlanePath::ToothpickTreeByCells->new (parts => 'octant');
my $path = Math::PlanePath::ToothpickTree->new (parts => 'octant');
for (my $depth = 2; $depth <= 36; $depth++) {
my $depth = $depth-2;
my $n = $cells->tree_depth_to_n($depth);
my $n2 = Math::BaseCnv::cnv($n,10,2);
my $f = oct_pow_by_formula($depth);
my $p = $path->tree_depth_to_n($depth);
my $diff = ($p == $n ? '' : '***');
my $d2 = $depth + 2;
if (is_pow2($d2)) { print "\n"; }
print "$depth $d2 n=$n $n2 p=$p$diff f=$f\n";
}
sub oct_pow_by_formula {
my ($depth) = @_;
$depth += 2; # parts=4 basis
if ($depth <= 2) { return 0; }
if ($depth == 4) { return 2; }
$depth /= 2;
return (4*oct_pow_by_formula($depth-2)
+ 4
- $depth
);
}
sub is_pow2 {
my ($n) = @_;
while ($n > 1) {
if ($n & 1) {
return 0;
}
lib/Math/PlanePath/OneOfEight.pm view on Meta::CPAN
my $powmult = $mult;
if ($rem <= 1) {
if ($rem == 0) {
### rem=0, oct(pow) only ...
} else { # $rem == 1
### rem=1, oct(pow)+1 ...
$n += $mult;
}
} else {
### formula ...
# oct(pow+rem) = oct(pow)
# + oct(rem+1)
# + 2*oct(rem)
# - floor(log2(rem+1))
# - rem - 3
my $rem1 = $rem + 1;
{
my ($lpow,$lexp) = round_down_pow ($rem1, 2);
$n -= ($lexp + $rem + 3)*$mult;
lib/Math/PlanePath/OneOfEight.pm view on Meta::CPAN
} else {
### rem=1, grow 3 ...
$added += 3 * $mult;
}
} else {
my $rem1 = $rem + 1;
if ($rem1 == $pow) {
### rem+1=pow, no lower part, 3/2 of pow ...
$added += ($pow/2) * (3*$mult);
} else {
### formula ...
# oadd(pow+rem) = oadd(rem+1) + 2*oadd(rem)
# + (is_pow2($rem+2) ? -2 : -1)
# upper/lower diagonal overlap, and no log2_extras in lower
$added -= (_is_pow2($rem+2) ? 2*$mult : $mult);
if (@new_depth && $new_depth[-1] == $rem1) {
### merge into previously pushed new_depth ...
# print "rem=$rem ",join(',',@new_depth),"\n";
$new_mult[-1] += $mult;
lib/Math/PlanePath/OneOfEight.pm view on Meta::CPAN
- floor(log2(rem+1)) # no pow2 points in lower
- rem # unduplicate diagonal upper/lower
- 3 # unduplicate centre points
oct(rem)+oct(rem+1) of upper and lower would count their common diagonal
twice, hence "-rem" being the length of that diagonal. The "centre" point
at X=pow,Y=pow is repeated by each of extend, upper, lower so "-2" to count
just once, and the X=pow+1,Y=pow point is repeated by extend and upper, so
"-1" to count it just once.
The 2*total(rem)+total(rem+1) in the formula is the same recurrence as the
toothpick pattern and the approach there can calculate it as a set of
pending depths and pow subtractions. See
L<Math::PlanePath::ToothpickTree/Depth to N>.
The other patterns can be expressed as combinations of octants,
parts=4 total = 8*oct(n) - 4*n - 7
parts=1 total = 2*oct(n) - n
3mid V2 total = 2*oct(n+1) + 4*oct(n)
- 3n - 2*floor(log(n+1)) - 6
3side V1 total = oct(n+1) + 3*oct(n) + 2*oct(n-1)
- 3n - floor(log(n+1)) - floor(log(n)) - 4
The depth offsets n,n+1,etc in these formulas become initial pending depth
for the toothpick style depth to N algorithm (and with respective initial
multipliers).
From the V1,V2 formulas it can be seen that V2(n)+V2(n+1) gives the same
combination of 1,3,2 times oct n-1,n,n+1 which is in V1, and that therefore
as noted in the Ndepth part of L</Three Side> above
V1(n) = (V2(n) + V2(n-1) + 1) / 2
=head1 OEIS
This cellular automaton is in Sloane's Online Encyclopedia of Integer
Sequences as
lib/Math/PlanePath/ToothpickTree.pm view on Meta::CPAN
# Do a binary search for the bits of depth which give Ndepth <= N.
#
# Ndepth grows as roughly depth*depth, so this is about log4(N) many
# compares. For large N wouldn't want to a table to search through to
# sqrt(N).
sub _n0_to_depth_and_rem {
my ($self, $n) = @_;
### _n0_to_depth_and_rem(): "n=$n parts=$self->{'parts'}"
# For parts=4 have depth=2^exp formula
# T[2^exp] = parts*(4^exp-1)*2/3 + 3
# parts*(4^exp-1)*2/3 + 3 = N
# 4^exp = (N-3)*3/2parts, round down
# but must be bigger ... (WHY-IS-IT-SO?)
#
my ($pow,$exp) = round_down_pow (12*$n, # /$self->{'parts'}
4);
if (is_infinite($exp)) {
return ($exp,0);
}
lib/Math/PlanePath/ToothpickTree.pm view on Meta::CPAN
my $basemult = $mult; # multiple of oct(2^k) base part
if ($rem == 0) {
### rem==0, so just the oct(2^k) part ...
} elsif ($rem == 1) {
### rem==1 "A" ...
$n += $mult;
} else {
### rem >= 2, formula ...
# formula oct(pow+rem) = oct(pow) + oct(rem+1) + 2*oct(rem) - rem + 4
$n += (4-$rem)*$mult;
$rem += 1; # to give rem+1
if ($rem == $pow) {
### rem+1==pow so oct(2^k) by increasing basemult ...
$basemult += $mult;
} elsif (@new_pending && $new_pending[-1] == $rem) {
### combine rem+1 here with rem of previous ...
$new_mult[-1] += $mult;
} else {
lib/Math/PlanePath/ToothpickTree.pm view on Meta::CPAN
### $mult
### $rem
### assert: $rem >= 0 && $rem <= $pow
if ($rem == 0 || $rem == $pow-1) {
### rem==0, A of each, add=1 ...
### or depth=2*pow-1, add=1 ...
$add += $mult;
} else {
### rem >= 2, formula ...
# A(pow+rem) = A(rem+1) + 2A(rem) - 1
$add -= $mult;
$rem += 1; # to make rem+1
if (@new_depth_list && $new_depth_list[-1] == $rem) {
# add to previously pushed pending depth
# print "rem=$rem ",join(',',@new_depth_list),"\n";
$new_mult_list[-1] += $mult;
} else {
push @new_depth_list, $rem;
lib/Math/PlanePath/ToothpickTree.pm view on Meta::CPAN
| | | | | |
For example the vertical depth "4" toothpick which is above the X=Y diagonal
folds down to become the horizontal "4" in the lower octant. Similarly the
block of five 8,10,12,12,12 above the diagonal fold down to make five
horizontals. And the final 12 above becomes the horizontal 12.
However the horizontals which are on the central diagonal spine don't have
corresponding verticals above. These are marked "....." in the octant shown
above. This means 1 toothpick missing at every second row and therefore the
floor(depth/2) in the oct() formula above.
The key is that a quadrant has the upper octant running 1 growth row behind
the lower. So the horizontals in the lower correspond to the verticals in
the upper (and vice-versa).
The correspondence can be seen algebraically in the formula for a quadrant,
quad(d) = oct(d) + oct(d-1) + d
reversed to
oct(d) = quad(d) - oct(d-1)
and the oct(d-1) term repeatedly expanded
oct(d) = quad(d) - (quad(d-1) - oct(d-2) + d-1) + d
lib/Math/PlanePath/ToothpickTree.pm view on Meta::CPAN
= ...
= quad(d)-quad(d-1)+1 + quad(d-2)-quad(d-3)+1 + ...
= quad(d)-quad(d-1) + quad(d-2)-quad(d-3) + ... + floor(d/2)
The difference quad(d)-quad(d-1) is the number of toothpicks added to make
depth d, and similarly quad(d-2)-quad(d-3) the number added to make depth
d-2. This means the number added at every second row, so if d is even then
this counts only the vertical toothpicks added. Or if d is odd then only
the horizontals.
The +d, +(d-1), +(d-2) additions from the quad(d) formula have alternating
signs and so cancel out to be +1 per pair, giving floor(d/2).
The parts=wedge pattern is two octants and therefore the wedge corresponds
to the horizontals or verticals of parts=2 which is two quadrants. But
there's an adjustment to make there though since parts=2 doesn't have a
toothpick at the origin the way the wedge does.
=head2 Quadrant and 2^k-1 Sums
In OEIS A168002 (L<http://oeis.org/A168002>) Omar Pol observed that the
lib/Math/PlanePath/ToothpickTree.pm view on Meta::CPAN
for 2 <= rem < 2^k
Taking these modulo 2
Q(2^k) == 0 mod 2 since (4^k-4)/6 always even
Q(2^k+1) == 1 mod 2
Q(2^k + rem) == Q(rem+1) mod 2 2 <= rem < 2^k
Q(2^k-1 + rem) == Q(rem) mod 2 1 <= rem < 2^k-1
The last formula is the key, being the same as the ways(2^k-1 + rem)
recurrence above. Then Q(2^k)=0 corresponds to ways(2^k-2)=0, and
Q(2^k+1)=1 corresponding to ways(2^k-1)=1. And therefore
Q(d-2) == ways(d) mod 2
=head1 FUNCTIONS
See L<Math::PlanePath/FUNCTIONS> for behaviour common to all path classes.
=over 4
lib/Math/PlanePath/ToothpickTree.pm view on Meta::CPAN
pattern. For example parts=two_horiz completes levels of the eighths
nearest the X axis (the "oct 2 behind" shown in L</Two Horizontal> above).
=back
=head1 FORMULAS
=head2 Depth to N
The first N at given depth is the total count of toothpicks in the preceding
rows. The paper by Applegate, Pol and Sloane above gives formulas for
parts=4 and parts=1. A similar formula can be made for parts=octant.
The depth in all the following is per the full pattern, which means the
octant starts at depth=2. So oct(2)=0 then oct(3)=1. This reckoning keeps
the replications on 2^k boundaries and is convenient for relating an octant
to the full pattern. Note though that C<tree_depth_to_n()> always counts
from C<$depth = 0> so an adjustment +1 or +2 is applied there.
for depth >= 2
depth = pow + rem where pow=2^k is the high bit of depth
so 0 <= rem < 2^k
t/ToothpickTree.t view on Meta::CPAN
foreach my $level (0 .. $#$depth_aref) {
my $depth = $depth_aref->[$level];
my $n_end = $path->tree_depth_to_n_end($depth);
my ($n_lo,$n_hi) = $path->level_to_n_range($level);
ok ($n_hi, $n_end);
}
}
#------------------------------------------------------------------------------
# oct() formulas
{
my $path = Math::PlanePath::ToothpickTree->new (parts => 'octant');
sub octant {
my ($d) = @_;
die "oops octant($d)" if $d < 2;
return $path->tree_depth_to_n($d-2);
}
}
{
xt/oeis/OneOfEight-oeis.t view on Meta::CPAN
# eg V2 (44+63+1)/2=54
require Math::NumSeq::OEIS::File;
my $V1seq = Math::NumSeq::OEIS::File->new(anum=>'A170879');
my $V2seq = Math::NumSeq::OEIS::File->new(anum=>'A170880');
my $V1_count = 0; while ($V1seq->next) { $V1_count++ }
my $V2_count = 0; while ($V2seq->next) { $V2_count++ }
my $max_count = min($V1_count,$V2_count);
MyOEIS::compare_values
(anum => 'A170879',
name => 'V2=A170880 formula',
func => sub {
my ($count) = @_;
my @got = (0);
my $prev = 0;
for (my $n = 1; @got < $count; $n++) {
my ($i, $value) = $V2seq->next;
push @got, ($value + $prev + 1) / 2;
$prev = $value;
}
return \@got;
});
}
sub my_V1_from_V2 {
my ($n) = @_;
if ($n == 0) { return 0; }
my $path = make_path('3mid'); # V2=3mid
return ($path->tree_depth_to_n($n)
+ $path->tree_depth_to_n($n-1)
+ 1) / 2;
# return (formula_V2($n) + formula_V2($n-1) + 5)/2;
}
{
# my_V1_from_V2() vs 3side
my $path = make_path('3side');
for (my $depth = 0; $depth < 1024; $depth++) {
my $f = my_V1_from_V2($depth);
my $p = $path->tree_depth_to_n($depth);
if ($f != $p) {
warn "depth=$depth f=$f 3side=$p";
}
xt/oeis/OneOfEight-oeis.t view on Meta::CPAN
my $t = my_total($depth);
my $ot = my_total_from_octant($depth);
if ($ot != $t) {
die "depth=$depth t=$t ot=$ot";
}
}
ok (1,1,'my_octant() vs my_total()');
}
#------------------------------------------------------------------------------
# "v2" added by toothpick paper formula
# 0 1 5 5 11 7 15 19 23 7
# n=0 v2=0 so first depth to depth+1 at n=1 v2=1
sub v2_formula {
my ($n) = @_;
### v2_formula(): $n
if ($n <= 0) {
return 0;
}
if ($n == 1) {
return 1;
}
my ($pow,$k) = round_down_pow($n,2);
my $i = $n - $pow;
if ($i == 0) {
return 3*$pow - 1;
}
if ($i == $pow-1) {
return 2*v2_formula($i) + v2_formula($i+1) - 2;
}
return 2*v2_formula($i) + v2_formula($i+1);
}
BEGIN {
use Memoize;
Memoize::memoize('v2_formula');
}
MyOEIS::compare_values
(anum => 'A151728',
name => 'v2_formula()',
func => sub {
my ($count) = @_;
my @got;
for (my $n = 1; @got < $count; $n++) {
push @got, v2_formula($n);
}
return \@got;
});
{
# v2_formula() vs path parts=3mid width
my $limit = ($class eq 'Math::PlanePath::OneOfEightByCells'
? 64 # small when ByCells
: 32768);
my $path = make_path('3mid');
for (my $depth = 0; $depth < $limit; $depth++) {
my $n = $depth+1;
my $v2 = v2_formula($n);
my $added = $path->tree_depth_to_width($depth);
if ($added != $v2) {
die "depth=$depth n=$n added=$added v2=$v2";
}
}
ok (1,1,'v2 against path parts=3mid');
}
{
# v2_formula() vs path parts=1 width at 2^k offset
my $limit = ($class eq 'Math::PlanePath::OneOfEightByCells'
? 64 # small when ByCells
: 32768);
my $offset = $limit;
my $path = make_path('1');
for (my $n = 0; $n < $limit; $n++) {
my $v2 = v2_formula($n+1);
my $added = $path->tree_depth_to_width($n+$offset);
if ($added != $v2) {
die "n=$n added=$added v2=$v2";
}
}
ok (1,1,'v2 against path parts=1');
}
#------------------------------------------------------------------------------
# A170880 "V2", 3mid total
xt/oeis/OneOfEight-oeis.t view on Meta::CPAN
# (anum => 'A151747',
# func => sub {
# my ($count) = @_;
# my @got = (0,1,3,5);
# for (my $depth = scalar(@got); @got < $count; $depth++) {
# push @got, v1_from_sides($depth);
# }
# return \@got;
# });
# foreach my $n (4 .. 30000) {
# my $sides = v1_formula($n);
# my $formula = v1_from_sides($n);
# if ($sides != $formula) {
# die "n=$n sides=$sides formula=$formula";
# }
# }
#------------------------------------------------------------------------------
# A170879 V1 total, 3side total
# cumulative A151747 "v1"
MyOEIS::compare_values
(anum => 'A170879',
name => 'by 3side tree_depth_to_n()',
xt/oeis/OneOfEight-oeis.t view on Meta::CPAN
#
# A151747
# except n<=3 a(n) = 2n-1
# j=0
# a(2^k + 0) = (3*k+1)*2^(k-2) + 1
# a(2^k + 1) = 3*2^(k-1) + 3 = (3*2^k + 6)/2
# a(2^k + j) = 2*a(j) + a(j+1) 2 <= j <= 2^k-1
# a(2^k + j) = 2*a(j) + a(j+1) - 1 if j=2^k-1
# n=0 v1=0 so first depth to depth+1 at n=1 v1=1
sub v1_formula {
my ($n) = @_;
if ($n <= 0) {
return 0;
}
if ($n <= 3) {
return 2*$n-1;
}
my ($pow,$k) = round_down_pow($n,2);
my $j = $n - $pow;
if ($j == 0) {
return 1+(3*$k+1)*2**($k-2);
}
if ($j == 1) {
return 3+ 3*2**($k-1);
}
if ($j == $pow-1) {
return 2*v1_formula($j) + v1_formula($j+1) - 1;
}
return 2*v1_formula($j) + v1_formula($j+1)
}
BEGIN {
use Memoize;
Memoize::memoize('v1_formula');
}
MyOEIS::compare_values
(anum => 'A151747',
name => 'v1_formula()',
func => sub {
my ($count) = @_;
my @got;
for (my $n = 0; @got < $count; $n++) {
push @got, v1_formula($n);
}
return \@got;
});
#------------------------------------------------------------------------------
# "v" total added by formula
# n=0 v=0 so first depth to depth+1 at n=1 v=1
sub v_formula {
my ($n) = @_;
if ($n == 0) {
return 0;
}
if ($n == 1) {
return 1;
}
my ($pow,$k) = round_down_pow($n,2);
my $i = $n - $pow;
if ($i == 0) {
return 6*$pow - 4;
}
return 4*v2_formula($i);
}
MyOEIS::compare_values
(anum => 'A151726',
name => 'v_formula()',
func => sub {
my ($count) = @_;
my @got;
for (my $n = 0; @got < $count; $n++) {
push @got, v_formula($n);
}
return \@got;
});
{
# v_formula() vs path parts=4 width at 2^k offset
my $limit = ($class eq 'Math::PlanePath::OneOfEightByCells'
? 64 # small when ByCells
: 32768);
my $offset = $limit;
my $path = make_path('4');
for (my $n = 1; $n < $limit; $n++) {
my $v = v_formula($n+1);
my $added = $path->tree_depth_to_width($n);
if ($added != $v) {
die "n=$n added=$added v=$v";
}
}
ok (1,1,'v against ByCells');
}
#------------------------------------------------------------------------------
( run in 0.656 second using v1.01-cache-2.11-cpan-26ccb49234f )