/lib/phpspreadsheet/vendor/phpoffice/phpspreadsheet/src/PhpSpreadsheet/Shared/JAMA/Matrix.php

https://github.com/markn86/moodle · PHP · 1202 lines · 776 code · 150 blank · 276 comment · 204 complexity · d1bd67a6abc6751ae6d7c8f8a3193133 MD5 · raw file

  1. <?php
  2. namespace PhpOffice\PhpSpreadsheet\Shared\JAMA;
  3. use PhpOffice\PhpSpreadsheet\Calculation\Exception as CalculationException;
  4. use PhpOffice\PhpSpreadsheet\Calculation\Functions;
  5. use PhpOffice\PhpSpreadsheet\Shared\StringHelper;
  6. /**
  7. * Matrix class.
  8. *
  9. * @author Paul Meagher
  10. * @author Michael Bommarito
  11. * @author Lukasz Karapuda
  12. * @author Bartek Matosiuk
  13. *
  14. * @version 1.8
  15. *
  16. * @see https://math.nist.gov/javanumerics/jama/
  17. */
  18. class Matrix
  19. {
  20. const POLYMORPHIC_ARGUMENT_EXCEPTION = 'Invalid argument pattern for polymorphic function.';
  21. const ARGUMENT_TYPE_EXCEPTION = 'Invalid argument type.';
  22. const ARGUMENT_BOUNDS_EXCEPTION = 'Invalid argument range.';
  23. const MATRIX_DIMENSION_EXCEPTION = 'Matrix dimensions are not equal.';
  24. const ARRAY_LENGTH_EXCEPTION = 'Array length must be a multiple of m.';
  25. const MATRIX_SPD_EXCEPTION = 'Can only perform operation on symmetric positive definite matrix.';
  26. /**
  27. * Matrix storage.
  28. *
  29. * @var array
  30. */
  31. public $A = [];
  32. /**
  33. * Matrix row dimension.
  34. *
  35. * @var int
  36. */
  37. private $m;
  38. /**
  39. * Matrix column dimension.
  40. *
  41. * @var int
  42. */
  43. private $n;
  44. /**
  45. * Polymorphic constructor.
  46. *
  47. * As PHP has no support for polymorphic constructors, we use tricks to make our own sort of polymorphism using func_num_args, func_get_arg, and gettype. In essence, we're just implementing a simple RTTI filter and calling the appropriate constructor.
  48. */
  49. public function __construct(...$args)
  50. {
  51. if (count($args) > 0) {
  52. $match = implode(',', array_map('gettype', $args));
  53. switch ($match) {
  54. //Rectangular matrix - m x n initialized from 2D array
  55. case 'array':
  56. $this->m = count($args[0]);
  57. $this->n = count($args[0][0]);
  58. $this->A = $args[0];
  59. break;
  60. //Square matrix - n x n
  61. case 'integer':
  62. $this->m = $args[0];
  63. $this->n = $args[0];
  64. $this->A = array_fill(0, $this->m, array_fill(0, $this->n, 0));
  65. break;
  66. //Rectangular matrix - m x n
  67. case 'integer,integer':
  68. $this->m = $args[0];
  69. $this->n = $args[1];
  70. $this->A = array_fill(0, $this->m, array_fill(0, $this->n, 0));
  71. break;
  72. //Rectangular matrix - m x n initialized from packed array
  73. case 'array,integer':
  74. $this->m = $args[1];
  75. if ($this->m != 0) {
  76. $this->n = count($args[0]) / $this->m;
  77. } else {
  78. $this->n = 0;
  79. }
  80. if (($this->m * $this->n) == count($args[0])) {
  81. for ($i = 0; $i < $this->m; ++$i) {
  82. for ($j = 0; $j < $this->n; ++$j) {
  83. $this->A[$i][$j] = $args[0][$i + $j * $this->m];
  84. }
  85. }
  86. } else {
  87. throw new CalculationException(self::ARRAY_LENGTH_EXCEPTION);
  88. }
  89. break;
  90. default:
  91. throw new CalculationException(self::POLYMORPHIC_ARGUMENT_EXCEPTION);
  92. break;
  93. }
  94. } else {
  95. throw new CalculationException(self::POLYMORPHIC_ARGUMENT_EXCEPTION);
  96. }
  97. }
  98. /**
  99. * getArray.
  100. *
  101. * @return array Matrix array
  102. */
  103. public function getArray()
  104. {
  105. return $this->A;
  106. }
  107. /**
  108. * getRowDimension.
  109. *
  110. * @return int Row dimension
  111. */
  112. public function getRowDimension()
  113. {
  114. return $this->m;
  115. }
  116. /**
  117. * getColumnDimension.
  118. *
  119. * @return int Column dimension
  120. */
  121. public function getColumnDimension()
  122. {
  123. return $this->n;
  124. }
  125. /**
  126. * get.
  127. *
  128. * Get the i,j-th element of the matrix.
  129. *
  130. * @param int $i Row position
  131. * @param int $j Column position
  132. *
  133. * @return mixed Element (int/float/double)
  134. */
  135. public function get($i = null, $j = null)
  136. {
  137. return $this->A[$i][$j];
  138. }
  139. /**
  140. * getMatrix.
  141. *
  142. * Get a submatrix
  143. *
  144. * @return Matrix Submatrix
  145. */
  146. public function getMatrix(...$args)
  147. {
  148. if (count($args) > 0) {
  149. $match = implode(',', array_map('gettype', $args));
  150. switch ($match) {
  151. //A($i0...; $j0...)
  152. case 'integer,integer':
  153. [$i0, $j0] = $args;
  154. if ($i0 >= 0) {
  155. $m = $this->m - $i0;
  156. } else {
  157. throw new CalculationException(self::ARGUMENT_BOUNDS_EXCEPTION);
  158. }
  159. if ($j0 >= 0) {
  160. $n = $this->n - $j0;
  161. } else {
  162. throw new CalculationException(self::ARGUMENT_BOUNDS_EXCEPTION);
  163. }
  164. $R = new self($m, $n);
  165. for ($i = $i0; $i < $this->m; ++$i) {
  166. for ($j = $j0; $j < $this->n; ++$j) {
  167. $R->set($i, $j, $this->A[$i][$j]);
  168. }
  169. }
  170. return $R;
  171. break;
  172. //A($i0...$iF; $j0...$jF)
  173. case 'integer,integer,integer,integer':
  174. [$i0, $iF, $j0, $jF] = $args;
  175. if (($iF > $i0) && ($this->m >= $iF) && ($i0 >= 0)) {
  176. $m = $iF - $i0;
  177. } else {
  178. throw new CalculationException(self::ARGUMENT_BOUNDS_EXCEPTION);
  179. }
  180. if (($jF > $j0) && ($this->n >= $jF) && ($j0 >= 0)) {
  181. $n = $jF - $j0;
  182. } else {
  183. throw new CalculationException(self::ARGUMENT_BOUNDS_EXCEPTION);
  184. }
  185. $R = new self($m + 1, $n + 1);
  186. for ($i = $i0; $i <= $iF; ++$i) {
  187. for ($j = $j0; $j <= $jF; ++$j) {
  188. $R->set($i - $i0, $j - $j0, $this->A[$i][$j]);
  189. }
  190. }
  191. return $R;
  192. break;
  193. //$R = array of row indices; $C = array of column indices
  194. case 'array,array':
  195. [$RL, $CL] = $args;
  196. if (count($RL) > 0) {
  197. $m = count($RL);
  198. } else {
  199. throw new CalculationException(self::ARGUMENT_BOUNDS_EXCEPTION);
  200. }
  201. if (count($CL) > 0) {
  202. $n = count($CL);
  203. } else {
  204. throw new CalculationException(self::ARGUMENT_BOUNDS_EXCEPTION);
  205. }
  206. $R = new self($m, $n);
  207. for ($i = 0; $i < $m; ++$i) {
  208. for ($j = 0; $j < $n; ++$j) {
  209. $R->set($i, $j, $this->A[$RL[$i]][$CL[$j]]);
  210. }
  211. }
  212. return $R;
  213. break;
  214. //A($i0...$iF); $CL = array of column indices
  215. case 'integer,integer,array':
  216. [$i0, $iF, $CL] = $args;
  217. if (($iF > $i0) && ($this->m >= $iF) && ($i0 >= 0)) {
  218. $m = $iF - $i0;
  219. } else {
  220. throw new CalculationException(self::ARGUMENT_BOUNDS_EXCEPTION);
  221. }
  222. if (count($CL) > 0) {
  223. $n = count($CL);
  224. } else {
  225. throw new CalculationException(self::ARGUMENT_BOUNDS_EXCEPTION);
  226. }
  227. $R = new self($m, $n);
  228. for ($i = $i0; $i < $iF; ++$i) {
  229. for ($j = 0; $j < $n; ++$j) {
  230. $R->set($i - $i0, $j, $this->A[$i][$CL[$j]]);
  231. }
  232. }
  233. return $R;
  234. break;
  235. //$RL = array of row indices
  236. case 'array,integer,integer':
  237. [$RL, $j0, $jF] = $args;
  238. if (count($RL) > 0) {
  239. $m = count($RL);
  240. } else {
  241. throw new CalculationException(self::ARGUMENT_BOUNDS_EXCEPTION);
  242. }
  243. if (($jF >= $j0) && ($this->n >= $jF) && ($j0 >= 0)) {
  244. $n = $jF - $j0;
  245. } else {
  246. throw new CalculationException(self::ARGUMENT_BOUNDS_EXCEPTION);
  247. }
  248. $R = new self($m, $n + 1);
  249. for ($i = 0; $i < $m; ++$i) {
  250. for ($j = $j0; $j <= $jF; ++$j) {
  251. $R->set($i, $j - $j0, $this->A[$RL[$i]][$j]);
  252. }
  253. }
  254. return $R;
  255. break;
  256. default:
  257. throw new CalculationException(self::POLYMORPHIC_ARGUMENT_EXCEPTION);
  258. break;
  259. }
  260. } else {
  261. throw new CalculationException(self::POLYMORPHIC_ARGUMENT_EXCEPTION);
  262. }
  263. }
  264. /**
  265. * checkMatrixDimensions.
  266. *
  267. * Is matrix B the same size?
  268. *
  269. * @param Matrix $B Matrix B
  270. *
  271. * @return bool
  272. */
  273. public function checkMatrixDimensions($B = null)
  274. {
  275. if ($B instanceof self) {
  276. if (($this->m == $B->getRowDimension()) && ($this->n == $B->getColumnDimension())) {
  277. return true;
  278. }
  279. throw new CalculationException(self::MATRIX_DIMENSION_EXCEPTION);
  280. }
  281. throw new CalculationException(self::ARGUMENT_TYPE_EXCEPTION);
  282. }
  283. // function checkMatrixDimensions()
  284. /**
  285. * set.
  286. *
  287. * Set the i,j-th element of the matrix.
  288. *
  289. * @param int $i Row position
  290. * @param int $j Column position
  291. * @param mixed $c Int/float/double value
  292. *
  293. * @return mixed Element (int/float/double)
  294. */
  295. public function set($i = null, $j = null, $c = null)
  296. {
  297. // Optimized set version just has this
  298. $this->A[$i][$j] = $c;
  299. }
  300. // function set()
  301. /**
  302. * identity.
  303. *
  304. * Generate an identity matrix.
  305. *
  306. * @param int $m Row dimension
  307. * @param int $n Column dimension
  308. *
  309. * @return Matrix Identity matrix
  310. */
  311. public function identity($m = null, $n = null)
  312. {
  313. return $this->diagonal($m, $n, 1);
  314. }
  315. /**
  316. * diagonal.
  317. *
  318. * Generate a diagonal matrix
  319. *
  320. * @param int $m Row dimension
  321. * @param int $n Column dimension
  322. * @param mixed $c Diagonal value
  323. *
  324. * @return Matrix Diagonal matrix
  325. */
  326. public function diagonal($m = null, $n = null, $c = 1)
  327. {
  328. $R = new self($m, $n);
  329. for ($i = 0; $i < $m; ++$i) {
  330. $R->set($i, $i, $c);
  331. }
  332. return $R;
  333. }
  334. /**
  335. * getMatrixByRow.
  336. *
  337. * Get a submatrix by row index/range
  338. *
  339. * @param int $i0 Initial row index
  340. * @param int $iF Final row index
  341. *
  342. * @return Matrix Submatrix
  343. */
  344. public function getMatrixByRow($i0 = null, $iF = null)
  345. {
  346. if (is_int($i0)) {
  347. if (is_int($iF)) {
  348. return $this->getMatrix($i0, 0, $iF + 1, $this->n);
  349. }
  350. return $this->getMatrix($i0, 0, $i0 + 1, $this->n);
  351. }
  352. throw new CalculationException(self::ARGUMENT_TYPE_EXCEPTION);
  353. }
  354. /**
  355. * getMatrixByCol.
  356. *
  357. * Get a submatrix by column index/range
  358. *
  359. * @param int $j0 Initial column index
  360. * @param int $jF Final column index
  361. *
  362. * @return Matrix Submatrix
  363. */
  364. public function getMatrixByCol($j0 = null, $jF = null)
  365. {
  366. if (is_int($j0)) {
  367. if (is_int($jF)) {
  368. return $this->getMatrix(0, $j0, $this->m, $jF + 1);
  369. }
  370. return $this->getMatrix(0, $j0, $this->m, $j0 + 1);
  371. }
  372. throw new CalculationException(self::ARGUMENT_TYPE_EXCEPTION);
  373. }
  374. /**
  375. * transpose.
  376. *
  377. * Tranpose matrix
  378. *
  379. * @return Matrix Transposed matrix
  380. */
  381. public function transpose()
  382. {
  383. $R = new self($this->n, $this->m);
  384. for ($i = 0; $i < $this->m; ++$i) {
  385. for ($j = 0; $j < $this->n; ++$j) {
  386. $R->set($j, $i, $this->A[$i][$j]);
  387. }
  388. }
  389. return $R;
  390. }
  391. // function transpose()
  392. /**
  393. * trace.
  394. *
  395. * Sum of diagonal elements
  396. *
  397. * @return float Sum of diagonal elements
  398. */
  399. public function trace()
  400. {
  401. $s = 0;
  402. $n = min($this->m, $this->n);
  403. for ($i = 0; $i < $n; ++$i) {
  404. $s += $this->A[$i][$i];
  405. }
  406. return $s;
  407. }
  408. /**
  409. * uminus.
  410. *
  411. * Unary minus matrix -A
  412. *
  413. * @return Matrix Unary minus matrix
  414. */
  415. public function uminus()
  416. {
  417. }
  418. /**
  419. * plus.
  420. *
  421. * A + B
  422. *
  423. * @return Matrix Sum
  424. */
  425. public function plus(...$args)
  426. {
  427. if (count($args) > 0) {
  428. $match = implode(',', array_map('gettype', $args));
  429. switch ($match) {
  430. case 'object':
  431. if ($args[0] instanceof self) {
  432. $M = $args[0];
  433. } else {
  434. throw new CalculationException(self::ARGUMENT_TYPE_EXCEPTION);
  435. }
  436. break;
  437. case 'array':
  438. $M = new self($args[0]);
  439. break;
  440. default:
  441. throw new CalculationException(self::POLYMORPHIC_ARGUMENT_EXCEPTION);
  442. break;
  443. }
  444. $this->checkMatrixDimensions($M);
  445. for ($i = 0; $i < $this->m; ++$i) {
  446. for ($j = 0; $j < $this->n; ++$j) {
  447. $M->set($i, $j, $M->get($i, $j) + $this->A[$i][$j]);
  448. }
  449. }
  450. return $M;
  451. }
  452. throw new CalculationException(self::POLYMORPHIC_ARGUMENT_EXCEPTION);
  453. }
  454. /**
  455. * plusEquals.
  456. *
  457. * A = A + B
  458. *
  459. * @return $this
  460. */
  461. public function plusEquals(...$args)
  462. {
  463. if (count($args) > 0) {
  464. $match = implode(',', array_map('gettype', $args));
  465. switch ($match) {
  466. case 'object':
  467. if ($args[0] instanceof self) {
  468. $M = $args[0];
  469. } else {
  470. throw new CalculationException(self::ARGUMENT_TYPE_EXCEPTION);
  471. }
  472. break;
  473. case 'array':
  474. $M = new self($args[0]);
  475. break;
  476. default:
  477. throw new CalculationException(self::POLYMORPHIC_ARGUMENT_EXCEPTION);
  478. break;
  479. }
  480. $this->checkMatrixDimensions($M);
  481. for ($i = 0; $i < $this->m; ++$i) {
  482. for ($j = 0; $j < $this->n; ++$j) {
  483. $validValues = true;
  484. $value = $M->get($i, $j);
  485. if ((is_string($this->A[$i][$j])) && (strlen($this->A[$i][$j]) > 0) && (!is_numeric($this->A[$i][$j]))) {
  486. $this->A[$i][$j] = trim($this->A[$i][$j], '"');
  487. $validValues &= StringHelper::convertToNumberIfFraction($this->A[$i][$j]);
  488. }
  489. if ((is_string($value)) && (strlen($value) > 0) && (!is_numeric($value))) {
  490. $value = trim($value, '"');
  491. $validValues &= StringHelper::convertToNumberIfFraction($value);
  492. }
  493. if ($validValues) {
  494. $this->A[$i][$j] += $value;
  495. } else {
  496. $this->A[$i][$j] = Functions::NAN();
  497. }
  498. }
  499. }
  500. return $this;
  501. }
  502. throw new CalculationException(self::POLYMORPHIC_ARGUMENT_EXCEPTION);
  503. }
  504. /**
  505. * minus.
  506. *
  507. * A - B
  508. *
  509. * @return Matrix Sum
  510. */
  511. public function minus(...$args)
  512. {
  513. if (count($args) > 0) {
  514. $match = implode(',', array_map('gettype', $args));
  515. switch ($match) {
  516. case 'object':
  517. if ($args[0] instanceof self) {
  518. $M = $args[0];
  519. } else {
  520. throw new CalculationException(self::ARGUMENT_TYPE_EXCEPTION);
  521. }
  522. break;
  523. case 'array':
  524. $M = new self($args[0]);
  525. break;
  526. default:
  527. throw new CalculationException(self::POLYMORPHIC_ARGUMENT_EXCEPTION);
  528. break;
  529. }
  530. $this->checkMatrixDimensions($M);
  531. for ($i = 0; $i < $this->m; ++$i) {
  532. for ($j = 0; $j < $this->n; ++$j) {
  533. $M->set($i, $j, $M->get($i, $j) - $this->A[$i][$j]);
  534. }
  535. }
  536. return $M;
  537. }
  538. throw new CalculationException(self::POLYMORPHIC_ARGUMENT_EXCEPTION);
  539. }
  540. /**
  541. * minusEquals.
  542. *
  543. * A = A - B
  544. *
  545. * @return $this
  546. */
  547. public function minusEquals(...$args)
  548. {
  549. if (count($args) > 0) {
  550. $match = implode(',', array_map('gettype', $args));
  551. switch ($match) {
  552. case 'object':
  553. if ($args[0] instanceof self) {
  554. $M = $args[0];
  555. } else {
  556. throw new CalculationException(self::ARGUMENT_TYPE_EXCEPTION);
  557. }
  558. break;
  559. case 'array':
  560. $M = new self($args[0]);
  561. break;
  562. default:
  563. throw new CalculationException(self::POLYMORPHIC_ARGUMENT_EXCEPTION);
  564. break;
  565. }
  566. $this->checkMatrixDimensions($M);
  567. for ($i = 0; $i < $this->m; ++$i) {
  568. for ($j = 0; $j < $this->n; ++$j) {
  569. $validValues = true;
  570. $value = $M->get($i, $j);
  571. if ((is_string($this->A[$i][$j])) && (strlen($this->A[$i][$j]) > 0) && (!is_numeric($this->A[$i][$j]))) {
  572. $this->A[$i][$j] = trim($this->A[$i][$j], '"');
  573. $validValues &= StringHelper::convertToNumberIfFraction($this->A[$i][$j]);
  574. }
  575. if ((is_string($value)) && (strlen($value) > 0) && (!is_numeric($value))) {
  576. $value = trim($value, '"');
  577. $validValues &= StringHelper::convertToNumberIfFraction($value);
  578. }
  579. if ($validValues) {
  580. $this->A[$i][$j] -= $value;
  581. } else {
  582. $this->A[$i][$j] = Functions::NAN();
  583. }
  584. }
  585. }
  586. return $this;
  587. }
  588. throw new CalculationException(self::POLYMORPHIC_ARGUMENT_EXCEPTION);
  589. }
  590. /**
  591. * arrayTimes.
  592. *
  593. * Element-by-element multiplication
  594. * Cij = Aij * Bij
  595. *
  596. * @return Matrix Matrix Cij
  597. */
  598. public function arrayTimes(...$args)
  599. {
  600. if (count($args) > 0) {
  601. $match = implode(',', array_map('gettype', $args));
  602. switch ($match) {
  603. case 'object':
  604. if ($args[0] instanceof self) {
  605. $M = $args[0];
  606. } else {
  607. throw new CalculationException(self::ARGUMENT_TYPE_EXCEPTION);
  608. }
  609. break;
  610. case 'array':
  611. $M = new self($args[0]);
  612. break;
  613. default:
  614. throw new CalculationException(self::POLYMORPHIC_ARGUMENT_EXCEPTION);
  615. break;
  616. }
  617. $this->checkMatrixDimensions($M);
  618. for ($i = 0; $i < $this->m; ++$i) {
  619. for ($j = 0; $j < $this->n; ++$j) {
  620. $M->set($i, $j, $M->get($i, $j) * $this->A[$i][$j]);
  621. }
  622. }
  623. return $M;
  624. }
  625. throw new CalculationException(self::POLYMORPHIC_ARGUMENT_EXCEPTION);
  626. }
  627. /**
  628. * arrayTimesEquals.
  629. *
  630. * Element-by-element multiplication
  631. * Aij = Aij * Bij
  632. *
  633. * @return $this
  634. */
  635. public function arrayTimesEquals(...$args)
  636. {
  637. if (count($args) > 0) {
  638. $match = implode(',', array_map('gettype', $args));
  639. switch ($match) {
  640. case 'object':
  641. if ($args[0] instanceof self) {
  642. $M = $args[0];
  643. } else {
  644. throw new CalculationException(self::ARGUMENT_TYPE_EXCEPTION);
  645. }
  646. break;
  647. case 'array':
  648. $M = new self($args[0]);
  649. break;
  650. default:
  651. throw new CalculationException(self::POLYMORPHIC_ARGUMENT_EXCEPTION);
  652. break;
  653. }
  654. $this->checkMatrixDimensions($M);
  655. for ($i = 0; $i < $this->m; ++$i) {
  656. for ($j = 0; $j < $this->n; ++$j) {
  657. $validValues = true;
  658. $value = $M->get($i, $j);
  659. if ((is_string($this->A[$i][$j])) && (strlen($this->A[$i][$j]) > 0) && (!is_numeric($this->A[$i][$j]))) {
  660. $this->A[$i][$j] = trim($this->A[$i][$j], '"');
  661. $validValues &= StringHelper::convertToNumberIfFraction($this->A[$i][$j]);
  662. }
  663. if ((is_string($value)) && (strlen($value) > 0) && (!is_numeric($value))) {
  664. $value = trim($value, '"');
  665. $validValues &= StringHelper::convertToNumberIfFraction($value);
  666. }
  667. if ($validValues) {
  668. $this->A[$i][$j] *= $value;
  669. } else {
  670. $this->A[$i][$j] = Functions::NAN();
  671. }
  672. }
  673. }
  674. return $this;
  675. }
  676. throw new CalculationException(self::POLYMORPHIC_ARGUMENT_EXCEPTION);
  677. }
  678. /**
  679. * arrayRightDivide.
  680. *
  681. * Element-by-element right division
  682. * A / B
  683. *
  684. * @return Matrix Division result
  685. */
  686. public function arrayRightDivide(...$args)
  687. {
  688. if (count($args) > 0) {
  689. $match = implode(',', array_map('gettype', $args));
  690. switch ($match) {
  691. case 'object':
  692. if ($args[0] instanceof self) {
  693. $M = $args[0];
  694. } else {
  695. throw new CalculationException(self::ARGUMENT_TYPE_EXCEPTION);
  696. }
  697. break;
  698. case 'array':
  699. $M = new self($args[0]);
  700. break;
  701. default:
  702. throw new CalculationException(self::POLYMORPHIC_ARGUMENT_EXCEPTION);
  703. break;
  704. }
  705. $this->checkMatrixDimensions($M);
  706. for ($i = 0; $i < $this->m; ++$i) {
  707. for ($j = 0; $j < $this->n; ++$j) {
  708. $validValues = true;
  709. $value = $M->get($i, $j);
  710. if ((is_string($this->A[$i][$j])) && (strlen($this->A[$i][$j]) > 0) && (!is_numeric($this->A[$i][$j]))) {
  711. $this->A[$i][$j] = trim($this->A[$i][$j], '"');
  712. $validValues &= StringHelper::convertToNumberIfFraction($this->A[$i][$j]);
  713. }
  714. if ((is_string($value)) && (strlen($value) > 0) && (!is_numeric($value))) {
  715. $value = trim($value, '"');
  716. $validValues &= StringHelper::convertToNumberIfFraction($value);
  717. }
  718. if ($validValues) {
  719. if ($value == 0) {
  720. // Trap for Divide by Zero error
  721. $M->set($i, $j, '#DIV/0!');
  722. } else {
  723. $M->set($i, $j, $this->A[$i][$j] / $value);
  724. }
  725. } else {
  726. $M->set($i, $j, Functions::NAN());
  727. }
  728. }
  729. }
  730. return $M;
  731. }
  732. throw new CalculationException(self::POLYMORPHIC_ARGUMENT_EXCEPTION);
  733. }
  734. /**
  735. * arrayRightDivideEquals.
  736. *
  737. * Element-by-element right division
  738. * Aij = Aij / Bij
  739. *
  740. * @return Matrix Matrix Aij
  741. */
  742. public function arrayRightDivideEquals(...$args)
  743. {
  744. if (count($args) > 0) {
  745. $match = implode(',', array_map('gettype', $args));
  746. switch ($match) {
  747. case 'object':
  748. if ($args[0] instanceof self) {
  749. $M = $args[0];
  750. } else {
  751. throw new CalculationException(self::ARGUMENT_TYPE_EXCEPTION);
  752. }
  753. break;
  754. case 'array':
  755. $M = new self($args[0]);
  756. break;
  757. default:
  758. throw new CalculationException(self::POLYMORPHIC_ARGUMENT_EXCEPTION);
  759. break;
  760. }
  761. $this->checkMatrixDimensions($M);
  762. for ($i = 0; $i < $this->m; ++$i) {
  763. for ($j = 0; $j < $this->n; ++$j) {
  764. $this->A[$i][$j] = $this->A[$i][$j] / $M->get($i, $j);
  765. }
  766. }
  767. return $M;
  768. }
  769. throw new CalculationException(self::POLYMORPHIC_ARGUMENT_EXCEPTION);
  770. }
  771. /**
  772. * arrayLeftDivide.
  773. *
  774. * Element-by-element Left division
  775. * A / B
  776. *
  777. * @return Matrix Division result
  778. */
  779. public function arrayLeftDivide(...$args)
  780. {
  781. if (count($args) > 0) {
  782. $match = implode(',', array_map('gettype', $args));
  783. switch ($match) {
  784. case 'object':
  785. if ($args[0] instanceof self) {
  786. $M = $args[0];
  787. } else {
  788. throw new CalculationException(self::ARGUMENT_TYPE_EXCEPTION);
  789. }
  790. break;
  791. case 'array':
  792. $M = new self($args[0]);
  793. break;
  794. default:
  795. throw new CalculationException(self::POLYMORPHIC_ARGUMENT_EXCEPTION);
  796. break;
  797. }
  798. $this->checkMatrixDimensions($M);
  799. for ($i = 0; $i < $this->m; ++$i) {
  800. for ($j = 0; $j < $this->n; ++$j) {
  801. $M->set($i, $j, $M->get($i, $j) / $this->A[$i][$j]);
  802. }
  803. }
  804. return $M;
  805. }
  806. throw new CalculationException(self::POLYMORPHIC_ARGUMENT_EXCEPTION);
  807. }
  808. /**
  809. * arrayLeftDivideEquals.
  810. *
  811. * Element-by-element Left division
  812. * Aij = Aij / Bij
  813. *
  814. * @return Matrix Matrix Aij
  815. */
  816. public function arrayLeftDivideEquals(...$args)
  817. {
  818. if (count($args) > 0) {
  819. $match = implode(',', array_map('gettype', $args));
  820. switch ($match) {
  821. case 'object':
  822. if ($args[0] instanceof self) {
  823. $M = $args[0];
  824. } else {
  825. throw new CalculationException(self::ARGUMENT_TYPE_EXCEPTION);
  826. }
  827. break;
  828. case 'array':
  829. $M = new self($args[0]);
  830. break;
  831. default:
  832. throw new CalculationException(self::POLYMORPHIC_ARGUMENT_EXCEPTION);
  833. break;
  834. }
  835. $this->checkMatrixDimensions($M);
  836. for ($i = 0; $i < $this->m; ++$i) {
  837. for ($j = 0; $j < $this->n; ++$j) {
  838. $this->A[$i][$j] = $M->get($i, $j) / $this->A[$i][$j];
  839. }
  840. }
  841. return $M;
  842. }
  843. throw new CalculationException(self::POLYMORPHIC_ARGUMENT_EXCEPTION);
  844. }
  845. /**
  846. * times.
  847. *
  848. * Matrix multiplication
  849. *
  850. * @return Matrix Product
  851. */
  852. public function times(...$args)
  853. {
  854. if (count($args) > 0) {
  855. $match = implode(',', array_map('gettype', $args));
  856. switch ($match) {
  857. case 'object':
  858. if ($args[0] instanceof self) {
  859. $B = $args[0];
  860. } else {
  861. throw new CalculationException(self::ARGUMENT_TYPE_EXCEPTION);
  862. }
  863. if ($this->n == $B->m) {
  864. $C = new self($this->m, $B->n);
  865. for ($j = 0; $j < $B->n; ++$j) {
  866. $Bcolj = [];
  867. for ($k = 0; $k < $this->n; ++$k) {
  868. $Bcolj[$k] = $B->A[$k][$j];
  869. }
  870. for ($i = 0; $i < $this->m; ++$i) {
  871. $Arowi = $this->A[$i];
  872. $s = 0;
  873. for ($k = 0; $k < $this->n; ++$k) {
  874. $s += $Arowi[$k] * $Bcolj[$k];
  875. }
  876. $C->A[$i][$j] = $s;
  877. }
  878. }
  879. return $C;
  880. }
  881. throw new CalculationException(self::MATRIX_DIMENSION_EXCEPTION);
  882. case 'array':
  883. $B = new self($args[0]);
  884. if ($this->n == $B->m) {
  885. $C = new self($this->m, $B->n);
  886. for ($i = 0; $i < $C->m; ++$i) {
  887. for ($j = 0; $j < $C->n; ++$j) {
  888. $s = '0';
  889. for ($k = 0; $k < $C->n; ++$k) {
  890. $s += $this->A[$i][$k] * $B->A[$k][$j];
  891. }
  892. $C->A[$i][$j] = $s;
  893. }
  894. }
  895. return $C;
  896. }
  897. throw new CalculationException(self::MATRIX_DIMENSION_EXCEPTION);
  898. case 'integer':
  899. $C = new self($this->A);
  900. for ($i = 0; $i < $C->m; ++$i) {
  901. for ($j = 0; $j < $C->n; ++$j) {
  902. $C->A[$i][$j] *= $args[0];
  903. }
  904. }
  905. return $C;
  906. case 'double':
  907. $C = new self($this->m, $this->n);
  908. for ($i = 0; $i < $C->m; ++$i) {
  909. for ($j = 0; $j < $C->n; ++$j) {
  910. $C->A[$i][$j] = $args[0] * $this->A[$i][$j];
  911. }
  912. }
  913. return $C;
  914. case 'float':
  915. $C = new self($this->A);
  916. for ($i = 0; $i < $C->m; ++$i) {
  917. for ($j = 0; $j < $C->n; ++$j) {
  918. $C->A[$i][$j] *= $args[0];
  919. }
  920. }
  921. return $C;
  922. default:
  923. throw new CalculationException(self::POLYMORPHIC_ARGUMENT_EXCEPTION);
  924. }
  925. } else {
  926. throw new CalculationException(self::POLYMORPHIC_ARGUMENT_EXCEPTION);
  927. }
  928. }
  929. /**
  930. * power.
  931. *
  932. * A = A ^ B
  933. *
  934. * @return $this
  935. */
  936. public function power(...$args)
  937. {
  938. if (count($args) > 0) {
  939. $match = implode(',', array_map('gettype', $args));
  940. switch ($match) {
  941. case 'object':
  942. if ($args[0] instanceof self) {
  943. $M = $args[0];
  944. } else {
  945. throw new CalculationException(self::ARGUMENT_TYPE_EXCEPTION);
  946. }
  947. break;
  948. case 'array':
  949. $M = new self($args[0]);
  950. break;
  951. default:
  952. throw new CalculationException(self::POLYMORPHIC_ARGUMENT_EXCEPTION);
  953. break;
  954. }
  955. $this->checkMatrixDimensions($M);
  956. for ($i = 0; $i < $this->m; ++$i) {
  957. for ($j = 0; $j < $this->n; ++$j) {
  958. $validValues = true;
  959. $value = $M->get($i, $j);
  960. if ((is_string($this->A[$i][$j])) && (strlen($this->A[$i][$j]) > 0) && (!is_numeric($this->A[$i][$j]))) {
  961. $this->A[$i][$j] = trim($this->A[$i][$j], '"');
  962. $validValues &= StringHelper::convertToNumberIfFraction($this->A[$i][$j]);
  963. }
  964. if ((is_string($value)) && (strlen($value) > 0) && (!is_numeric($value))) {
  965. $value = trim($value, '"');
  966. $validValues &= StringHelper::convertToNumberIfFraction($value);
  967. }
  968. if ($validValues) {
  969. $this->A[$i][$j] = $this->A[$i][$j] ** $value;
  970. } else {
  971. $this->A[$i][$j] = Functions::NAN();
  972. }
  973. }
  974. }
  975. return $this;
  976. }
  977. throw new CalculationException(self::POLYMORPHIC_ARGUMENT_EXCEPTION);
  978. }
  979. /**
  980. * concat.
  981. *
  982. * A = A & B
  983. *
  984. * @return $this
  985. */
  986. public function concat(...$args)
  987. {
  988. if (count($args) > 0) {
  989. $match = implode(',', array_map('gettype', $args));
  990. switch ($match) {
  991. case 'object':
  992. if ($args[0] instanceof self) {
  993. $M = $args[0];
  994. } else {
  995. throw new CalculationException(self::ARGUMENT_TYPE_EXCEPTION);
  996. }
  997. break;
  998. case 'array':
  999. $M = new self($args[0]);
  1000. break;
  1001. default:
  1002. throw new CalculationException(self::POLYMORPHIC_ARGUMENT_EXCEPTION);
  1003. break;
  1004. }
  1005. $this->checkMatrixDimensions($M);
  1006. for ($i = 0; $i < $this->m; ++$i) {
  1007. for ($j = 0; $j < $this->n; ++$j) {
  1008. $this->A[$i][$j] = trim($this->A[$i][$j], '"') . trim($M->get($i, $j), '"');
  1009. }
  1010. }
  1011. return $this;
  1012. }
  1013. throw new CalculationException(self::POLYMORPHIC_ARGUMENT_EXCEPTION);
  1014. }
  1015. /**
  1016. * Solve A*X = B.
  1017. *
  1018. * @param Matrix $B Right hand side
  1019. *
  1020. * @return Matrix ... Solution if A is square, least squares solution otherwise
  1021. */
  1022. public function solve($B)
  1023. {
  1024. if ($this->m == $this->n) {
  1025. $LU = new LUDecomposition($this);
  1026. return $LU->solve($B);
  1027. }
  1028. $QR = new QRDecomposition($this);
  1029. return $QR->solve($B);
  1030. }
  1031. /**
  1032. * Matrix inverse or pseudoinverse.
  1033. *
  1034. * @return Matrix ... Inverse(A) if A is square, pseudoinverse otherwise.
  1035. */
  1036. public function inverse()
  1037. {
  1038. return $this->solve($this->identity($this->m, $this->m));
  1039. }
  1040. /**
  1041. * det.
  1042. *
  1043. * Calculate determinant
  1044. *
  1045. * @return float Determinant
  1046. */
  1047. public function det()
  1048. {
  1049. $L = new LUDecomposition($this);
  1050. return $L->det();
  1051. }
  1052. }