87        $this->d = $this->V[$this->n - 1];
 
   90        for (
$i = $this->n - 1; 
$i > 0; --
$i) {
 
   94            $scale += array_sum(array_map(
'abs', $this->d));
 
   96                $this->e[
$i] = $this->d[$i_];
 
   97                $this->d = array_slice($this->V[$i_], 0, $i_);
 
   98                for ($j = 0; $j < 
$i; ++$j) {
 
   99                    $this->V[$j][
$i] = $this->V[
$i][$j] = 0.0;
 
  103                for ($k = 0; $k < 
$i; ++$k) {
 
  104                    $this->d[$k] /= $scale;
 
  105                    $h += $this->d[$k] ** 2;
 
  112                $this->e[
$i] = $scale * $g;
 
  114                $this->d[$i_] = 
$f - $g;
 
  115                for ($j = 0; $j < 
$i; ++$j) {
 
  119                for ($j = 0; $j < 
$i; ++$j) {
 
  121                    $this->V[$j][
$i] = 
$f;
 
  122                    $g = $this->e[$j] + $this->V[$j][$j] * 
$f;
 
  123                    for ($k = $j + 1; $k <= $i_; ++$k) {
 
  124                        $g += $this->V[$k][$j] * $this->d[$k];
 
  125                        $this->e[$k] += $this->V[$k][$j] * 
$f;
 
  130                for ($j = 0; $j < 
$i; ++$j) {
 
  132                    $f += $this->e[$j] * $this->d[$j];
 
  135                for ($j = 0; $j < 
$i; ++$j) {
 
  136                    $this->e[$j] -= $hh * $this->d[$j];
 
  138                for ($j = 0; $j < 
$i; ++$j) {
 
  141                    for ($k = $j; $k <= $i_; ++$k) {
 
  142                        $this->V[$k][$j] -= (
$f * $this->e[$k] + $g * $this->d[$k]);
 
  144                    $this->d[$j] = $this->V[
$i - 1][$j];
 
  145                    $this->V[
$i][$j] = 0.0;
 
  152        for (
$i = 0; 
$i < $this->n - 1; ++
$i) {
 
  153            $this->V[$this->n - 1][
$i] = $this->V[
$i][
$i];
 
  154            $this->V[
$i][
$i] = 1.0;
 
  155            $h = $this->d[
$i + 1];
 
  157                for ($k = 0; $k <= 
$i; ++$k) {
 
  158                    $this->d[$k] = $this->V[$k][
$i + 1] / 
$h;
 
  160                for ($j = 0; $j <= 
$i; ++$j) {
 
  162                    for ($k = 0; $k <= 
$i; ++$k) {
 
  163                        $g += $this->V[$k][
$i + 1] * $this->V[$k][$j];
 
  165                    for ($k = 0; $k <= 
$i; ++$k) {
 
  166                        $this->V[$k][$j] -= $g * $this->d[$k];
 
  170            for ($k = 0; $k <= 
$i; ++$k) {
 
  171                $this->V[$k][
$i + 1] = 0.0;
 
  175        $this->d = $this->V[$this->n - 1];
 
  176        $this->V[$this->n - 1] = array_fill(0, $j, 0.0);
 
  177        $this->V[$this->n - 1][$this->n - 1] = 1.0;
 
  192            $this->e[
$i - 1] = $this->e[
$i];
 
  194        $this->e[$this->n - 1] = 0.0;
 
  197        $eps = 2.0 ** (-52.0);
 
  201            $tst1 = max($tst1, abs($this->d[
$l]) + abs($this->e[
$l]));
 
  203            while ($m < $this->n) {
 
  204                if (abs($this->e[
$m]) <= 
$eps * $tst1) {
 
  218                    $p = ($this->d[
$l + 1] - $g) / (2.0 * $this->e[
$l]);
 
  223                    $this->d[
$l] = $this->e[
$l] / ($p + 
$r);
 
  224                    $this->d[
$l + 1] = $this->e[
$l] * ($p + 
$r);
 
  225                    $dl1 = $this->d[
$l + 1];
 
  226                    $h = $g - $this->d[
$l];
 
  235                    $el1 = $this->e[
$l + 1];
 
  241                        $g = 
$c * $this->e[
$i];
 
  244                        $this->e[
$i + 1] = 
$s * 
$r;
 
  247                        $p = 
$c * $this->d[
$i] - 
$s * $g;
 
  248                        $this->d[
$i + 1] = 
$h + 
$s * (
$c * $g + 
$s * $this->d[
$i]);
 
  251                            $h = $this->V[$k][
$i + 1];
 
  252                            $this->V[$k][
$i + 1] = 
$s * $this->V[$k][
$i] + 
$c * 
$h;
 
  253                            $this->V[$k][
$i] = 
$c * $this->V[$k][
$i] - 
$s * 
$h;
 
  256                    $p = -
$s * $s2 * $c3 * $el1 * $this->e[
$l] / $dl1;
 
  257                    $this->e[
$l] = 
$s * $p;
 
  258                    $this->d[
$l] = 
$c * $p;
 
  260                } 
while (abs($this->e[
$l]) > 
$eps * $tst1);
 
  262            $this->d[
$l] = $this->d[
$l] + 
$f;
 
  267        for (
$i = 0; 
$i < $this->n - 1; ++
$i) {
 
  271                if ($this->d[$j] < $p) {
 
  277                $this->d[$k] = $this->d[
$i];
 
  280                    $p = $this->V[$j][
$i];
 
  281                    $this->V[$j][
$i] = $this->V[$j][$k];
 
  282                    $this->V[$j][$k] = $p;
 
  299        $high = $this->n - 1;
 
  301        for (
$m = $low + 1; 
$m <= $high - 1; ++
$m) {
 
  305                $scale = $scale + abs($this->H[
$i][
$m - 1]);
 
  311                    $this->ort[
$i] = $this->H[
$i][
$m - 1] / $scale;
 
  312                    $h += $this->ort[
$i] * $this->ort[
$i];
 
  315                if ($this->ort[
$m] > 0) {
 
  318                $h -= $this->ort[
$m] * $g;
 
  319                $this->ort[
$m] -= $g;
 
  325                        $f += $this->ort[
$i] * $this->H[
$i][$j];
 
  329                        $this->H[
$i][$j] -= 
$f * $this->ort[
$i];
 
  332                for (
$i = 0; 
$i <= $high; ++
$i) {
 
  334                    for ($j = $high; $j >= 
$m; --$j) {
 
  335                        $f += $this->ort[$j] * $this->H[
$i][$j];
 
  338                    for ($j = 
$m; $j <= $high; ++$j) {
 
  339                        $this->H[
$i][$j] -= 
$f * $this->ort[$j];
 
  342                $this->ort[
$m] = $scale * $this->ort[
$m];
 
  343                $this->H[
$m][
$m - 1] = $scale * $g;
 
  350                $this->V[
$i][$j] = (
$i == $j ? 1.0 : 0.0);
 
  353        for (
$m = $high - 1; 
$m >= $low + 1; --
$m) {
 
  354            if ($this->H[
$m][
$m - 1] != 0.0) {
 
  355                for (
$i = 
$m + 1; 
$i <= $high; ++
$i) {
 
  356                    $this->ort[
$i] = $this->H[
$i][
$m - 1];
 
  358                for ($j = 
$m; $j <= $high; ++$j) {
 
  361                        $g += $this->ort[
$i] * $this->V[
$i][$j];
 
  364                    $g = ($g / $this->ort[
$m]) / $this->H[
$m][
$m - 1];
 
  366                        $this->V[
$i][$j] += $g * $this->ort[
$i];
 
  381    private function cdiv($xr, $xi, $yr, $yi): void
 
  383        if (abs($yr) > abs($yi)) {
 
  386            $this->cdivr = ($xr + 
$r * $xi) / 
$d;
 
  387            $this->cdivi = ($xi - 
$r * $xr) / 
$d;
 
  391            $this->cdivr = (
$r * $xr + $xi) / 
$d;
 
  392            $this->cdivi = (
$r * $xi - $xr) / 
$d;
 
  411        $eps = 2.0 ** (-52.0);
 
  413        $p = $q = 
$r = 
$s = $z = 0;
 
  417        for (
$i = 0; 
$i < $nn; ++
$i) {
 
  418            if ((
$i < $low) || (
$i > $high)) {
 
  419                $this->d[
$i] = $this->H[
$i][
$i];
 
  422            for ($j = max(
$i - 1, 0); $j < $nn; ++$j) {
 
  423                $norm = $norm + abs($this->H[
$i][$j]);
 
  433                $s = abs($this->H[
$l - 1][
$l - 1]) + abs($this->H[
$l][
$l]);
 
  445                $this->H[
$n][
$n] = $this->H[
$n][
$n] + $exshift;
 
  446                $this->d[
$n] = $this->H[
$n][
$n];
 
  451            } elseif (
$l == 
$n - 1) {
 
  452                $w = $this->H[
$n][
$n - 1] * $this->H[
$n - 1][
$n];
 
  453                $p = ($this->H[
$n - 1][
$n - 1] - $this->H[
$n][
$n]) / 2.0;
 
  456                $this->H[
$n][
$n] = $this->H[
$n][
$n] + $exshift;
 
  457                $this->H[
$n - 1][
$n - 1] = $this->H[
$n - 1][
$n - 1] + $exshift;
 
  466                    $this->d[
$n - 1] = 
$x + $z;
 
  467                    $this->d[
$n] = $this->d[
$n - 1];
 
  469                        $this->d[
$n] = 
$x - 
$w / $z;
 
  471                    $this->e[
$n - 1] = 0.0;
 
  473                    $x = $this->H[
$n][
$n - 1];
 
  474                    $s = abs(
$x) + abs($z);
 
  477                    $r = sqrt($p * $p + $q * $q);
 
  481                    for ($j = 
$n - 1; $j < $nn; ++$j) {
 
  482                        $z = $this->H[
$n - 1][$j];
 
  483                        $this->H[
$n - 1][$j] = $q * $z + $p * $this->H[
$n][$j];
 
  484                        $this->H[
$n][$j] = $q * $this->H[
$n][$j] - $p * $z;
 
  488                        $z = $this->H[
$i][
$n - 1];
 
  489                        $this->H[
$i][
$n - 1] = $q * $z + $p * $this->H[
$i][
$n];
 
  490                        $this->H[
$i][
$n] = $q * $this->H[
$i][
$n] - $p * $z;
 
  493                    for (
$i = $low; 
$i <= $high; ++
$i) {
 
  494                        $z = $this->V[
$i][
$n - 1];
 
  495                        $this->V[
$i][
$n - 1] = $q * $z + $p * $this->V[
$i][
$n];
 
  496                        $this->V[
$i][
$n] = $q * $this->V[
$i][
$n] - $p * $z;
 
  500                    $this->d[
$n - 1] = 
$x + $p;
 
  501                    $this->d[
$n] = 
$x + $p;
 
  502                    $this->e[
$n - 1] = $z;
 
  514                    $y = $this->H[
$n - 1][
$n - 1];
 
  515                    $w = $this->H[
$n][
$n - 1] * $this->H[
$n - 1][
$n];
 
  523                    $s = abs($this->H[
$n][
$n - 1]) + abs($this->H[
$n - 1][
$n - 2]);
 
  549                    $z = $this->H[
$m][
$m];
 
  552                    $p = (
$r * 
$s - 
$w) / $this->H[
$m + 1][
$m] + $this->H[
$m][
$m + 1];
 
  553                    $q = $this->H[
$m + 1][
$m + 1] - $z - 
$r - 
$s;
 
  554                    $r = $this->H[
$m + 2][
$m + 1];
 
  555                    $s = abs($p) + abs($q) + abs(
$r);
 
  563                        abs($this->H[
$m][
$m - 1]) * (abs($q) + abs(
$r)) <
 
  564                        $eps * (abs($p) * (abs($this->H[
$m - 1][
$m - 1]) + abs($z) + abs($this->H[
$m + 1][
$m + 1])))
 
  571                    $this->H[
$i][
$i - 2] = 0.0;
 
  573                        $this->H[
$i][
$i - 3] = 0.0;
 
  577                for ($k = 
$m; $k <= 
$n - 1; ++$k) {
 
  578                    $notlast = ($k != 
$n - 1);
 
  580                        $p = $this->H[$k][$k - 1];
 
  581                        $q = $this->H[$k + 1][$k - 1];
 
  582                        $r = ($notlast ? $this->H[$k + 2][$k - 1] : 0.0);
 
  583                        $x = abs($p) + abs($q) + abs(
$r);
 
  593                    $s = sqrt($p * $p + $q * $q + 
$r * 
$r);
 
  599                            $this->H[$k][$k - 1] = -
$s * 
$x;
 
  600                        } elseif (
$l != 
$m) {
 
  601                            $this->H[$k][$k - 1] = -$this->H[$k][$k - 1];
 
  610                        for ($j = $k; $j < $nn; ++$j) {
 
  611                            $p = $this->H[$k][$j] + $q * $this->H[$k + 1][$j];
 
  613                                $p = $p + 
$r * $this->H[$k + 2][$j];
 
  614                                $this->H[$k + 2][$j] = $this->H[$k + 2][$j] - $p * $z;
 
  616                            $this->H[$k][$j] = $this->H[$k][$j] - $p * 
$x;
 
  617                            $this->H[$k + 1][$j] = $this->H[$k + 1][$j] - $p * 
$y;
 
  620                        $iMax = min(
$n, $k + 3);
 
  621                        for (
$i = 0; 
$i <= $iMax; ++
$i) {
 
  622                            $p = 
$x * $this->H[
$i][$k] + 
$y * $this->H[
$i][$k + 1];
 
  624                                $p = $p + $z * $this->H[
$i][$k + 2];
 
  625                                $this->H[
$i][$k + 2] = $this->H[
$i][$k + 2] - $p * 
$r;
 
  627                            $this->H[
$i][$k] = $this->H[
$i][$k] - $p;
 
  628                            $this->H[
$i][$k + 1] = $this->H[
$i][$k + 1] - $p * $q;
 
  631                        for (
$i = $low; 
$i <= $high; ++
$i) {
 
  632                            $p = 
$x * $this->V[
$i][$k] + 
$y * $this->V[
$i][$k + 1];
 
  634                                $p = $p + $z * $this->V[
$i][$k + 2];
 
  635                                $this->V[
$i][$k + 2] = $this->V[
$i][$k + 2] - $p * 
$r;
 
  637                            $this->V[
$i][$k] = $this->V[
$i][$k] - $p;
 
  638                            $this->V[
$i][$k + 1] = $this->V[
$i][$k + 1] - $p * $q;
 
  650        for (
$n = $nn - 1; 
$n >= 0; --
$n) {
 
  656                $this->H[
$n][
$n] = 1.0;
 
  658                    $w = $this->H[
$i][
$i] - $p;
 
  660                    for ($j = 
$l; $j <= 
$n; ++$j) {
 
  661                        $r = 
$r + $this->H[
$i][$j] * $this->H[$j][
$n];
 
  663                    if ($this->e[
$i] < 0.0) {
 
  668                        if ($this->e[
$i] == 0.0) {
 
  676                            $x = $this->H[
$i][
$i + 1];
 
  677                            $y = $this->H[
$i + 1][
$i];
 
  678                            $q = ($this->d[
$i] - $p) * ($this->d[
$i] - $p) + $this->e[
$i] * $this->e[
$i];
 
  681                            if (abs(
$x) > abs($z)) {
 
  688                        $t = abs($this->H[
$i][
$n]);
 
  690                            for ($j = 
$i; $j <= 
$n; ++$j) {
 
  691                                $this->H[$j][
$n] = $this->H[$j][
$n] / 
$t;
 
  700                if (abs($this->H[
$n][
$n - 1]) > abs($this->H[
$n - 1][
$n])) {
 
  701                    $this->H[
$n - 1][
$n - 1] = $q / $this->H[
$n][
$n - 1];
 
  702                    $this->H[
$n - 1][
$n] = -($this->H[
$n][
$n] - $p) / $this->H[
$n][
$n - 1];
 
  704                    $this->
cdiv(0.0, -$this->H[
$n - 1][
$n], $this->H[
$n - 1][
$n - 1] - $p, $q);
 
  708                $this->H[
$n][
$n - 1] = 0.0;
 
  709                $this->H[
$n][
$n] = 1.0;
 
  714                    for ($j = 
$l; $j <= 
$n; ++$j) {
 
  715                        $ra = $ra + $this->H[
$i][$j] * $this->H[$j][
$n - 1];
 
  716                        $sa = $sa + $this->H[
$i][$j] * $this->H[$j][
$n];
 
  718                    $w = $this->H[
$i][
$i] - $p;
 
  719                    if ($this->e[
$i] < 0.0) {
 
  725                        if ($this->e[
$i] == 0) {
 
  726                            $this->
cdiv(-$ra, -$sa, 
$w, $q);
 
  731                            $x = $this->H[
$i][
$i + 1];
 
  732                            $y = $this->H[
$i + 1][
$i];
 
  733                            $vr = ($this->d[
$i] - $p) * ($this->d[
$i] - $p) + $this->e[
$i] * $this->e[
$i] - $q * $q;
 
  734                            $vi = ($this->d[
$i] - $p) * 2.0 * $q;
 
  735                            if ($vr == 0.0 & $vi == 0.0) {
 
  736                                $vr = 
$eps * $norm * (abs(
$w) + abs($q) + abs(
$x) + abs(
$y) + abs($z));
 
  738                            $this->
cdiv(
$x * 
$r - $z * $ra + $q * $sa, 
$x * 
$s - $z * $sa - $q * $ra, $vr, $vi);
 
  741                            if (abs(
$x) > (abs($z) + abs($q))) {
 
  742                                $this->H[
$i + 1][
$n - 1] = (-$ra - 
$w * $this->H[
$i][
$n - 1] + $q * $this->H[
$i][
$n]) / 
$x;
 
  743                                $this->H[
$i + 1][
$n] = (-$sa - 
$w * $this->H[
$i][
$n] - $q * $this->H[
$i][
$n - 1]) / 
$x;
 
  751                        $t = max(abs($this->H[
$i][
$n - 1]), abs($this->H[
$i][
$n]));
 
  753                            for ($j = 
$i; $j <= 
$n; ++$j) {
 
  754                                $this->H[$j][
$n - 1] = $this->H[$j][
$n - 1] / 
$t;
 
  755                                $this->H[$j][
$n] = $this->H[$j][
$n] / 
$t;
 
  764        for (
$i = 0; 
$i < $nn; ++
$i) {
 
  765            if ($i < $low | $i > $high) {
 
  766                for ($j = 
$i; $j < $nn; ++$j) {
 
  767                    $this->V[
$i][$j] = $this->H[
$i][$j];
 
  773        for ($j = $nn - 1; $j >= $low; --$j) {
 
  774            for (
$i = $low; 
$i <= $high; ++
$i) {
 
  776                $kMax = min($j, $high);
 
  777                for ($k = $low; $k <= $kMax; ++$k) {
 
  778                    $z = $z + $this->V[
$i][$k] * $this->H[$k][$j];
 
  780                $this->V[
$i][$j] = $z;
 
  798        for ($j = 0; ($j < 
$this->n) & $issymmetric; ++$j) {
 
  800                $issymmetric = ($this->A[
$i][$j] == $this->A[$j][
$i]);
 
  827        return new Matrix($this->V, $this->n, $this->n);
 
  859            $D[
$i] = array_fill(0, $this->n, 0.0);
 
  860            $D[
$i][
$i] = $this->d[
$i];
 
  861            if ($this->e[
$i] == 0) {
 
  864            $o = ($this->e[
$i] > 0) ? 
$i + 1 : 
$i - 1;
 
  865            $D[
$i][$o] = $this->e[
$i];
 
hypo($a, $b)
Pythagorean Theorem:.
An exception for terminatinating execution or to throw for unit testing.
Class to obtain eigenvalues and eigenvectors of a real matrix.
tred2()
Symmetric Householder reduction to tridiagonal form.
getImagEigenvalues()
Return the imaginary parts of the eigenvalues.
orthes()
Nonsymmetric reduction to Hessenberg form.
cdiv($xr, $xi, $yr, $yi)
Performs complex division.
hqr2()
Nonsymmetric reduction from Hessenberg to real Schur form.
getV()
Return the eigenvector matrix.
getRealEigenvalues()
Return the real parts of the eigenvalues.
tql2()
Symmetric tridiagonal QL algorithm.
getD()
Return the block diagonal eigenvalue matrix.
__construct(Matrix $Arg)
Constructor: Check for symmetry, then construct the eigenvalue decomposition.
getColumnDimension()
getColumnDimension.
Class for the creating "special" Matrices.