Vector: added subtraction and corrected addition interface in ILinearAlgebra

Signed-off-by: Marcus Cuda <marcus@cuda.net>

src/Numerics/Algorithms/LinearAlgebra/ILinearAlgebra.cs

@@ -30,10 +30,20 @@ namespace MathNet.Numerics.Algorithms.LinearAlgebra
     public interface ILinearAlgebra
     {
         /// <summary>
-        /// Adds the two arrays together: <c>a += c</c>.
+        /// Adds a scaled vector to another: <c>y += alpha*x</c>.
         /// </summary>
-        /// <param name="a">One of the arrays to add.</param>
-        /// <param name="b">The other array to add.</param>
-        void AddArrays(double[] a, double[] b);
+        /// <param name="y">The vector to update.</param>
+        /// <param name="alpha">The value to scale <param name="x"/> by.</param>
+        /// <param name="x">The vector to add to <paramref name="y"/>.</param>
+        /// <remarks>This is equivalent to the AXPY BLAS routine.</remarks>
+        void AddVectorToScaledVector(double[] y, double alpha, double[] x);
+
+        /// <summary>
+        /// Scales an array. Can be used to scale a vector and a matrix.
+        /// </summary>
+        /// <param name="alpha">The scalar.</param>
+        /// <param name="x">The values to scale.</param>
+        /// <remarks>This is equivalent to the SCAL BLAS routine.</remarks>
+        void ScaleArray(double alpha, double[] x);
     }
 }

src/Numerics/Algorithms/LinearAlgebra/ManagedLinearAlgebra.cs

@@ -31,37 +31,52 @@ namespace MathNet.Numerics.Algorithms.LinearAlgebra
     /// </summary>
     internal class ManagedLinearAlgebra : ILinearAlgebra
     {
-        #region ILinearAlgebra Members
-
         /// <summary>
-        /// Adds the two arrays together: <c>a += c</c>.
+        /// Adds a scaled vector to another: <c>y += alpha*x</c>.
         /// </summary>
-        /// <param name="a">
-        /// One of the arrays to add.
-        /// </param>
-        /// <param name="b">
-        /// The other array to add.
-        /// </param>
-        public void AddArrays(double[] a, double[] b)
+        /// <param name="y">The vector to update.</param>
+        /// <param name="alpha">The value to scale <param name="x"/> by.</param>
+        /// <param name="x">The vector to add to <paramref name="y"/>.</param>
+        /// <remarks>This equivalent to the AXPY BLAS routine.</remarks>
+        public void AddVectorToScaledVector(double[] y, double alpha, double[] x)
         {
-            if (a == null)
+            if (y == null)
             {
-                throw new ArgumentNullException("a");
+                throw new ArgumentNullException("y");
             }
 
-            if (b == null)
+            if (x == null)
             {
-                throw new ArgumentNullException("b");
+                throw new ArgumentNullException("x");
             }
 
-            if (a.Length != b.Length)
+            if (y.Length != x.Length)
             {
                 throw new ArgumentException(Properties.Resources.ArgumentVectorsSameLength);
             }
 
-            Parallel.For(0, a.Length, i => a[i] += b[i]);
+            if (alpha.AlmostZero())
+            {
+                return;
+            }
+
+            if (alpha.AlmostEqual(1.0))
+            {
+                Parallel.For(0, y.Length, i => y[i] += x[i]);
+            }
+            else
+            {
+                Parallel.For(0, y.Length, i => y[i] += alpha * x[i]);
+            }
         }
 
-        #endregion
+        public void ScaleArray(double alpha, double[] x)
+        {
+            if (alpha.AlmostEqual(1.0))
+            {
+                return;
+            }
+            Parallel.For(0, x.Length, i => x[i] = alpha * x[i]);
+        }
     }
 }

src/Numerics/LinearAlgebra/Double/DenseVector.cs

@@ -265,7 +265,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double
             }
             else
             {
-                _linearAlgebra.AddArrays(Data, denseVector.Data);
+                _linearAlgebra.AddVectorToScaledVector(Data, 1.0, denseVector.Data);
             }
         }
 
@@ -354,5 +354,169 @@ namespace MathNet.Numerics.LinearAlgebra.Double
             ret.Add(rightSide);
             return ret;
         }
