[Summary view]

Differences Between: [Versions 39 and 311] [Versions 39 and 400] [Versions 39 and 401]

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