File Coverage

blib/lib/Math/NumSeq/FractionDigits.pm

Criterion	Covered	Total	%
statement	105	113	92.9
branch	17	24	70.8
condition	8	36	22.2
subroutine	17	17	100.0
pod	5	5	100.0
total	152	195	77.9

# Copyright 2010, 2011, 2012, 2013, 2014 Kevin Ryde

# This file is part of Math-NumSeq.
#
# Math-NumSeq 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-NumSeq 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-NumSeq.  If not, see .

package Math::NumSeq::FractionDigits;
use 5.004;
use strict;
use List::Util 'max';

use Math::NumSeq;
*_is_infinite = \&Math::NumSeq::_is_infinite;
*_to_bigint = \&Math::NumSeq::_to_bigint;

use vars '$VERSION', '@ISA';
$VERSION = 72;
use Math::NumSeq::Base::Digits;
@ISA = ('Math::NumSeq::Base::Digits');

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


# use constant name => Math::NumSeq::__('Fraction Digits');
use constant description => Math::NumSeq::__('A given fraction number written out in binary.');
use constant i_start => 0;

use constant parameter_info_array =>
  [ Math::NumSeq::Base::Digits->parameter_info_list,
    { name       => 'fraction',
      display    => Math::NumSeq::__('Fraction'),
      type       => 'string',
      type_hint  => 'fraction',
      width      => 10,
      default    => '5/29', # an arbitrary choice
      description => Math::NumSeq::__('The fraction to show, for example 5/29.  Press Return when ready to display the expression.'),
    },
  ];

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

my @oeis_anum;

$oeis_anum[10] =
  {
   # Any fixed-length repeating sequence is a fraction of some sort in
   # some radix.  There's many more not expressed here, and constant
   # digits sequences can be done by more than one fraction, etc.

   '1/7'   => 'A020806',   # 1/7 decimal
   # OEIS-Catalogue: A020806 fraction=1/7

   # OFFSET=1, unlike other fractions which are OFFSET=0
   # # '22/7'  => 'A068028',   # 22/7 decimal
   # # # OEIS-Catalogue: A068028 fraction=22/7

   '1/9'   => 'A000012',   # 1/9 decimal, is just 1,1,1,1
   # pending something better for a constant sequence
   # OEIS-Catalogue: A000012 fraction=1/9

   '1/11'  => 'A010680',   # 1/11 decimal
   # OEIS-Catalogue: A010680 fraction=1/11

   # OEIS-Catalogue: A021015 fraction=1/11 # duplicate of A010680
   # OEIS-Catalogue: A021016 fraction=1/12
   # OEIS-Catalogue: A021017 fraction=1/13
   # OEIS-Catalogue: A021018 fraction=1/14
   # OEIS-Catalogue: A021019 fraction=1/15
   # OEIS-Catalogue: A021020 fraction=1/16

   '1/17'  => 'A007450',   # 1/17 decimal
   # OEIS-Catalogue: A007450 fraction=1/17

   # Math::NumSeq::OEIS::Catalogue::Plugin::FractionDigits has A021022
   # through A021999, being 1/18 to 1/995.
   # A022000 is not 1/996, that fraction missing apparently.
   #
   # OEIS-Catalogue: A022001 fraction=1/997
   # OEIS-Catalogue: A022002 fraction=1/998
   # OEIS-Catalogue: A022003 fraction=1/999

   # OFFSET ?
   # '1/999999'  => 'A172051',
   # # OEIS-Catalogue: A172051 fraction=1/999999

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

   # extra 10 in the denominator to give the leading 0
   '13717421/1111111110' => 'A010888',  # .012345678912...
   # OEIS-Catalogue: A010888 fraction=13717421/1111111110

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

   # constant digits 3,3,3,...
   '10/3' => 'A010701',
   # ENHANCE-ME: of course can generate 3s more efficiently just as a
   # constant sequence, in which case would prefer that over this for the
   # catalogue.
   # OEIS-Catalogue: A010701 fraction=10/3
  };
sub oeis_anum {
  my ($self) = @_;
  ### oeis_anum() ...
  my $radix = $self->{'radix'};
  my $fraction = $self->{'fraction'};
  if (my $anum = $oeis_anum[$radix]->{$fraction}) {
    return $anum;
  }
  if ($radix == 10
      && $fraction =~ m{(\d+)/(\d+)}
      && $1 == 1
      && $2 >= 12 && $2 <= 999 && $2 != 996
      && ($2 % 10) != 0
      && $2 != 25) {
    return 'A0'.($2 + 21016-12);
  }
  ### $fraction
  return undef;
}

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