+
+        /// <summary>
+        /// Subtracts a scalar from each element of the vector.
+        /// </summary>
+        /// <param name="scalar">The scalar to subtract.</param>
+        public override void Subtract(double scalar)
+        {
+            if (scalar.AlmostZero())
+            {
+                return;
+            }
+
+            Parallel.For(0, Count, i => Data[i] -= scalar);
+        }
+
+        /// <summary>
+        /// Subtracts a scalar from each element of the vector and stores the result in the result vector.
+        /// </summary>
+        /// <param name="scalar">The scalar to subtract.</param>
+        /// <param name="result">The vector to store the result of the subtraction.</param>
+        /// <exception cref="ArgumentNullException">If the result vector is <see langword="null"/>.</exception>
+        /// <exception cref="ArgumentException">If this vector and <paramref name="result"/> are not the same size.</exception>
+        public override void Subtract(double scalar, Vector result)
+        {
+            if (result == null)
+            {
+                throw new ArgumentNullException("result");
+            }
+
+            if (Count != result.Count)
+            {
+                throw new ArgumentException("result", Resources.ArgumentVectorsSameLength);
+            }
+
+            CopyTo(result);
+            result.Subtract(scalar);
+        }
+
+        /// <summary>
+        /// Subtracts another vector from this vector.
+        /// </summary>
+        /// <param name="other">The vector to subtract from this one.</param>
+        /// <exception cref="ArgumentNullException">If the other vector is <see langword="null"/>.</exception>
+        /// <exception cref="ArgumentException">If this vector and <paramref name="other"/> are not the same size.</exception>
+         public override void Subtract(Vector other)
+        {
+            if (other == null)
+            {
+                throw new ArgumentNullException("other");
+            }
+
+            if (Count != other.Count)
+            {
+                throw new ArgumentException("other", Resources.ArgumentVectorsSameLength);
+            }
+
+            var denseVector = other as DenseVector;
+
+            if (denseVector == null)
+            {
+                base.Subtract(other);
+            }
+            else
+            {
+                _linearAlgebra.AddVectorToScaledVector(Data, -1.0, denseVector.Data);
+            }
+        }
+
+         /// <summary>
+         /// Subtracts another vector to this vector and stores the result into the result vector.
+         /// </summary>
+         /// <param name="other">The vector to subtract from this one.</param>
+         /// <param name="result">The vector to store the result of the subtraction.</param>
+         /// <exception cref="ArgumentNullException">If the other vector is <see langword="null"/>.</exception>
+         /// <exception cref="ArgumentNullException">If the result vector is <see langword="null"/>.</exception>
+         /// <exception cref="ArgumentException">If this vector and <paramref name="other"/> are not the same size.</exception>
+         /// <exception cref="ArgumentException">If this vector and <paramref name="result"/> are not the same size.</exception>
+        public override void Subtract(Vector other, Vector result)
+        {
+            if (result == null)
+            {
+                throw new ArgumentNullException("result");
+            }
+
+            if (Count != other.Count)
+            {
+                throw new ArgumentException("other", Resources.ArgumentVectorsSameLength);
+            }
+
+            if (Count != result.Count)
+            {
+                throw new ArgumentException("result", Resources.ArgumentVectorsSameLength);
+            }
+
+            if (ReferenceEquals(this, result) || ReferenceEquals(other, result))
+            {
+                var tmp = result.CreateVector(result.Count);
+                Subtract(other, tmp);
+                tmp.CopyTo(result);
+            }
+            else
+            {
+                CopyTo(result);
+                result.Subtract(other);
+            }
+        }
+
+        /// <summary>
+        /// Returns a <strong>Vector</strong> containing the negated values of rightSide. 
+        /// </summary>
+        /// <param name="rightSide">The vector to get the values from.</param>
+        /// <returns>A vector containing the negated values as <paramref name="rightSide"/>.</returns>
+        /// <exception cref="ArgumentNullException">If <paramref name="rightSide"/> is <see langword="null" />.</exception>
+        public static Vector operator -(DenseVector rightSide)
+        {
+            if (rightSide == null)
+            {
+                throw new ArgumentNullException("rightSide");
+            }
+
+            return rightSide.Negate();
+        }
+
+        /// <summary>
+        /// Subtracts two <strong>Vectors</strong> and returns the results.
+        /// </summary>
+        /// <param name="leftSide">The vector to subtract from.</param>
+        /// <param name="rightSide">The vector to subtract.</param>
+        /// <returns>The result of the subtraction.</returns>
+        /// <exception cref="ArgumentException">If <paramref name="leftSide"/> and <paramref name="rightSide"/> are not the same size.</exception>
+        /// <exception cref="ArgumentNullException">If <paramref name="leftSide"/> or <paramref name="rightSide"/> is <see langword="null" />.</exception>
+        public static Vector operator -(DenseVector leftSide, DenseVector rightSide)
+        {
+            if (rightSide == null)
+            {
+                throw new ArgumentNullException("rightSide");
+            }
+
+            if (leftSide == null)
+            {
+                throw new ArgumentNullException("leftSide");
+            }
+
+            if (leftSide.Count != rightSide.Count)
+            {
+                throw new ArgumentException("rightSide", Resources.ArgumentVectorsSameLength);
+            }
+
+            var ret = leftSide.Clone();
+            ret.Subtract(rightSide);
+            return ret;
+        }
+
+        /// <summary>
+        /// Returns a negated vector.
+        /// </summary>
+        /// <returns>The negated vector.</returns>
+        /// <remarks>Added as an alternative to the unary negation operator.</remarks>
+        public virtual Vector Negate()
+        {
+            var result = new DenseVector(Count);
+            Parallel.For(0, Count, i => result[i] = -Data[i]);
+            return result;
+        }
     }
 }

