1
|
|
|
<?php |
2
|
|
|
|
3
|
|
|
namespace Samsara\Fermat\Provider; |
4
|
|
|
|
5
|
|
|
use Samsara\Exceptions\UsageError\IntegrityConstraint; |
6
|
|
|
use Samsara\Fermat\Numbers; |
7
|
|
|
use Samsara\Fermat\Types\Base\DecimalInterface; |
8
|
|
|
use Samsara\Fermat\Types\Base\NumberInterface; |
9
|
|
|
use Samsara\Fermat\Values\ImmutableNumber; |
10
|
|
|
|
11
|
|
|
class SequenceProvider |
12
|
|
|
{ |
13
|
|
|
|
14
|
|
|
const EULER_ZIGZAG = [ |
15
|
|
|
'1', // 0 |
16
|
|
|
'1', // 1 |
17
|
|
|
'1', // 2 |
18
|
|
|
'2', // 3 |
19
|
|
|
'5', // 4 |
20
|
|
|
'16', // 5 |
21
|
|
|
'61', // 6 |
22
|
|
|
'272', // 7 |
23
|
|
|
'1385', // 8 |
24
|
|
|
'7936', // 9 |
25
|
|
|
'50521', // 10 |
26
|
|
|
'353792', // 11 |
27
|
|
|
'2702765', // 12 |
28
|
|
|
'22368256', // 13 |
29
|
|
|
'199360981', // 14 |
30
|
|
|
'1903757312', // 15 |
31
|
|
|
'19391512145', // 16 |
32
|
|
|
'209865342976', // 17 |
33
|
|
|
'2404879675441', // 18 |
34
|
|
|
'29088885112832', // 19 |
35
|
|
|
'370371188237525', // 20 |
36
|
|
|
'4951498053124096', // 21 |
37
|
|
|
'69348874393137901', // 22 |
38
|
|
|
'1015423886506852352', // 23 |
39
|
|
|
'15514534163557086905', // 24 |
40
|
|
|
'246921480190207983616', // 25 |
41
|
|
|
'4087072509293123892361', // 26 |
42
|
|
|
'70251601603943959887872', // 27 |
43
|
|
|
'1252259641403629865468285', // 28 |
44
|
|
|
'23119184187809597841473536', // 29 |
45
|
|
|
'441543893249023104553682821', // 30 |
46
|
|
|
'8713962757125169296170811392', // 31 |
47
|
|
|
'177519391579539289436664789665', // 32 |
48
|
|
|
'3729407703720529571097509625856', // 33 |
49
|
|
|
'80723299235887898062168247453281', // 34 |
50
|
|
|
'1798651693450888780071750349094912', // 35 |
51
|
|
|
'41222060339517702122347079671259045', // 36 |
52
|
|
|
'970982810785059112379399707952152576', // 37 |
53
|
|
|
'23489580527043108252017828576198947741', // 38 |
54
|
|
|
'583203324917310043943191641625494290432', // 39 |
55
|
|
|
'14851150718114980017877156781405826684425', // 40 |
56
|
|
|
'387635983772083031828014624002175135645696', // 41 |
57
|
|
|
'10364622733519612119397957304745185976310201', // 42 |
58
|
|
|
'283727921907431909304183316295787837183229952', // 43 |
59
|
|
|
'7947579422597592703608040510088070619519273805', // 44 |
60
|
|
|
'227681379129930886488600284336316164603920777216', // 45 |
61
|
|
|
'6667537516685544977435028474773748197524107684661', // 46 |
62
|
|
|
'199500252157859031027160499643195658166340757225472', // 47 |
63
|
|
|
'6096278645568542158691685742876843153976539044435185', // 48 |
64
|
|
|
'190169564657928428175235445073924928592047775873499136', // 49 |
65
|
|
|
'6053285248188621896314383785111649088103498225146815121', // 50 |
66
|
|
|
]; |
67
|
|
|
|
68
|
|
|
/** |
69
|
|
|
* OEIS: A005408 |
70
|
|
|
* |
71
|
|
|
* @param $n |
72
|
|
|
* |
73
|
|
|
* @return DecimalInterface|NumberInterface |
74
|
|
|
*/ |
75
|
15 |
|
public static function nthOddNumber($n) |
76
|
|
|
{ |
77
|
|
|
|
78
|
15 |
|
$n = Numbers::makeOrDont(Numbers::IMMUTABLE, $n, 100); |
79
|
|
|
|
80
|
15 |
|
return $n->multiply(2)->add(1); |
|
|
|
|
81
|
|
|
|
82
|
|
|
} |
83
|
|
|
|
84
|
|
|
/** |
85
|
|
|
* OEIS: A005843 |
86
|
|
|
* |
87
|
|
|
* @param $n |
88
|
|
|
* |
89
|
|
|
* @return DecimalInterface|NumberInterface |
90
|
|
|
*/ |
91
|
2 |
|
public static function nthEvenNumber($n) |
92
|
|
|
{ |
93
|
|
|
|
94
|
2 |
|
$n = Numbers::makeOrDont(Numbers::IMMUTABLE, $n, 100); |
95
|
|
|
|
96
|
2 |
|
return $n->multiply(2); |
|
|
|
|
97
|
|
|
|
98
|
|
|
} |
99
|
|
|
|
100
|
|
|
/** |
101
|
|
|
* OEIS: A033999 |
102
|
|
|
* |
103
|
|
|
* @param $n |
104
|
|
|
* |
105
|
|
|
* @return DecimalInterface|NumberInterface |
106
|
|
|
*/ |
107
|
8 |
|
public static function nthPowerNegativeOne($n) |
108
|
|
|
{ |
109
|
|
|
|
110
|
8 |
|
$negOne = Numbers::makeOrDont(Numbers::IMMUTABLE, -1, 100); |
111
|
|
|
|
112
|
8 |
|
return $negOne->pow($n); |
|
|
|
|
113
|
|
|
|
114
|
|
|
} |
115
|
|
|
|
116
|
|
|
/** |
117
|
|
|
* OEIS: A000111 |
118
|
|
|
* |
119
|
|
|
* @param $n |
120
|
|
|
* |
121
|
|
|
* @return DecimalInterface|NumberInterface |
122
|
|
|
* @throws IntegrityConstraint |
123
|
|
|
*/ |
124
|
4 |
|
public static function nthEulerZigzag($n) |
125
|
|
|
{ |
126
|
|
|
|
127
|
4 |
|
$n = Numbers::makeOrDont(Numbers::IMMUTABLE, $n, 100); |
128
|
|
|
|
129
|
4 |
|
if ($n->isGreaterThan(50)) { |
130
|
1 |
|
throw new IntegrityConstraint( |
131
|
1 |
|
'$n cannot be larger than 50', |
132
|
1 |
|
'Limit your use of the Euler Zigzag Sequence to the 50th index', |
133
|
1 |
|
'This library does not support the Euler Zigzag Sequence (OEIS: A000111) beyond E(50)' |
134
|
|
|
); |
135
|
|
|
} |
136
|
|
|
|
137
|
4 |
|
return Numbers::make(Numbers::IMMUTABLE, static::EULER_ZIGZAG[$n->asInt()], 100); |
|
|
|
|
138
|
|
|
|
139
|
|
|
} |
140
|
|
|
|
141
|
|
|
/** |
142
|
|
|
* WARNING: This function is VERY unoptimized. Be careful of large m values. |
143
|
|
|
* |
144
|
|
|
* @param $n |
145
|
|
|
* |
146
|
|
|
* @return DecimalInterface|NumberInterface |
147
|
|
|
*/ |
148
|
1 |
|
public static function nthBernoulliNumber($n) |
149
|
|
|
{ |
150
|
|
|
|
151
|
1 |
|
$n = Numbers::makeOrDont(Numbers::IMMUTABLE, $n, 100); |
152
|
|
|
|
153
|
1 |
|
if ($n->isEqual(0)) { |
154
|
1 |
|
return Numbers::makeOne(); |
155
|
1 |
|
} elseif ($n->isEqual(1)) { |
156
|
1 |
|
return Numbers::make(Numbers::IMMUTABLE, '0.5', 100); |
|
|
|
|
157
|
|
|
} |
158
|
|
|
|
159
|
1 |
|
if ($n->isGreaterThan(1) && $n->modulo(2)->isEqual(1)) { |
|
|
|
|
160
|
1 |
|
return Numbers::makeZero(); |
161
|
|
|
} |
162
|
|
|
|
163
|
1 |
|
$b = Numbers::make(Numbers::IMMUTABLE, -1, 100); |
164
|
1 |
|
$two = Numbers::make(Numbers::IMMUTABLE, 2, 100); |
165
|
1 |
|
$four = Numbers::make(Numbers::IMMUTABLE, 4, 100); |
166
|
|
|
|
167
|
1 |
|
$b = $b->pow($n->divide(2)->floor()) |
|
|
|
|
168
|
1 |
|
->multiply($n->divide($two->pow($n)->subtract($four->pow($n)))) |
|
|
|
|
169
|
1 |
|
->multiply(static::nthEulerZigzag($n)); |
170
|
|
|
|
171
|
1 |
|
return $b; |
|
|
|
|
172
|
|
|
|
173
|
|
|
} |
174
|
|
|
|
175
|
|
|
/** |
176
|
|
|
* OEIS: A000045 |
177
|
|
|
* |
178
|
|
|
* This uses an implementation of the fast-doubling Karatsuba multiplication algorithm as described by 'Nayuki': |
179
|
|
|
* |
180
|
|
|
* https://www.nayuki.io/page/fast-fibonacci-algorithms |
181
|
|
|
* |
182
|
|
|
* @param $n |
183
|
|
|
* @return ImmutableNumber |
184
|
|
|
* @throws IntegrityConstraint |
185
|
|
|
*/ |
186
|
3 |
|
public static function nthFibonacciNumber($n): ImmutableNumber |
187
|
|
|
{ |
188
|
3 |
|
$n = Numbers::makeOrDont(Numbers::IMMUTABLE, $n); |
189
|
3 |
|
if (!$n->isInt()) { |
|
|
|
|
190
|
1 |
|
throw new IntegrityConstraint( |
191
|
1 |
|
'Sequences can only have integer term numbers', |
192
|
1 |
|
'Provide a valid term number', |
193
|
1 |
|
'The nthFibonacciNumber function takes the term number as its argument; provide an integer term number' |
194
|
|
|
); |
195
|
|
|
} |
196
|
|
|
|
197
|
2 |
|
if ($n->isLessThan(0)) { |
198
|
1 |
|
throw new IntegrityConstraint( |
199
|
1 |
|
'Negative term numbers not valid for Fibonacci Sequence', |
200
|
1 |
|
'Provide a positive term number', |
201
|
1 |
|
'A negative term number for the Fibonacci sequence was requested; provide a positive term number' |
202
|
|
|
); |
203
|
|
|
} |
204
|
|
|
|
205
|
1 |
|
$fastFib = static::_fib($n); |
|
|
|
|
206
|
|
|
|
207
|
1 |
|
return $fastFib[0]; |
208
|
|
|
} |
209
|
|
|
|
210
|
1 |
|
private static function _fib(ImmutableNumber $number): array |
211
|
|
|
{ |
212
|
1 |
|
if ($number->isEqual(0)) { |
213
|
1 |
|
return [Numbers::makeZero(), Numbers::makeOne()]; |
214
|
|
|
} |
215
|
|
|
|
216
|
|
|
/** |
217
|
|
|
* @var ImmutableNumber $a |
218
|
|
|
* @var ImmutableNumber $b |
219
|
|
|
* @var ImmutableNumber $prevCall |
220
|
|
|
*/ |
221
|
1 |
|
$prevCall = $number->divide(2)->floor(); |
222
|
1 |
|
list($a, $b) = static::_fib($prevCall); |
223
|
1 |
|
$c = $a->multiply($b->multiply(2)->subtract($a)); |
224
|
1 |
|
$d = $a->multiply($a)->add($b->multiply($b)); |
225
|
1 |
|
if ($number->modulo(2)->isEqual(0)) { |
226
|
1 |
|
return [$c, $d]; |
227
|
|
|
} else { |
228
|
1 |
|
return [$d, $c->add($d)]; |
229
|
|
|
} |
230
|
|
|
} |
231
|
|
|
|
232
|
|
|
} |