From 75160b12821f7f4299cce7f0b69c83c1502ae071 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Anton=20Luka=20=C5=A0ijanec?= Date: Mon, 27 May 2024 13:08:29 +0200 Subject: 2024-02-19 upstream --- .../Shared/JAMA/CholeskyDecomposition.php | 147 +++++++++++++++++++++ 1 file changed, 147 insertions(+) create mode 100644 vendor/phpoffice/phpspreadsheet/src/PhpSpreadsheet/Shared/JAMA/CholeskyDecomposition.php (limited to 'vendor/phpoffice/phpspreadsheet/src/PhpSpreadsheet/Shared/JAMA/CholeskyDecomposition.php') diff --git a/vendor/phpoffice/phpspreadsheet/src/PhpSpreadsheet/Shared/JAMA/CholeskyDecomposition.php b/vendor/phpoffice/phpspreadsheet/src/PhpSpreadsheet/Shared/JAMA/CholeskyDecomposition.php new file mode 100644 index 0000000..4f7c666 --- /dev/null +++ b/vendor/phpoffice/phpspreadsheet/src/PhpSpreadsheet/Shared/JAMA/CholeskyDecomposition.php @@ -0,0 +1,147 @@ +L = $A->getArray(); + $this->m = $A->getRowDimension(); + + for ($i = 0; $i < $this->m; ++$i) { + for ($j = $i; $j < $this->m; ++$j) { + for ($sum = $this->L[$i][$j], $k = $i - 1; $k >= 0; --$k) { + $sum -= $this->L[$i][$k] * $this->L[$j][$k]; + } + if ($i == $j) { + if ($sum >= 0) { + $this->L[$i][$i] = sqrt($sum); + } else { + $this->isspd = false; + } + } else { + if ($this->L[$i][$i] != 0) { + $this->L[$j][$i] = $sum / $this->L[$i][$i]; + } + } + } + + for ($k = $i + 1; $k < $this->m; ++$k) { + $this->L[$i][$k] = 0.0; + } + } + } + + /** + * Is the matrix symmetric and positive definite? + * + * @return bool + */ + public function isSPD() + { + return $this->isspd; + } + + /** + * getL. + * + * Return triangular factor. + * + * @return Matrix Lower triangular matrix + */ + public function getL() + { + return new Matrix($this->L); + } + + /** + * Solve A*X = B. + * + * @param $B Row-equal matrix + * + * @return Matrix L * L' * X = B + */ + public function solve(Matrix $B) + { + if ($B->getRowDimension() == $this->m) { + if ($this->isspd) { + $X = $B->getArrayCopy(); + $nx = $B->getColumnDimension(); + + for ($k = 0; $k < $this->m; ++$k) { + for ($i = $k + 1; $i < $this->m; ++$i) { + for ($j = 0; $j < $nx; ++$j) { + $X[$i][$j] -= $X[$k][$j] * $this->L[$i][$k]; + } + } + for ($j = 0; $j < $nx; ++$j) { + $X[$k][$j] /= $this->L[$k][$k]; + } + } + + for ($k = $this->m - 1; $k >= 0; --$k) { + for ($j = 0; $j < $nx; ++$j) { + $X[$k][$j] /= $this->L[$k][$k]; + } + for ($i = 0; $i < $k; ++$i) { + for ($j = 0; $j < $nx; ++$j) { + $X[$i][$j] -= $X[$k][$j] * $this->L[$k][$i]; + } + } + } + + return new Matrix($X, $this->m, $nx); + } + + throw new CalculationException(Matrix::MATRIX_SPD_EXCEPTION); + } + + throw new CalculationException(Matrix::MATRIX_DIMENSION_EXCEPTION); + } +} -- cgit v1.2.3