diff options
Diffstat (limited to 'admin/survey/excel/PHPExcel/Shared/JAMA/CholeskyDecomposition.php')
-rw-r--r-- | admin/survey/excel/PHPExcel/Shared/JAMA/CholeskyDecomposition.php | 149 |
1 files changed, 149 insertions, 0 deletions
diff --git a/admin/survey/excel/PHPExcel/Shared/JAMA/CholeskyDecomposition.php b/admin/survey/excel/PHPExcel/Shared/JAMA/CholeskyDecomposition.php new file mode 100644 index 0000000..2d1fdfb --- /dev/null +++ b/admin/survey/excel/PHPExcel/Shared/JAMA/CholeskyDecomposition.php @@ -0,0 +1,149 @@ +<?php
+/**
+ * @package JAMA
+ *
+ * Cholesky decomposition class
+ *
+ * For a symmetric, positive definite matrix A, the Cholesky decomposition
+ * is an lower triangular matrix L so that A = L*L'.
+ *
+ * If the matrix is not symmetric or positive definite, the constructor
+ * returns a partial decomposition and sets an internal flag that may
+ * be queried by the isSPD() method.
+ *
+ * @author Paul Meagher
+ * @author Michael Bommarito
+ * @version 1.2
+ */
+class CholeskyDecomposition {
+
+ /**
+ * Decomposition storage
+ * @var array
+ * @access private
+ */
+ private $L = array();
+
+ /**
+ * Matrix row and column dimension
+ * @var int
+ * @access private
+ */
+ private $m;
+
+ /**
+ * Symmetric positive definite flag
+ * @var boolean
+ * @access private
+ */
+ private $isspd = true;
+
+
+ /**
+ * CholeskyDecomposition
+ *
+ * Class constructor - decomposes symmetric positive definite matrix
+ * @param mixed Matrix square symmetric positive definite matrix
+ */
+ public function __construct($A = null) {
+ if ($A instanceof Matrix) {
+ $this->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;
+ }
+ }
+ } else {
+ throw new Exception(JAMAError(ArgumentTypeException));
+ }
+ } // function __construct()
+
+
+ /**
+ * Is the matrix symmetric and positive definite?
+ *
+ * @return boolean
+ */
+ public function isSPD() {
+ return $this->isspd;
+ } // function isSPD()
+
+
+ /**
+ * getL
+ *
+ * Return triangular factor.
+ * @return Matrix Lower triangular matrix
+ */
+ public function getL() {
+ return new Matrix($this->L);
+ } // function getL()
+
+
+ /**
+ * Solve A*X = B
+ *
+ * @param $B Row-equal matrix
+ * @return Matrix L * L' * X = B
+ */
+ public function solve($B = null) {
+ if ($B instanceof Matrix) {
+ 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);
+ } else {
+ throw new Exception(JAMAError(MatrixSPDException));
+ }
+ } else {
+ throw new Exception(JAMAError(MatrixDimensionException));
+ }
+ } else {
+ throw new Exception(JAMAError(ArgumentTypeException));
+ }
+ } // function solve()
+
+} // class CholeskyDecomposition
|