sub new {
  my $self = shift->SUPER::new(@_);

  my $radix = $self->{'radix'};
  my $fraction = $self->{'fraction'};

  my $num = 0;  # 0/0 if unrecognised
  my $den = 0;
  ($num, $den) = ($fraction =~ m{^\s*
                                 ([.[:digit:]]+)?
                                 \s*
                                 (?:/\s*
                                   ([.[:digit:]]+)?
                                 )?
                                 \s*$}x);
  if (! defined $num) { $num = 1; }
  if (! defined $den) { $den = 1; }
  ### $num
  ### $den
  $fraction = "$num/$den";

  # decimals like 1.5/2.75 become 150/275
  {
    ($num, my $num_decimals) = _to_int_and_decimals ($num);
    ($den, my $den_decimals) = _to_int_and_decimals ($den);
    $num .= '0' x max(0, $den_decimals - $num_decimals);
    $den .= '0' x max(0, $num_decimals - $den_decimals);
  }

  if (max(length($num),length($den)) >= length(int (~0 / $radix))) {
    $num = _to_bigint($num);
    $den = _to_bigint($den);
  }

  # increase den so first digit is 0 to radix-1
  while ($den != 0 && $num >= $den) {
    $den *= $radix;
  }

  ### create
  ### $num
  ### $den
  $self->{'fraction'} = $fraction;
  $self->{'initial_num'} = $num;
  $self->{'den'} = $den;

  $self->rewind;
  return $self;
}

sub _to_int_and_decimals {
  my ($n) = @_;
  if ($n =~ m{^(\d*)\.(\d*?)0*$}) {
    return ($1 . $2,
            length($2));
  } else {
    return ($n, 0);
  }
}

sub rewind {
  my ($self) = @_;
  $self->{'i'} = $self->i_start;
  $self->{'num'} = $self->{'initial_num'};
}

sub next {
  my ($self) = @_;

  my $num   = $self->{'num'} || return;  # num==0 exact radix frac
  my $den   = $self->{'den'} || return;  # den==0 invalid
  my $radix = $self->{'radix'};
  ### FractionDigits next(): "$self->{'i'}  $num/$den"

  $num *= $radix;
  my $quot = int ($num / $den);
  $self->{'num'} = $num - $quot * $den;

  ### $quot
  ### rem: $self->{'num'}

  return ($self->{'i'}++, $quot);
}

sub ith {
  my ($self, $i) = @_;

  if ($i < 0 || _is_infinite($i)) {
    return undef;
  }

  my $radix = $self->{'radix'};
  my $den = $self->{'den'};
  my $num = (($self->{'initial_num'} * _modpow ($self->{'radix'}, $i, $den))
             % $den);
  return int (($num * $radix) / $den);
}

use constant 1.02 _UV_MAX_SQRT => do {
  my $uv_max = ~0;
  my $bit = 1;
  my $shift = 0;
  for (;;) {
    $shift++;
    my $try_bit = $bit << 1;
    if ($uv_max >> $shift < $try_bit) {
      last;
    }
    $bit = $try_bit;
  }
  ### $bit

  my $uv_max_sqrt = $bit;
  for (;;) {
    $bit >>= 1;
    last if $bit == 0;
    my $try_sqrt = $uv_max_sqrt + $bit;
    if (int($uv_max / $try_sqrt) <= $try_sqrt) {
      $uv_max_sqrt = $try_sqrt;
    }
  }

  ### $uv_max_sqrt
  ### uv_max_sqrt: sprintf '%#X', $uv_max_sqrt
  ### squared: $uv_max_sqrt*$uv_max_sqrt
  ### squared: sprintf '%#X', $uv_max_sqrt*$uv_max_sqrt

  $uv_max_sqrt
};

sub _modpow {
  my ($base, $exp, $mod) = @_;
  ### _modpow(): "$base $exp $mod"

  my $ret = 1;
  if (ref $mod || $mod > _UV_MAX_SQRT) {
    return _to_bigint($base)->bmodpow($exp,$mod);
  }

  # only if base and mod have no common factor ...
  # $exp %= $mod-1;

  my $power = $base;
  for (;;) {
    ### $exp
    if ($exp % 2) {
      ### step: "power=$power"
      $ret = ($ret * $power) % $mod;
    }
    $exp = int($exp/2) || last;
    $power = ($power*$power) % $mod;
  }
  return $ret;
}


# ENHANCE-ME: only some digits occur, being the modulo den residue class
# containing num.
# sub pred {
#   my ($self, $value) = @_;
# }
#
# =item C<$bool = $seq-Epred($value)>
#
# Return true if C<$value> occurs as a digit in the fraction.

1;
__END__