private function EigenvalueDecomposition::tql2 in Loft Data Grids 7.2

Same name and namespace in other branches

1. 6.2 vendor/phpoffice/phpexcel/Classes/PHPExcel/Shared/JAMA/EigenvalueDecomposition.php \EigenvalueDecomposition::tql2()

* Symmetric tridiagonal QL algorithm. * * This is derived from the Algol procedures tql2, by * Bowdler, Martin, Reinsch, and Wilkinson, Handbook for * Auto. Comp., Vol.ii-Linear Algebra, and the corresponding * Fortran subroutine in EISPACK. * * @access private

1 call to EigenvalueDecomposition::tql2()

EigenvalueDecomposition::__construct in vendor/phpoffice/phpexcel/Classes/PHPExcel/Shared/JAMA/EigenvalueDecomposition.php

* Constructor: Check for symmetry, then construct the eigenvalue decomposition * * @access public *

File

vendor/phpoffice/phpexcel/Classes/PHPExcel/Shared/JAMA/EigenvalueDecomposition.php, line 185

Class

EigenvalueDecomposition

@package JAMA

Code

private function tql2() {
  for ($i = 1; $i < $this->n; ++$i) {
    $this->e[$i - 1] = $this->e[$i];
  }
  $this->e[$this->n - 1] = 0.0;
  $f = 0.0;
  $tst1 = 0.0;
  $eps = pow(2.0, -52.0);
  for ($l = 0; $l < $this->n; ++$l) {

    // Find small subdiagonal element
    $tst1 = max($tst1, abs($this->d[$l]) + abs($this->e[$l]));
    $m = $l;
    while ($m < $this->n) {
      if (abs($this->e[$m]) <= $eps * $tst1) {
        break;
      }
      ++$m;
    }

    // If m == l, $this->d[l] is an eigenvalue,
    // otherwise, iterate.
    if ($m > $l) {
      $iter = 0;
      do {

        // Could check iteration count here.
        $iter += 1;

        // Compute implicit shift
        $g = $this->d[$l];
        $p = ($this->d[$l + 1] - $g) / (2.0 * $this->e[$l]);
        $r = hypo($p, 1.0);
        if ($p < 0) {
          $r *= -1;
        }
        $this->d[$l] = $this->e[$l] / ($p + $r);
        $this->d[$l + 1] = $this->e[$l] * ($p + $r);
        $dl1 = $this->d[$l + 1];
        $h = $g - $this->d[$l];
        for ($i = $l + 2; $i < $this->n; ++$i) {
          $this->d[$i] -= $h;
        }
        $f += $h;

        // Implicit QL transformation.
        $p = $this->d[$m];
        $c = 1.0;
        $c2 = $c3 = $c;
        $el1 = $this->e[$l + 1];
        $s = $s2 = 0.0;
        for ($i = $m - 1; $i >= $l; --$i) {
          $c3 = $c2;
          $c2 = $c;
          $s2 = $s;
          $g = $c * $this->e[$i];
          $h = $c * $p;
          $r = hypo($p, $this->e[$i]);
          $this->e[$i + 1] = $s * $r;
          $s = $this->e[$i] / $r;
          $c = $p / $r;
          $p = $c * $this->d[$i] - $s * $g;
          $this->d[$i + 1] = $h + $s * ($c * $g + $s * $this->d[$i]);

          // Accumulate transformation.
          for ($k = 0; $k < $this->n; ++$k) {
            $h = $this->V[$k][$i + 1];
            $this->V[$k][$i + 1] = $s * $this->V[$k][$i] + $c * $h;
            $this->V[$k][$i] = $c * $this->V[$k][$i] - $s * $h;
          }
        }
        $p = -$s * $s2 * $c3 * $el1 * $this->e[$l] / $dl1;
        $this->e[$l] = $s * $p;
        $this->d[$l] = $c * $p;

        // Check for convergence.
      } while (abs($this->e[$l]) > $eps * $tst1);
    }
    $this->d[$l] = $this->d[$l] + $f;
    $this->e[$l] = 0.0;
  }

  // Sort eigenvalues and corresponding vectors.
  for ($i = 0; $i < $this->n - 1; ++$i) {
    $k = $i;
    $p = $this->d[$i];
    for ($j = $i + 1; $j < $this->n; ++$j) {
      if ($this->d[$j] < $p) {
        $k = $j;
        $p = $this->d[$j];
      }
    }
    if ($k != $i) {
      $this->d[$k] = $this->d[$i];
      $this->d[$i] = $p;
      for ($j = 0; $j < $this->n; ++$j) {
        $p = $this->V[$j][$i];
        $this->V[$j][$i] = $this->V[$j][$k];
        $this->V[$j][$k] = $p;
      }
    }
  }
}