src/Numerics/LinearAlgebra/Double/Vector.cs

@@ -329,9 +329,9 @@ namespace MathNet.Numerics.LinearAlgebra.Double
         }
 
         /// <summary>
-        /// Returns a clone of this vector.
+        /// Returns this vector.
         /// </summary>
-        /// <returns>A clone of this vector.</returns>
+        /// <returns>This vector.</returns>
         /// <remarks>Added as an alternative to the unary addition operator.</remarks>
         public virtual Vector Plus()
         {
@@ -398,7 +398,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double
         /// </summary>
         /// <remarks>This method is included for completeness.</remarks>
         /// <param name="rightSide">The vector to get the values from.</param>
-        /// <returns>A vector containing a the same values as <paramref name="rightSide"/>.</returns>
+        /// <returns>A vector containing the same values as <paramref name="rightSide"/>.</returns>
         /// <exception cref="ArgumentNullException">If <paramref name="rightSide"/> is <see langword="null" />.</exception>
         public static Vector operator +(Vector rightSide)
         {
@@ -440,6 +440,157 @@ namespace MathNet.Numerics.LinearAlgebra.Double
             return ret;
         }
 
+        /// <summary>
+        /// Subtracts a scalar from each element of the vector.
+        /// </summary>
+        /// <param name="scalar">The scalar to subtract.</param>
+        public virtual void Subtract(double scalar)
+        {
+            if (scalar.AlmostZero())
+            {
+                return;
+            }
+
+            Parallel.For(0, Count, i => this[i] -= scalar);
+        }
+
+        /// <summary>
+        ///  Subtracts a scalar from each element of the vector and stores the result in the result vector.
+        /// </summary>
+        /// <param name="scalar">The scalar to subtract.</param>
+        /// <param name="result">The vector to store the result of the subtraction.</param>
+        /// <exception cref="ArgumentNullException">If the result vector is <see langword="null" />.</exception>
+        /// <exception cref="ArgumentException">If this vector and <paramref name="result"/> are not the same size.</exception>
+        public virtual void Subtract(double scalar, Vector result)
+        {
+            if (result == null)
+            {
+                throw new ArgumentNullException("result");
+            }
+
+            if (Count != result.Count)
+            {
+                throw new ArgumentException(Resources.ArgumentVectorsSameLength, "result");
+            }
+
+            CopyTo(result);
+            result.Subtract(scalar);
+        }
+
+        /// <summary>
+        /// Returns a negated vector.
+        /// </summary>
+        /// <returns>The negated vector.</returns>
+        /// <remarks>Added as an alternative to the unary negation operator.</remarks>
+        public virtual Vector Negate()
+        {
+            var result = CreateVector(Count);
+            Parallel.For(0, Count, i => result[i] = -this[i]);
+            return result;
+        }
+
+        /// <summary>
+        /// Subtracts another vector from this vector.
+        /// </summary>
+        /// <param name="other">The vector to subtract from this one.</param>
+        /// <exception cref="ArgumentNullException">If the other vector is <see langword="null" />.</exception> 
+        /// <exception cref="ArgumentException">If this vector and <paramref name="other"/> are not the same size.</exception>
+        public virtual void Subtract(Vector other)
+        {
+            if (other == null)
+            {
+                throw new ArgumentNullException("other");
+            }
+
+            if (Count != other.Count)
+            {
+                throw new ArgumentException(Resources.ArgumentVectorsSameLength, "other");
+            }
+
+            Parallel.For(0, Count, i => this[i] -= other[i]);
+        }
+
+        /// <summary>
+        /// Subtracts another vector to this vector and stores the result into the result vector.
+        /// </summary>
+        /// <param name="other">The vector to subtract from this one.</param>
+        /// <param name="result">The vector to store the result of the subtraction.</param>
+        /// <exception cref="ArgumentNullException">If the other vector is <see langword="null" />.</exception> 
+        /// <exception cref="ArgumentNullException">If the result vector is <see langword="null" />.</exception> 
+        /// <exception cref="ArgumentException">If this vector and <paramref name="other"/> are not the same size.</exception>
+        /// <exception cref="ArgumentException">If this vector and <paramref name="result"/> are not the same size.</exception>
+        public virtual void Subtract(Vector other, Vector result)
+        {
+            if (result == null)
+            {
+                throw new ArgumentNullException("result");
+            }
+
+            if (Count != result.Count)
+            {
+                throw new ArgumentException(Resources.ArgumentVectorsSameLength, "result");
+            }
+
+            if (ReferenceEquals(this, result) || ReferenceEquals(other, result))
+            {
+                var tmp = result.CreateVector(result.Count);
+                Subtract(other, tmp);
+                tmp.CopyTo(result);
+            }
+            else
+            {
+                CopyTo(result);
+                result.Subtract(other);
+            }
+        }
+
+        /// <summary>
+        /// Returns a <strong>Vector</strong> containing the negated values of rightSide. 
+        /// </summary>
+        /// <param name="rightSide">The vector to get the values from.</param>
+        /// <returns>A vector containing the negated values as <paramref name="rightSide"/>.</returns>
+        /// <exception cref="ArgumentNullException">If <paramref name="rightSide"/> is <see langword="null" />.</exception>
+        public static Vector operator -(Vector rightSide)
+        {
+            if (rightSide == null)
+            {
+                throw new ArgumentNullException("rightSide");
+            }
+
+            return rightSide.Negate();
+        }
+
+        /// <summary>
+        /// Subtracts two <strong>Vectors</strong> and returns the results.
+        /// </summary>
+        /// <param name="leftSide">The vector to subtract from.</param>
+        /// <param name="rightSide">The vector to subtract.</param>
+        /// <returns>The result of the subtraction.</returns>
+        /// <exception cref="ArgumentException">If <paramref name="leftSide"/> and <paramref name="rightSide"/> are not the same size.</exception>
+        /// <exception cref="ArgumentNullException">If <paramref name="leftSide"/> or <paramref name="rightSide"/> is <see langword="null" />.</exception>
+        public static Vector operator -(Vector leftSide, Vector rightSide)
+        {
+            if (rightSide == null)
+            {
+                throw new ArgumentNullException("rightSide");
+            }
+
+            if (leftSide == null)
+            {
+                throw new ArgumentNullException("leftSide");
+            }
+
+            if (leftSide.Count != rightSide.Count)
+            {
+                throw new ArgumentException(Resources.ArgumentVectorsSameLength, "rightSide");
+            }
+
+            var ret = leftSide.Clone();
+            ret.Subtract(rightSide);
+            return ret;
+        }
+
+        
         #region Implemented Interfaces
 
         #region ICloneable

src/UnitTests/LinearAlgebraTests/Double/DenseVectorTests.cs

@@ -6,7 +6,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Double
 
     public class DenseVectorTests : VectorTests
     {
-        protected override MathNet.Numerics.LinearAlgebra.Double.Vector CreateVector(int size)
+        protected override Vector CreateVector(int size)
         {
             return new DenseVector(size);
         }

src/UnitTests/LinearAlgebraTests/Double/VectorTests.cs

@@ -346,7 +346,6 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Double
             }
         }
 
