Changeset 10077

Timestamp:

10/21/13 23:43:43 (11 years ago)

Author:

ascheibe

Message:

#1886 fixed a bug in the point calculation from distance matrices

File:

• branches/HeuristicLab.Analysis.AlgorithmBehavior/HeuristicLab.Analysis.AlgorithmBehavior.Analyzers/3.3/DistanceMatrixToPoints.cs (modified) (3 diffs)

branches/HeuristicLab.Analysis.AlgorithmBehavior/HeuristicLab.Analysis.AlgorithmBehavior.Analyzers/3.3/DistanceMatrixToPoints.cs

@@ -21 +21 @@
 
 using System;
+using System.Linq;
+using HeuristicLab.Common;
 
 namespace HeuristicLab.Analysis.AlgorithmBehavior.Analyzers {
@@ -33 +35 @@
     public static double[][] ConvertDistanceMatrixToPoints(double[][] dm) {
       double[][] points = new double[dm.Length][];
-      for (int i = 0; i < dm.Length; i++) {
-        points[i] = new double[dm.Length];
-      }
       double[,] m = new double[dm.Length, dm.Length];
       double[] q = new double[dm.Length]; //eigenvalues
-      double[,] v = new double[dm.Length, dm.Length];  //eigenvectors
+      double[,] v = new double[dm.Length, dm.Length]; //eigenvectors
 
       for (int i = 0; i < dm.Length; i++) {
@@ -46 +45 @@
       }
 
-      bool res = alglib.smatrixevd(m, dm.Length, 1, true, out q, out v);
+      //QR decomposition to get the upper part for smatrixevd
+      double[] tau;
+      double[,] r;
+      alglib.rmatrixqr(ref m, dm.Length, dm.Length, out tau);
+      alglib.rmatrixqrunpackr(m, dm.Length, dm.Length, out r);
+
+      bool res = alglib.smatrixevd(r, dm.Length, 1, true, out q, out v);
       if (!res) throw new Exception("Eigenvalue computation did not converge!");
 
+      int zeroCnt = q.Count(x => x.IsAlmost(0));
       for (int i = 0; i < dm.Length; i++) {
-        for (int j = 0; j < dm.Length; j++) {
-          points[j][i] = Math.Sqrt(q[i]) * v[j, i];
+        points[i] = new double[dm.Length - zeroCnt];
+      }
+
+      int pi = 0;
+      for (int i = 0; i < dm.Length; i++) {
+        if (q[i].IsAlmost(0)) {
+          for (int j = 0; j < dm.Length; j++) {
+            points[j][pi] = Math.Sqrt(q[i]) * v[j, i];
+          }
+          pi++;
         }
       }