30 private $LU = array();
66 $this->LU = $A->getArray();
67 $this->m = $A->getRowDimension();
68 $this->n = $A->getColumnDimension();
73 $LUrowi = $LUcolj = array();
79 $LUcolj[
$i] = &$this->LU[
$i][$j];
83 $LUrowi = $this->LU[
$i];
87 for ($k = 0; $k < $kmax; ++$k) {
88 $s += $LUrowi[$k] * $LUcolj[$k];
90 $LUrowi[$j] = $LUcolj[
$i] -=
$s;
95 if (abs($LUcolj[
$i]) > abs($LUcolj[$p])) {
101 $t = $this->LU[$p][$k];
102 $this->LU[$p][$k] = $this->LU[$j][$k];
103 $this->LU[$j][$k] =
$t;
106 $this->piv[$p] = $this->piv[$j];
108 $this->pivsign = $this->pivsign * -1;
111 if (($j < $this->m) && ($this->LU[$j][$j] != 0.0)) {
113 $this->LU[
$i][$j] /= $this->LU[$j][$j];
132 $L[
$i][$j] = $this->LU[
$i][$j];
133 } elseif (
$i == $j) {
153 $U[
$i][$j] = $this->LU[
$i][$j];
190 if ($this->LU[$j][$j] == 0) {
204 if ($this->m == $this->n) {
207 $d *= $this->LU[$j][$j];
225 if ($B->getRowDimension() == $this->m) {
228 $nx = $B->getColumnDimension();
229 $X = $B->getMatrix($this->piv, 0, $nx-1);
233 for ($j = 0; $j < $nx; ++$j) {
234 $X->A[
$i][$j] -= $X->A[$k][$j] * $this->LU[
$i][$k];
239 for ($k = $this->n-1; $k >= 0; --$k) {
240 for ($j = 0; $j < $nx; ++$j) {
241 $X->A[$k][$j] /= $this->LU[$k][$k];
243 for (
$i = 0;
$i < $k; ++
$i) {
244 for ($j = 0; $j < $nx; ++$j) {
245 $X->A[
$i][$j] -= $X->A[$k][$j] * $this->LU[
$i][$k];
An exception for terminatinating execution or to throw for unit testing.
__construct($A)
LU Decomposition constructor.
getPivot()
Return pivot permutation vector.
getL()
Get lower triangular factor.
getDoublePivot()
Alias for getPivot.
const MatrixSingularException
const MatrixSquareException
getU()
Get upper triangular factor.
isNonsingular()
Is the matrix nonsingular?
const ArgumentTypeException
const MatrixDimensionException
for( $i=6;$i< 13;$i++) for($i=1; $i< 13; $i++) $d
const MatrixSingularException