-
         [Test]
         [MultipleAsserts]
         public void CanAddVectorToItselfUsingResultVector()
@@ -377,6 +376,231 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Double
             }
         }
 
+        [Test]
+        public void CanCallNegate()
+        {
+            var vector = CreateVector(_data);
+            var other = vector.Negate();
+            for (var i = 0; i < _data.Length; i++)
+            {
+                Assert.AreEqual(-_data[i], other[i]);
+            }
+        }
+
+        [Test]
+        public void OperatorNegateThrowsArgumentNullExceptionWhenCallOnNullVector()
+        {
+            Vector vector = null;
+            Vector other = null;
+            Assert.Throws<ArgumentNullException>(() => other = -vector);
+        }
+
+        [Test]
+        public void CanCallUnaryNegationOperator()
+        {
+            var vector = CreateVector(_data);
+            var other = -vector;
+            for (var i = 0; i < _data.Length; i++)
+            {
+                Assert.AreEqual(-_data[i], other[i]);
+            }
+        }
+
+        [Test]
+        [MultipleAsserts]
+        public void CanSubtractScalarFromVector()
+        {
+            var vector = CreateVector(_data);
+            vector.Subtract(2.0);
+
+            for (var i = 0; i < _data.Length; i++)
+            {
+                Assert.AreEqual(_data[i] - 2.0, vector[i]);
+            }
+
+            vector.Subtract(0.0);
+            for (var i = 0; i < _data.Length; i++)
+            {
+                Assert.AreEqual(_data[i] - 2.0, vector[i]);
+            }
+        }
+
+        [Test]
+        [MultipleAsserts]
+        public void CanSubtractScalarFromVectorUsingResultVector()
+        {
+            var vector = CreateVector(_data);
+            var result = CreateVector(_data.Length);
+            vector.Subtract(2.0, result);
+
+            for (var i = 0; i < _data.Length; i++)
+            {
+                Assert.AreEqual(_data[i], vector[i], "Making sure the original vector wasn't modified.");
+                Assert.AreEqual(_data[i] - 2.0, result[i]);
+            }
+
+            vector.Subtract(0.0, result);
+            for (var i = 0; i < _data.Length; i++)
+            {
+                Assert.AreEqual(_data[i], result[i]);
+            }
+        }
+
+        [Test]
+        public void ThrowsArgumentNullExceptionWhenSubtractingScalarWithNullResultVector()
+        {
+            var vector = CreateVector(_data.Length);
+            Assert.Throws<ArgumentNullException>(() => vector.Subtract(0.0, null));
+        }
+
+        [Test]
+        public void ThrowsArgumentExceptionWhenSubtractingScalarWithWrongSizeResultVector()
+        {
+            var vector = CreateVector(_data.Length);
+            var result = CreateVector(_data.Length + 1);
+            Assert.Throws<ArgumentException>(() => vector.Subtract(0.0, result));
+        }
+
+        [Test]
+        public void ThrowsArgumentNullExceptionWhenSubtractingTwoVectorsAndOneIsNull()
+        {
+            var vector = CreateVector(_data);
+            Assert.Throws<ArgumentNullException>(() => vector.Subtract(null));
+        }
+
+        [Test]
+        public void ThrowsArgumentExceptionWhenSubtractingTwoVectorsOfDifferingSize()
+        {
+            var vector = CreateVector(_data.Length);
+            var other = CreateVector(_data.Length + 1);
+            Assert.Throws<ArgumentException>(() => vector.Subtract(other));
+        }
+
+        [Test]
+        public void ThrowsArgumentNullExceptionWhenSubtractingTwoVectorsAndResultIsNull()
+        {
+            var vector = CreateVector(_data.Length);
+            var other = CreateVector(_data.Length + 1);
+            Assert.Throws<ArgumentNullException>(() => vector.Subtract(other, null));
+        }
+
+        [Test]
+        public void ThrowsArgumentExceptionWhenSubtractingTwoVectorsAndResultIsDifferentSize()
+        {
+            var vector = CreateVector(_data.Length);
+            var other = CreateVector(_data.Length);
+            var result = CreateVector(_data.Length + 1);
+            Assert.Throws<ArgumentException>(() => vector.Subtract(other, result));
+        }
+
+        [Test]
+        public void SubtractionOperatorThrowsArgumentNullExpectionIfAVectorIsNull()
+        {
+            Vector a = null;
+            var b = CreateVector(_data.Length);
+            Assert.Throws<ArgumentNullException>(() => a -= b);
+
+            a = b;
+            b = null;
+            Assert.Throws<ArgumentNullException>(() => a -= b);
+        }
+
+        [Test]
+        public void SubtractionOperatorThrowsArgumentExpectionIfVectorsAreDifferentSize()
+        {
+            var a = CreateVector(_data.Length);
+            var b = CreateVector(_data.Length + 1);
+            Assert.Throws<ArgumentException>(() => a -= b);
+        }
+
+        [Test]
+        public void CanSubtractTwoVectors()
+        {
+            var vector = CreateVector(_data);
+            var other = CreateVector(_data);
+            vector.Subtract(other);
+
+            for (var i = 0; i < _data.Length; i++)
+            {
+                Assert.AreEqual(0.0, vector[i]);
+            }
+        }
+
+        [Test]
+        [MultipleAsserts]
+        public void CanSubtractTwoVectorsUsingResultVector()
+        {
+            var vector = CreateVector(_data);
+            var other = CreateVector(_data);
+            var result = CreateVector(_data.Length);
+            vector.Subtract(other, result);
+
+            for (var i = 0; i < _data.Length; i++)
+            {
+                Assert.AreEqual(_data[i], vector[i], "Making sure the original vector wasn't modified.");
+                Assert.AreEqual(_data[i], other[i], "Making sure the original vector wasn't modified.");
+                Assert.AreEqual(0.0, result[i]);
+            }
+        }
+
+        [Test]
+        [MultipleAsserts]
+        public void CanSubtractTwoVectorsUsingOperator()
+        {
+            var vector = CreateVector(_data);
+            var other = CreateVector(_data);
+            var result = vector - other;
+
+            for (var i = 0; i < _data.Length; i++)
+            {
+                Assert.AreEqual(_data[i], vector[i], "Making sure the original vector wasn't modified.");
+                Assert.AreEqual(_data[i], other[i], "Making sure the original vector wasn't modified.");
+                Assert.AreEqual(0.0, result[i]);
+            }
+        }
+
+        [Test]
+        public void CanSubtractVectorFromItself()
+        {
+            var vector = CreateVector(_data);
+            vector.Subtract(vector);
+
+            for (var i = 0; i < _data.Length; i++)
+            {
+                Assert.AreEqual(0.0, vector[i]);
+            }
+        }
+
+
+        [Test]
+        [MultipleAsserts]
+        public void CanSubtractVectorFromItselfUsingResultVector()
+        {
+            var vector = CreateVector(_data);
+            var result = CreateVector(_data.Length);
+            vector.Subtract(vector, result);
+
+            for (var i = 0; i < _data.Length; i++)
+            {
+                Assert.AreEqual(_data[i], vector[i], "Making sure the original vector wasn't modified.");
+                Assert.AreEqual(0.0, result[i]);
+            }
+        }
+
+        [Test]
+        [MultipleAsserts]
+        public void CanSubtractTwoVectorsUsingItselfAsResultVector()
+        {
+            var vector = CreateVector(_data);
+            var other = CreateVector(_data);
+            vector.Subtract(other, vector);
+
+            for (var i = 0; i < _data.Length; i++)
+            {
+                Assert.AreEqual(_data[i], other[i], "Making sure the original vector wasn't modified.");
+                Assert.AreEqual(0.0, vector[i]);
+            }
+        }
 
         protected abstract Vector CreateVector(int size);