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

   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 http://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                     list($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                     list($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                     list($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 - $i0, $j - $j0, $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                     list($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[$RL[$i]][$j]);
 261                         }
 262                     }
 263
 264                     return $R;
 265
 266                     break;
 267                 //$RL = array of row indices
 268                 case 'array,integer,integer':
 269                     list($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                             for ($k = 0; $k < $this->n; ++$k) {
1017                                 $Bcolj[$k] = $B->A[$k][$j];
1018                             }
1019                             for ($i = 0; $i < $this->m; ++$i) {
1020                                 $Arowi = $this->A[$i];
1021                                 $s = 0;
1022                                 for ($k = 0; $k < $this->n; ++$k) {
1023                                     $s += $Arowi[$k] * $Bcolj[$k];
1024                                 }
1025                                 $C->A[$i][$j] = $s;
1026                             }
1027                         }
1028
1029                         return $C;
1030                     }
1031
1032                     throw new CalculationException(self::MATRIX_DIMENSION_EXCEPTION);
1033                 case 'array':
1034                     $B = new self($args[0]);
1035                     if ($this->n == $B->m) {
1036                         $C = new self($this->m, $B->n);
1037                         for ($i = 0; $i < $C->m; ++$i) {
1038                             for ($j = 0; $j < $C->n; ++$j) {
1039                                 $s = '0';
1040                                 for ($k = 0; $k < $C->n; ++$k) {
1041                                     $s += $this->A[$i][$k] * $B->A[$k][$j];
1042                                 }
1043                                 $C->A[$i][$j] = $s;
1044                             }
1045                         }
1046
1047                         return $C;
1048                     }
1049
1050                     throw new CalculationException(self::MATRIX_DIMENSION_EXCEPTION);
1051                 case 'integer':
1052                     $C = new self($this->A);
1053                     for ($i = 0; $i < $C->m; ++$i) {
1054                         for ($j = 0; $j < $C->n; ++$j) {
1055                             $C->A[$i][$j] *= $args[0];
1056                         }
1057                     }
1058
1059                     return $C;
1060                 case 'double':
1061                     $C = new self($this->m, $this->n);
1062                     for ($i = 0; $i < $C->m; ++$i) {
1063                         for ($j = 0; $j < $C->n; ++$j) {
1064                             $C->A[$i][$j] = $args[0] * $this->A[$i][$j];
1065                         }
1066                     }
1067
1068                     return $C;
1069                 case 'float':
1070                     $C = new self($this->A);
1071                     for ($i = 0; $i < $C->m; ++$i) {
1072                         for ($j = 0; $j < $C->n; ++$j) {
1073                             $C->A[$i][$j] *= $args[0];
1074                         }
1075                     }
1076
1077                     return $C;
1078                 default:
1079                     throw new CalculationException(self::POLYMORPHIC_ARGUMENT_EXCEPTION);
1080             }
1081         } else {
1082             throw new CalculationException(self::POLYMORPHIC_ARGUMENT_EXCEPTION);
1083         }
1084     }
1085
1086     /**
1087      * power.
1088      *
1089      *    A = A ^ B
1090      *
1091      * @param mixed $B Matrix/Array
1092      *
1093      * @return Matrix Sum
1094      */
1095     public function power(...$args)
1096     {
1097         if (count($args) > 0) {
1098             $match = implode(',', array_map('gettype', $args));
1099
1100             switch ($match) {
1101                 case 'object':
1102                     if ($args[0] instanceof self) {
1103                         $M = $args[0];
1104                     } else {
1105                         throw new CalculationException(self::ARGUMENT_TYPE_EXCEPTION);
1106                     }
1107
1108                     break;
1109                 case 'array':
1110                     $M = new self($args[0]);
1111
1112                     break;
1113                 default:
1114                     throw new CalculationException(self::POLYMORPHIC_ARGUMENT_EXCEPTION);
1115
1116                     break;
1117             }
1118             $this->checkMatrixDimensions($M);
1119             for ($i = 0; $i < $this->m; ++$i) {
1120                 for ($j = 0; $j < $this->n; ++$j) {
1121                     $validValues = true;
1122                     $value = $M->get($i, $j);
1123                     if ((is_string($this->A[$i][$j])) && (strlen($this->A[$i][$j]) > 0) && (!is_numeric($this->A[$i][$j]))) {
1124                         $this->A[$i][$j] = trim($this->A[$i][$j], '"');
1125                         $validValues &= StringHelper::convertToNumberIfFraction($this->A[$i][$j]);
1126                     }
1127                     if ((is_string($value)) && (strlen($value) > 0) && (!is_numeric($value))) {
1128                         $value = trim($value, '"');
1129                         $validValues &= StringHelper::convertToNumberIfFraction($value);
1130                     }
1131                     if ($validValues) {
1132                         $this->A[$i][$j] = pow($this->A[$i][$j], $value);
1133                     } else {
1134                         $this->A[$i][$j] = Functions::NAN();
1135                     }
1136                 }
1137             }
1138
1139             return $this;
1140         }
1141
1142         throw new CalculationException(self::POLYMORPHIC_ARGUMENT_EXCEPTION);
1143     }
1144
1145     /**
1146      * concat.
1147      *
1148      *    A = A & B
1149      *
1150      * @param mixed $B Matrix/Array
1151      *
1152      * @return Matrix Sum
1153      */
1154     public function concat(...$args)
1155     {
1156         if (count($args) > 0) {
1157             $match = implode(',', array_map('gettype', $args));
1158
1159             switch ($match) {
1160                 case 'object':
1161                     if ($args[0] instanceof self) {
1162                         $M = $args[0];
1163                     } else {
1164                         throw new CalculationException(self::ARGUMENT_TYPE_EXCEPTION);
1165                     }
1166
1167                     break;
1168                 case 'array':
1169                     $M = new self($args[0]);
1170
1171                     break;
1172                 default:
1173                     throw new CalculationException(self::POLYMORPHIC_ARGUMENT_EXCEPTION);
1174
1175                     break;
1176             }
1177             $this->checkMatrixDimensions($M);
1178             for ($i = 0; $i < $this->m; ++$i) {
1179                 for ($j = 0; $j < $this->n; ++$j) {
1180                     $this->A[$i][$j] = trim($this->A[$i][$j], '"') . trim($M->get($i, $j), '"');
1181                 }
1182             }
1183
1184             return $this;
1185         }
1186
1187         throw new CalculationException(self::POLYMORPHIC_ARGUMENT_EXCEPTION);
1188     }
1189
1190     /**
1191      * Solve A*X = B.
1192      *
1193      * @param Matrix $B Right hand side
1194      *
1195      * @return Matrix ... Solution if A is square, least squares solution otherwise
1196      */
1197     public function solve($B)
1198     {
1199         if ($this->m == $this->n) {
1200             $LU = new LUDecomposition($this);
1201
1202             return $LU->solve($B);
1203         }
1204         $QR = new QRDecomposition($this);
1205
1206         return $QR->solve($B);
1207     }
1208
1209     /**
1210      * Matrix inverse or pseudoinverse.
1211      *
1212      * @return Matrix ... Inverse(A) if A is square, pseudoinverse otherwise.
1213      */
1214     public function inverse()
1215     {
1216         return $this->solve($this->identity($this->m, $this->m));
1217     }
1218
1219     /**
1220      * det.
1221      *
1222      *    Calculate determinant
1223      *
1224      * @return float Determinant
1225      */
1226     public function det()
1227     {
1228         $L = new LUDecomposition($this);
1229
1230         return $L->det();
1231     }
1232 }