File Coverage

blib/lib/Math/NumSeq/SqrtContinuedPeriod.pm
Criterion Covered Total %
statement 52 54 96.3
branch 6 8 75.0
condition n/a
subroutine 13 13 100.0
pod 1 1 100.0
total 72 76 94.7


line stmt bran cond sub pod time code
1             # Copyright 2011, 2012, 2013, 2014 Kevin Ryde
2              
3             # This file is part of Math-NumSeq.
4             #
5             # Math-NumSeq 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-NumSeq 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-NumSeq. If not, see .
17              
18             package Math::NumSeq::SqrtContinuedPeriod;
19 2     2   55 use 5.004;
  2         7  
  2         105  
20 2     2   14 use strict;
  2         5  
  2         87  
21              
22 2     2   11 use vars '$VERSION', '@ISA';
  2         5  
  2         157  
23             $VERSION = 71;
24              
25 2     2   11 use Math::NumSeq;
  2         3  
  2         53  
26 2     2   598 use Math::NumSeq::Base::IterateIth;
  2         4  
  2         123  
27             @ISA = ('Math::NumSeq::Base::IterateIth',
28             'Math::NumSeq');
29             *_is_infinite = \&Math::NumSeq::_is_infinite;
30              
31             # uncomment this to run the ### lines
32             #use Smart::Comments;
33              
34              
35             # use constant name => Math::NumSeq::__('Sqrt Continued Fraction Period');
36 2     2   12 use constant description => Math::NumSeq::__('Period of the continued fraction expansion of sqrt(i), or 0 for perfect squares (where the expansion is finite)');
  2         3  
  2         12  
37 2     2   10 use constant characteristic_count => 1;
  2         4  
  2         74  
38 2     2   8 use constant characteristic_smaller => 1;
  2         4  
  2         81  
39 2     2   10 use constant characteristic_increasing => 0;
  2         3  
  2         117  
40 2     2   9 use constant default_i_start => 1;
  2         4  
  2         73  
41 2     2   9 use constant values_min => 0;
  2         3  
  2         87  
42              
43             # cf A097853 contfrac sqrt(n) period, or 1 if square
44             # A054269 contfrac sqrt(prime(i)) period
45             #
46 2     2   8 use constant oeis_anum => 'A003285'; # sqrt period, or 0 if square
  2         4  
  2         386  
47              
48             sub ith {
49 200     200 1 306 my ($self, $i) = @_;
50             ### SqrtContinuedPeriod ith(): $i
51              
52 200 50       494 if (_is_infinite($i)) {
53 0         0 return $i;
54             }
55 200 50       469 if ($i <= 0) {
56 0         0 return 0;
57             }
58              
59             # initial a[1] = floor(sqrt(i)) = $root
60             # then root + 1/x = sqrt(i)
61             # 1/x = sqrt(i)-root
62             # x = 1/(sqrt(i)-root)
63             # x = (sqrt(i)+root)/(i-root*root)
64             # so P = root
65             # Q = i - root*root
66             #
67 200         353 my $p = my $root = int(sqrt($i));
68 200         287 my $q = $i - $root*$root;
69 200 100       419 if ($q <= 0) {
70             # perfect square
71 13         54 return 0;
72             }
73              
74 187         203 my %seen;
75 187         267 my $count = 0;
76 187         186 for (;;) {
77 1273 100       7306 if ($seen{"$p,$q"}++) {
78 187         1580 return $count;
79             }
80 1086         1099 $count++;
81              
82 1086         1887 my $value = int (($root + $p) / $q);
83 1086         1429 $p -= $value*$q;
84 1086         11122 ($p, $q) = (-$p,
85             ($i - $p*$p) / $q);
86              
87             ### assert: $p >= 0
88             ### assert: $q > 0
89             ### assert: $p <= $root
90             ### assert: $q <= 2*$root+1
91             ### assert: (($p*$p - $i) % $q) == 0
92             }
93             }
94              
95             1;
96             __END__