commons-dev mailing list archives

Site index · List index
Message view « Date » · « Thread »
Top « Date » · « Thread »
From Phil Steitz <phil.ste...@gmail.com>
Subject Re: svn commit: r832577 - in /commons/proper/math/trunk/src: main/java/org/apache/commons/math/linear/EigenDecompositionImpl.java site/xdoc/changes.xml test/java/org/apache/commons/math/linear/EigenDecompositionImplTest.java
Date Tue, 03 Nov 2009 23:56:11 GMT
Kudos indeed. Nice work!

Phil

luc@apache.org wrote:
> Author: luc
> Date: Tue Nov  3 22:04:08 2009
> New Revision: 832577
> 
> URL: http://svn.apache.org/viewvc?rev=832577&view=rev
> Log:
> fixed an ArrayIndexOutOfBoundsException
> Kudos to Dimitri who debugged this mess of fortran/java array indices translation
> JIRA: MATH-308
> 
> Modified:
>     commons/proper/math/trunk/src/main/java/org/apache/commons/math/linear/EigenDecompositionImpl.java
>     commons/proper/math/trunk/src/site/xdoc/changes.xml
>     commons/proper/math/trunk/src/test/java/org/apache/commons/math/linear/EigenDecompositionImplTest.java
> 
> Modified: commons/proper/math/trunk/src/main/java/org/apache/commons/math/linear/EigenDecompositionImpl.java
> URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math/linear/EigenDecompositionImpl.java?rev=832577&r1=832576&r2=832577&view=diff
> ==============================================================================
> --- commons/proper/math/trunk/src/main/java/org/apache/commons/math/linear/EigenDecompositionImpl.java
(original)
> +++ commons/proper/math/trunk/src/main/java/org/apache/commons/math/linear/EigenDecompositionImpl.java
Tue Nov  3 22:04:08 2009
> @@ -595,8 +595,12 @@
>          }
>  
>          final double dCurrent = main[m - 1];
> -        work[lowerStart + m - 1] = dCurrent - eCurrent;
> -        work[upperStart + m - 1] = dCurrent + eCurrent;
> +        final double lower = dCurrent - eCurrent;
> +        work[lowerStart + m - 1] = lower;
> +        lowerSpectra = Math.min(lowerSpectra, lower);
> +        final double upper = dCurrent + eCurrent;
> +        work[upperStart + m - 1] = upper;
> +        upperSpectra = Math.max(upperSpectra, upper);
>          minPivot = MathUtils.SAFE_MIN * Math.max(1.0, eMax * eMax);
>  
>      }
> @@ -647,7 +651,7 @@
>  
>                  tau = (range[1] - range[0]) * MathUtils.EPSILON * n + 2 * minPivot;
>  
> -                // decompose T&lambda;I as LDL<sup>T</sup>
> +                // decompose T-&lambda;I as LDL<sup>T</sup>
>                  ldlTDecomposition(lambda, begin, n);
>  
>                  // apply general dqd/dqds method
> @@ -899,8 +903,8 @@
>                      diagMax    = work[4 * i0];
>                      offDiagMin = work[4 * i0 + 2];
>                      double previousEMin = work[4 * i0 + 3];
> -                    for (int i = 4 * i0; i < 4 * n0 - 11; i += 4) {
> -                        if ((work[i + 3] <= TOLERANCE_2 * work[i]) &&
> +                    for (int i = 4 * i0; i < 4 * n0 - 16; i += 4) {
> +                        if ((work[i + 3] <= TOLERANCE_2 * work[i]) ||
>                              (work[i + 2] <= TOLERANCE_2 * sigma)) {
>                              // insert a split
>                              work[i + 2]  = -sigma;
> @@ -1537,7 +1541,7 @@
>                  double a2 = (work[np - 8] / b2) * (1 + work[np - 4] / b1);
>  
>                  // approximate contribution to norm squared from i < nn-2.
> -                if (end - start > 2) {
> +                if (end - start > 3) {
>                      b2 = work[nn - 13] / work[nn - 15];
>                      a2 = a2 + b2;
>                      for (int i4 = nn - 17; i4 >= 4 * start + 2 + pingPong; i4 -=
4) {
> 
> Modified: commons/proper/math/trunk/src/site/xdoc/changes.xml
> URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/site/xdoc/changes.xml?rev=832577&r1=832576&r2=832577&view=diff
> ==============================================================================
> --- commons/proper/math/trunk/src/site/xdoc/changes.xml (original)
> +++ commons/proper/math/trunk/src/site/xdoc/changes.xml Tue Nov  3 22:04:08 2009
> @@ -39,6 +39,10 @@
>    </properties>
>    <body>
>      <release version="2.1" date="TBD" description="TBD">
> +      <action dev="luc" type="fix" issue="MATH-308" due-to="Dimitri Pourbaix">
> +        Fixed an ArrayIndexOutOfBoundsException in eigen decomposition. Kudos to Dimitri
> +        for debugging this, it was really difficult.
> +      </action>
>        <action dev="psteitz" type="fix" issue="MATH-309" due-to="Mikkel Meyer Andersen">
>          Fixed parameter test in RandomDataImpl#nextExponential. The method now throws
>          IllegalArgumentException for mean = 0.
> 
> Modified: commons/proper/math/trunk/src/test/java/org/apache/commons/math/linear/EigenDecompositionImplTest.java
> URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/test/java/org/apache/commons/math/linear/EigenDecompositionImplTest.java?rev=832577&r1=832576&r2=832577&view=diff
> ==============================================================================
> --- commons/proper/math/trunk/src/test/java/org/apache/commons/math/linear/EigenDecompositionImplTest.java
(original)
> +++ commons/proper/math/trunk/src/test/java/org/apache/commons/math/linear/EigenDecompositionImplTest.java
Tue Nov  3 22:04:08 2009
> @@ -108,6 +108,40 @@
>          assertEquals(0.1, ed.getRealEigenvalue(3), 1.0e-15);
>      }
>  
> +    // the following test triggered an ArrayIndexOutOfBoundsException in commons-math
2.0
> +    public void testMath308() {
> +
> +        double[] mainTridiagonal = {
> +            22.330154644539597, 46.65485522478641, 17.393672330044705, 54.46687435351116,
80.17800767709437
> +        };
> +        double[] secondaryTridiagonal = {
> +            13.04450406501361, -5.977590941539671, 2.9040909856707517, 7.1570352792841225
> +        };
> +
> +        // the reference values have been computed using routine DSTEMR
> +        // from the fortran library LAPACK version 3.2.1
> +        double[] refEigenValues = {
> +            82.044413207204002, 53.456697699894512, 52.536278520113882, 18.847969733754262,
14.138204224043099
> +        };
> +        RealVector[] refEigenVectors = {
> +            new ArrayRealVector(new double[] { -0.000462690386766, -0.002118073109055,
 0.011530080757413,  0.252322434584915,  0.967572088232592 }),
> +            new ArrayRealVector(new double[] {  0.314647769490148,  0.750806415553905,
-0.167700312025760, -0.537092972407375,  0.143854968127780 }),
> +            new ArrayRealVector(new double[] {  0.222368839324646,  0.514921891363332,
-0.021377019336614,  0.801196801016305, -0.207446991247740 }),
> +            new ArrayRealVector(new double[] {  0.713933751051495, -0.190582113553930,
 0.671410443368332, -0.056056055955050,  0.006541576993581 }),
> +            new ArrayRealVector(new double[] {  0.584677060845929, -0.367177264979103,
-0.721453187784497,  0.052971054621812, -0.005740715188257 })
> +        };
> +
> +        EigenDecomposition decomposition =
> +            new EigenDecompositionImpl(mainTridiagonal, secondaryTridiagonal, MathUtils.SAFE_MIN);
> +
> +        double[] eigenValues = decomposition.getRealEigenvalues();
> +        for (int i = 0; i < refEigenValues.length; ++i) {
> +            assertEquals(refEigenValues[i], eigenValues[i], 1.0e-5);
> +            assertEquals(0, refEigenVectors[i].subtract(decomposition.getEigenvector(i)).getNorm(),
2.0e-7);
> +        }
> +
> +    }
> +
>      /** test a matrix already in tridiagonal form. */
>      public void testTridiagonal() {
>          Random r = new Random(4366663527842l);
> 
> 


---------------------------------------------------------------------
To unsubscribe, e-mail: dev-unsubscribe@commons.apache.org
For additional commands, e-mail: dev-help@commons.apache.org


Mime
View raw message