wpdatatables / lib / phpoffice / phpspreadsheet / src / PhpSpreadsheet / Calculation / Engineering / BesselJ.php

BesselJ.php in wpDataTables – WordPress Data Table, Dynamic Tables & Table Charts Plugin 6.5.1.7, at lib/phpoffice/phpspreadsheet/src/PhpSpreadsheet/Calculation/Engineering/BesselJ.php

181 lines 5.8 KB

Raw Download Zip

1	<?php
2
3	namespace PhpOffice\PhpSpreadsheet\Calculation\Engineering;
4
5	use PhpOffice\PhpSpreadsheet\Calculation\ArrayEnabled;
6	use PhpOffice\PhpSpreadsheet\Calculation\Exception;
7	use PhpOffice\PhpSpreadsheet\Calculation\Information\ExcelError;
8
9	class BesselJ
10	{
11	use ArrayEnabled;
12
13	/**
14	* BESSELJ.
15	*
16	* Returns the Bessel function
17	*
18	* Excel Function:
19	* BESSELJ(x,ord)
20	*
21	* NOTE: The MS Excel implementation of the BESSELJ function is still not accurate, particularly for higher order
22	* values with x < -8 and x > 8. This code provides a more accurate calculation
23	*
24	* @param mixed $x A float value at which to evaluate the function.
25	* If x is nonnumeric, BESSELJ returns the #VALUE! error value.
26	* Or can be an array of values
27	* @param mixed $ord The integer order of the Bessel function.
28	* If ord is not an integer, it is truncated.
29	* If $ord is nonnumeric, BESSELJ returns the #VALUE! error value.
30	* If $ord < 0, BESSELJ returns the #NUM! error value.
31	* Or can be an array of values
32	*
33	* @return array\|float\|string Result, or a string containing an error
34	* If an array of numbers is passed as an argument, then the returned result will also be an array
35	* with the same dimensions
36	*/
37	public static function BESSELJ($x, $ord)
38	{
39	if (is_array($x) \|\| is_array($ord)) {
40	return self::evaluateArrayArguments([self::class, __FUNCTION__], $x, $ord);
41	}
42
43	try {
44	$x = EngineeringValidations::validateFloat($x);
45	$ord = EngineeringValidations::validateInt($ord);
46	} catch (Exception $e) {
47	return $e->getMessage();
48	}
49
50	if ($ord < 0) {
51	return ExcelError::NAN();
52	}
53
54	$fResult = self::calculate($x, $ord);
55
56	return (is_nan($fResult)) ? ExcelError::NAN() : $fResult;
57	}
58
59	private static function calculate(float $x, int $ord): float
60	{
61	// special cases
62	switch ($ord) {
63	case 0:
64	return self::besselJ0($x);
65	case 1:
66	return self::besselJ1($x);
67	}
68
69	return self::besselJ2($x, $ord);
70	}
71
72	private static function besselJ0(float $x): float
73	{
74	$ax = abs($x);
75
76	if ($ax < 8.0) {
77	$y = $x * $x;
78	$ans1 = 57568490574.0 + $y * (-13362590354.0 + $y * (651619640.7 + $y * (-11214424.18 + $y *
79	(77392.33017 + $y * (-184.9052456)))));
80	$ans2 = 57568490411.0 + $y * (1029532985.0 + $y * (9494680.718 + $y * (59272.64853 + $y *
81	(267.8532712 + $y * 1.0))));
82
83	return $ans1 / $ans2;
84	}
85
86	$z = 8.0 / $ax;
87	$y = $z * $z;
88	$xx = $ax - 0.785398164;
89	$ans1 = 1.0 + $y * (-0.1098628627e-2 + $y * (0.2734510407e-4 + $y * (-0.2073370639e-5 + $y * 0.2093887211e-6)));
90	$ans2 = -0.1562499995e-1 + $y * (0.1430488765e-3 + $y * (-0.6911147651e-5 + $y *
91	(0.7621095161e-6 - $y * 0.934935152e-7)));
92
93	return sqrt(0.636619772 / $ax) * (cos($xx) * $ans1 - $z * sin($xx) * $ans2);
94	}
95
96	private static function besselJ1(float $x): float
97	{
98	$ax = abs($x);
99
100	if ($ax < 8.0) {
101	$y = $x * $x;
102	$ans1 = $x * (72362614232.0 + $y * (-7895059235.0 + $y * (242396853.1 + $y *
103	(-2972611.439 + $y * (15704.48260 + $y * (-30.16036606))))));
104	$ans2 = 144725228442.0 + $y * (2300535178.0 + $y * (18583304.74 + $y * (99447.43394 + $y *
105	(376.9991397 + $y * 1.0))));
106
107	return $ans1 / $ans2;
108	}
109
110	$z = 8.0 / $ax;
111	$y = $z * $z;
112	$xx = $ax - 2.356194491;
113
114	$ans1 = 1.0 + $y * (0.183105e-2 + $y * (-0.3516396496e-4 + $y * (0.2457520174e-5 + $y * (-0.240337019e-6))));
115	$ans2 = 0.04687499995 + $y * (-0.2002690873e-3 + $y * (0.8449199096e-5 + $y *
116	(-0.88228987e-6 + $y * 0.105787412e-6)));
117	$ans = sqrt(0.636619772 / $ax) * (cos($xx) * $ans1 - $z * sin($xx) * $ans2);
118
119	return ($x < 0.0) ? -$ans : $ans;
120	}
121
122	private static function besselJ2(float $x, int $ord): float
123	{
124	$ax = abs($x);
125	if ($ax === 0.0) {
126	return 0.0;
127	}
128
129	if ($ax > $ord) {
130	return self::besselj2a($ax, $ord, $x);
131	}
132
133	return self::besselj2b($ax, $ord, $x);
134	}
135
136	private static function besselj2a(float $ax, int $ord, float $x): float
137	{
138	$tox = 2.0 / $ax;
139	$bjm = self::besselJ0($ax);
140	$bj = self::besselJ1($ax);
141	for ($j = 1; $j < $ord; ++$j) {
142	$bjp = $j * $tox * $bj - $bjm;
143	$bjm = $bj;
144	$bj = $bjp;
145	}
146	$ans = $bj;
147
148	return ($x < 0.0 && ($ord % 2) == 1) ? -$ans : $ans;
149	}
150
151	private static function besselj2b(float $ax, int $ord, float $x): float
152	{
153	$tox = 2.0 / $ax;
154	$jsum = false;
155	$bjp = $ans = $sum = 0.0;
156	$bj = 1.0;
157	for ($j = 2 * ($ord + (int) sqrt(40.0 * $ord)); $j > 0; --$j) {
158	$bjm = $j * $tox * $bj - $bjp;
159	$bjp = $bj;
160	$bj = $bjm;
161	if (abs($bj) > 1.0e+10) {
162	$bj *= 1.0e-10;
163	$bjp *= 1.0e-10;
164	$ans *= 1.0e-10;
165	$sum *= 1.0e-10;
166	}
167	if ($jsum === true) {
168	$sum += $bj;
169	}
170	$jsum = $jsum === false;
171	if ($j === $ord) {
172	$ans = $bjp;
173	}
174	}
175	$sum = 2.0 * $sum - $bj;
176	$ans /= $sum;
177
178	return ($x < 0.0 && ($ord % 2) === 1) ? -$ans : $ans;
179	}
180	}
181