1
|
|
|
<?php |
2
|
|
|
|
3
|
|
|
namespace Tdn\PhpTypes\Math\Library; |
4
|
|
|
|
5
|
|
|
use Tdn\PhpTypes\Exception\InvalidNumberException; |
6
|
|
|
use Tdn\PhpTypes\Math\DefaultMathAdapter; |
7
|
|
|
use Tdn\PhpTypes\Type\StringType; |
8
|
|
|
|
9
|
|
|
/** |
10
|
|
|
* Class Spl. |
11
|
|
|
*/ |
12
|
|
|
class Spl implements MathLibraryInterface |
13
|
|
|
{ |
14
|
|
|
/** |
15
|
|
|
* @var int |
16
|
|
|
*/ |
17
|
|
|
private $roundingStrategy; |
18
|
|
|
|
19
|
|
|
/** |
20
|
|
|
* @param int $roundingStrategy a PHP_ROUND_HALF_* integer |
21
|
|
|
*/ |
22
|
106 |
|
public function __construct(int $roundingStrategy) |
23
|
|
|
{ |
24
|
106 |
|
$this->roundingStrategy = $roundingStrategy; |
25
|
106 |
|
} |
26
|
|
|
|
27
|
|
|
/** |
28
|
|
|
* Add two arbitrary precision numbers. |
29
|
|
|
* |
30
|
|
|
* @param string $leftOperand |
31
|
|
|
* @param string $rightOperand |
32
|
|
|
* @param int $precision |
33
|
|
|
* |
34
|
|
|
* @return string |
35
|
|
|
*/ |
36
|
1 |
|
public function add(string $leftOperand, string $rightOperand, int $precision = 0): string |
37
|
|
|
{ |
38
|
1 |
|
return (string) ($this->isIntOperation($precision) ? (intval($leftOperand) + intval($rightOperand)) : |
39
|
1 |
|
round(floatval($leftOperand) + floatval($rightOperand), $precision, $this->roundingStrategy)); |
40
|
|
|
} |
41
|
|
|
|
42
|
|
|
/** |
43
|
|
|
* Subtract two arbitrary precision numbers. |
44
|
|
|
* |
45
|
|
|
* @param string $leftOperand |
46
|
|
|
* @param string $rightOperand |
47
|
|
|
* @param int $precision |
48
|
|
|
* |
49
|
|
|
* @return string |
50
|
|
|
*/ |
51
|
1 |
|
public function subtract(string $leftOperand, string $rightOperand, int $precision = 0): string |
52
|
|
|
{ |
53
|
1 |
|
return (string) ($this->isIntOperation($precision) ? (intval($leftOperand) - intval($rightOperand)) : |
54
|
1 |
|
round($leftOperand - $rightOperand, $precision, $this->roundingStrategy)); |
55
|
|
|
} |
56
|
|
|
|
57
|
|
|
/** |
58
|
|
|
* Multiply two arbitrary precision numbers. |
59
|
|
|
* |
60
|
|
|
* @param string $leftOperand |
61
|
|
|
* @param string $rightOperand |
62
|
|
|
* @param int $precision |
63
|
|
|
* |
64
|
|
|
* @return string |
65
|
|
|
*/ |
66
|
1 |
|
public function multiply(string $leftOperand, string $rightOperand, int $precision = 0): string |
67
|
|
|
{ |
68
|
1 |
|
return (string) ($this->isIntOperation($precision) ? (intval($leftOperand) * intval($rightOperand)) : |
69
|
1 |
|
round($leftOperand * $rightOperand, ($precision ?? 0), $this->roundingStrategy)); |
70
|
|
|
} |
71
|
|
|
|
72
|
|
|
/** |
73
|
|
|
* Divide two arbitrary precision numbers. |
74
|
|
|
* |
75
|
|
|
* @param string $leftOperand |
76
|
|
|
* @param string $rightOperand |
77
|
|
|
* @param int $precision |
78
|
|
|
* |
79
|
|
|
* @return string |
80
|
|
|
*/ |
81
|
1 |
|
public function divide(string $leftOperand, string $rightOperand, int $precision = 0): string |
82
|
|
|
{ |
83
|
1 |
|
return (string) ($this->isIntOperation($precision) ? (intval($leftOperand) / intval($rightOperand)) : |
84
|
1 |
|
round($leftOperand / $rightOperand, $precision, $this->roundingStrategy)); |
85
|
|
|
} |
86
|
|
|
|
87
|
|
|
/** |
88
|
|
|
* Compare two arbitrary precision numbers. |
89
|
|
|
* |
90
|
|
|
* @param string $leftOperand |
91
|
|
|
* @param string $rightOperand |
92
|
|
|
* @param int $precision |
93
|
|
|
* |
94
|
|
|
* @return string |
95
|
|
|
*/ |
96
|
2 |
|
public function compare(string $leftOperand, string $rightOperand, int $precision = 0): string |
97
|
|
|
{ |
98
|
2 |
|
return strval($leftOperand <=> $rightOperand); |
99
|
|
|
} |
100
|
|
|
|
101
|
|
|
/** |
102
|
|
|
* Get modulus of an arbitrary precision number. |
103
|
|
|
* |
104
|
|
|
* @param string $operand |
105
|
|
|
* @param string $modulus |
106
|
|
|
* @param int $precision |
107
|
|
|
* |
108
|
|
|
* @return string |
109
|
|
|
*/ |
110
|
3 |
|
public function modulus(string $operand, string $modulus, int $precision = 0): string |
111
|
|
|
{ |
112
|
3 |
|
return (string) round( |
113
|
3 |
|
fmod( |
114
|
3 |
|
floatval($operand), |
115
|
3 |
|
floatval($modulus) |
116
|
|
|
), |
117
|
3 |
|
($precision ?? 0), |
118
|
3 |
|
$this->roundingStrategy |
119
|
|
|
); |
120
|
|
|
} |
121
|
|
|
|
122
|
|
|
/** |
123
|
|
|
* Raise an arbitrary precision number to another. |
124
|
|
|
* |
125
|
|
|
* @param string $leftOperand |
126
|
|
|
* @param string $rightOperand |
127
|
|
|
* @param int $precision |
128
|
|
|
* |
129
|
|
|
* @return string |
130
|
|
|
*/ |
131
|
1 |
|
public function power(string $leftOperand, string $rightOperand, int $precision = 0): string |
132
|
|
|
{ |
133
|
1 |
|
return (string) round( |
134
|
1 |
|
pow( |
135
|
1 |
|
floatval($leftOperand), |
136
|
1 |
|
floatval($rightOperand) |
137
|
|
|
), |
138
|
1 |
|
($precision ?? 0), |
139
|
1 |
|
$this->roundingStrategy |
140
|
|
|
); |
141
|
|
|
} |
142
|
|
|
|
143
|
|
|
/** |
144
|
|
|
* Get the square root of an arbitrary precision number. |
145
|
|
|
* |
146
|
|
|
* @param string $operand |
147
|
|
|
* @param int $precision |
148
|
|
|
* |
149
|
|
|
* @return string |
150
|
|
|
*/ |
151
|
3 |
|
public function squareRoot(string $operand, int $precision = 0): string |
152
|
|
|
{ |
153
|
3 |
|
return (string) round(sqrt(floatval($operand)), ($precision ?? 0), $this->roundingStrategy); |
154
|
|
|
} |
155
|
|
|
|
156
|
|
|
/** |
157
|
|
|
* Returns absolute value of operand. |
158
|
|
|
* |
159
|
|
|
* @param string $operand |
160
|
|
|
* |
161
|
|
|
* @return string |
162
|
|
|
*/ |
163
|
3 |
|
public function absolute(string $operand): string |
164
|
|
|
{ |
165
|
3 |
|
return (string) abs($operand); |
166
|
|
|
} |
167
|
|
|
|
168
|
|
|
/** |
169
|
|
|
* Negates a number. Opposite of absolute/abs. |
170
|
|
|
* |
171
|
|
|
* @param string $operand |
172
|
|
|
* |
173
|
|
|
* @return string |
174
|
|
|
*/ |
175
|
3 |
|
public function negate(string $operand): string |
176
|
|
|
{ |
177
|
3 |
|
return strval($operand * -1); |
178
|
|
|
} |
179
|
|
|
|
180
|
|
|
/** |
181
|
|
|
* Returns the factorial of operand. |
182
|
|
|
* |
183
|
|
|
* @param string $operand |
184
|
|
|
* |
185
|
|
|
* @return string |
186
|
|
|
*/ |
187
|
1 |
|
public function factorial(string $operand): string |
188
|
|
|
{ |
189
|
1 |
|
if (StringType::create($operand)->contains('.')) { |
190
|
1 |
|
++$operand; |
191
|
|
|
|
192
|
1 |
|
return $this->gamma((string) $operand); |
193
|
|
|
} |
194
|
|
|
|
195
|
|
|
$factorial = function (string $num) use (&$factorial) { |
196
|
1 |
|
if ($num < 2) { |
197
|
1 |
|
return 1; |
198
|
|
|
} |
199
|
|
|
|
200
|
1 |
|
return $factorial(strval($num - 1)) * $num; |
201
|
1 |
|
}; |
202
|
|
|
|
203
|
1 |
|
return (string) $factorial($operand); |
204
|
|
|
} |
205
|
|
|
|
206
|
|
|
/** |
207
|
|
|
* Greatest common divisor. |
208
|
|
|
* |
209
|
|
|
* @param string $leftOperand |
210
|
|
|
* @param string $rightOperand |
211
|
|
|
* |
212
|
|
|
* @return string |
213
|
|
|
*/ |
214
|
|
|
public function gcd(string $leftOperand, string $rightOperand): string |
215
|
|
|
{ |
216
|
1 |
|
$gcd = function (string $a, string $b) use (&$gcd) { |
217
|
1 |
|
return $b ? $gcd($b, strval($a % $b)) : $a; |
218
|
1 |
|
}; |
219
|
|
|
|
220
|
1 |
|
$exponent = $this->getSmallestDecimalPlaceCount($leftOperand, $rightOperand); |
221
|
|
|
|
222
|
|
|
return (string) ( |
223
|
1 |
|
$gcd( |
224
|
1 |
|
strval($leftOperand * (pow(10, $exponent))), |
225
|
1 |
|
strval($rightOperand * (pow(10, $exponent))) |
226
|
1 |
|
) / pow(10, $exponent) |
227
|
|
|
); |
228
|
|
|
} |
229
|
|
|
|
230
|
|
|
/** |
231
|
|
|
* Calculates to the nth root. |
232
|
|
|
* |
233
|
|
|
* @param string $operand |
234
|
|
|
* @param int $nth |
235
|
|
|
* |
236
|
|
|
* @return string |
237
|
|
|
*/ |
238
|
1 |
|
public function root(string $operand, int $nth): string |
239
|
|
|
{ |
240
|
|
|
//Implement something later. |
241
|
1 |
|
throw new \RuntimeException('Not a valid library for root^n.'); |
242
|
|
|
} |
243
|
|
|
|
244
|
|
|
/** |
245
|
|
|
* Gets the next prime after operand. |
246
|
|
|
* |
247
|
|
|
* @param string $operand |
248
|
|
|
* |
249
|
|
|
* @return string |
250
|
|
|
*/ |
251
|
2 |
|
public function nextPrime(string $operand): string |
252
|
|
|
{ |
253
|
2 |
|
$operand = (intval($operand) + 1); |
254
|
2 |
|
for ($i = $operand;; ++$i) { |
255
|
2 |
|
if ($this->isPrime(strval($i))) { |
256
|
2 |
|
break; |
257
|
|
|
} |
258
|
|
|
} |
259
|
|
|
|
260
|
2 |
|
return (string) $i; |
261
|
|
|
} |
262
|
|
|
|
263
|
|
|
/** |
264
|
|
|
* @param string $operand |
265
|
|
|
* @param int $reps |
266
|
|
|
* |
267
|
|
|
* @return bool |
268
|
|
|
*/ |
269
|
3 |
|
public function isPrime(string $operand, int $reps = 10): bool |
270
|
|
|
{ |
271
|
3 |
|
$x = floor(sqrt(floatval($operand))); |
272
|
3 |
|
for ($i = 2; $i <= $x; ++$i) { |
273
|
3 |
|
if (($operand % $i) == 0) { |
274
|
1 |
|
break; |
275
|
|
|
} |
276
|
|
|
} |
277
|
|
|
|
278
|
3 |
|
return ($x == ($i - 1)) ? true : false; |
279
|
|
|
} |
280
|
|
|
|
281
|
|
|
/** |
282
|
|
|
* Checks if operand is perfect square. |
283
|
|
|
* |
284
|
|
|
* @param string $operand |
285
|
|
|
* @param int $precision |
286
|
|
|
* |
287
|
|
|
* @return bool |
288
|
|
|
*/ |
289
|
2 |
|
public function isPerfectSquare(string $operand, int $precision = 0): bool |
290
|
|
|
{ |
291
|
2 |
|
$candidate = $this->squareRoot($operand, $precision + 1); |
292
|
|
|
|
293
|
2 |
|
return $candidate == intval($candidate); |
294
|
|
|
} |
295
|
|
|
|
296
|
|
|
/** |
297
|
|
|
* @return bool |
298
|
|
|
*/ |
299
|
9 |
|
public function isEnabled(): bool |
300
|
|
|
{ |
301
|
9 |
|
return true; |
302
|
|
|
} |
303
|
|
|
|
304
|
|
|
/** |
305
|
|
|
* @param string $type |
306
|
|
|
* |
307
|
|
|
* @return bool |
308
|
|
|
*/ |
309
|
9 |
|
public function supportsOperationType(string $type): bool |
310
|
|
|
{ |
311
|
|
|
// Supports both float and int. |
312
|
9 |
|
return true; |
313
|
|
|
} |
314
|
|
|
|
315
|
|
|
/** |
316
|
|
|
* @param string $operand |
317
|
|
|
* |
318
|
|
|
* @return string |
319
|
|
|
*/ |
320
|
5 |
|
public function gamma(string $operand): string |
321
|
|
|
{ |
322
|
5 |
|
if ($operand <= 0.0) { |
323
|
1 |
|
throw new InvalidNumberException('Operand must be a positive number.'); |
324
|
|
|
} |
325
|
|
|
|
326
|
|
|
// Euler's gamma constant |
327
|
4 |
|
$gamma = 0.577215664901532860606512090; |
328
|
|
|
|
329
|
4 |
|
if ($operand < 0.001) { |
330
|
1 |
|
return strval(1.0 / ($operand * (1.0 + $gamma * $operand))); |
331
|
|
|
} |
332
|
|
|
|
333
|
4 |
|
if ($operand < 12.0) { |
334
|
3 |
|
$y = $operand; |
335
|
3 |
|
$n = 0; |
336
|
3 |
|
$lessThanOne = ($y < 1.0); |
337
|
3 |
|
if ($lessThanOne) { |
338
|
2 |
|
$y += 1.0; |
339
|
|
|
} else { |
340
|
|
|
// will use n later |
341
|
3 |
|
$n = floor($y) - 1; |
342
|
3 |
|
$y -= $n; |
343
|
|
|
} |
344
|
|
|
|
345
|
|
|
$p = [ |
346
|
3 |
|
-1.71618513886549492533811E+0, |
347
|
|
|
2.47656508055759199108314E+1, |
348
|
|
|
-3.79804256470945635097577E+2, |
349
|
|
|
6.29331155312818442661052E+2, |
350
|
|
|
8.66966202790413211295064E+2, |
351
|
|
|
-3.14512729688483675254357E+4, |
352
|
|
|
-3.61444134186911729807069E+4, |
353
|
|
|
6.64561438202405440627855E+4, |
354
|
|
|
]; |
355
|
|
|
|
356
|
|
|
$q = [ |
357
|
3 |
|
-3.08402300119738975254353E+1, |
358
|
|
|
3.15350626979604161529144E+2, |
359
|
|
|
-1.01515636749021914166146E+3, |
360
|
|
|
-3.10777167157231109440444E+3, |
361
|
|
|
2.25381184209801510330112E+4, |
362
|
|
|
4.75584627752788110767815E+3, |
363
|
|
|
-1.34659959864969306392456E+5, |
364
|
|
|
-1.15132259675553483497211E+5, |
365
|
|
|
]; |
366
|
|
|
|
367
|
3 |
|
$num = 0.0; |
368
|
3 |
|
$den = 1.0; |
369
|
3 |
|
$z = $y - 1; |
370
|
|
|
|
371
|
3 |
|
for ($i = 0; $i < 8; ++$i) { |
372
|
3 |
|
$num = ($num + $p[$i]) * $z; |
373
|
3 |
|
$den = $den * $z + $q[$i]; |
374
|
|
|
} |
375
|
|
|
|
376
|
3 |
|
$result = $num / $den + 1.0; |
377
|
|
|
|
378
|
3 |
|
if ($lessThanOne) { |
379
|
2 |
|
$result /= ($y - 1.0); |
380
|
|
|
} else { |
381
|
3 |
|
for ($i = 0; $i < $n; ++$i) { |
382
|
3 |
|
$result *= $y++; |
383
|
|
|
} |
384
|
|
|
} |
385
|
|
|
|
386
|
3 |
|
return (string) $result; |
387
|
|
|
} |
388
|
|
|
|
389
|
2 |
|
if ($operand > 171.624) { |
390
|
1 |
|
throw new \RuntimeException('Number too large.'); |
391
|
|
|
} |
392
|
|
|
|
393
|
1 |
|
return (string) exp(floatval($this->logGamma((string) $operand))); |
394
|
|
|
} |
395
|
|
|
|
396
|
|
|
/** |
397
|
|
|
* @param string $operand |
398
|
|
|
* |
399
|
|
|
* @return string |
400
|
|
|
*/ |
401
|
3 |
|
public function logGamma(string $operand): string |
402
|
|
|
{ |
403
|
3 |
|
if ($operand <= 0.0) { |
404
|
1 |
|
throw new InvalidNumberException('Operand must be a positive number.'); |
405
|
|
|
} |
406
|
|
|
|
407
|
2 |
|
if ($operand < 12.0) { |
408
|
1 |
|
return (string) log(abs($this->gamma((string) $operand))); |
409
|
|
|
} |
410
|
|
|
|
411
|
|
|
$c = [ |
412
|
2 |
|
1.0 / 12.0, |
413
|
|
|
-1.0 / 360.0, |
414
|
|
|
1.0 / 1260.0, |
415
|
|
|
-1.0 / 1680.0, |
416
|
|
|
1.0 / 1188.0, |
417
|
|
|
-691.0 / 360360.0, |
418
|
|
|
1.0 / 156.0, |
419
|
|
|
-3617.0 / 122400.0, |
420
|
|
|
]; |
421
|
|
|
|
422
|
2 |
|
$z = 1.0 / ($operand * $operand); |
423
|
2 |
|
$sum = $c[7]; |
424
|
2 |
|
for ($i = 6; $i >= 0; --$i) { |
425
|
2 |
|
$sum *= $z; |
426
|
2 |
|
$sum += $c[$i]; |
427
|
|
|
} |
428
|
|
|
|
429
|
2 |
|
$series = $sum / $operand; |
430
|
2 |
|
$halfLogTwoPi = 0.91893853320467274178032973640562; |
431
|
2 |
|
$logGamma = (floatval($operand) - 0.5) * log(floatval($operand)) - $operand + $halfLogTwoPi + $series; |
432
|
|
|
|
433
|
2 |
|
return (string) $logGamma; |
434
|
|
|
} |
435
|
|
|
|
436
|
|
|
/** |
437
|
|
|
* Figures out the smallest number of decimal places between the two numbers and returns that count. |
438
|
|
|
* Eg. (1.005, 2.4) => 1, (1.005, 2.5399) => 3. |
439
|
|
|
* |
440
|
|
|
* @param string $leftOperand |
441
|
|
|
* @param string $rightOperand |
442
|
|
|
* |
443
|
|
|
* @return int |
444
|
|
|
*/ |
445
|
1 |
|
private function getSmallestDecimalPlaceCount(string $leftOperand, string $rightOperand): int |
446
|
|
|
{ |
447
|
1 |
|
$leftPrecision = DefaultMathAdapter::getNumberPrecision($leftOperand); |
448
|
1 |
|
$rightPrecision = DefaultMathAdapter::getNumberPrecision($rightOperand); |
449
|
|
|
|
450
|
1 |
|
return $leftPrecision < $rightPrecision ? $leftPrecision : $rightPrecision; |
451
|
|
|
} |
452
|
|
|
|
453
|
|
|
/** |
454
|
|
|
* Ensures that an operation is meant to be an integer operation, float operation otherwise. |
455
|
|
|
* |
456
|
|
|
* @param int $precision |
457
|
|
|
* |
458
|
|
|
* @return bool |
459
|
|
|
*/ |
460
|
4 |
|
private function isIntOperation(int $precision = 0): bool |
461
|
|
|
{ |
462
|
4 |
|
return $precision === null || $precision === 0; |
463
|
|
|
} |
464
|
|
|
} |
465
|
|
|
|