From 35f12a62e6f0b7a5e81e43f7234522ecf5524f62 Mon Sep 17 00:00:00 2001 From: Christoph Ruegg Date: Mon, 28 Dec 2015 08:58:51 +0100 Subject: [PATCH] Linear Algebra: Matrix Rank should use effective epsilon #334 Previously it was using the normalized epsilon at 1.0. See http://discuss.mathdotnet.com/t/wrong-compute-of-the-matrix-rank/120 for discussion. --- .../LinearAlgebra/Complex/DenseVector.cs | 18 ------------------ .../LinearAlgebra/Complex/Factorization/Svd.cs | 2 +- src/Numerics/LinearAlgebra/Complex/Vector.cs | 2 +- .../LinearAlgebra/Complex32/DenseVector.cs | 18 ------------------ .../Complex32/Factorization/Svd.cs | 2 +- src/Numerics/LinearAlgebra/Complex32/Vector.cs | 2 +- .../LinearAlgebra/Double/DenseVector.cs | 18 ------------------ .../LinearAlgebra/Double/Factorization/Svd.cs | 2 +- .../LinearAlgebra/Single/DenseVector.cs | 18 ------------------ .../LinearAlgebra/Single/Factorization/Svd.cs | 2 +- .../Double/Factorization/SvdTests.cs | 14 ++++++++++++++ .../Single/Factorization/SvdTests.cs | 14 ++++++++++++++ 12 files changed, 34 insertions(+), 78 deletions(-) diff --git a/src/Numerics/LinearAlgebra/Complex/DenseVector.cs b/src/Numerics/LinearAlgebra/Complex/DenseVector.cs index c5bdc6bd..d3d39d09 100644 --- a/src/Numerics/LinearAlgebra/Complex/DenseVector.cs +++ b/src/Numerics/LinearAlgebra/Complex/DenseVector.cs @@ -519,24 +519,6 @@ namespace MathNet.Numerics.LinearAlgebra.Complex return index; } - /// - /// Returns the value of the absolute minimum element. - /// - /// The value of the absolute minimum element. - public override Complex AbsoluteMinimum() - { - return _values[AbsoluteMinimumIndex()].Magnitude; - } - - /// - /// Returns the value of the absolute maximum element. - /// - /// The value of the absolute maximum element. - public override Complex AbsoluteMaximum() - { - return _values[AbsoluteMaximumIndex()].Magnitude; - } - /// /// Returns the index of the absolute maximum element. /// diff --git a/src/Numerics/LinearAlgebra/Complex/Factorization/Svd.cs b/src/Numerics/LinearAlgebra/Complex/Factorization/Svd.cs index a2f778c3..45348c9f 100644 --- a/src/Numerics/LinearAlgebra/Complex/Factorization/Svd.cs +++ b/src/Numerics/LinearAlgebra/Complex/Factorization/Svd.cs @@ -71,7 +71,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization { get { - double tolerance = Precision.DoublePrecision*Math.Max(U.RowCount, VT.RowCount); + double tolerance = Precision.EpsilonOf(S.AbsoluteMaximum().Magnitude)*Math.Max(U.RowCount, VT.RowCount); return S.Count(t => t.Magnitude > tolerance); } } diff --git a/src/Numerics/LinearAlgebra/Complex/Vector.cs b/src/Numerics/LinearAlgebra/Complex/Vector.cs index a057e15f..068634cb 100644 --- a/src/Numerics/LinearAlgebra/Complex/Vector.cs +++ b/src/Numerics/LinearAlgebra/Complex/Vector.cs @@ -323,7 +323,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex /// Returns the value of the absolute minimum element. /// /// The value of the absolute minimum element. - public override Complex AbsoluteMinimum() + public sealed override Complex AbsoluteMinimum() { return At(AbsoluteMinimumIndex()).Magnitude; } diff --git a/src/Numerics/LinearAlgebra/Complex32/DenseVector.cs b/src/Numerics/LinearAlgebra/Complex32/DenseVector.cs index d41c99bb..178aa5a5 100644 --- a/src/Numerics/LinearAlgebra/Complex32/DenseVector.cs +++ b/src/Numerics/LinearAlgebra/Complex32/DenseVector.cs @@ -514,24 +514,6 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 return index; } - /// - /// Returns the value of the absolute minimum element. - /// - /// The value of the absolute minimum element. - public override Complex32 AbsoluteMinimum() - { - return _values[AbsoluteMinimumIndex()].Magnitude; - } - - /// - /// Returns the value of the absolute maximum element. - /// - /// The value of the absolute maximum element. - public override Complex32 AbsoluteMaximum() - { - return _values[AbsoluteMaximumIndex()].Magnitude; - } - /// /// Returns the index of the absolute maximum element. /// diff --git a/src/Numerics/LinearAlgebra/Complex32/Factorization/Svd.cs b/src/Numerics/LinearAlgebra/Complex32/Factorization/Svd.cs index c96d2b4f..1ac1f3de 100644 --- a/src/Numerics/LinearAlgebra/Complex32/Factorization/Svd.cs +++ b/src/Numerics/LinearAlgebra/Complex32/Factorization/Svd.cs @@ -66,7 +66,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization { get { - double tolerance = Precision.SinglePrecision*Math.Max(U.RowCount, VT.RowCount); + double tolerance = Precision.EpsilonOf(S.AbsoluteMaximum().Magnitude)*Math.Max(U.RowCount, VT.RowCount); return S.Count(t => t.Magnitude > tolerance); } } diff --git a/src/Numerics/LinearAlgebra/Complex32/Vector.cs b/src/Numerics/LinearAlgebra/Complex32/Vector.cs index 41d93489..392a1d08 100644 --- a/src/Numerics/LinearAlgebra/Complex32/Vector.cs +++ b/src/Numerics/LinearAlgebra/Complex32/Vector.cs @@ -318,7 +318,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 /// Returns the value of the absolute minimum element. /// /// The value of the absolute minimum element. - public override Complex32 AbsoluteMinimum() + public sealed override Complex32 AbsoluteMinimum() { return At(AbsoluteMinimumIndex()).Magnitude; } diff --git a/src/Numerics/LinearAlgebra/Double/DenseVector.cs b/src/Numerics/LinearAlgebra/Double/DenseVector.cs index 7d7dc876..7714c050 100644 --- a/src/Numerics/LinearAlgebra/Double/DenseVector.cs +++ b/src/Numerics/LinearAlgebra/Double/DenseVector.cs @@ -556,24 +556,6 @@ namespace MathNet.Numerics.LinearAlgebra.Double return index; } - /// - /// Returns the value of the absolute minimum element. - /// - /// The value of the absolute minimum element. - public override double AbsoluteMinimum() - { - return Math.Abs(_values[AbsoluteMinimumIndex()]); - } - - /// - /// Returns the value of the absolute maximum element. - /// - /// The value of the absolute maximum element. - public override double AbsoluteMaximum() - { - return Math.Abs(_values[AbsoluteMaximumIndex()]); - } - /// /// Returns the index of the absolute maximum element. /// diff --git a/src/Numerics/LinearAlgebra/Double/Factorization/Svd.cs b/src/Numerics/LinearAlgebra/Double/Factorization/Svd.cs index ba39ffed..baa3db1a 100644 --- a/src/Numerics/LinearAlgebra/Double/Factorization/Svd.cs +++ b/src/Numerics/LinearAlgebra/Double/Factorization/Svd.cs @@ -64,7 +64,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization { get { - double tolerance = Precision.DoublePrecision*Math.Max(U.RowCount, VT.RowCount); + double tolerance = Precision.EpsilonOf(S.Maximum())*Math.Max(U.RowCount, VT.RowCount); return S.Count(t => Math.Abs(t) > tolerance); } } diff --git a/src/Numerics/LinearAlgebra/Single/DenseVector.cs b/src/Numerics/LinearAlgebra/Single/DenseVector.cs index be93a48a..b98cf348 100644 --- a/src/Numerics/LinearAlgebra/Single/DenseVector.cs +++ b/src/Numerics/LinearAlgebra/Single/DenseVector.cs @@ -546,24 +546,6 @@ namespace MathNet.Numerics.LinearAlgebra.Single return index; } - /// - /// Returns the value of the absolute minimum element. - /// - /// The value of the absolute minimum element. - public override float AbsoluteMinimum() - { - return Math.Abs(_values[AbsoluteMinimumIndex()]); - } - - /// - /// Returns the value of the absolute maximum element. - /// - /// The value of the absolute maximum element. - public override float AbsoluteMaximum() - { - return Math.Abs(_values[AbsoluteMaximumIndex()]); - } - /// /// Returns the index of the absolute maximum element. /// diff --git a/src/Numerics/LinearAlgebra/Single/Factorization/Svd.cs b/src/Numerics/LinearAlgebra/Single/Factorization/Svd.cs index d83ba9a1..17f4001d 100644 --- a/src/Numerics/LinearAlgebra/Single/Factorization/Svd.cs +++ b/src/Numerics/LinearAlgebra/Single/Factorization/Svd.cs @@ -64,7 +64,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization { get { - double tolerance = Precision.SinglePrecision*Math.Max(U.RowCount, VT.RowCount); + double tolerance = Precision.EpsilonOf(S.Maximum())*Math.Max(U.RowCount, VT.RowCount); return S.Count(t => Math.Abs(t) > tolerance); } } diff --git a/src/UnitTests/LinearAlgebraTests/Double/Factorization/SvdTests.cs b/src/UnitTests/LinearAlgebraTests/Double/Factorization/SvdTests.cs index b18a3b36..7ac2e126 100644 --- a/src/UnitTests/LinearAlgebraTests/Double/Factorization/SvdTests.cs +++ b/src/UnitTests/LinearAlgebraTests/Double/Factorization/SvdTests.cs @@ -180,6 +180,20 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Double.Factorization Assert.AreEqual(factorSvd.Rank, order - 1); } + [Test] + public void RankAcceptance() + { + // http://discuss.mathdotnet.com/t/wrong-compute-of-the-matrix-rank/120 + Matrix m = DenseMatrix.OfArray(new double[,] { + { 4, 4, 1, 3 }, + { 1,-2, 1, 0 }, + { 4, 0, 2, 2 }, + { 7, 6, 2, 5 } }); + + Assert.That(m.Svd(true).Rank, Is.EqualTo(2)); + Assert.That(m.Svd(false).Rank, Is.EqualTo(2)); + } + /// /// Solve for matrix if vectors are not computed throws InvalidOperationException. /// diff --git a/src/UnitTests/LinearAlgebraTests/Single/Factorization/SvdTests.cs b/src/UnitTests/LinearAlgebraTests/Single/Factorization/SvdTests.cs index 27b613af..31ad9eb7 100644 --- a/src/UnitTests/LinearAlgebraTests/Single/Factorization/SvdTests.cs +++ b/src/UnitTests/LinearAlgebraTests/Single/Factorization/SvdTests.cs @@ -180,6 +180,20 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single.Factorization Assert.AreEqual(factorSvd.Rank, order - 1); } + [Test] + public void RankAcceptance() + { + // http://discuss.mathdotnet.com/t/wrong-compute-of-the-matrix-rank/120 + Matrix m = DenseMatrix.OfArray(new float[,] { + { 4, 4, 1, 3 }, + { 1,-2, 1, 0 }, + { 4, 0, 2, 2 }, + { 7, 6, 2, 5 } }); + + Assert.That(m.Svd(true).Rank, Is.EqualTo(2)); + Assert.That(m.Svd(false).Rank, Is.EqualTo(2)); + } + /// /// Solve for matrix if vectors are not computed throws InvalidOperationException. ///