Browse Source

refactored the Vector classes to use an intermediate, type specific class

la-knuth
Marcus Cuda 16 years ago
parent
commit
d9bbaa4a90
  1. 6
      src/MathNet.Numerics.5.1.ReSharper
  2. 186
      src/Numerics/LinearAlgebra/Complex/DenseVector.cs
  3. 4
      src/Numerics/LinearAlgebra/Complex/Solvers/Iterative/TFQMR.cs
  4. 2
      src/Numerics/LinearAlgebra/Complex/Solvers/Preconditioners/Ilutp.cs
  5. 2
      src/Numerics/LinearAlgebra/Complex/Solvers/StopCriterium/DivergenceStopCriterium.cs
  6. 2
      src/Numerics/LinearAlgebra/Complex/Solvers/StopCriterium/FailureStopCriterium.cs
  7. 6
      src/Numerics/LinearAlgebra/Complex/Solvers/StopCriterium/ResidualStopCriterium.cs
  8. 170
      src/Numerics/LinearAlgebra/Complex/SparseVector.cs
  9. 446
      src/Numerics/LinearAlgebra/Complex/Vector.cs
  10. 138
      src/Numerics/LinearAlgebra/Complex32/DenseVector.cs
  11. 2
      src/Numerics/LinearAlgebra/Complex32/Factorization/UserGramSchmidt.cs
  12. 2
      src/Numerics/LinearAlgebra/Complex32/Solvers/Iterative/MlkBiCgStab.cs
  13. 4
      src/Numerics/LinearAlgebra/Complex32/Solvers/Iterative/TFQMR.cs
  14. 2
      src/Numerics/LinearAlgebra/Complex32/Solvers/Preconditioners/Ilutp.cs
  15. 2
      src/Numerics/LinearAlgebra/Complex32/Solvers/StopCriterium/DivergenceStopCriterium.cs
  16. 2
      src/Numerics/LinearAlgebra/Complex32/Solvers/StopCriterium/FailureStopCriterium.cs
  17. 20
      src/Numerics/LinearAlgebra/Complex32/Solvers/StopCriterium/ResidualStopCriterium.cs
  18. 178
      src/Numerics/LinearAlgebra/Complex32/SparseVector.cs
  19. 446
      src/Numerics/LinearAlgebra/Complex32/Vector.cs
  20. 161
      src/Numerics/LinearAlgebra/Double/DenseVector.cs
  21. 148
      src/Numerics/LinearAlgebra/Double/SparseVector.cs
  22. 470
      src/Numerics/LinearAlgebra/Double/Vector.cs
  23. 21
      src/Numerics/LinearAlgebra/Generic/Common.cs
  24. 6
      src/Numerics/LinearAlgebra/Generic/Matrix.cs
  25. 416
      src/Numerics/LinearAlgebra/Generic/Vector.cs
  26. 171
      src/Numerics/LinearAlgebra/Single/DenseVector.cs
  27. 18
      src/Numerics/LinearAlgebra/Single/Solvers/StopCriterium/ResidualStopCriterium.cs
  28. 162
      src/Numerics/LinearAlgebra/Single/SparseVector.cs
  29. 473
      src/Numerics/LinearAlgebra/Single/Vector.cs
  30. 5
      src/Numerics/Numerics.csproj
  31. 64
      src/Numerics/Precision.cs
  32. 52
      src/Numerics/Threading/CommonParallel.cs
  33. 15
      src/Silverlight/Silverlight.csproj
  34. 126
      src/UnitTests/LinearAlgebraTests/Complex/UserDefinedVectorTests.cs
  35. 4
      src/UnitTests/LinearAlgebraTests/Complex32/MatrixTests.Arithmetic.cs
  36. 34
      src/UnitTests/LinearAlgebraTests/Complex32/Solvers/StopCriterium/ResidualStopCriteriumTest.cs
  37. 136
      src/UnitTests/LinearAlgebraTests/Complex32/UserDefinedVectorTests.cs
  38. 10
      src/UnitTests/LinearAlgebraTests/Complex32/VectorTests.Norm.cs
  39. 132
      src/UnitTests/LinearAlgebraTests/Double/UserDefinedVectorTests.cs
  40. 34
      src/UnitTests/LinearAlgebraTests/Single/Solvers/StopCriterium/ResidualStopCriteriumTest.cs
  41. 159
      src/UnitTests/LinearAlgebraTests/Single/UserDefinedVectorTests.cs
  42. 10
      src/UnitTests/LinearAlgebraTests/Single/VectorTests.Norm.cs

6
src/MathNet.Numerics.5.1.ReSharper

@ -19,7 +19,11 @@ Matlab
Matlab
Matlab
Endian
indices</UserWords>
indices
&amp;lt
&amp;gt
Frobenius
Pointwise</UserWords>
</CustomDictionary>
</Dictionaries>
</CustomDictionaries>

186
src/Numerics/LinearAlgebra/Complex/DenseVector.cs

@ -39,7 +39,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
/// <summary>
/// A vector using dense storage.
/// </summary>
public class DenseVector : Vector<Complex>
public class DenseVector : Vector
{
/// <summary>
/// Initializes a new instance of the <see cref="DenseVector"/> class with a given size.
@ -767,7 +767,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
/// Returns the value of the absolute minimum element.
/// </summary>
/// <returns>The value of the absolute minimum element.</returns>
public override double AbsoluteMinimum()
public override Complex AbsoluteMinimum()
{
return Data[AbsoluteMinimumIndex()].Magnitude;
}
@ -776,7 +776,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
/// Returns the value of the absolute maximum element.
/// </summary>
/// <returns>The value of the absolute maximum element.</returns>
public override double AbsoluteMaximum()
public override Complex AbsoluteMaximum()
{
return Data[AbsoluteMaximumIndex()].Magnitude;
}
@ -869,28 +869,22 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
/// <returns>The sum of the vector's elements.</returns>
public override Complex Sum()
{
var result = Complex.Zero;
for (var i = 0; i < Count; i++)
{
result += Data[i];
}
return result;
return CommonParallel.Aggregate(
0,
Count,
i => Data[i]);
}
/// <summary>
/// Computes the sum of the absolute value of the vector's elements.
/// </summary>
/// <returns>The sum of the absolute value of the vector's elements.</returns>
public override double SumMagnitudes()
public override Complex SumMagnitudes()
{
double result = 0;
for (var i = 0; i < Count; i++)
{
result += Data[i].Magnitude;
}
return result;
return CommonParallel.Aggregate(
0,
Count,
i => Data[i].Magnitude);
}
/// <summary>
@ -1101,58 +1095,6 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
return matrix;
}
/// <summary>
/// Generates a vector with random elements
/// </summary>
/// <param name="length">Number of elements in the vector.</param>
/// <param name="randomDistribution">Continuous Random Distribution or Source</param>
/// <returns>
/// A vector with n-random elements distributed according
/// to the specified random distribution.
/// </returns>
/// <exception cref="ArgumentNullException">If the n vector is non positive<see langword="null" />.</exception>
public override Vector<Complex> Random(int length, IContinuousDistribution randomDistribution)
{
if (length < 1)
{
throw new ArgumentException(Resources.ArgumentMustBePositive, "length");
}
var v = (DenseVector)CreateVector(length);
for (var index = 0; index < v.Data.Length; index++)
{
v.Data[index] = new Complex(randomDistribution.Sample(), randomDistribution.Sample());
}
return v;
}
/// <summary>
/// Generates a vector with random elements
/// </summary>
/// <param name="length">Number of elements in the vector.</param>
/// <param name="randomDistribution">Continuous Random Distribution or Source</param>
/// <returns>
/// A vector with n-random elements distributed according
/// to the specified random distribution.
/// </returns>
/// <exception cref="ArgumentNullException">If the n vector is non positive<see langword="null" />.</exception>
public override Vector<Complex> Random(int length, IDiscreteDistribution randomDistribution)
{
if (length < 1)
{
throw new ArgumentException(Resources.ArgumentMustBePositive, "length");
}
var v = (DenseVector)CreateVector(length);
for (var index = 0; index < v.Data.Length; index++)
{
v.Data[index] = new Complex(randomDistribution.Sample(), randomDistribution.Sample());
}
return v;
}
/// <summary>
/// Outer product of this and another vector.
/// </summary>
@ -1171,7 +1113,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
/// </summary>
/// <param name="p">The p value.</param>
/// <returns>Scalar <c>ret = (sum(abs(this[i])^p))^(1/p)</c></returns>
public override double Norm(double p)
public override Complex Norm(double p)
{
if (p < 0.0)
{
@ -1208,34 +1150,6 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
return Math.Pow(sum, 1.0 / p);
}
/// <summary>
/// Normalizes this vector to a unit vector with respect to the p-norm.
/// </summary>
/// <param name="p">
/// The p value.
/// </param>
/// <returns>
/// This vector normalized to a unit vector with respect to the p-norm.
/// </returns>
public override Vector<Complex> Normalize(double p)
{
if (p < 0.0)
{
throw new ArgumentOutOfRangeException("p");
}
var norm = Norm(p);
var clone = Clone();
if (norm == 0.0)
{
return clone;
}
clone.Multiply(1.0 / norm, clone);
return clone;
}
#region Parse Functions
/// <summary>
@ -1395,24 +1309,6 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
#endregion
/// <summary>
/// Returns the index of the absolute maximum element.
/// </summary>
/// <returns>The index of absolute maximum element.</returns>
public override int MaximumIndex()
{
throw new NotSupportedException();
}
/// <summary>
/// Returns the index of the minimum element.
/// </summary>
/// <returns>The index of minimum element.</returns>
public override int MinimumIndex()
{
throw new NotSupportedException();
}
/// <summary>
/// Resets all values to zero.
/// </summary>
@ -1457,61 +1353,5 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
index => otherVector.Data[index] = Data[index].Conjugate());
}
}
#region Simple arithmetic of type T
/// <summary>
/// Add two values T+T
/// </summary>
/// <param name="val1">Left operand value</param>
/// <param name="val2">Right operand value</param>
/// <returns>Result of addition</returns>
protected sealed override Complex AddT(Complex val1, Complex val2)
{
return val1 + val2;
}
/// <summary>
/// Subtract two values T-T
/// </summary>
/// <param name="val1">Left operand value</param>
/// <param name="val2">Right operand value</param>
/// <returns>Result of subtract</returns>
protected sealed override Complex SubtractT(Complex val1, Complex val2)
{
return val1 - val2;
}
/// <summary>
/// Multiply two values T*T
/// </summary>
/// <param name="val1">Left operand value</param>
/// <param name="val2">Right operand value</param>
/// <returns>Result of multiplication</returns>
protected sealed override Complex MultiplyT(Complex val1, Complex val2)
{
return val1 * val2;
}
/// <summary>
/// Divide two values T/T
/// </summary>
/// <param name="val1">Left operand value</param>
/// <param name="val2">Right operand value</param>
/// <returns>Result of divide</returns>
protected sealed override Complex DivideT(Complex val1, Complex val2)
{
return val1 / val2;
}
/// <summary>
/// Take absolute value
/// </summary>
/// <param name="val1">Source alue</param>
/// <returns>True if one; otherwise false</returns>
protected sealed override double AbsoluteT(Complex val1)
{
return val1.Magnitude;
}
#endregion
}
}

4
src/Numerics/LinearAlgebra/Complex/Solvers/Iterative/TFQMR.cs

@ -286,7 +286,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Solvers.Iterative
Complex eta = 0;
double theta = 0;
var tau = startNorm;
var tau = startNorm.Real;
Complex rho = tau * tau;
// Calculate the initial values for v
@ -342,7 +342,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Solvers.Iterative
yinternal.Add(temp, d);
// theta = ||pseudoResiduals||_2 / tau
theta = pseudoResiduals.Norm(2) / tau;
theta = pseudoResiduals.Norm(2).Real / tau;
var c = 1 / Math.Sqrt(1 + (theta * theta));
// tau = tau * theta * c

2
src/Numerics/LinearAlgebra/Complex/Solvers/Preconditioners/Ilutp.cs

@ -442,7 +442,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Solvers.Preconditioners
// {
// w(j) = 0
// }
if (workVector[j].Magnitude <= _dropTolerance * vectorNorm)
if (workVector[j].Magnitude <= _dropTolerance * vectorNorm.Real)
{
workVector[j] = 0.0;
}

2
src/Numerics/LinearAlgebra/Complex/Solvers/StopCriterium/DivergenceStopCriterium.cs

@ -252,7 +252,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Solvers.StopCriterium
// Store the infinity norms of both the solution and residual vectors
// These values will be used to calculate the relative drop in residuals later on.
_residualHistory[_residualHistory.Length - 1] = residualVector.Norm(Double.PositiveInfinity);
_residualHistory[_residualHistory.Length - 1] = residualVector.Norm(Double.PositiveInfinity).Real;
// Check if we have NaN's. If so we've gone way beyond normal divergence.
// Stop the iteration.

2
src/Numerics/LinearAlgebra/Complex/Solvers/StopCriterium/FailureStopCriterium.cs

@ -110,7 +110,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Solvers.StopCriterium
var residualNorm = residualVector.Norm(Double.PositiveInfinity);
var solutionNorm = solutionVector.Norm(Double.PositiveInfinity);
if (Double.IsNaN(solutionNorm) || Double.IsNaN(residualNorm))
if (Double.IsNaN(solutionNorm.Real) || Double.IsNaN(residualNorm.Real))
{
SetStatusToFailed();
}

6
src/Numerics/LinearAlgebra/Complex/Solvers/StopCriterium/ResidualStopCriterium.cs

@ -263,12 +263,12 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Solvers.StopCriterium
// Check the residuals by calculating:
// ||r_i|| <= stop_tol * ||b||
var stopCriterium = ComputeStopCriterium(sourceVector.Norm(Double.PositiveInfinity));
var stopCriterium = ComputeStopCriterium(sourceVector.Norm(Double.PositiveInfinity).Real);
// First check that we have real numbers not NaN's.
// NaN's can occur when the iterative process diverges so we
// stop if that is the case.
if (double.IsNaN(stopCriterium) || double.IsNaN(residualNorm))
if (double.IsNaN(stopCriterium) || double.IsNaN(residualNorm.Real))
{
_iterationCount = 0;
SetStatusToDiverged();
@ -278,7 +278,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Solvers.StopCriterium
// ||r_i|| <= stop_tol * ||b||
// Stop the calculation if it's clearly smaller than the tolerance
var decimalMagnitude = Math.Abs(stopCriterium.Magnitude()) + 1;
if (residualNorm.IsSmallerWithDecimalPlaces(stopCriterium, decimalMagnitude))
if (residualNorm.Real.IsSmallerWithDecimalPlaces(stopCriterium, decimalMagnitude))
{
if (_lastIteration <= iterationNumber)
{

170
src/Numerics/LinearAlgebra/Complex/SparseVector.cs

@ -30,7 +30,6 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
using System.Collections.Generic;
using System.Linq;
using System.Numerics;
using Distributions;
using Generic;
using NumberTheory;
using Properties;
@ -39,7 +38,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
/// <summary>
/// A vector with sparse storage.
/// </summary>
public class SparseVector : Vector<Complex>
public class SparseVector : Vector
{
/// <summary>
/// Lock object for the indexer.
@ -581,11 +580,12 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
var sparseother = other as SparseVector;
if (sparseother == null)
{
sparse.AddScaledSparseVector(Complex.One, sparseother);
base.Add(other, result);
}
else
{
base.Add(other, result);
CopyTo(result);
sparse.AddScaledSparseVector(Complex.One, sparseother);
}
}
}
@ -753,11 +753,12 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
var sparseother = other as SparseVector;
if (sparseother == null)
{
sparse.AddScaledSparseVector(-Complex.One, sparseother);
base.Subtract(other, result);
}
else
{
base.Subtract(other, result);
CopyTo(result);
sparse.AddScaledSparseVector(-Complex.One, sparseother);
}
}
}
@ -1075,7 +1076,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
/// Computes the sum of the absolute value of the vector's elements.
/// </summary>
/// <returns>The sum of the absolute value of the vector's elements.</returns>
public override double SumMagnitudes()
public override Complex SumMagnitudes()
{
double result = 0;
for (var i = 0; i < NonZerosCount; i++)
@ -1153,58 +1154,6 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
return matrix;
}
/// <summary>
/// Generates a vector with random elements
/// </summary>
/// <param name="length">Number of elements in the vector.</param>
/// <param name="randomDistribution">Continuous Random Distribution or Source</param>
/// <returns>
/// A vector with n-random elements distributed according
/// to the specified random distribution.
/// </returns>
/// <exception cref="ArgumentNullException">If the length vector is non positive<see langword="null" />.</exception>
public override Vector<Complex> Random(int length, IContinuousDistribution randomDistribution)
{
if (length < 1)
{
throw new ArgumentException(Resources.ArgumentMustBePositive, "length");
}
var v = (SparseVector)CreateVector(length);
for (var index = 0; index < v.Count; index++)
{
v[index] = new Complex(randomDistribution.Sample(), randomDistribution.Sample());
}
return v;
}
/// <summary>
/// Generates a vector with random elements
/// </summary>
/// <param name="length">Number of elements in the vector.</param>
/// <param name="randomDistribution">Continuous Random Distribution or Source</param>
/// <returns>
/// A vector with n-random elements distributed according
/// to the specified random distribution.
/// </returns>
/// <exception cref="ArgumentNullException">If the n vector is non positive<see langword="null" />.</exception>
public override Vector<Complex> Random(int length, IDiscreteDistribution randomDistribution)
{
if (length < 1)
{
throw new ArgumentException(Resources.ArgumentMustBePositive, "length");
}
var v = (SparseVector)CreateVector(length);
for (var index = 0; index < v.Count; index++)
{
v[index] = new Complex(randomDistribution.Sample(), randomDistribution.Sample());
}
return v;
}
/// <summary>
/// Outer product of this and another vector.
/// </summary>
@ -1227,7 +1176,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
/// </summary>
/// <param name="p">The p value.</param>
/// <returns>Scalar <c>ret = (sum(abs(this[i])^p))^(1/p)</c></returns>
public override double Norm(double p)
public override Complex Norm(double p)
{
if (1 > p)
{
@ -1257,33 +1206,6 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
return Math.Pow(sum, 1.0 / p);
}
/// <summary>
/// Normalizes this vector to a unit vector with respect to the p-norm.
/// </summary>
/// <param name="p">
/// The p value.
/// </param>
/// <returns>
/// This vector normalized to a unit vector with respect to the p-norm.
/// </returns>
public override Vector<Complex> Normalize(double p)
{
if (p < 0.0)
{
throw new ArgumentOutOfRangeException("p");
}
var norm = Norm(p);
var clone = Clone();
if (norm == 0.0)
{
return clone;
}
clone.Multiply(1.0 / norm, clone);
return clone;
}
#endregion
#region Parse Functions
@ -1610,79 +1532,5 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
}
#endregion
/// <summary>
/// Returns the index of the absolute maximum element.
/// </summary>
/// <returns>The index of absolute maximum element.</returns>
public override int MaximumIndex()
{
throw new NotSupportedException();
}
/// <summary>
/// Returns the index of the minimum element.
/// </summary>
/// <returns>The index of minimum element.</returns>
public override int MinimumIndex()
{
throw new NotSupportedException();
}
#region Simple arithmetic of type T
/// <summary>
/// Add two values T+T
/// </summary>
/// <param name="val1">Left operand value</param>
/// <param name="val2">Right operand value</param>
/// <returns>Result of addition</returns>
protected sealed override Complex AddT(Complex val1, Complex val2)
{
return val1 + val2;
}
/// <summary>
/// Subtract two values T-T
/// </summary>
/// <param name="val1">Left operand value</param>
/// <param name="val2">Right operand value</param>
/// <returns>Result of subtract</returns>
protected sealed override Complex SubtractT(Complex val1, Complex val2)
{
return val1 - val2;
}
/// <summary>
/// Multiply two values T*T
/// </summary>
/// <param name="val1">Left operand value</param>
/// <param name="val2">Right operand value</param>
/// <returns>Result of multiplication</returns>
protected sealed override Complex MultiplyT(Complex val1, Complex val2)
{
return val1 * val2;
}
/// <summary>
/// Divide two values T/T
/// </summary>
/// <param name="val1">Left operand value</param>
/// <param name="val2">Right operand value</param>
/// <returns>Result of divide</returns>
protected sealed override Complex DivideT(Complex val1, Complex val2)
{
return val1 / val2;
}
/// <summary>
/// Take absolute value
/// </summary>
/// <param name="val1">Source alue</param>
/// <returns>True if one; otherwise false</returns>
protected sealed override double AbsoluteT(Complex val1)
{
return val1.Magnitude;
}
#endregion
}
}

446
src/Numerics/LinearAlgebra/Complex/Vector.cs

@ -0,0 +1,446 @@
// <copyright file="Vector.cs" company="Math.NET">
// Math.NET Numerics, part of the Math.NET Project
// http://numerics.mathdotnet.com
// http://github.com/mathnet/mathnet-numerics
// http://mathnetnumerics.codeplex.com
// Copyright (c) 2009-2010 Math.NET
// Permission is hereby granted, free of charge, to any person
// obtaining a copy of this software and associated documentation
// files (the "Software"), to deal in the Software without
// restriction, including without limitation the rights to use,
// copy, modify, merge, publish, distribute, sublicense, and/or sell
// copies of the Software, and to permit persons to whom the
// Software is furnished to do so, subject to the following
// conditions:
// The above copyright notice and this permission notice shall be
// included in all copies or substantial portions of the Software.
// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES
// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT
// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY,
// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
// OTHER DEALINGS IN THE SOFTWARE.
// </copyright>
namespace MathNet.Numerics.LinearAlgebra.Complex
{
using System;
using System.Numerics;
using Distributions;
using Generic;
using Properties;
using Threading;
/// <summary>
/// <c>Complex</c> version of the <see cref="Vector{T}"/> class.
/// </summary>
public abstract class Vector : Vector<Complex>
{
/// <summary>
/// Initializes a new instance of the Vector class.
/// Constructs a <strong>Vector</strong> with the given size.
/// </summary>
/// <param name="size">
/// The size of the <strong>Vector</strong> to construct.
/// </param>
/// <exception cref="ArgumentException">
/// If <paramref name="size"/> is less than one.
/// </exception>
protected Vector(int size) : base(size)
{
}
/// <summary>
/// Adds a scalar to each element of the vector and stores the result in the result vector.
/// </summary>
/// <param name="scalar">
/// The scalar to add.
/// </param>
/// <param name="result">
/// The vector to store the result of the addition.
/// </param>
protected override void DoAdd(Complex scalar, Vector<Complex> result)
{
CommonParallel.For(
0,
Count,
index => result[index] = result[index] + scalar);
}
/// <summary>
/// Adds another vector to this vector and stores the result into the result vector.
/// </summary>
/// <param name="other">
/// The vector to add to this one.
/// </param>
/// <param name="result">
/// The vector to store the result of the addition.
/// </param>
protected override void DoAdd(Vector<Complex> other, Vector<Complex> result)
{
CommonParallel.For(
0,
Count,
index => result[index] = this[index] + other[index]);
}
/// <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>
protected override void DoSubtract(Complex scalar, Vector<Complex> result)
{
DoAdd(-scalar, result);
}
/// <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>
protected override void DoSubtract(Vector<Complex> other, Vector<Complex> result)
{
CopyTo(result);
CommonParallel.For(
0,
Count,
index => result[index] = this[index] - other[index]);
}
/// <summary>
/// Multiplies a scalar to each element of the vector and stores the result in the result vector.
/// </summary>
/// <param name="scalar">
/// The scalar to multiply.
/// </param>
/// <param name="result">
/// The vector to store the result of the multiplication.
/// </param>
protected override void DoMultiply(Complex scalar, Vector<Complex> result)
{
CommonParallel.For(
0,
Count,
index => result[index] = result[index] * scalar);
}
/// <summary>
/// Divides each element of the vector by a scalar and stores the result in the result vector.
/// </summary>
/// <param name="scalar">
/// The scalar to divide with.
/// </param>
/// <param name="result">
/// The vector to store the result of the division.
/// </param>
protected override void DoDivide(Complex scalar, Vector<Complex> result)
{
CommonParallel.For(
0,
Count,
index => result[index] = result[index] / scalar);
}
/// <summary>
/// Pointwise multiplies this vector with another vector and stores the result into the result vector.
/// </summary>
/// <param name="other">The vector to pointwise multiply with this one.</param>
/// <param name="result">The vector to store the result of the pointwise multiplication.</param>
protected override void DoPointwiseMultiply(Vector<Complex> other, Vector<Complex> result)
{
CommonParallel.For(
0,
Count,
index => result[index] = this[index] * other[index]);
}
/// <summary>
/// Pointwise divide this vector with another vector and stores the result into the result vector.
/// </summary>
/// <param name="other">The vector to pointwise divide this one by.</param>
/// <param name="result">The vector to store the result of the pointwise division.</param>
protected override void DoPointwiseDivide(Vector<Complex> other, Vector<Complex> result)
{
CommonParallel.For(
0,
Count,
index => result[index] = this[index] / other[index]);
}
/// <summary>
/// Computes the dot product between this vector and another vector.
/// </summary>
/// <param name="other">
/// The other vector to add.
/// </param>
/// <returns>s
/// The result of the addition.
/// </returns>
protected override Complex DoDotProduct(Vector<Complex> other)
{
return CommonParallel.Aggregate(
0,
Count,
i => this[i] * other[i]);
}
/// <summary>
/// Returns the value of the absolute minimum element.
/// </summary>
/// <returns>The value of the absolute minimum element.</returns>
public override Complex AbsoluteMinimum()
{
return this[AbsoluteMinimumIndex()].Magnitude;
}
/// <summary>
/// Returns the index of the absolute minimum element.
/// </summary>
/// <returns>The index of absolute minimum element.</returns>
public override int AbsoluteMinimumIndex()
{
var index = 0;
var min = this[index].Magnitude;
for (var i = 1; i < Count; i++)
{
var test = this[i].Magnitude;
if (test < min)
{
index = i;
min = test;
}
}
return index;
}
/// <summary>
/// Returns the value of the absolute maximum element.
/// </summary>
/// <returns>The value of the absolute maximum element.</returns>
public override Complex AbsoluteMaximum()
{
return this[AbsoluteMaximumIndex()].Magnitude;
}
/// <summary>
/// Returns the index of the absolute maximum element.
/// </summary>
/// <returns>The index of absolute maximum element.</returns>
public override int AbsoluteMaximumIndex()
{
var index = 0;
var max = this[index].Magnitude;
for (var i = 1; i < Count; i++)
{
var test = this[i].Magnitude;
if (test > max)
{
index = i;
max = test;
}
}
return index;
}
/// <summary>
/// Computes the sum of the vector's elements.
/// </summary>
/// <returns>The sum of the vector's elements.</returns>
public override Complex Sum()
{
return CommonParallel.Aggregate(
0,
Count,
i => this[i]);
}
/// <summary>
/// Computes the sum of the absolute value of the vector's elements.
/// </summary>
/// <returns>The sum of the absolute value of the vector's elements.</returns>
public override Complex SumMagnitudes()
{
return CommonParallel.Aggregate(
0,
Count,
i => this[i].Magnitude);
}
/// <summary>
/// Computes the p-Norm.
/// </summary>
/// <param name="p">
/// The p value.
/// </param>
/// <returns>
/// <c>Scalar ret = (sum(abs(this[i])^p))^(1/p)</c>
/// </returns>
public override Complex Norm(double p)
{
if (p < 0.0)
{
throw new ArgumentOutOfRangeException("p");
}
if (double.IsPositiveInfinity(p))
{
return CommonParallel.Select(
0,
Count,
(index, localData) => Math.Max(localData, this[index].Magnitude),
Math.Max);
}
var sum = CommonParallel.Aggregate(
0,
Count,
index => Math.Pow(this[index].Magnitude, p));
return Math.Pow(sum, 1.0 / p);
}
/// <summary>
/// Conjugates vector and save result to <paramref name="target"/>
/// </summary>
/// <param name="target">Target vector</param>
protected override void DoConjugate(Vector<Complex> target)
{
CopyTo(target);
CommonParallel.For(
0,
Count,
index => target[index] = this[index].Conjugate());
}
/// <summary>
/// Returns a negated vector.
/// </summary>
/// <returns>
/// The negated vector.
/// </returns>
/// <remarks>
/// Added as an alternative to the unary negation operator.
/// </remarks>
public override Vector<Complex> Negate()
{
var result = CreateVector(Count);
CommonParallel.For(
0,
Count,
index => result[index] = -this[index]);
return result;
}
/// <summary>
/// Generates a vector with random elements
/// </summary>
/// <param name="length">Number of elements in the vector.</param>
/// <param name="randomDistribution">Continuous Random Distribution or Source</param>
/// <returns>
/// A vector with n-random elements distributed according
/// to the specified random distribution.
/// </returns>
/// <exception cref="ArgumentException">If the n vector is non-positive.</exception>
public override Vector<Complex> Random(int length, IContinuousDistribution randomDistribution)
{
if (length < 1)
{
throw new ArgumentException(Resources.ArgumentMustBePositive, "length");
}
var vector = CreateVector(length);
for (var index = 0; index < length; index++)
{
vector[index] = new Complex(randomDistribution.Sample(), randomDistribution.Sample());
}
return vector;
}
/// <summary>
/// Generates a vector with random elements
/// </summary>
/// <param name="length">Number of elements in the vector.</param>
/// <param name="randomDistribution">Continuous Random Distribution or Source</param>
/// <returns>
/// A vector with n-random elements distributed according
/// to the specified random distribution.
/// </returns>
/// <exception cref="ArgumentException">If the n vector is not positive.</exception>
public override Vector<Complex> Random(int length, IDiscreteDistribution randomDistribution)
{
if (length < 1)
{
throw new ArgumentException(Resources.ArgumentMustBePositive, "length");
}
var vector = CreateVector(length);
for (var index = 0; index < length; index++)
{
vector[index] = new Complex(randomDistribution.Sample(), randomDistribution.Sample());
}
return vector;
}
/// <summary>
/// Returns the index of the absolute maximum element.
/// </summary>
/// <returns>The index of absolute maximum element.</returns>
public override int MaximumIndex()
{
throw new NotSupportedException();
}
/// <summary>
/// Returns the index of the minimum element.
/// </summary>
/// <returns>The index of minimum element.</returns>
public override int MinimumIndex()
{
throw new NotSupportedException();
}
/// <summary>
/// Normalizes this vector to a unit vector with respect to the p-norm.
/// </summary>
/// <param name="p">
/// The p value.
/// </param>
/// <returns>
/// This vector normalized to a unit vector with respect to the p-norm.
/// </returns>
public override Vector<Complex> Normalize(double p)
{
if (p < 0.0)
{
throw new ArgumentOutOfRangeException("p");
}
var norm = Norm(p);
var clone = Clone();
if (norm.Real == 0.0)
{
return clone;
}
clone.Multiply(1.0 / norm, clone);
return clone;
}
}
}

138
src/Numerics/LinearAlgebra/Complex32/DenseVector.cs

@ -39,7 +39,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
/// <summary>
/// A vector using dense storage.
/// </summary>
public class DenseVector : Vector<Complex32>
public class DenseVector : Vector
{
/// <summary>
/// Initializes a new instance of the <see cref="DenseVector"/> class with a given size.
@ -767,7 +767,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
/// Returns the value of the absolute minimum element.
/// </summary>
/// <returns>The value of the absolute minimum element.</returns>
public override double AbsoluteMinimum()
public override Complex32 AbsoluteMinimum()
{
return Data[AbsoluteMinimumIndex()].Magnitude;
}
@ -776,7 +776,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
/// Returns the value of the absolute maximum element.
/// </summary>
/// <returns>The value of the absolute maximum element.</returns>
public override double AbsoluteMaximum()
public override Complex32 AbsoluteMaximum()
{
return Data[AbsoluteMaximumIndex()].Magnitude;
}
@ -869,28 +869,22 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
/// <returns>The sum of the vector's elements.</returns>
public override Complex32 Sum()
{
var result = Complex32.Zero;
for (var i = 0; i < Count; i++)
{
result += Data[i];
}
return result;
return CommonParallel.Aggregate(
0,
Count,
i => Data[i]);
}
/// <summary>
/// Computes the sum of the absolute value of the vector's elements.
/// </summary>
/// <returns>The sum of the absolute value of the vector's elements.</returns>
public override double SumMagnitudes()
public override Complex32 SumMagnitudes()
{
double result = 0;
for (var i = 0; i < Count; i++)
{
result += Data[i].Magnitude;
}
return result;
return CommonParallel.Aggregate(
0,
Count,
i => Data[i].Magnitude);
}
/// <summary>
@ -1171,7 +1165,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
/// </summary>
/// <param name="p">The p value.</param>
/// <returns>Scalar <c>ret = (sum(abs(this[i])^p))^(1/p)</c></returns>
public override double Norm(double p)
public override Complex32 Norm(double p)
{
if (p < 0.0)
{
@ -1197,7 +1191,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
0,
Count,
(index, localData) => Math.Max(localData, Data[index].Magnitude),
Math.Max);
Common.Max);
}
var sum = CommonParallel.Aggregate(
@ -1205,35 +1199,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
Count,
index => Math.Pow(Data[index].Magnitude, p));
return Math.Pow(sum, 1.0 / p);
}
/// <summary>
/// Normalizes this vector to a unit vector with respect to the p-norm.
/// </summary>
/// <param name="p">
/// The p value.
/// </param>
/// <returns>
/// This vector normalized to a unit vector with respect to the p-norm.
/// </returns>
public override Vector<Complex32> Normalize(double p)
{
if (p < 0.0)
{
throw new ArgumentOutOfRangeException("p");
}
var norm = Norm(p);
var clone = Clone();
if (norm == 0.0)
{
return clone;
}
clone.Multiply(Complex32.One / (float)norm, clone);
return clone;
return (float)Math.Pow(sum, 1.0 / p);
}
#region Parse Functions
@ -1395,24 +1361,6 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
#endregion
/// <summary>
/// Returns the index of the absolute maximum element.
/// </summary>
/// <returns>The index of absolute maximum element.</returns>
public override int MaximumIndex()
{
throw new NotSupportedException();
}
/// <summary>
/// Returns the index of the minimum element.
/// </summary>
/// <returns>The index of minimum element.</returns>
public override int MinimumIndex()
{
throw new NotSupportedException();
}
/// <summary>
/// Resets all values to zero.
/// </summary>
@ -1457,61 +1405,5 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
index => otherVector.Data[index] = Data[index].Conjugate());
}
}
#region Simple arithmetic of type T
/// <summary>
/// Add two values T+T
/// </summary>
/// <param name="val1">Left operand value</param>
/// <param name="val2">Right operand value</param>
/// <returns>Result of addition</returns>
protected sealed override Complex32 AddT(Complex32 val1, Complex32 val2)
{
return val1 + val2;
}
/// <summary>
/// Subtract two values T-T
/// </summary>
/// <param name="val1">Left operand value</param>
/// <param name="val2">Right operand value</param>
/// <returns>Result of subtract</returns>
protected sealed override Complex32 SubtractT(Complex32 val1, Complex32 val2)
{
return val1 - val2;
}
/// <summary>
/// Multiply two values T*T
/// </summary>
/// <param name="val1">Left operand value</param>
/// <param name="val2">Right operand value</param>
/// <returns>Result of multiplication</returns>
protected sealed override Complex32 MultiplyT(Complex32 val1, Complex32 val2)
{
return val1 * val2;
}
/// <summary>
/// Divide two values T/T
/// </summary>
/// <param name="val1">Left operand value</param>
/// <param name="val2">Right operand value</param>
/// <returns>Result of divide</returns>
protected sealed override Complex32 DivideT(Complex32 val1, Complex32 val2)
{
return val1 / val2;
}
/// <summary>
/// Take absolute value
/// </summary>
/// <param name="val1">Source alue</param>
/// <returns>True if one; otherwise false</returns>
protected sealed override double AbsoluteT(Complex32 val1)
{
return val1.Magnitude;
}
#endregion
}
}

2
src/Numerics/LinearAlgebra/Complex32/Factorization/UserGramSchmidt.cs

@ -70,7 +70,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization
for (var k = 0; k < MatrixQ.ColumnCount; k++)
{
var norm = (float)MatrixQ.Column(k).Norm(2);
var norm = MatrixQ.Column(k).Norm(2).Real;
if (norm == 0.0f)
{
throw new ArgumentException(Resources.ArgumentMatrixNotRankDeficient);

2
src/Numerics/LinearAlgebra/Complex32/Solvers/Iterative/MlkBiCgStab.cs

@ -652,7 +652,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Solvers.Iterative
result.Add(orthogonalMatrix.Column(i));
// Normalize the result vector
result[i].Multiply(1 / (float)result[i].Norm(2), result[i]);
result[i].Multiply(1 / result[i].Norm(2).Real, result[i]);
}
return result;

4
src/Numerics/LinearAlgebra/Complex32/Solvers/Iterative/TFQMR.cs

@ -286,7 +286,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Solvers.Iterative
Complex32 eta = 0;
float theta = 0;
var tau = (float)startNorm;
var tau = startNorm.Real;
Complex32 rho = tau * tau;
// Calculate the initial values for v
@ -342,7 +342,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Solvers.Iterative
yinternal.Add(temp, d);
// theta = ||pseudoResiduals||_2 / tau
theta = (float)pseudoResiduals.Norm(2) / tau;
theta = pseudoResiduals.Norm(2).Real / tau;
var c = 1 / (float)Math.Sqrt(1 + (theta * theta));
// tau = tau * theta * c

2
src/Numerics/LinearAlgebra/Complex32/Solvers/Preconditioners/Ilutp.cs

@ -442,7 +442,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Solvers.Preconditioners
// {
// w(j) = 0
// }
if (workVector[j].Magnitude <= _dropTolerance * vectorNorm)
if (workVector[j].Magnitude <= _dropTolerance * vectorNorm.Real)
{
workVector[j] = 0.0f;
}

2
src/Numerics/LinearAlgebra/Complex32/Solvers/StopCriterium/DivergenceStopCriterium.cs

@ -252,7 +252,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Solvers.StopCriterium
// Store the infinity norms of both the solution and residual vectors
// These values will be used to calculate the relative drop in residuals later on.
_residualHistory[_residualHistory.Length - 1] = residualVector.Norm(Double.PositiveInfinity);
_residualHistory[_residualHistory.Length - 1] = residualVector.Norm(Double.PositiveInfinity).Real;
// Check if we have NaN's. If so we've gone way beyond normal divergence.
// Stop the iteration.

2
src/Numerics/LinearAlgebra/Complex32/Solvers/StopCriterium/FailureStopCriterium.cs

@ -110,7 +110,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Solvers.StopCriterium
var residualNorm = residualVector.Norm(Double.PositiveInfinity);
var solutionNorm = solutionVector.Norm(Double.PositiveInfinity);
if (Double.IsNaN(solutionNorm) || Double.IsNaN(residualNorm))
if (Single.IsNaN(solutionNorm.Real) || Single.IsNaN(residualNorm.Real))
{
SetStatusToFailed();
}

20
src/Numerics/LinearAlgebra/Complex32/Solvers/StopCriterium/ResidualStopCriterium.cs

@ -46,7 +46,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Solvers.StopCriterium
/// <summary>
/// The default value for the maximum value of the residual.
/// </summary>
public const double DefaultMaximumResidual = 1e-6;
public const float DefaultMaximumResidual = 1e-6f;
/// <summary>
/// The default value for the minimum number of iterations.
@ -66,7 +66,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Solvers.StopCriterium
/// <summary>
/// The maximum value for the residual below which the calculation is considered converged.
/// </summary>
private double _maximum;
private float _maximum;
/// <summary>
/// The minimum number of iterations for which the residual has to be below the maximum before
@ -102,7 +102,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Solvers.StopCriterium
/// maximum residual and the default minimum number of iterations.
/// </summary>
/// <param name="maximum">The maximum value for the residual below which the calculation is considered converged.</param>
public ResidualStopCriterium(double maximum) : this(maximum, DefaultMinimumIterationsBelowMaximum)
public ResidualStopCriterium(float maximum) : this(maximum, DefaultMinimumIterationsBelowMaximum)
{
}
@ -129,7 +129,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Solvers.StopCriterium
/// The minimum number of iterations for which the residual has to be below the maximum before
/// the calculation is considered converged.
/// </param>
public ResidualStopCriterium(double maximum, int minimumIterationsBelowMaximum)
public ResidualStopCriterium(float maximum, int minimumIterationsBelowMaximum)
{
if (maximum < 0)
{
@ -150,7 +150,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Solvers.StopCriterium
/// converged.
/// </summary>
/// <exception cref="ArgumentOutOfRangeException">Thrown if the <c>Maximum</c> is set to a negative value.</exception>
public double Maximum
public float Maximum
{
[DebuggerStepThrough]
get
@ -259,16 +259,16 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Solvers.StopCriterium
// Store the infinity norms of both the solution and residual vectors
// These values will be used to calculate the relative drop in residuals
// later on.
var residualNorm = residualVector.Norm(Double.PositiveInfinity);
var residualNorm = residualVector.Norm(float.PositiveInfinity);
// Check the residuals by calculating:
// ||r_i|| <= stop_tol * ||b||
var stopCriterium = ComputeStopCriterium(sourceVector.Norm(Double.PositiveInfinity));
var stopCriterium = ComputeStopCriterium(sourceVector.Norm(float.PositiveInfinity).Real);
// First check that we have real numbers not NaN's.
// NaN's can occur when the iterative process diverges so we
// stop if that is the case.
if (double.IsNaN(stopCriterium) || double.IsNaN(residualNorm))
if (float.IsNaN(stopCriterium) || float.IsNaN(residualNorm.Real))
{
_iterationCount = 0;
SetStatusToDiverged();
@ -278,7 +278,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Solvers.StopCriterium
// ||r_i|| <= stop_tol * ||b||
// Stop the calculation if it's clearly smaller than the tolerance
var decimalMagnitude = Math.Abs(stopCriterium.Magnitude()) + 1;
if (residualNorm.IsSmallerWithDecimalPlaces(stopCriterium, decimalMagnitude))
if (residualNorm.Real.IsSmallerWithDecimalPlaces(stopCriterium, decimalMagnitude))
{
if (_lastIteration <= iterationNumber)
{
@ -307,7 +307,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Solvers.StopCriterium
/// </summary>
/// <param name="solutionNorm">Solution vector norm</param>
/// <returns>Criterium value</returns>
private double ComputeStopCriterium(double solutionNorm)
private float ComputeStopCriterium(float solutionNorm)
{
// This is criterium 1 from Templates for the solution of linear systems.
// The problem with this criterium is that it's not limiting enough. For now

178
src/Numerics/LinearAlgebra/Complex32/SparseVector.cs

@ -29,7 +29,6 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
using System;
using System.Collections.Generic;
using System.Linq;
using Distributions;
using Generic;
using NumberTheory;
using Numerics;
@ -39,7 +38,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
/// <summary>
/// A vector with sparse storage.
/// </summary>
public class SparseVector : Vector<Complex32>
public class SparseVector : Vector
{
/// <summary>
/// Lock object for the indexer.
@ -581,11 +580,12 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
var sparseother = other as SparseVector;
if (sparseother == null)
{
sparse.AddScaledSparseVector(Complex32.One, sparseother);
base.Add(other, result);
}
else
{
base.Add(other, result);
CopyTo(result);
sparse.AddScaledSparseVector(Complex32.One, sparseother);
}
}
}
@ -753,11 +753,12 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
var sparseother = other as SparseVector;
if (sparseother == null)
{
sparse.AddScaledSparseVector(-Complex32.One, sparseother);
base.Subtract(other, result);
}
else
{
base.Subtract(other, result);
CopyTo(result);
sparse.AddScaledSparseVector(-Complex32.One, sparseother);
}
}
}
@ -1075,9 +1076,9 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
/// Computes the sum of the absolute value of the vector's elements.
/// </summary>
/// <returns>The sum of the absolute value of the vector's elements.</returns>
public override double SumMagnitudes()
public override Complex32 SumMagnitudes()
{
double result = 0;
var result = 0.0f;
for (var i = 0; i < NonZerosCount; i++)
{
result += _nonZeroValues[i].Magnitude;
@ -1153,58 +1154,6 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
return matrix;
}
/// <summary>
/// Generates a vector with random elements
/// </summary>
/// <param name="length">Number of elements in the vector.</param>
/// <param name="randomDistribution">Continuous Random Distribution or Source</param>
/// <returns>
/// A vector with n-random elements distributed according
/// to the specified random distribution.
/// </returns>
/// <exception cref="ArgumentNullException">If the length vector is non positive<see langword="null" />.</exception>
public override Vector<Complex32> Random(int length, IContinuousDistribution randomDistribution)
{
if (length < 1)
{
throw new ArgumentException(Resources.ArgumentMustBePositive, "length");
}
var v = (SparseVector)CreateVector(length);
for (var index = 0; index < v.Count; index++)
{
v[index] = new Complex32((float)randomDistribution.Sample(), (float)randomDistribution.Sample());
}
return v;
}
/// <summary>
/// Generates a vector with random elements
/// </summary>
/// <param name="length">Number of elements in the vector.</param>
/// <param name="randomDistribution">Continuous Random Distribution or Source</param>
/// <returns>
/// A vector with n-random elements distributed according
/// to the specified random distribution.
/// </returns>
/// <exception cref="ArgumentNullException">If the n vector is non positive<see langword="null" />.</exception>
public override Vector<Complex32> Random(int length, IDiscreteDistribution randomDistribution)
{
if (length < 1)
{
throw new ArgumentException(Resources.ArgumentMustBePositive, "length");
}
var v = (SparseVector)CreateVector(length);
for (var index = 0; index < v.Count; index++)
{
v[index] = new Complex32(randomDistribution.Sample(), randomDistribution.Sample());
}
return v;
}
/// <summary>
/// Outer product of this and another vector.
/// </summary>
@ -1227,7 +1176,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
/// </summary>
/// <param name="p">The p value.</param>
/// <returns>Scalar <c>ret = (sum(abs(this[i])^p))^(1/p)</c></returns>
public override double Norm(double p)
public override Complex32 Norm(double p)
{
if (1 > p)
{
@ -1236,7 +1185,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
if (NonZerosCount == 0)
{
return 0.0;
return Complex32.Zero;
}
if (2.0 == p)
@ -1246,7 +1195,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
if (Double.IsPositiveInfinity(p))
{
return CommonParallel.Select(0, NonZerosCount, (index, localData) => Math.Max(localData, _nonZeroValues[index].Magnitude), Math.Max);
return CommonParallel.Select(0, NonZerosCount, (index, localData) => Math.Max(localData, _nonZeroValues[index].Magnitude), Common.Max);
}
var sum = CommonParallel.Aggregate(
@ -1254,36 +1203,9 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
NonZerosCount,
index => Math.Pow(_nonZeroValues[index].Magnitude, p));
return Math.Pow(sum, 1.0 / p);
return (float)Math.Pow(sum, 1.0 / p);
}
/// <summary>
/// Normalizes this vector to a unit vector with respect to the p-norm.
/// </summary>
/// <param name="p">
/// The p value.
/// </param>
/// <returns>
/// This vector normalized to a unit vector with respect to the p-norm.
/// </returns>
public override Vector<Complex32> Normalize(double p)
{
if (p < 0.0)
{
throw new ArgumentOutOfRangeException("p");
}
var norm = Norm(p);
var clone = Clone();
if (norm == 0.0)
{
return clone;
}
clone.Multiply(Complex32.One / (float)norm, clone);
return clone;
}
#endregion
#region Parse Functions
@ -1610,79 +1532,5 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
}
#endregion
/// <summary>
/// Returns the index of the absolute maximum element.
/// </summary>
/// <returns>The index of absolute maximum element.</returns>
public override int MaximumIndex()
{
throw new NotSupportedException();
}
/// <summary>
/// Returns the index of the minimum element.
/// </summary>
/// <returns>The index of minimum element.</returns>
public override int MinimumIndex()
{
throw new NotSupportedException();
}
#region Simple arithmetic of type T
/// <summary>
/// Add two values T+T
/// </summary>
/// <param name="val1">Left operand value</param>
/// <param name="val2">Right operand value</param>
/// <returns>Result of addition</returns>
protected sealed override Complex32 AddT(Complex32 val1, Complex32 val2)
{
return val1 + val2;
}
/// <summary>
/// Subtract two values T-T
/// </summary>
/// <param name="val1">Left operand value</param>
/// <param name="val2">Right operand value</param>
/// <returns>Result of subtract</returns>
protected sealed override Complex32 SubtractT(Complex32 val1, Complex32 val2)
{
return val1 - val2;
}
/// <summary>
/// Multiply two values T*T
/// </summary>
/// <param name="val1">Left operand value</param>
/// <param name="val2">Right operand value</param>
/// <returns>Result of multiplication</returns>
protected sealed override Complex32 MultiplyT(Complex32 val1, Complex32 val2)
{
return val1 * val2;
}
/// <summary>
/// Divide two values T/T
/// </summary>
/// <param name="val1">Left operand value</param>
/// <param name="val2">Right operand value</param>
/// <returns>Result of divide</returns>
protected sealed override Complex32 DivideT(Complex32 val1, Complex32 val2)
{
return val1 / val2;
}
/// <summary>
/// Take absolute value
/// </summary>
/// <param name="val1">Source alue</param>
/// <returns>True if one; otherwise false</returns>
protected sealed override double AbsoluteT(Complex32 val1)
{
return val1.Magnitude;
}
#endregion
}
}

446
src/Numerics/LinearAlgebra/Complex32/Vector.cs

@ -0,0 +1,446 @@
// <copyright file="Vector.cs" company="Math.NET">
// Math.NET Numerics, part of the Math.NET Project
// http://numerics.mathdotnet.com
// http://github.com/mathnet/mathnet-numerics
// http://mathnetnumerics.codeplex.com
// Copyright (c) 2009-2010 Math.NET
// Permission is hereby granted, free of charge, to any person
// obtaining a copy of this software and associated documentation
// files (the "Software"), to deal in the Software without
// restriction, including without limitation the rights to use,
// copy, modify, merge, publish, distribute, sublicense, and/or sell
// copies of the Software, and to permit persons to whom the
// Software is furnished to do so, subject to the following
// conditions:
// The above copyright notice and this permission notice shall be
// included in all copies or substantial portions of the Software.
// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES
// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT
// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY,
// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
// OTHER DEALINGS IN THE SOFTWARE.
// </copyright>
namespace MathNet.Numerics.LinearAlgebra.Complex32
{
using System;
using Distributions;
using Generic;
using Properties;
using Threading;
using Complex32 = Numerics.Complex32;
/// <summary>
/// <c>Complex32</c> version of the <see cref="Vector{T}"/> class.
/// </summary>
public abstract class Vector : Vector<Complex32>
{
/// <summary>
/// Initializes a new instance of the Vector class.
/// Constructs a <strong>Vector</strong> with the given size.
/// </summary>
/// <param name="size">
/// The size of the <strong>Vector</strong> to construct.
/// </param>
/// <exception cref="ArgumentException">
/// If <paramref name="size"/> is less than one.
/// </exception>
protected Vector(int size) : base(size)
{
}
/// <summary>
/// Adds a scalar to each element of the vector and stores the result in the result vector.
/// </summary>
/// <param name="scalar">
/// The scalar to add.
/// </param>
/// <param name="result">
/// The vector to store the result of the addition.
/// </param>
protected override void DoAdd(Complex32 scalar, Vector<Complex32> result)
{
CommonParallel.For(
0,
Count,
index => result[index] = result[index] + scalar);
}
/// <summary>
/// Adds another vector to this vector and stores the result into the result vector.
/// </summary>
/// <param name="other">
/// The vector to add to this one.
/// </param>
/// <param name="result">
/// The vector to store the result of the addition.
/// </param>
protected override void DoAdd(Vector<Complex32> other, Vector<Complex32> result)
{
CommonParallel.For(
0,
Count,
index => result[index] = this[index] + other[index]);
}
/// <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>
protected override void DoSubtract(Complex32 scalar, Vector<Complex32> result)
{
DoAdd(-scalar, result);
}
/// <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>
protected override void DoSubtract(Vector<Complex32> other, Vector<Complex32> result)
{
CopyTo(result);
CommonParallel.For(
0,
Count,
index => result[index] = this[index] - other[index]);
}
/// <summary>
/// Multiplies a scalar to each element of the vector and stores the result in the result vector.
/// </summary>
/// <param name="scalar">
/// The scalar to multiply.
/// </param>
/// <param name="result">
/// The vector to store the result of the multiplication.
/// </param>
protected override void DoMultiply(Complex32 scalar, Vector<Complex32> result)
{
CommonParallel.For(
0,
Count,
index => result[index] = result[index] * scalar);
}
/// <summary>
/// Divides each element of the vector by a scalar and stores the result in the result vector.
/// </summary>
/// <param name="scalar">
/// The scalar to divide with.
/// </param>
/// <param name="result">
/// The vector to store the result of the division.
/// </param>
protected override void DoDivide(Complex32 scalar, Vector<Complex32> result)
{
CommonParallel.For(
0,
Count,
index => result[index] = result[index] / scalar);
}
/// <summary>
/// Pointwise multiplies this vector with another vector and stores the result into the result vector.
/// </summary>
/// <param name="other">The vector to pointwise multiply with this one.</param>
/// <param name="result">The vector to store the result of the pointwise multiplication.</param>
protected override void DoPointwiseMultiply(Vector<Complex32> other, Vector<Complex32> result)
{
CommonParallel.For(
0,
Count,
index => result[index] = this[index] * other[index]);
}
/// <summary>
/// Pointwise divide this vector with another vector and stores the result into the result vector.
/// </summary>
/// <param name="other">The vector to pointwise divide this one by.</param>
/// <param name="result">The vector to store the result of the pointwise division.</param>
protected override void DoPointwiseDivide(Vector<Complex32> other, Vector<Complex32> result)
{
CommonParallel.For(
0,
Count,
index => result[index] = this[index] / other[index]);
}
/// <summary>
/// Computes the dot product between this vector and another vector.
/// </summary>
/// <param name="other">
/// The other vector to add.
/// </param>
/// <returns>s
/// The result of the addition.
/// </returns>
protected override Complex32 DoDotProduct(Vector<Complex32> other)
{
return CommonParallel.Aggregate(
0,
Count,
i => this[i] * other[i]);
}
/// <summary>
/// Returns the value of the absolute minimum element.
/// </summary>
/// <returns>The value of the absolute minimum element.</returns>
public override Complex32 AbsoluteMinimum()
{
return this[AbsoluteMinimumIndex()].Magnitude;
}
/// <summary>
/// Returns the index of the absolute minimum element.
/// </summary>
/// <returns>The index of absolute minimum element.</returns>
public override int AbsoluteMinimumIndex()
{
var index = 0;
var min = this[index].Magnitude;
for (var i = 1; i < Count; i++)
{
var test = this[i].Magnitude;
if (test < min)
{
index = i;
min = test;
}
}
return index;
}
/// <summary>
/// Returns the value of the absolute maximum element.
/// </summary>
/// <returns>The value of the absolute maximum element.</returns>
public override Complex32 AbsoluteMaximum()
{
return this[AbsoluteMaximumIndex()].Magnitude;
}
/// <summary>
/// Returns the index of the absolute maximum element.
/// </summary>
/// <returns>The index of absolute maximum element.</returns>
public override int AbsoluteMaximumIndex()
{
var index = 0;
var max = this[index].Magnitude;
for (var i = 1; i < Count; i++)
{
var test = this[i].Magnitude;
if (test > max)
{
index = i;
max = test;
}
}
return index;
}
/// <summary>
/// Computes the sum of the vector's elements.
/// </summary>
/// <returns>The sum of the vector's elements.</returns>
public override Complex32 Sum()
{
return CommonParallel.Aggregate(
0,
Count,
i => this[i]);
}
/// <summary>
/// Computes the sum of the absolute value of the vector's elements.
/// </summary>
/// <returns>The sum of the absolute value of the vector's elements.</returns>
public override Complex32 SumMagnitudes()
{
return CommonParallel.Aggregate(
0,
Count,
i => this[i].Magnitude);
}
/// <summary>
/// Computes the p-Norm.
/// </summary>
/// <param name="p">
/// The p value.
/// </param>
/// <returns>
/// <c>Scalar ret = (sum(abs(this[i])^p))^(1/p)</c>
/// </returns>
public override Complex32 Norm(double p)
{
if (p < 0.0)
{
throw new ArgumentOutOfRangeException("p");
}
if (double.IsPositiveInfinity(p))
{
return CommonParallel.Select(
0,
Count,
(index, localData) => Math.Max(localData, this[index].Magnitude),
Common.Max);
}
var sum = CommonParallel.Aggregate(
0,
Count,
index => Math.Pow(this[index].Magnitude, p));
return (float)Math.Pow(sum, 1.0 / p);
}
/// <summary>
/// Conjugates vector and save result to <paramref name="target"/>
/// </summary>
/// <param name="target">Target vector</param>
protected override void DoConjugate(Vector<Complex32> target)
{
CopyTo(target);
CommonParallel.For(
0,
Count,
index => target[index] = this[index].Conjugate());
}
/// <summary>
/// Returns a negated vector.
/// </summary>
/// <returns>
/// The negated vector.
/// </returns>
/// <remarks>
/// Added as an alternative to the unary negation operator.
/// </remarks>
public override Vector<Complex32> Negate()
{
var result = CreateVector(Count);
CommonParallel.For(
0,
Count,
index => result[index] = -this[index]);
return result;
}
/// <summary>
/// Returns the index of the absolute maximum element.
/// </summary>
/// <returns>The index of absolute maximum element.</returns>
public override int MaximumIndex()
{
throw new NotSupportedException();
}
/// <summary>
/// Returns the index of the minimum element.
/// </summary>
/// <returns>The index of minimum element.</returns>
public override int MinimumIndex()
{
throw new NotSupportedException();
}
/// <summary>
/// Normalizes this vector to a unit vector with respect to the p-norm.
/// </summary>
/// <param name="p">
/// The p value.
/// </param>
/// <returns>
/// This vector normalized to a unit vector with respect to the p-norm.
/// </returns>
public override Vector<Complex32> Normalize(double p)
{
if (p < 0.0)
{
throw new ArgumentOutOfRangeException("p");
}
var norm = Norm(p);
var clone = Clone();
if (norm.Real == 0.0f)
{
return clone;
}
clone.Multiply(1.0f / norm, clone);
return clone;
}
/// <summary>
/// Generates a vector with random elements
/// </summary>
/// <param name="length">Number of elements in the vector.</param>
/// <param name="randomDistribution">Continuous Random Distribution or Source</param>
/// <returns>
/// A vector with n-random elements distributed according
/// to the specified random distribution.
/// </returns>
/// <exception cref="ArgumentException">If the n vector is non-positive.</exception>
public override Vector<Complex32> Random(int length, IContinuousDistribution randomDistribution)
{
if (length < 1)
{
throw new ArgumentException(Resources.ArgumentMustBePositive, "length");
}
var vector = CreateVector(length);
for (var index = 0; index < length; index++)
{
vector[index] = new Complex32(Convert.ToSingle(randomDistribution.Sample()), Convert.ToSingle(randomDistribution.Sample()));
}
return vector;
}
/// <summary>
/// Generates a vector with random elements
/// </summary>
/// <param name="length">Number of elements in the vector.</param>
/// <param name="randomDistribution">Continuous Random Distribution or Source</param>
/// <returns>
/// A vector with n-random elements distributed according
/// to the specified random distribution.
/// </returns>
/// <exception cref="ArgumentException">If the n vector is not positive.</exception>
public override Vector<Complex32> Random(int length, IDiscreteDistribution randomDistribution)
{
if (length < 1)
{
throw new ArgumentException(Resources.ArgumentMustBePositive, "length");
}
var vector = CreateVector(length);
for (var index = 0; index < length; index++)
{
vector[index] = new Complex32(Convert.ToSingle(randomDistribution.Sample()), Convert.ToSingle(randomDistribution.Sample()));
}
return vector;
}
}
}

161
src/Numerics/LinearAlgebra/Double/DenseVector.cs

@ -39,7 +39,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double
/// <summary>
/// A vector using dense storage.
/// </summary>
public class DenseVector : Vector<double>
public class DenseVector : Vector
{
/// <summary>
/// Initializes a new instance of the <see cref="DenseVector"/> class with a given size.
@ -961,13 +961,10 @@ namespace MathNet.Numerics.LinearAlgebra.Double
/// <returns>The sum of the vector's elements.</returns>
public override double Sum()
{
double result = 0;
for (var i = 0; i < Count; i++)
{
result += Data[i];
}
return result;
return CommonParallel.Aggregate(
0,
Count,
i => Data[i]);
}
/// <summary>
@ -976,13 +973,10 @@ namespace MathNet.Numerics.LinearAlgebra.Double
/// <returns>The sum of the absolute value of the vector's elements.</returns>
public override double SumMagnitudes()
{
double result = 0;
for (var i = 0; i < Count; i++)
{
result += Math.Abs(Data[i]);
}
return result;
return CommonParallel.Aggregate(
0,
Count,
i => Math.Abs(Data[i]));
}
/// <summary>
@ -1193,58 +1187,6 @@ namespace MathNet.Numerics.LinearAlgebra.Double
return matrix;
}
/// <summary>
/// Generates a vector with random elements
/// </summary>
/// <param name="length">Number of elements in the vector.</param>
/// <param name="randomDistribution">Continuous Random Distribution or Source</param>
/// <returns>
/// A vector with n-random elements distributed according
/// to the specified random distribution.
/// </returns>
/// <exception cref="ArgumentNullException">If the n vector is non positive<see langword="null" />.</exception>
public override Vector<double> Random(int length, IContinuousDistribution randomDistribution)
{
if (length < 1)
{
throw new ArgumentException(Resources.ArgumentMustBePositive, "length");
}
var v = (DenseVector)CreateVector(length);
for (var index = 0; index < v.Data.Length; index++)
{
v.Data[index] = randomDistribution.Sample();
}
return v;
}
/// <summary>
/// Generates a vector with random elements
/// </summary>
/// <param name="length">Number of elements in the vector.</param>
/// <param name="randomDistribution">Continuous Random Distribution or Source</param>
/// <returns>
/// A vector with n-random elements distributed according
/// to the specified random distribution.
/// </returns>
/// <exception cref="ArgumentNullException">If the n vector is non positive<see langword="null" />.</exception>
public override Vector<double> Random(int length, IDiscreteDistribution randomDistribution)
{
if (length < 1)
{
throw new ArgumentException(Resources.ArgumentMustBePositive, "length");
}
var v = (DenseVector)CreateVector(length);
for (var index = 0; index < v.Data.Length; index++)
{
v.Data[index] = randomDistribution.Sample();
}
return v;
}
/// <summary>
/// Outer product of this and another vector.
/// </summary>
@ -1301,34 +1243,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double
return Math.Pow(sum, 1.0 / p);
}
/// <summary>
/// Normalizes this vector to a unit vector with respect to the p-norm.
/// </summary>
/// <param name="p">
/// The p value.
/// </param>
/// <returns>
/// This vector normalized to a unit vector with respect to the p-norm.
/// </returns>
public override Vector<double> Normalize(double p)
{
if (p < 0.0)
{
throw new ArgumentOutOfRangeException("p");
}
var norm = Norm(p);
var clone = Clone();
if (norm == 0.0)
{
return clone;
}
clone.Multiply(1.0 / norm, clone);
return clone;
}
#endregion
#region Parse Functions
@ -1497,61 +1412,5 @@ namespace MathNet.Numerics.LinearAlgebra.Double
{
Array.Clear(Data, 0, Data.Length);
}
#region Simple arithmetic of type T
/// <summary>
/// Add two values T+T
/// </summary>
/// <param name="val1">Left operand value</param>
/// <param name="val2">Right operand value</param>
/// <returns>Result of addition</returns>
protected sealed override double AddT(double val1, double val2)
{
return val1 + val2;
}
/// <summary>
/// Subtract two values T-T
/// </summary>
/// <param name="val1">Left operand value</param>
/// <param name="val2">Right operand value</param>
/// <returns>Result of subtract</returns>
protected sealed override double SubtractT(double val1, double val2)
{
return val1 - val2;
}
/// <summary>
/// Multiply two values T*T
/// </summary>
/// <param name="val1">Left operand value</param>
/// <param name="val2">Right operand value</param>
/// <returns>Result of multiplication</returns>
protected sealed override double MultiplyT(double val1, double val2)
{
return val1 * val2;
}
/// <summary>
/// Divide two values T/T
/// </summary>
/// <param name="val1">Left operand value</param>
/// <param name="val2">Right operand value</param>
/// <returns>Result of divide</returns>
protected sealed override double DivideT(double val1, double val2)
{
return val1 / val2;
}
/// <summary>
/// Take absolute value
/// </summary>
/// <param name="val1">Source alue</param>
/// <returns>True if one; otherwise false</returns>
protected sealed override double AbsoluteT(double val1)
{
return Math.Abs(val1);
}
#endregion
}
}

148
src/Numerics/LinearAlgebra/Double/SparseVector.cs

@ -30,7 +30,6 @@ namespace MathNet.Numerics.LinearAlgebra.Double
using System.Collections.Generic;
using System.Globalization;
using System.Linq;
using Distributions;
using Generic;
using NumberTheory;
using Properties;
@ -39,7 +38,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double
/// <summary>
/// A vector with sparse storage.
/// </summary>
public class SparseVector : Vector<double>
public class SparseVector : Vector
{
/// <summary>
/// Lock object for the indexer.
@ -530,11 +529,12 @@ namespace MathNet.Numerics.LinearAlgebra.Double
var sparseother = other as SparseVector;
if (sparseother == null)
{
sparse.AddScaledSparseVector(1.0, sparseother);
base.Add(other, result);
}
else
{
base.Add(other, result);
CopyTo(result);
sparse.AddScaledSparseVector(1.0, sparseother);
}
}
}
@ -702,11 +702,12 @@ namespace MathNet.Numerics.LinearAlgebra.Double
var sparseother = other as SparseVector;
if (sparseother == null)
{
sparse.AddScaledSparseVector(-1.0, sparseother);
base.Subtract(other, result);
}
else
{
base.Subtract(other, result);
CopyTo(result);
sparse.AddScaledSparseVector(-1.0, sparseother);
}
}
}
@ -1152,58 +1153,6 @@ namespace MathNet.Numerics.LinearAlgebra.Double
return matrix;
}
/// <summary>
/// Generates a vector with random elements
/// </summary>
/// <param name="length">Number of elements in the vector.</param>
/// <param name="randomDistribution">Continuous Random Distribution or Source</param>
/// <returns>
/// A vector with n-random elements distributed according
/// to the specified random distribution.
/// </returns>
/// <exception cref="ArgumentNullException">If the length vector is non positive<see langword="null" />.</exception>
public override Vector<double> Random(int length, IContinuousDistribution randomDistribution)
{
if (length < 1)
{
throw new ArgumentException(Resources.ArgumentMustBePositive, "length");
}
var v = (SparseVector)CreateVector(length);
for (var index = 0; index < v.Count; index++)
{
v[index] = randomDistribution.Sample();
}
return v;
}
/// <summary>
/// Generates a vector with random elements
/// </summary>
/// <param name="length">Number of elements in the vector.</param>
/// <param name="randomDistribution">Continuous Random Distribution or Source</param>
/// <returns>
/// A vector with n-random elements distributed according
/// to the specified random distribution.
/// </returns>
/// <exception cref="ArgumentNullException">If the n vector is non positive<see langword="null" />.</exception>
public override Vector<double> Random(int length, IDiscreteDistribution randomDistribution)
{
if (length < 1)
{
throw new ArgumentException(Resources.ArgumentMustBePositive, "length");
}
var v = (SparseVector)CreateVector(length);
for (var index = 0; index < v.Count; index++)
{
v[index] = randomDistribution.Sample();
}
return v;
}
/// <summary>
/// Outer product of this and another vector.
/// </summary>
@ -1256,33 +1205,6 @@ namespace MathNet.Numerics.LinearAlgebra.Double
return Math.Pow(sum, 1.0 / p);
}
/// <summary>
/// Normalizes this vector to a unit vector with respect to the p-norm.
/// </summary>
/// <param name="p">
/// The p value.
/// </param>
/// <returns>
/// This vector normalized to a unit vector with respect to the p-norm.
/// </returns>
public override Vector<double> Normalize(double p)
{
if (p < 0.0)
{
throw new ArgumentOutOfRangeException("p");
}
var norm = Norm(p);
var clone = Clone();
if (norm == 0.0)
{
return clone;
}
clone.Multiply(1.0 / norm, clone);
return clone;
}
#endregion
#region Parse Functions
@ -1623,61 +1545,5 @@ namespace MathNet.Numerics.LinearAlgebra.Double
return true;
}
#region Simple arithmetic of type T
/// <summary>
/// Add two values T+T
/// </summary>
/// <param name="val1">Left operand value</param>
/// <param name="val2">Right operand value</param>
/// <returns>Result of addition</returns>
protected sealed override double AddT(double val1, double val2)
{
return val1 + val2;
}
/// <summary>
/// Subtract two values T-T
/// </summary>
/// <param name="val1">Left operand value</param>
/// <param name="val2">Right operand value</param>
/// <returns>Result of subtract</returns>
protected sealed override double SubtractT(double val1, double val2)
{
return val1 - val2;
}
/// <summary>
/// Multiply two values T*T
/// </summary>
/// <param name="val1">Left operand value</param>
/// <param name="val2">Right operand value</param>
/// <returns>Result of multiplication</returns>
protected sealed override double MultiplyT(double val1, double val2)
{
return val1 * val2;
}
/// <summary>
/// Divide two values T/T
/// </summary>
/// <param name="val1">Left operand value</param>
/// <param name="val2">Right operand value</param>
/// <returns>Result of divide</returns>
protected sealed override double DivideT(double val1, double val2)
{
return val1 / val2;
}
/// <summary>
/// Take absolute value
/// </summary>
/// <param name="val1">Source alue</param>
/// <returns>True if one; otherwise false</returns>
protected sealed override double AbsoluteT(double val1)
{
return Math.Abs(val1);
}
#endregion
}
}

470
src/Numerics/LinearAlgebra/Double/Vector.cs

@ -0,0 +1,470 @@
// <copyright file="Vector.cs" company="Math.NET">
// Math.NET Numerics, part of the Math.NET Project
// http://numerics.mathdotnet.com
// http://github.com/mathnet/mathnet-numerics
// http://mathnetnumerics.codeplex.com
// Copyright (c) 2009-2010 Math.NET
// Permission is hereby granted, free of charge, to any person
// obtaining a copy of this software and associated documentation
// files (the "Software"), to deal in the Software without
// restriction, including without limitation the rights to use,
// copy, modify, merge, publish, distribute, sublicense, and/or sell
// copies of the Software, and to permit persons to whom the
// Software is furnished to do so, subject to the following
// conditions:
// The above copyright notice and this permission notice shall be
// included in all copies or substantial portions of the Software.
// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES
// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT
// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY,
// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
// OTHER DEALINGS IN THE SOFTWARE.
// </copyright>
namespace MathNet.Numerics.LinearAlgebra.Double
{
using System;
using Distributions;
using Generic;
using Properties;
using Threading;
/// <summary>
/// <c>double</c> version of the <see cref="Vector{T}"/> class.
/// </summary>
public abstract class Vector : Vector<double>
{
/// <summary>
/// Initializes a new instance of the Vector class.
/// Constructs a <strong>Vector</strong> with the given size.
/// </summary>
/// <param name="size">
/// The size of the <strong>Vector</strong> to construct.
/// </param>
/// <exception cref="ArgumentException">
/// If <paramref name="size"/> is less than one.
/// </exception>
protected Vector(int size) : base(size)
{
}
/// <summary>
/// Adds a scalar to each element of the vector and stores the result in the result vector.
/// </summary>
/// <param name="scalar">
/// The scalar to add.
/// </param>
/// <param name="result">
/// The vector to store the result of the addition.
/// </param>
protected override void DoAdd(double scalar, Vector<double> result)
{
CommonParallel.For(
0,
Count,
index => result[index] = result[index] + scalar);
}
/// <summary>
/// Adds another vector to this vector and stores the result into the result vector.
/// </summary>
/// <param name="other">
/// The vector to add to this one.
/// </param>
/// <param name="result">
/// The vector to store the result of the addition.
/// </param>
protected override void DoAdd(Vector<double> other, Vector<double> result)
{
CommonParallel.For(
0,
Count,
index => result[index] = this[index] + other[index]);
}
/// <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>
protected override void DoSubtract(double scalar, Vector<double> result)
{
DoAdd(-scalar, result);
}
/// <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>
protected override void DoSubtract(Vector<double> other, Vector<double> result)
{
CopyTo(result);
CommonParallel.For(
0,
Count,
index => result[index] = this[index] - other[index]);
}
/// <summary>
/// Multiplies a scalar to each element of the vector and stores the result in the result vector.
/// </summary>
/// <param name="scalar">
/// The scalar to multiply.
/// </param>
/// <param name="result">
/// The vector to store the result of the multiplication.
/// </param>
protected override void DoMultiply(double scalar, Vector<double> result)
{
CommonParallel.For(
0,
Count,
index => result[index] = result[index] * scalar);
}
/// <summary>
/// Divides each element of the vector by a scalar and stores the result in the result vector.
/// </summary>
/// <param name="scalar">
/// The scalar to divide with.
/// </param>
/// <param name="result">
/// The vector to store the result of the division.
/// </param>
protected override void DoDivide(double scalar, Vector<double> result)
{
CommonParallel.For(
0,
Count,
index => result[index] = result[index] / scalar);
}
/// <summary>
/// Pointwise multiplies this vector with another vector and stores the result into the result vector.
/// </summary>
/// <param name="other">The vector to pointwise multiply with this one.</param>
/// <param name="result">The vector to store the result of the pointwise multiplication.</param>
protected override void DoPointwiseMultiply(Vector<double> other, Vector<double> result)
{
CommonParallel.For(
0,
Count,
index => result[index] = this[index] * other[index]);
}
/// <summary>
/// Pointwise divide this vector with another vector and stores the result into the result vector.
/// </summary>
/// <param name="other">The vector to pointwise divide this one by.</param>
/// <param name="result">The vector to store the result of the pointwise division.</param>
protected override void DoPointwiseDivide(Vector<double> other, Vector<double> result)
{
CommonParallel.For(
0,
Count,
index => result[index] = this[index] / other[index]);
}
/// <summary>
/// Computes the dot product between this vector and another vector.
/// </summary>
/// <param name="other">
/// The other vector to add.
/// </param>
/// <returns>s
/// The result of the addition.
/// </returns>
protected override double DoDotProduct(Vector<double> other)
{
return CommonParallel.Aggregate(
0,
Count,
i => this[i] * other[i]);
}
/// <summary>
/// Returns the value of the absolute minimum element.
/// </summary>
/// <returns>The value of the absolute minimum element.</returns>
public override double AbsoluteMinimum()
{
return Math.Abs(this[AbsoluteMinimumIndex()]);
}
/// <summary>
/// Returns the index of the absolute minimum element.
/// </summary>
/// <returns>The index of absolute minimum element.</returns>
public override int AbsoluteMinimumIndex()
{
var index = 0;
var min = Math.Abs(this[index]);
for (var i = 1; i < Count; i++)
{
var test = Math.Abs(this[i]);
if (test < min)
{
index = i;
min = test;
}
}
return index;
}
/// <summary>
/// Returns the value of the absolute maximum element.
/// </summary>
/// <returns>The value of the absolute maximum element.</returns>
public override double AbsoluteMaximum()
{
return Math.Abs(this[AbsoluteMaximumIndex()]);
}
/// <summary>
/// Returns the index of the absolute maximum element.
/// </summary>
/// <returns>The index of absolute maximum element.</returns>
public override int AbsoluteMaximumIndex()
{
var index = 0;
var max = Math.Abs(this[index]);
for (var i = 1; i < Count; i++)
{
var test = Math.Abs(this[i]);
if (test > max)
{
index = i;
max = test;
}
}
return index;
}
/// <summary>
/// Computes the sum of the vector's elements.
/// </summary>
/// <returns>The sum of the vector's elements.</returns>
public override double Sum()
{
return CommonParallel.Aggregate(
0,
Count,
i => this[i]);
}
/// <summary>
/// Computes the sum of the absolute value of the vector's elements.
/// </summary>
/// <returns>The sum of the absolute value of the vector's elements.</returns>
public override double SumMagnitudes()
{
return CommonParallel.Aggregate(
0,
Count,
i => Math.Abs(this[i]));
}
/// <summary>
/// Computes the p-Norm.
/// </summary>
/// <param name="p">
/// The p value.
/// </param>
/// <returns>
/// <c>Scalar ret = (sum(abs(this[i])^p))^(1/p)</c>
/// </returns>
public override double Norm(double p)
{
if (p < 0.0)
{
throw new ArgumentOutOfRangeException("p");
}
if (Double.IsPositiveInfinity(p))
{
return CommonParallel.Select(
0,
Count,
(index, localData) => Math.Max(localData, Math.Abs(this[index])),
Math.Max);
}
var sum = CommonParallel.Aggregate(
0,
Count,
index => Math.Pow(Math.Abs(this[index]), p));
return Math.Pow(sum, 1.0 / p);
}
/// <summary>
/// Conjugates vector and save result to <paramref name="target"/>
/// </summary>
/// <param name="target">Target vector</param>
protected override void DoConjugate(Vector<double> target)
{
if (ReferenceEquals(this, target))
{
return;
}
CopyTo(target);
}
/// <summary>
/// Returns a negated vector.
/// </summary>
/// <returns>
/// The negated vector.
/// </returns>
/// <remarks>
/// Added as an alternative to the unary negation operator.
/// </remarks>
public override Vector<double> Negate()
{
var result = CreateVector(Count);
CommonParallel.For(
0,
Count,
index => result[index] = -this[index]);
return result;
}
/// <summary>
/// Returns the index of the absolute maximum element.
/// </summary>
/// <returns>The index of absolute maximum element.</returns>
public override int MaximumIndex()
{
var index = 0;
var max = this[index];
for (var i = 1; i < Count; i++)
{
var test = this[i];
if (test > max)
{
index = i;
max = test;
}
}
return index;
}
/// <summary>
/// Returns the index of the minimum element.
/// </summary>
/// <returns>The index of minimum element.</returns>
public override int MinimumIndex()
{
var index = 0;
var min = this[index];
for (var i = 1; i < Count; i++)
{
var test = this[i];
if (test < min)
{
index = i;
min = test;
}
}
return index;
}
/// <summary>
/// Normalizes this vector to a unit vector with respect to the p-norm.
/// </summary>
/// <param name="p">
/// The p value.
/// </param>
/// <returns>
/// This vector normalized to a unit vector with respect to the p-norm.
/// </returns>
public override Vector<double> Normalize(double p)
{
if (p < 0.0)
{
throw new ArgumentOutOfRangeException("p");
}
var norm = Norm(p);
var clone = Clone();
if (norm == 0.0)
{
return clone;
}
clone.Multiply(1.0 / norm, clone);
return clone;
}
/// <summary>
/// Generates a vector with random elements
/// </summary>
/// <param name="length">Number of elements in the vector.</param>
/// <param name="randomDistribution">Continuous Random Distribution or Source</param>
/// <returns>
/// A vector with n-random elements distributed according
/// to the specified random distribution.
/// </returns>
/// <exception cref="ArgumentNullException">If the n vector is non positive<see langword="null" />.</exception>
public override Vector<double> Random(int length, IContinuousDistribution randomDistribution)
{
if (length < 1)
{
throw new ArgumentException(Resources.ArgumentMustBePositive, "length");
}
var v = CreateVector(length);
for (var index = 0; index < Count; index++)
{
v[index] = randomDistribution.Sample();
}
return v;
}
/// <summary>
/// Generates a vector with random elements
/// </summary>
/// <param name="length">Number of elements in the vector.</param>
/// <param name="randomDistribution">Continuous Random Distribution or Source</param>
/// <returns>
/// A vector with n-random elements distributed according
/// to the specified random distribution.
/// </returns>
/// <exception cref="ArgumentNullException">If the n vector is non positive<see langword="null" />.</exception>
public override Vector<double> Random(int length, IDiscreteDistribution randomDistribution)
{
if (length < 1)
{
throw new ArgumentException(Resources.ArgumentMustBePositive, "length");
}
var v = CreateVector(length);
for (var index = 0; index < Count; index++)
{
this[index] = randomDistribution.Sample();
}
return v;
}
}
}

21
src/Numerics/LinearAlgebra/Generic/Common.cs

@ -0,0 +1,21 @@
using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;
namespace MathNet.Numerics.LinearAlgebra.Generic
{
internal static class Common
{
/// <summary>
/// Returns the maximum value.
/// </summary>
/// <param name="a">The first value.</param>
/// <param name="b">The second value.</param>
/// <returns>The maximum value.</returns>
public static float Max(float a, float b)
{
return Math.Max(a, b);
}
}
}

6
src/Numerics/LinearAlgebra/Generic/Matrix.cs

@ -1888,15 +1888,15 @@ namespace MathNet.Numerics.LinearAlgebra.Generic
/// <summary>
/// Take absolute value
/// </summary>
/// <param name="val1">Source alue</param>
/// <returns>True if one; otherwise false</returns>
/// <param name="val1">Source value</param>
/// <returns>True if one; otherwise <c>false</c></returns>
protected abstract double AbsoluteT(T val1);
/// <summary>
/// Is equal to one?
/// </summary>
/// <param name="val1">Value to check</param>
/// <returns>True if one; otherwise false</returns>
/// <returns>True if one; otherwise <c>false</c></returns>
private static bool IsOneT(T val1)
{
if (typeof(T) == typeof(Complex))

416
src/Numerics/LinearAlgebra/Generic/Vector.cs

@ -49,6 +49,16 @@ namespace MathNet.Numerics.LinearAlgebra.Generic
#endif
where T : struct, IEquatable<T>, IFormattable
{
/// <summary>
/// The zero value for type T.
/// </summary>
private static readonly T Zero = default(T);
/// <summary>
/// The value on 1.0 for type T.
/// </summary>
private static readonly T One = SetOne();
/// <summary>
/// Initializes a new instance of the Vector class.
/// Constructs a <strong>Vector</strong> with the given size.
@ -127,14 +137,14 @@ namespace MathNet.Numerics.LinearAlgebra.Generic
/// <returns>A copy of the vector with the scalar added.</returns>
public virtual Vector<T> Add(T scalar)
{
if (scalar.Equals(default(T)))
if (scalar.Equals(Zero))
{
return Clone();
}
var copy = Clone();
Add(scalar, copy);
return copy;
var result = CreateVector(Count);
Add(scalar, result);
return result;
}
/// <summary>
@ -169,12 +179,20 @@ namespace MathNet.Numerics.LinearAlgebra.Generic
CopyTo(result);
}
CommonParallel.For(
0,
Count,
index => result[index] = AddT(result[index], scalar));
DoAdd(scalar, result);
}
/// <summary>
/// Adds a scalar to each element of the vector and stores the result in the result vector.
/// </summary>
/// <param name="scalar">
/// The scalar to add.
/// </param>
/// <param name="result">
/// The vector to store the result of the addition.
/// </param>
protected abstract void DoAdd(T scalar, Vector<T> result);
/// <summary>
/// Returns a copy of this vector.
/// </summary>
@ -214,9 +232,9 @@ namespace MathNet.Numerics.LinearAlgebra.Generic
throw new ArgumentException(Resources.ArgumentVectorsSameLength, "other");
}
var copy = Clone();
Add(other, copy);
return copy;
var result = CreateVector(Count);
Add(other, result);
return result;
}
/// <summary>
@ -259,13 +277,21 @@ namespace MathNet.Numerics.LinearAlgebra.Generic
}
else
{
CommonParallel.For(
0,
Count,
index => result[index] = AddT(this[index], other[index]));
DoAdd(other, result);
}
}
/// <summary>
/// Adds another vector to this vector and stores the result into the result vector.
/// </summary>
/// <param name="other">
/// The vector to add to this one.
/// </param>
/// <param name="result">
/// The vector to store the result of the addition.
/// </param>
protected abstract void DoAdd(Vector<T> other, Vector<T> result);
/// <summary>
/// Subtracts a scalar from each element of the vector.
/// </summary>
@ -280,9 +306,9 @@ namespace MathNet.Numerics.LinearAlgebra.Generic
return Clone();
}
var copy = Clone();
Subtract(scalar, copy);
return copy;
var result = CreateVector(Count);
Subtract(scalar, result);
return result;
}
/// <summary>
@ -317,12 +343,20 @@ namespace MathNet.Numerics.LinearAlgebra.Generic
CopyTo(result);
}
CommonParallel.For(
0,
Count,
index => result[index] = SubtractT(result[index], scalar));
DoSubtract(scalar, result);
}
/// <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>
protected abstract void DoSubtract(T scalar, Vector<T> result);
/// <summary>
/// Returns a negated vector.
/// </summary>
@ -359,9 +393,9 @@ namespace MathNet.Numerics.LinearAlgebra.Generic
throw new ArgumentException(Resources.ArgumentVectorsSameLength, "other");
}
var copy = Clone();
Subtract(other, copy);
return copy;
var result = CreateVector(Count);
Subtract(other, result);
return result;
}
/// <summary>
@ -404,14 +438,21 @@ namespace MathNet.Numerics.LinearAlgebra.Generic
}
else
{
CopyTo(result);
CommonParallel.For(
0,
Count,
index => result[index] = SubtractT(this[index], other[index]));
DoSubtract(other, result);
}
}
/// <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>
protected abstract void DoSubtract(Vector<T> other, Vector<T> result);
/// <summary>
/// Multiplies a scalar to each element of the vector.
/// </summary>
@ -421,14 +462,14 @@ namespace MathNet.Numerics.LinearAlgebra.Generic
/// <returns>A new vector that is the multiplication of the vector and the scalar.</returns>
public virtual Vector<T> Multiply(T scalar)
{
if (IsOneT(scalar))
if (scalar.Equals(One))
{
return Clone();
}
var copy = Clone();
Multiply(scalar, copy);
return copy;
var result = CreateVector(Count);
Multiply(scalar, result);
return result;
}
/// <summary>
@ -463,19 +504,27 @@ namespace MathNet.Numerics.LinearAlgebra.Generic
CopyTo(result);
}
CommonParallel.For(
0,
Count,
index => result[index] = MultiplyT(result[index], scalar));
DoMultiply(scalar, result);
}
/// <summary>
/// Multiplies a scalar to each element of the vector and stores the result in the result vector.
/// </summary>
/// <param name="scalar">
/// The scalar to multiply.
/// </param>
/// <param name="result">
/// The vector to store the result of the multiplication.
/// </param>
protected abstract void DoMultiply(T scalar, Vector<T> result);
/// <summary>
/// Computes the dot product between this vector and another vector.
/// </summary>
/// <param name="other">
/// The other vector to add.
/// </param>
/// <returns>
/// <returns>s
/// The result of the addition.
/// </returns>
/// <exception cref="ArgumentException">
@ -496,15 +545,20 @@ namespace MathNet.Numerics.LinearAlgebra.Generic
throw new ArgumentException(Resources.ArgumentVectorsSameLength, "other");
}
var dot = default(T);
for (var i = 0; i < Count; i++)
{
dot = AddT(dot, MultiplyT(this[i], other[i]));
}
return dot;
return DoDotProduct(other);
}
/// <summary>
/// Computes the dot product between this vector and another vector.
/// </summary>
/// <param name="other">
/// The other vector to add.
/// </param>
/// <returns>s
/// The result of the addition.
/// </returns>
protected abstract T DoDotProduct(Vector<T> other);
/// <summary>
/// Divides each element of the vector by a scalar.
/// </summary>
@ -514,14 +568,14 @@ namespace MathNet.Numerics.LinearAlgebra.Generic
/// <returns>A new vector that is the division of the vector and the scalar.</returns>
public virtual Vector<T> Divide(T scalar)
{
if (IsOneT(scalar))
if (scalar.Equals(One))
{
return Clone();
}
var copy = Clone();
Divide(scalar, copy);
return copy;
var result = CreateVector(Count);
Divide(scalar, result);
return result;
}
/// <summary>
@ -556,12 +610,20 @@ namespace MathNet.Numerics.LinearAlgebra.Generic
CopyTo(result);
}
CommonParallel.For(
0,
Count,
index => result[index] = DivideT(result[index], scalar));
DoDivide(scalar, result);
}
/// <summary>
/// Divides each element of the vector by a scalar and stores the result in the result vector.
/// </summary>
/// <param name="scalar">
/// The scalar to divide with.
/// </param>
/// <param name="result">
/// The vector to store the result of the division.
/// </param>
protected abstract void DoDivide(T scalar, Vector<T> result);
/// <summary>
/// Pointwise multiplies this vector with another vector.
/// </summary>
@ -581,9 +643,9 @@ namespace MathNet.Numerics.LinearAlgebra.Generic
throw new ArgumentException(Resources.ArgumentVectorsSameLength, "other");
}
var copy = Clone();
PointwiseMultiply(other, copy);
return copy;
var result = CreateVector(Count);
PointwiseMultiply(other, result);
return result;
}
/// <summary>
@ -624,13 +686,17 @@ namespace MathNet.Numerics.LinearAlgebra.Generic
}
else
{
CommonParallel.For(
0,
Count,
index => result[index] = MultiplyT(this[index], other[index]));
DoPointwiseMultiply(other, result);
}
}
/// <summary>
/// Pointwise multiplies this vector with another vector and stores the result into the result vector.
/// </summary>
/// <param name="other">The vector to pointwise multiply with this one.</param>
/// <param name="result">The vector to store the result of the pointwise multiplication.</param>
protected abstract void DoPointwiseMultiply(Vector<T> other, Vector<T> result);
/// <summary>
/// Pointwise divide this vector with another vector.
/// </summary>
@ -650,9 +716,9 @@ namespace MathNet.Numerics.LinearAlgebra.Generic
throw new ArgumentException(Resources.ArgumentVectorsSameLength, "other");
}
var copy = Clone();
PointwiseDivide(other, copy);
return copy;
var result = CreateVector(Count);
PointwiseDivide(other, result);
return result;
}
/// <summary>
@ -693,13 +759,17 @@ namespace MathNet.Numerics.LinearAlgebra.Generic
}
else
{
CommonParallel.For(
0,
Count,
index => result[index] = DivideT(this[index], other[index]));
DoPointwiseDivide(other, result);
}
}
/// <summary>
/// Pointwise divide this vector with another vector and stores the result into the result vector.
/// </summary>
/// <param name="other">The vector to pointwise divide this one by.</param>
/// <param name="result">The vector to store the result of the pointwise division.</param>
protected abstract void DoPointwiseDivide(Vector<T> other, Vector<T> result);
/// <summary>
/// Outer product of two vectors
/// </summary>
@ -769,61 +839,25 @@ namespace MathNet.Numerics.LinearAlgebra.Generic
/// Returns the value of the absolute minimum element.
/// </summary>
/// <returns>The value of the absolute minimum element.</returns>
public virtual double AbsoluteMinimum()
{
return AbsoluteT(this[AbsoluteMinimumIndex()]);
}
public abstract T AbsoluteMinimum();
/// <summary>
/// Returns the index of the absolute minimum element.
/// </summary>
/// <returns>The index of absolute minimum element.</returns>
public virtual int AbsoluteMinimumIndex()
{
var index = 0;
var min = AbsoluteT(this[index]);
for (var i = 1; i < Count; i++)
{
var test = AbsoluteT(this[i]);
if (test < min)
{
index = i;
min = test;
}
}
public abstract int AbsoluteMinimumIndex();
return index;
}
/// <summary>
/// Returns the value of the absolute maximum element.
/// </summary>
/// <returns>The value of the absolute maximum element.</returns>
public virtual double AbsoluteMaximum()
{
return AbsoluteT(this[AbsoluteMaximumIndex()]);
}
public abstract T AbsoluteMaximum();
/// <summary>
/// Returns the index of the absolute maximum element.
/// </summary>
/// <returns>The index of absolute maximum element.</returns>
public virtual int AbsoluteMaximumIndex()
{
var index = 0;
var max = AbsoluteT(this[index]);
for (var i = 1; i < Count; i++)
{
var test = AbsoluteT(this[i]);
if (test > max)
{
index = i;
max = test;
}
}
return index;
}
public abstract int AbsoluteMaximumIndex();
/// <summary>
/// Returns the value of maximum element.
@ -859,31 +893,13 @@ namespace MathNet.Numerics.LinearAlgebra.Generic
/// Computes the sum of the vector's elements.
/// </summary>
/// <returns>The sum of the vector's elements.</returns>
public virtual T Sum()
{
var result = default(T);
for (var i = 0; i < Count; i++)
{
result = AddT(result, this[i]);
}
return result;
}
public abstract T Sum();
/// <summary>
/// Computes the sum of the absolute value of the vector's elements.
/// </summary>
/// <returns>The sum of the absolute value of the vector's elements.</returns>
public virtual double SumMagnitudes()
{
double result = 0;
for (var i = 0; i < Count; i++)
{
result += AbsoluteT(this[i]);
}
return result;
}
public abstract T SumMagnitudes();
#endregion
@ -1070,29 +1086,7 @@ namespace MathNet.Numerics.LinearAlgebra.Generic
/// <returns>
/// <c>Scalar ret = (sum(abs(this[i])^p))^(1/p)</c>
/// </returns>
public virtual double Norm(double p)
{
if (p < 0.0)
{
throw new ArgumentOutOfRangeException("p");
}
if (Double.IsPositiveInfinity(p))
{
return CommonParallel.Select(
0,
Count,
(index, localData) => Math.Max(localData, AbsoluteT(this[index])),
Math.Max);
}
var sum = CommonParallel.Aggregate(
0,
Count,
index => Math.Pow(AbsoluteT(this[index]), p));
return Math.Pow(sum, 1.0 / p);
}
public abstract T Norm(double p);
/// <summary>
/// Normalizes this vector to a unit vector with respect to the p-norm.
@ -1124,13 +1118,6 @@ namespace MathNet.Numerics.LinearAlgebra.Generic
/// <param name="target">Target vector</param>
public virtual void Conjugate(Vector<T> target)
{
// In case of real return copy of vector
if (typeof(T) == typeof(double) || (typeof(T) == typeof(float)))
{
CopyTo(target);
return;
}
if (target == null)
{
throw new ArgumentNullException("target");
@ -1141,19 +1128,15 @@ namespace MathNet.Numerics.LinearAlgebra.Generic
throw new ArgumentException(Resources.ArgumentVectorsSameLength, "target");
}
if (ReferenceEquals(this, target))
{
var tmp = CreateVector(Count);
Conjugate(tmp);
tmp.CopyTo(target);
}
CommonParallel.For(
0,
Count,
index => target[index] = ConjugateT(this[index]));
DoConjugate(target);
}
/// <summary>
/// Conjugates vector and save result to <paramref name="target"/>
/// </summary>
/// <param name="target">Target vector</param>
protected abstract void DoConjugate(Vector<T> target);
#region Copying and Conversion
/// <summary>
@ -1604,114 +1587,33 @@ namespace MathNet.Numerics.LinearAlgebra.Generic
CommonParallel.For(0, Count, index => this[index] = default(T));
}
#region Simple arithmetic of type T
/// <summary>
/// Add two values T+T
/// </summary>
/// <param name="val1">Left operand value</param>
/// <param name="val2">Right operand value</param>
/// <returns>Result of addition</returns>
protected abstract T AddT(T val1, T val2);
/// <summary>
/// Subtract two values T-T
/// </summary>
/// <param name="val1">Left operand value</param>
/// <param name="val2">Right operand value</param>
/// <returns>Result of subtract</returns>
protected abstract T SubtractT(T val1, T val2);
/// <summary>
/// Multiply two values T*T
/// </summary>
/// <param name="val1">Left operand value</param>
/// <param name="val2">Right operand value</param>
/// <returns>Result of multiplication</returns>
protected abstract T MultiplyT(T val1, T val2);
/// <summary>
/// Divide two values T/T
/// </summary>
/// <param name="val1">Left operand value</param>
/// <param name="val2">Right operand value</param>
/// <returns>Result of divide</returns>
protected abstract T DivideT(T val1, T val2);
/// <summary>
/// Take absolute value
/// </summary>
/// <param name="val1">Source alue</param>
/// <returns>True if one; otherwise false</returns>
protected abstract double AbsoluteT(T val1);
/// <summary>
/// Is equal to one?
/// </summary>
/// <param name="val1">Value to check</param>
/// <returns>True if one; otherwise false</returns>
private static bool IsOneT(T val1)
{
if (typeof(T) == typeof(Complex))
{
object obj1 = val1;
return Complex.One.AlmostEqual((Complex)obj1);
}
if (typeof(T) == typeof(Complex32))
{
object obj1 = val1;
return Complex32.One.AlmostEqual((Complex32)obj1);
}
if (typeof(T) == typeof(double))
{
object obj1 = val1;
return 1.0.AlmostEqualInDecimalPlaces((double)obj1, 15);
}
if (typeof(T) == typeof(float))
{
object obj1 = val1;
return 1.0f.AlmostEqualInDecimalPlaces((float)obj1, 7);
}
throw new NotSupportedException();
}
/// <summary>
/// Conjugate complex value. In real case the same value is returned
/// Sets the value of <c>1.0</c> for type T.
/// </summary>
/// <param name="val1">Value to conjugate</param>
/// <returns>Conjugated value (complex) or the same (real)</returns>
private static T ConjugateT(T val1)
/// <returns>The value of <c>1.0</c> for type T.</returns>
private static T SetOne()
{
if (typeof(T) == typeof(Complex))
{
object obj = val1;
object conj = Complex.Conjugate((Complex)obj);
return (T)conj;
return (T)(object)Complex.One;
}
if (typeof(T) == typeof(Complex32))
{
object obj = val1;
object conj = ((Complex32)obj).Conjugate();
return (T)conj;
return (T)(object)Complex32.One;
}
if (typeof(T) == typeof(double))
{
return val1;
return (T)(object)1.0;
}
if (typeof(T) == typeof(float))
{
return val1;
return (T)(object)1.0f;
}
throw new NotSupportedException();
}
#endregion
}
}

171
src/Numerics/LinearAlgebra/Single/DenseVector.cs

@ -39,7 +39,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single
/// <summary>
/// A vector using dense storage.
/// </summary>
public class DenseVector : Vector<float>
public class DenseVector : Vector
{
/// <summary>
/// Initializes a new instance of the <see cref="DenseVector"/> class with a given size.
@ -819,7 +819,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single
/// Returns the value of the absolute minimum element.
/// </summary>
/// <returns>The value of the absolute minimum element.</returns>
public override double AbsoluteMinimum()
public override float AbsoluteMinimum()
{
return Math.Abs(Data[AbsoluteMinimumIndex()]);
}
@ -828,7 +828,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single
/// Returns the value of the absolute maximum element.
/// </summary>
/// <returns>The value of the absolute maximum element.</returns>
public override double AbsoluteMaximum()
public override float AbsoluteMaximum()
{
return Math.Abs(Data[AbsoluteMaximumIndex()]);
}
@ -961,28 +961,22 @@ namespace MathNet.Numerics.LinearAlgebra.Single
/// <returns>The sum of the vector's elements.</returns>
public override float Sum()
{
float result = 0;
for (var i = 0; i < Count; i++)
{
result += Data[i];
}
return result;
return CommonParallel.Aggregate(
0,
Count,
i => Data[i]);
}
/// <summary>
/// Computes the sum of the absolute value of the vector's elements.
/// </summary>
/// <returns>The sum of the absolute value of the vector's elements.</returns>
public override double SumMagnitudes()
public override float SumMagnitudes()
{
float result = 0;
for (var i = 0; i < Count; i++)
{
result += Math.Abs(Data[i]);
}
return result;
return CommonParallel.Aggregate(
0,
Count,
i => Math.Abs(Data[i]));
}
/// <summary>
@ -1193,58 +1187,6 @@ namespace MathNet.Numerics.LinearAlgebra.Single
return matrix;
}
/// <summary>
/// Generates a vector with random elements
/// </summary>
/// <param name="length">Number of elements in the vector.</param>
/// <param name="randomDistribution">Continuous Random Distribution or Source</param>
/// <returns>
/// A vector with n-random elements distributed according
/// to the specified random distribution.
/// </returns>
/// <exception cref="ArgumentNullException">If the n vector is non positive<see langword="null" />.</exception>
public override Vector<float> Random(int length, IContinuousDistribution randomDistribution)
{
if (length < 1)
{
throw new ArgumentException(Resources.ArgumentMustBePositive, "length");
}
var v = (DenseVector)CreateVector(length);
for (var index = 0; index < v.Data.Length; index++)
{
v.Data[index] = (float)randomDistribution.Sample();
}
return v;
}
/// <summary>
/// Generates a vector with random elements
/// </summary>
/// <param name="length">Number of elements in the vector.</param>
/// <param name="randomDistribution">Continuous Random Distribution or Source</param>
/// <returns>
/// A vector with n-random elements distributed according
/// to the specified random distribution.
/// </returns>
/// <exception cref="ArgumentNullException">If the n vector is non positive<see langword="null" />.</exception>
public override Vector<float> Random(int length, IDiscreteDistribution randomDistribution)
{
if (length < 1)
{
throw new ArgumentException(Resources.ArgumentMustBePositive, "length");
}
var v = (DenseVector)CreateVector(length);
for (var index = 0; index < v.Data.Length; index++)
{
v.Data[index] = randomDistribution.Sample();
}
return v;
}
/// <summary>
/// Outer product of this and another vector.
/// </summary>
@ -1265,7 +1207,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single
/// </summary>
/// <param name="p">The p value.</param>
/// <returns>Scalar <c>ret = (sum(abs(this[i])^p))^(1/p)</c></returns>
public override double Norm(double p)
public override float Norm(double p)
{
if (p < 0.0)
{
@ -1291,7 +1233,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single
0,
Count,
(index, localData) => Math.Max(localData, Math.Abs(Data[index])),
Math.Max);
Common.Max);
}
var sum = CommonParallel.Aggregate(
@ -1299,36 +1241,9 @@ namespace MathNet.Numerics.LinearAlgebra.Single
Count,
index => Math.Pow(Math.Abs(Data[index]), p));
return Math.Pow(sum, 1.0 / p);
return (float)Math.Pow(sum, 1.0 / p);
}
/// <summary>
/// Normalizes this vector to a unit vector with respect to the p-norm.
/// </summary>
/// <param name="p">
/// The p value.
/// </param>
/// <returns>
/// This vector normalized to a unit vector with respect to the p-norm.
/// </returns>
public override Vector<float> Normalize(double p)
{
if (p < 0.0)
{
throw new ArgumentOutOfRangeException("p");
}
var norm = Norm(p);
var clone = Clone();
if (norm == 0.0)
{
return clone;
}
clone.Multiply(1.0f / (float)norm, clone);
return clone;
}
#endregion
#region Parse Functions
@ -1497,61 +1412,5 @@ namespace MathNet.Numerics.LinearAlgebra.Single
{
Array.Clear(Data, 0, Data.Length);
}
#region Simple arithmetic of type T
/// <summary>
/// Add two values T+T
/// </summary>
/// <param name="val1">Left operand value</param>
/// <param name="val2">Right operand value</param>
/// <returns>Result of addition</returns>
protected sealed override float AddT(float val1, float val2)
{
return val1 + val2;
}
/// <summary>
/// Subtract two values T-T
/// </summary>
/// <param name="val1">Left operand value</param>
/// <param name="val2">Right operand value</param>
/// <returns>Result of subtract</returns>
protected sealed override float SubtractT(float val1, float val2)
{
return val1 - val2;
}
/// <summary>
/// Multiply two values T*T
/// </summary>
/// <param name="val1">Left operand value</param>
/// <param name="val2">Right operand value</param>
/// <returns>Result of multiplication</returns>
protected sealed override float MultiplyT(float val1, float val2)
{
return val1 * val2;
}
/// <summary>
/// Divide two values T/T
/// </summary>
/// <param name="val1">Left operand value</param>
/// <param name="val2">Right operand value</param>
/// <returns>Result of divide</returns>
protected sealed override float DivideT(float val1, float val2)
{
return val1 / val2;
}
/// <summary>
/// Take absolute value
/// </summary>
/// <param name="val1">Source alue</param>
/// <returns>True if one; otherwise false</returns>
protected sealed override double AbsoluteT(float val1)
{
return Math.Abs(val1);
}
#endregion
}
}

18
src/Numerics/LinearAlgebra/Single/Solvers/StopCriterium/ResidualStopCriterium.cs

@ -45,7 +45,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Solvers.StopCriterium
/// <summary>
/// The default value for the maximum value of the residual.
/// </summary>
public const double DefaultMaximumResidual = 1e-6;
public const float DefaultMaximumResidual = 1e-6f;
/// <summary>
/// The default value for the minimum number of iterations.
@ -65,7 +65,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Solvers.StopCriterium
/// <summary>
/// The maximum value for the residual below which the calculation is considered converged.
/// </summary>
private double _maximum;
private float _maximum;
/// <summary>
/// The minimum number of iterations for which the residual has to be below the maximum before
@ -101,7 +101,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Solvers.StopCriterium
/// maximum residual and the default minimum number of iterations.
/// </summary>
/// <param name="maximum">The maximum value for the residual below which the calculation is considered converged.</param>
public ResidualStopCriterium(double maximum) : this(maximum, DefaultMinimumIterationsBelowMaximum)
public ResidualStopCriterium(float maximum) : this(maximum, DefaultMinimumIterationsBelowMaximum)
{
}
@ -128,7 +128,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Solvers.StopCriterium
/// The minimum number of iterations for which the residual has to be below the maximum before
/// the calculation is considered converged.
/// </param>
public ResidualStopCriterium(double maximum, int minimumIterationsBelowMaximum)
public ResidualStopCriterium(float maximum, int minimumIterationsBelowMaximum)
{
if (maximum < 0)
{
@ -149,7 +149,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Solvers.StopCriterium
/// converged.
/// </summary>
/// <exception cref="ArgumentOutOfRangeException">Thrown if the <c>Maximum</c> is set to a negative value.</exception>
public double Maximum
public float Maximum
{
[DebuggerStepThrough]
get
@ -258,16 +258,16 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Solvers.StopCriterium
// Store the infinity norms of both the solution and residual vectors
// These values will be used to calculate the relative drop in residuals
// later on.
var residualNorm = residualVector.Norm(Double.PositiveInfinity);
var residualNorm = residualVector.Norm(float.PositiveInfinity);
// Check the residuals by calculating:
// ||r_i|| <= stop_tol * ||b||
var stopCriterium = ComputeStopCriterium(sourceVector.Norm(Double.PositiveInfinity));
var stopCriterium = ComputeStopCriterium(sourceVector.Norm(float.PositiveInfinity));
// First check that we have real numbers not NaN's.
// NaN's can occur when the iterative process diverges so we
// stop if that is the case.
if (double.IsNaN(stopCriterium) || double.IsNaN(residualNorm))
if (float.IsNaN(stopCriterium) || float.IsNaN(residualNorm))
{
_iterationCount = 0;
SetStatusToDiverged();
@ -306,7 +306,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Solvers.StopCriterium
/// </summary>
/// <param name="solutionNorm">Solution vector norm</param>
/// <returns>Criterium value</returns>
private double ComputeStopCriterium(double solutionNorm)
private float ComputeStopCriterium(float solutionNorm)
{
// This is criterium 1 from Templates for the solution of linear systems.
// The problem with this criterium is that it's not limiting enough. For now

162
src/Numerics/LinearAlgebra/Single/SparseVector.cs

@ -30,7 +30,6 @@ namespace MathNet.Numerics.LinearAlgebra.Single
using System.Collections.Generic;
using System.Globalization;
using System.Linq;
using Distributions;
using Generic;
using NumberTheory;
using Properties;
@ -39,7 +38,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single
/// <summary>
/// A vector with sparse storage.
/// </summary>
public class SparseVector : Vector<float>
public class SparseVector : Vector
{
/// <summary>
/// Lock object for the indexer.
@ -530,11 +529,12 @@ namespace MathNet.Numerics.LinearAlgebra.Single
var sparseother = other as SparseVector;
if (sparseother == null)
{
sparse.AddScaledSparseVector(1.0f, sparseother);
base.Add(other, result);
}
else
{
base.Add(other, result);
CopyTo(result);
sparse.AddScaledSparseVector(1.0f, sparseother);
}
}
}
@ -702,11 +702,12 @@ namespace MathNet.Numerics.LinearAlgebra.Single
var sparseother = other as SparseVector;
if (sparseother == null)
{
sparse.AddScaledSparseVector(-1.0f, sparseother);
base.Subtract(other, result);
}
else
{
base.Subtract(other, result);
CopyTo(result);
sparse.AddScaledSparseVector(-1.0f, sparseother);
}
}
}
@ -1061,7 +1062,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single
/// <returns>The sum of the vector's elements.</returns>
public override float Sum()
{
float result = 0;
var result = 0.0f;
for (var i = 0; i < NonZerosCount; i++)
{
result += _nonZeroValues[i];
@ -1074,9 +1075,9 @@ namespace MathNet.Numerics.LinearAlgebra.Single
/// Computes the sum of the absolute value of the vector's elements.
/// </summary>
/// <returns>The sum of the absolute value of the vector's elements.</returns>
public override double SumMagnitudes()
public override float SumMagnitudes()
{
double result = 0;
var result = 0.0f;
for (var i = 0; i < NonZerosCount; i++)
{
result += Math.Abs(_nonZeroValues[i]);
@ -1152,58 +1153,6 @@ namespace MathNet.Numerics.LinearAlgebra.Single
return matrix;
}
/// <summary>
/// Generates a vector with random elements
/// </summary>
/// <param name="length">Number of elements in the vector.</param>
/// <param name="randomDistribution">Continuous Random Distribution or Source</param>
/// <returns>
/// A vector with n-random elements distributed according
/// to the specified random distribution.
/// </returns>
/// <exception cref="ArgumentNullException">If the length vector is non positive<see langword="null" />.</exception>
public override Vector<float> Random(int length, IContinuousDistribution randomDistribution)
{
if (length < 1)
{
throw new ArgumentException(Resources.ArgumentMustBePositive, "length");
}
var v = (SparseVector)CreateVector(length);
for (var index = 0; index < v.Count; index++)
{
v[index] = (float)randomDistribution.Sample();
}
return v;
}
/// <summary>
/// Generates a vector with random elements
/// </summary>
/// <param name="length">Number of elements in the vector.</param>
/// <param name="randomDistribution">Continuous Random Distribution or Source</param>
/// <returns>
/// A vector with n-random elements distributed according
/// to the specified random distribution.
/// </returns>
/// <exception cref="ArgumentNullException">If the n vector is non positive<see langword="null" />.</exception>
public override Vector<float> Random(int length, IDiscreteDistribution randomDistribution)
{
if (length < 1)
{
throw new ArgumentException(Resources.ArgumentMustBePositive, "length");
}
var v = (SparseVector)CreateVector(length);
for (var index = 0; index < v.Count; index++)
{
v[index] = randomDistribution.Sample();
}
return v;
}
/// <summary>
/// Outer product of this and another vector.
/// </summary>
@ -1226,7 +1175,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single
/// </summary>
/// <param name="p">The p value.</param>
/// <returns>Scalar <c>ret = (sum(abs(this[i])^p))^(1/p)</c></returns>
public override double Norm(double p)
public override float Norm(double p)
{
if (1 > p)
{
@ -1235,7 +1184,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single
if (NonZerosCount == 0)
{
return 0.0;
return 0.0f;
}
if (2.0 == p)
@ -1245,7 +1194,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single
if (Double.IsPositiveInfinity(p))
{
return CommonParallel.Select(0, NonZerosCount, (index, localData) => Math.Max(localData, Math.Abs(_nonZeroValues[index])), Math.Max);
return CommonParallel.Select(0, NonZerosCount, (index, localData) => Math.Max(localData, Math.Abs(_nonZeroValues[index])), Common.Max);
}
var sum = CommonParallel.Aggregate(
@ -1253,36 +1202,9 @@ namespace MathNet.Numerics.LinearAlgebra.Single
NonZerosCount,
index => Math.Pow(Math.Abs(_nonZeroValues[index]), p));
return Math.Pow(sum, 1.0 / p);
return (float)Math.Pow(sum, 1.0 / p);
}
/// <summary>
/// Normalizes this vector to a unit vector with respect to the p-norm.
/// </summary>
/// <param name="p">
/// The p value.
/// </param>
/// <returns>
/// This vector normalized to a unit vector with respect to the p-norm.
/// </returns>
public override Vector<float> Normalize(double p)
{
if (p < 0.0)
{
throw new ArgumentOutOfRangeException("p");
}
var norm = Norm(p);
var clone = Clone();
if (norm == 0.0)
{
return clone;
}
clone.Multiply(1.0f / (float)norm, clone);
return clone;
}
#endregion
#region Parse Functions
@ -1623,61 +1545,5 @@ namespace MathNet.Numerics.LinearAlgebra.Single
return true;
}
#region Simple arithmetic of type T
/// <summary>
/// Add two values T+T
/// </summary>
/// <param name="val1">Left operand value</param>
/// <param name="val2">Right operand value</param>
/// <returns>Result of addition</returns>
protected sealed override float AddT(float val1, float val2)
{
return val1 + val2;
}
/// <summary>
/// Subtract two values T-T
/// </summary>
/// <param name="val1">Left operand value</param>
/// <param name="val2">Right operand value</param>
/// <returns>Result of subtract</returns>
protected sealed override float SubtractT(float val1, float val2)
{
return val1 - val2;
}
/// <summary>
/// Multiply two values T*T
/// </summary>
/// <param name="val1">Left operand value</param>
/// <param name="val2">Right operand value</param>
/// <returns>Result of multiplication</returns>
protected sealed override float MultiplyT(float val1, float val2)
{
return val1 * val2;
}
/// <summary>
/// Divide two values T/T
/// </summary>
/// <param name="val1">Left operand value</param>
/// <param name="val2">Right operand value</param>
/// <returns>Result of divide</returns>
protected sealed override float DivideT(float val1, float val2)
{
return val1 / val2;
}
/// <summary>
/// Take absolute value
/// </summary>
/// <param name="val1">Source alue</param>
/// <returns>True if one; otherwise false</returns>
protected sealed override double AbsoluteT(float val1)
{
return Math.Abs(val1);
}
#endregion
}
}

473
src/Numerics/LinearAlgebra/Single/Vector.cs

@ -0,0 +1,473 @@
// <copyright file="Vector.cs" company="Math.NET">
// Math.NET Numerics, part of the Math.NET Project
// http://numerics.mathdotnet.com
// http://github.com/mathnet/mathnet-numerics
// http://mathnetnumerics.codeplex.com
// Copyright (c) 2009-2010 Math.NET
// Permission is hereby granted, free of charge, to any person
// obtaining a copy of this software and associated documentation
// files (the "Software"), to deal in the Software without
// restriction, including without limitation the rights to use,
// copy, modify, merge, publish, distribute, sublicense, and/or sell
// copies of the Software, and to permit persons to whom the
// Software is furnished to do so, subject to the following
// conditions:
// The above copyright notice and this permission notice shall be
// included in all copies or substantial portions of the Software.
// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES
// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT
// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY,
// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
// OTHER DEALINGS IN THE SOFTWARE.
// </copyright>
namespace MathNet.Numerics.LinearAlgebra.Single
{
using System;
using Distributions;
using Generic;
using Properties;
using Threading;
/// <summary>
/// <c>float</c> version of the <see cref="Vector{T}"/> class.
/// </summary>
public abstract class Vector : Vector<float>
{
/// <summary>
/// Initializes a new instance of the Vector class.
/// Constructs a <strong>Vector</strong> with the given size.
/// </summary>
/// <param name="size">
/// The size of the <strong>Vector</strong> to construct.
/// </param>
/// <exception cref="ArgumentException">
/// If <paramref name="size"/> is less than one.
/// </exception>
protected Vector(int size)
: base(size)
{
}
/// <summary>
/// Adds a scalar to each element of the vector and stores the result in the result vector.
/// </summary>
/// <param name="scalar">
/// The scalar to add.
/// </param>
/// <param name="result">
/// The vector to store the result of the addition.
/// </param>
protected override void DoAdd(float scalar, Vector<float> result)
{
CommonParallel.For(
0,
Count,
index => result[index] = result[index] + scalar);
}
/// <summary>
/// Adds another vector to this vector and stores the result into the result vector.
/// </summary>
/// <param name="other">
/// The vector to add to this one.
/// </param>
/// <param name="result">
/// The vector to store the result of the addition.
/// </param>
protected override void DoAdd(Vector<float> other, Vector<float> result)
{
CommonParallel.For(
0,
Count,
index => result[index] = this[index] + other[index]);
}
/// <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>
protected override void DoSubtract(float scalar, Vector<float> result)
{
DoAdd(-scalar, result);
}
/// <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>
protected override void DoSubtract(Vector<float> other, Vector<float> result)
{
CopyTo(result);
CommonParallel.For(
0,
Count,
index => result[index] = this[index] - other[index]);
}
/// <summary>
/// Multiplies a scalar to each element of the vector and stores the result in the result vector.
/// </summary>
/// <param name="scalar">
/// The scalar to multiply.
/// </param>
/// <param name="result">
/// The vector to store the result of the multiplication.
/// </param>
protected override void DoMultiply(float scalar, Vector<float> result)
{
CommonParallel.For(
0,
Count,
index => result[index] = result[index] * scalar);
}
/// <summary>
/// Divides each element of the vector by a scalar and stores the result in the result vector.
/// </summary>
/// <param name="scalar">
/// The scalar to divide with.
/// </param>
/// <param name="result">
/// The vector to store the result of the division.
/// </param>
protected override void DoDivide(float scalar, Vector<float> result)
{
CommonParallel.For(
0,
Count,
index => result[index] = result[index] / scalar);
}
/// <summary>
/// Pointwise multiplies this vector with another vector and stores the result into the result vector.
/// </summary>
/// <param name="other">The vector to pointwise multiply with this one.</param>
/// <param name="result">The vector to store the result of the pointwise multiplication.</param>
protected override void DoPointwiseMultiply(Vector<float> other, Vector<float> result)
{
CommonParallel.For(
0,
Count,
index => result[index] = this[index] * other[index]);
}
/// <summary>
/// Pointwise divide this vector with another vector and stores the result into the result vector.
/// </summary>
/// <param name="other">The vector to pointwise divide this one by.</param>
/// <param name="result">The vector to store the result of the pointwise division.</param>
protected override void DoPointwiseDivide(Vector<float> other, Vector<float> result)
{
CommonParallel.For(
0,
Count,
index => result[index] = this[index] / other[index]);
}
/// <summary>
/// Computes the dot product between this vector and another vector.
/// </summary>
/// <param name="other">
/// The other vector to add.
/// </param>
/// <returns>s
/// The result of the addition.
/// </returns>
protected override float DoDotProduct(Vector<float> other)
{
return CommonParallel.Aggregate(
0,
Count,
i => this[i] * other[i]);
}
/// <summary>
/// Returns the value of the absolute minimum element.
/// </summary>
/// <returns>The value of the absolute minimum element.</returns>
public override float AbsoluteMinimum()
{
return Math.Abs(this[AbsoluteMinimumIndex()]);
}
/// <summary>
/// Returns the index of the absolute minimum element.
/// </summary>
/// <returns>The index of absolute minimum element.</returns>
public override int AbsoluteMinimumIndex()
{
var index = 0;
var min = Math.Abs(this[index]);
for (var i = 1; i < Count; i++)
{
var test = Math.Abs(this[i]);
if (test < min)
{
index = i;
min = test;
}
}
return index;
}
/// <summary>
/// Returns the value of the absolute maximum element.
/// </summary>
/// <returns>The value of the absolute maximum element.</returns>
public override float AbsoluteMaximum()
{
return Math.Abs(this[AbsoluteMaximumIndex()]);
}
/// <summary>
/// Returns the index of the absolute maximum element.
/// </summary>
/// <returns>The index of absolute maximum element.</returns>
public override int AbsoluteMaximumIndex()
{
var index = 0;
var max = Math.Abs(this[index]);
for (var i = 1; i < Count; i++)
{
var test = Math.Abs(this[i]);
if (test > max)
{
index = i;
max = test;
}
}
return index;
}
/// <summary>
/// Computes the sum of the vector's elements.
/// </summary>
/// <returns>The sum of the vector's elements.</returns>
public override float Sum()
{
return CommonParallel.Aggregate(
0,
Count,
i => this[i]);
}
/// <summary>
/// Computes the sum of the absolute value of the vector's elements.
/// </summary>
/// <returns>The sum of the absolute value of the vector's elements.</returns>
public override float SumMagnitudes()
{
return CommonParallel.Aggregate(
0,
Count,
i => Math.Abs(this[i]));
}
/// <summary>
/// Computes the p-Norm.
/// </summary>
/// <param name="p">
/// The p value.
/// </param>
/// <returns>
/// <c>Scalar ret = (sum(abs(this[i])^p))^(1/p)</c>
/// </returns>
public override float Norm(double p)
{
if (p < 0.0)
{
throw new ArgumentOutOfRangeException("p");
}
if (float.IsPositiveInfinity((float)p))
{
return CommonParallel.Select(
0,
Count,
(index, localData) => Math.Max(localData, Math.Abs(this[index])),
Common.Max);
}
var sum = CommonParallel.Aggregate(
0,
Count,
index => Math.Pow(Math.Abs(this[index]), p));
return (float)Math.Pow(sum, 1.0 / p);
}
/// <summary>
/// Conjugates vector and save result to <paramref name="target"/>
/// </summary>
/// <param name="target">Target vector</param>
protected override void DoConjugate(Vector<float> target)
{
if (ReferenceEquals(this, target))
{
return;
}
CopyTo(target);
}
/// <summary>
/// Returns a negated vector.
/// </summary>
/// <returns>
/// The negated vector.
/// </returns>
/// <remarks>
/// Added as an alternative to the unary negation operator.
/// </remarks>
public override Vector<float> Negate()
{
var result = CreateVector(Count);
CommonParallel.For(
0,
Count,
index => result[index] = -this[index]);
return result;
}
/// <summary>
/// Returns the index of the absolute maximum element.
/// </summary>
/// <returns>The index of absolute maximum element.</returns>
public override int MaximumIndex()
{
var index = 0;
var max = this[index];
for (var i = 1; i < Count; i++)
{
var test = this[i];
if (test > max)
{
index = i;
max = test;
}
}
return index;
}
/// <summary>
/// Returns the index of the minimum element.
/// </summary>
/// <returns>The index of minimum element.</returns>
public override int MinimumIndex()
{
var index = 0;
var min = this[index];
for (var i = 1; i < Count; i++)
{
var test = this[i];
if (test < min)
{
index = i;
min = test;
}
}
return index;
}
/// <summary>
/// Normalizes this vector to a unit vector with respect to the p-norm.
/// </summary>
/// <param name="p">
/// The p value.
/// </param>
/// <returns>
/// This vector normalized to a unit vector with respect to the p-norm.
/// </returns>
public override Vector<float> Normalize(double p)
{
if (p < 0.0)
{
throw new ArgumentOutOfRangeException("p");
}
var norm = Norm(p);
var clone = Clone();
if (norm == 0.0)
{
return clone;
}
clone.Multiply(1.0f / norm, clone);
return clone;
}
/// <summary>
/// Generates a vector with random elements
/// </summary>
/// <param name="length">Number of elements in the vector.</param>
/// <param name="randomDistribution">Continuous Random Distribution or Source</param>
/// <returns>
/// A vector with n-random elements distributed according
/// to the specified random distribution.
/// </returns>
/// <exception cref="ArgumentNullException">If the n vector is non positive<see langword="null" />.</exception>
public override Vector<float> Random(int length, IContinuousDistribution randomDistribution)
{
if (length < 1)
{
throw new ArgumentException(Resources.ArgumentMustBePositive, "length");
}
var v = CreateVector(length);
for (var index = 0; index < Count; index++)
{
v[index] = Convert.ToSingle(randomDistribution.Sample());
}
return v;
}
/// <summary>
/// Generates a vector with random elements
/// </summary>
/// <param name="length">Number of elements in the vector.</param>
/// <param name="randomDistribution">Continuous Random Distribution or Source</param>
/// <returns>
/// A vector with n-random elements distributed according
/// to the specified random distribution.
/// </returns>
/// <exception cref="ArgumentNullException">If the n vector is non positive<see langword="null" />.</exception>
public override Vector<float> Random(int length, IDiscreteDistribution randomDistribution)
{
if (length < 1)
{
throw new ArgumentException(Resources.ArgumentMustBePositive, "length");
}
var v = CreateVector(length);
for (var index = 0; index < Count; index++)
{
v[index] = Convert.ToSingle(randomDistribution.Sample());
}
return v;
}
}
}

5
src/Numerics/Numerics.csproj

@ -126,9 +126,13 @@
<Compile Include="Distributions\Multivariate\Wishart.cs" />
<Compile Include="LinearAlgebra\Complex32\IO\DelimitedWriter.cs" />
<Compile Include="LinearAlgebra\Complex32\IO\MatlabReader.cs" />
<Compile Include="LinearAlgebra\Complex32\Vector.cs" />
<Compile Include="LinearAlgebra\Complex\IO\DelimitedWriter.cs" />
<Compile Include="LinearAlgebra\Complex\IO\MatlabReader.cs" />
<Compile Include="LinearAlgebra\Complex\Vector.cs" />
<Compile Include="LinearAlgebra\Double\IO\DelimitedWriter.cs" />
<Compile Include="LinearAlgebra\Double\Vector.cs" />
<Compile Include="LinearAlgebra\Generic\Common.cs" />
<Compile Include="LinearAlgebra\IO\MatlabReader.cs" />
<Compile Include="LinearAlgebra\IO\Matlab\ArrayClass.cs" />
<Compile Include="LinearAlgebra\IO\Matlab\ArrayFlags.cs" />
@ -305,6 +309,7 @@
<Compile Include="LinearAlgebra\Generic\Matrix.Arithmetic.cs" />
<Compile Include="LinearAlgebra\Generic\Matrix.cs" />
<Compile Include="LinearAlgebra\Generic\Vector.cs" />
<Compile Include="LinearAlgebra\Single\Vector.cs" />
<Compile Include="Permutation.cs" />
<Compile Include="Distributions\Continuous\Beta.cs" />
<Compile Include="Distributions\Continuous\ContinuousUniform.cs" />

64
src/Numerics/Precision.cs

@ -200,6 +200,42 @@ namespace MathNet.Numerics
#endif
}
/// <summary>
/// Returns the magnitude of the number.
/// </summary>
/// <param name="value">The value.</param>
/// <returns>The magnitude of the number.</returns>
public static int Magnitude(this float value)
{
// Can't do this with zero because the 10-log of zero doesn't exist.
if (value.Equals(0.0f))
{
return 0;
}
// Note that we need the absolute value of the input because Log10 doesn't
// work for negative numbers (obviously).
var magnitude = Convert.ToSingle(Math.Log10(Math.Abs(value)));
// To get the right number we need to know if the value is negative or positive
// truncating a positive number will always give use the correct magnitude
// truncating a negative number will give us a magnitude that is off by 1
if (magnitude < 0)
{
#if SILVERLIGHT
return (int)Truncate(magnitude - 1);
#else
return (int)Math.Truncate(magnitude - 1);
#endif
}
#if SILVERLIGHT
return (int)Truncate(magnitude);
#else
return (int)Math.Truncate(magnitude);
#endif
}
/// <summary>
/// Returns the number divided by it's magnitude, effectively returning a number between -10 and 10.
/// </summary>
@ -1582,6 +1618,34 @@ namespace MathNet.Numerics
return CompareToInDecimalPlaces(a, b, decimalPlaces) < 0;
}
///<summary>
/// Compares two floats and determines if the <c>first</c> value is smaller than the <c>second</c>
/// value to within the specified number of decimal places or not.
/// </summary>
/// <remarks>
/// <para>
/// The values are equal if the difference between the two numbers is smaller than 10^(-numberOfDecimalPlaces). We divide by
/// two so that we have half the range on each side of th<paramref name="decimalPlaces"/>g. if <paramref name="decimalPlaces"/> == 2, then 0.01 will equal between
/// 0.005 and 0.015, but not 0.02 and not 0.00
/// </para>
/// </remarks>
/// <param name="a">The first value.</param>
/// <param name="b">The second value.</param>
/// <param name="decimalPlaces">The number of decimal places.</param>
/// <returns><c>true</c> if the first value is smaller than the second value; otherwise <c>false</c>.</returns>
public static bool IsSmallerWithDecimalPlaces(this float a, float b, int decimalPlaces)
{
// If A or B are a NAN, return false. NANs are equal to nothing,
// not even themselves, and thus they're not bigger or
// smaller than anything either
if (double.IsNaN(a) || double.IsNaN(b))
{
return false;
}
return CompareToInDecimalPlaces(a, b, decimalPlaces) < 0;
}
/// <summary>
/// Compares two doubles and determines which double is bigger.

52
src/Numerics/Threading/CommonParallel.cs

@ -385,5 +385,57 @@ namespace MathNet.Numerics.Threading
return ret;
}
/// <summary>
/// Selects an item (such as Max or Min).
/// </summary>
/// <param name="fromInclusive">Starting index of the loop.</param>
/// <param name="toExclusive">Ending index of the loop</param>
/// <param name="body">The function to select items over a subset.</param>
/// <param name="localFinally">The function to select the item of selection from the subsets.</param>
/// <returns>The selected value.</returns>
public static float Select(int fromInclusive, int toExclusive, Func<int, float, float> body, Func<float, float, float> localFinally)
{
float ret = 0;
var syncLock = new object();
#if SILVERLIGHT
Parallel.For(
fromInclusive,
toExclusive,
() => 0.0f,
(i, localData) => localData += body(i, localData),
localResult =>
{
lock (syncLock)
{
ret = localFinally(ret, localResult);
}
});
#else
Parallel.ForEach(
Partitioner.Create(fromInclusive, toExclusive),
new ParallelOptions { MaxDegreeOfParallelism = Control.NumberOfParallelWorkerThreads },
() => 0.0f,
(range, loop, localData) =>
{
for (var i = range.Item1; i < range.Item2; i++)
{
localData = body(i, localData);
}
return localData;
},
localResult =>
{
lock (syncLock)
{
ret = localFinally(ret, localResult);
}
});
#endif
return ret;
}
}
}

15
src/Silverlight/Silverlight.csproj

@ -368,6 +368,9 @@
<Compile Include="..\Numerics\LinearAlgebra\Complex32\SparseVector.cs">
<Link>LinearAlgebra\Complex32\SparseVector.cs</Link>
</Compile>
<Compile Include="..\Numerics\LinearAlgebra\Complex32\Vector.cs">
<Link>LinearAlgebra\Complex32\Vector.cs</Link>
</Compile>
<Compile Include="..\Numerics\LinearAlgebra\Complex\DenseMatrix.cs">
<Link>LinearAlgebra\Complex\DenseMatrix.cs</Link>
</Compile>
@ -464,6 +467,9 @@
<Compile Include="..\Numerics\LinearAlgebra\Complex\SparseVector.cs">
<Link>LinearAlgebra\Complex\SparseVector.cs</Link>
</Compile>
<Compile Include="..\Numerics\LinearAlgebra\Complex\Vector.cs">
<Link>LinearAlgebra\Complex\Vector.cs</Link>
</Compile>
<Compile Include="..\Numerics\LinearAlgebra\Double\DenseMatrix.cs">
<Link>LinearAlgebra\Double\DenseMatrix.cs</Link>
</Compile>
@ -572,6 +578,12 @@
<Compile Include="..\Numerics\LinearAlgebra\Double\SparseVector.cs">
<Link>LinearAlgebra\Double\SparseVector.cs</Link>
</Compile>
<Compile Include="..\Numerics\LinearAlgebra\Double\Vector.cs">
<Link>LinearAlgebra\Double\Vector.cs</Link>
</Compile>
<Compile Include="..\Numerics\LinearAlgebra\Generic\Common.cs">
<Link>LinearAlgebra\Generic\Common.cs</Link>
</Compile>
<Compile Include="..\Numerics\LinearAlgebra\Generic\Factorization\Evd.cs">
<Link>LinearAlgebra\Generic\Factorization\Evd.cs</Link>
</Compile>
@ -743,6 +755,9 @@
<Compile Include="..\Numerics\LinearAlgebra\Generic\Vector.cs">
<Link>LinearAlgebra\Generic\Vector.cs</Link>
</Compile>
<Compile Include="..\Numerics\LinearAlgebra\Single\Vector.cs">
<Link>LinearAlgebra\Single\Vector.cs</Link>
</Compile>
<Compile Include="..\Numerics\NumberTheory\IntegerTheory.cs">
<Link>NumberTheory\IntegerTheory.cs</Link>
</Compile>

126
src/UnitTests/LinearAlgebraTests/Complex/UserDefinedVectorTests.cs

@ -34,10 +34,11 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex
using System.Numerics;
using System.Collections.Generic;
using Distributions;
using LinearAlgebra.Complex;
using LinearAlgebra.Generic;
using Threading;
internal class UserDefinedVector : Vector<Complex>
internal class UserDefinedVector : Vector
{
private readonly Complex[] _data;
@ -76,128 +77,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex
return new UserDefinedVector(size);
}
public override Vector<Complex> Negate()
{
var result = new UserDefinedVector(Count);
CommonParallel.For(
0,
_data.Length,
index => result[index] = -_data[index]);
return result;
}
public override Vector<Complex> Random(int length, IContinuousDistribution randomDistribution)
{
if (length < 1)
{
throw new ArgumentException("length");
}
var v = (UserDefinedVector)CreateVector(length);
for (var index = 0; index < v._data.Length; index++)
{
v._data[index] = new Complex(randomDistribution.Sample(), randomDistribution.Sample());
}
return v;
}
public override Vector<Complex> Random(int length, IDiscreteDistribution randomDistribution)
{
if (length < 1)
{
throw new ArgumentException("length");
}
var v = (UserDefinedVector)CreateVector(length);
for (var index = 0; index < v._data.Length; index++)
{
v._data[index] = new Complex(randomDistribution.Sample(), randomDistribution.Sample());
}
return v;
}
public override int MinimumIndex()
{
throw new NotSupportedException();
}
public override int MaximumIndex()
{
throw new NotSupportedException();
}
public override Vector<Complex> Normalize(double p)
{
if (p < 0.0)
{
throw new ArgumentOutOfRangeException("p");
}
var norm = Norm(p);
var clone = Clone();
if (norm == 0.0)
{
return clone;
}
clone.Multiply(1.0 / norm, clone);
return clone;
}
protected sealed override Complex AddT(Complex val1, Complex val2)
{
return val1 + val2;
}
protected sealed override Complex SubtractT(Complex val1, Complex val2)
{
return val1 - val2;
}
protected sealed override Complex MultiplyT(Complex val1, Complex val2)
{
return val1 * val2;
}
protected sealed override Complex DivideT(Complex val1, Complex val2)
{
return val1 / val2;
}
protected sealed override double AbsoluteT(Complex val1)
{
return val1.Magnitude;
}
public override void Conjugate(Vector<Complex> target)
{
if (target == null)
{
throw new ArgumentNullException("target");
}
if (Count != target.Count)
{
throw new ArgumentException("target");
}
if (ReferenceEquals(this, target))
{
var tmp = CreateVector(Count);
Conjugate(tmp);
tmp.CopyTo(target);
}
CommonParallel.For(
0,
Count,
index => target[index] = this[index].Conjugate());
}
}
}
public class UserDefinedVectorTests : VectorTests
{

4
src/UnitTests/LinearAlgebraTests/Complex32/MatrixTests.Arithmetic.cs

@ -707,7 +707,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32
for (var j = 0; j < result.ColumnCount; j++)
{
var col = result.Column(j);
AssertHelpers.AlmostEqual(Complex32.One, (float)col.Norm(pValue), 6);
AssertHelpers.AlmostEqual(Complex32.One, col.Norm(pValue).Real, 6);
}
}
@ -721,7 +721,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32
for (var i = 0; i < matrix.RowCount; i++)
{
var row = matrix.Row(i);
AssertHelpers.AlmostEqual(Complex32.One, (float)row.Norm(pValue), 6);
AssertHelpers.AlmostEqual(Complex32.One, row.Norm(pValue).Real, 6);
}
}

34
src/UnitTests/LinearAlgebraTests/Complex32/Solvers/StopCriterium/ResidualStopCriteriumTest.cs

@ -13,7 +13,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32.Solvers.StopCr
[ExpectedArgumentOutOfRangeException]
public void CreateWithNegativeMaximum()
{
new ResidualStopCriterium(-0.1);
new ResidualStopCriterium(-0.1f);
Assert.Fail();
}
@ -29,10 +29,10 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32.Solvers.StopCr
[MultipleAsserts]
public void Create()
{
var criterium = new ResidualStopCriterium(1e-8, 50);
var criterium = new ResidualStopCriterium(1e-6f, 50);
Assert.IsNotNull(criterium, "There should be a criterium");
Assert.AreEqual(1e-8, criterium.Maximum, "Incorrect maximum");
Assert.AreEqual(1e-6f, criterium.Maximum, "Incorrect maximum");
Assert.AreEqual(50, criterium.MinimumIterationsBelowMaximum, "Incorrect iteration count");
}
@ -40,7 +40,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32.Solvers.StopCr
[MultipleAsserts]
public void ResetMaximum()
{
var criterium = new ResidualStopCriterium(1e-8, 50);
var criterium = new ResidualStopCriterium(1e-6f, 50);
Assert.IsNotNull(criterium, "There should be a criterium");
criterium.ResetMaximumResidualToDefault();
@ -51,7 +51,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32.Solvers.StopCr
[MultipleAsserts]
public void ResetMinimumIterationsBelowMaximum()
{
var criterium = new ResidualStopCriterium(1e-8, 50);
var criterium = new ResidualStopCriterium(1e-6f, 50);
Assert.IsNotNull(criterium, "There should be a criterium");
criterium.ResetMinimumIterationsBelowMaximumToDefault();
@ -62,7 +62,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32.Solvers.StopCr
[ExpectedArgumentOutOfRangeException]
public void DetermineStatusWithIllegalIterationNumber()
{
var criterium = new ResidualStopCriterium(1e-8, 50);
var criterium = new ResidualStopCriterium(1e-6f, 50);
Assert.IsNotNull(criterium, "There should be a criterium");
criterium.DetermineStatus(-1,
@ -76,7 +76,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32.Solvers.StopCr
[ExpectedArgumentNullException]
public void DetermineStatusWithNullSolutionVector()
{
var criterium = new ResidualStopCriterium(1e-8, 50);
var criterium = new ResidualStopCriterium(1e-6f, 50);
Assert.IsNotNull(criterium, "There should be a criterium");
criterium.DetermineStatus(1,
@ -90,7 +90,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32.Solvers.StopCr
[ExpectedArgumentNullException]
public void DetermineStatusWithNullSourceVector()
{
var criterium = new ResidualStopCriterium(1e-8, 50);
var criterium = new ResidualStopCriterium(1e-6f, 50);
Assert.IsNotNull(criterium, "There should be a criterium");
criterium.DetermineStatus(1,
@ -104,7 +104,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32.Solvers.StopCr
[ExpectedArgumentNullException]
public void DetermineStatusWithNullResidualVector()
{
var criterium = new ResidualStopCriterium(1e-8, 50);
var criterium = new ResidualStopCriterium(1e-6f, 50);
Assert.IsNotNull(criterium, "There should be a criterium");
criterium.DetermineStatus(1,
@ -117,7 +117,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32.Solvers.StopCr
[ExpectedArgumentException]
public void DetermineStatusWithNonMatchingSolutionVector()
{
var criterium = new ResidualStopCriterium(1e-8, 50);
var criterium = new ResidualStopCriterium(1e-6f, 50);
Assert.IsNotNull(criterium, "There should be a criterium");
criterium.DetermineStatus(1,
@ -130,7 +130,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32.Solvers.StopCr
[ExpectedArgumentException]
public void DetermineStatusWithNonMatchingSourceVector()
{
var criterium = new ResidualStopCriterium(1e-8, 50);
var criterium = new ResidualStopCriterium(1e-6f, 50);
Assert.IsNotNull(criterium, "There should be a criterium");
criterium.DetermineStatus(1,
@ -143,7 +143,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32.Solvers.StopCr
[ExpectedArgumentException]
public void DetermineStatusWithNonMatchingResidualVector()
{
var criterium = new ResidualStopCriterium(1e-8, 50);
var criterium = new ResidualStopCriterium(1e-6f, 50);
Assert.IsNotNull(criterium, "There should be a criterium");
criterium.DetermineStatus(1,
@ -156,7 +156,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32.Solvers.StopCr
[MultipleAsserts]
public void DetermineStatusWithSourceNaN()
{
var criterium = new ResidualStopCriterium(1e-3, 10);
var criterium = new ResidualStopCriterium(1e-3f, 10);
Assert.IsNotNull(criterium, "There should be a criterium");
var solution = new DenseVector(new[] { new Complex32(1.0f, 1), new Complex32(1.0f, 1), new Complex32(2.0f, 1) });
@ -171,7 +171,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32.Solvers.StopCr
[MultipleAsserts]
public void DetermineStatusWithResidualNaN()
{
var criterium = new ResidualStopCriterium(1e-3, 10);
var criterium = new ResidualStopCriterium(1e-3f, 10);
Assert.IsNotNull(criterium, "There should be a criterium");
var solution = new DenseVector(new[] { new Complex32(1.0f, 1), new Complex32(1.0f, 1), new Complex32(2.0f, 1) });
@ -204,7 +204,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32.Solvers.StopCr
[MultipleAsserts]
public void DetermineStatus()
{
var criterium = new ResidualStopCriterium(1e-3, 10);
var criterium = new ResidualStopCriterium(1e-3f, 10);
Assert.IsNotNull(criterium, "There should be a criterium");
// Note that the solution vector isn't actually being used so ...
@ -227,7 +227,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32.Solvers.StopCr
[MultipleAsserts]
public void ResetCalculationState()
{
var criterium = new ResidualStopCriterium(1e-3, 10);
var criterium = new ResidualStopCriterium(1e-3f, 10);
Assert.IsNotNull(criterium, "There should be a criterium");
var solution = new DenseVector(new[] { new Complex32(0.001f, 1), new Complex32(0.001f, 1), new Complex32(0.002f, 1) });
@ -245,7 +245,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32.Solvers.StopCr
[MultipleAsserts]
public void Clone()
{
var criterium = new ResidualStopCriterium(1e-3, 10);
var criterium = new ResidualStopCriterium(1e-3f, 10);
Assert.IsNotNull(criterium, "There should be a criterium");
var clone = criterium.Clone();

136
src/UnitTests/LinearAlgebraTests/Complex32/UserDefinedVectorTests.cs

@ -3,9 +3,7 @@
// http://numerics.mathdotnet.com
// http://github.com/mathnet/mathnet-numerics
// http://mathnetnumerics.codeplex.com
//
// Copyright (c) 2009-2010 Math.NET
//
// Permission is hereby granted, free of charge, to any person
// obtaining a copy of this software and associated documentation
// files (the "Software"), to deal in the Software without
@ -14,10 +12,8 @@
// copies of the Software, and to permit persons to whom the
// Software is furnished to do so, subject to the following
// conditions:
//
// The above copyright notice and this permission notice shall be
// included in all copies or substantial portions of the Software.
//
// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES
// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
@ -30,14 +26,12 @@
namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32
{
using System;
using Numerics;
using System.Collections.Generic;
using Distributions;
using LinearAlgebra.Complex32;
using LinearAlgebra.Generic;
using Threading;
using Complex32 = Numerics.Complex32;
internal class UserDefinedVector : Vector<Complex32>
internal class UserDefinedVector : Vector
{
private readonly Complex32[] _data;
@ -75,128 +69,6 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32
{
return new UserDefinedVector(size);
}
public override Vector<Complex32> Negate()
{
var result = new UserDefinedVector(Count);
CommonParallel.For(
0,
_data.Length,
index => result[index] = -_data[index]);
return result;
}
public override Vector<Complex32> Random(int length, IContinuousDistribution randomDistribution)
{
if (length < 1)
{
throw new ArgumentException("length");
}
var v = (UserDefinedVector)CreateVector(length);
for (var index = 0; index < v._data.Length; index++)
{
v._data[index] = new Complex32((float)randomDistribution.Sample(), (float)randomDistribution.Sample());
}
return v;
}
public override Vector<Complex32> Random(int length, IDiscreteDistribution randomDistribution)
{
if (length < 1)
{
throw new ArgumentException("length");
}
var v = (UserDefinedVector)CreateVector(length);
for (var index = 0; index < v._data.Length; index++)
{
v._data[index] = new Complex32(randomDistribution.Sample(), randomDistribution.Sample());
}
return v;
}
public override int MinimumIndex()
{
throw new NotSupportedException();
}
public override int MaximumIndex()
{
throw new NotSupportedException();
}
public override Vector<Complex32> Normalize(double p)
{
if (p < 0.0)
{
throw new ArgumentOutOfRangeException("p");
}
var norm = Norm(p);
var clone = Clone();
if (norm == 0.0)
{
return clone;
}
clone.Multiply(1.0f / (float)norm, clone);
return clone;
}
protected sealed override Complex32 AddT(Complex32 val1, Complex32 val2)
{
return val1 + val2;
}
protected sealed override Complex32 SubtractT(Complex32 val1, Complex32 val2)
{
return val1 - val2;
}
protected sealed override Complex32 MultiplyT(Complex32 val1, Complex32 val2)
{
return val1 * val2;
}
protected sealed override Complex32 DivideT(Complex32 val1, Complex32 val2)
{
return val1 / val2;
}
protected sealed override double AbsoluteT(Complex32 val1)
{
return val1.Magnitude;
}
public override void Conjugate(Vector<Complex32> target)
{
if (target == null)
{
throw new ArgumentNullException("target");
}
if (Count != target.Count)
{
throw new ArgumentException("target");
}
if (ReferenceEquals(this, target))
{
var tmp = CreateVector(Count);
Conjugate(tmp);
tmp.CopyTo(target);
}
CommonParallel.For(
0,
Count,
index => target[index] = this[index].Conjugate());
}
}
public class UserDefinedVectorTests : VectorTests
@ -217,4 +89,4 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32
return vector;
}
}
}
}

10
src/UnitTests/LinearAlgebraTests/Complex32/VectorTests.Norm.cs

@ -40,21 +40,21 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32
public void CanComputeNorm()
{
var vector = CreateVector(Data);
AssertHelpers.AlmostEqual(7.7459666f, (float)vector.Norm(2), 7);
AssertHelpers.AlmostEqual(7.7459666f, vector.Norm(2).Real, 7);
}
[Test]
public void CanComputeNorm1()
{
var vector = CreateVector(Data);
AssertHelpers.AlmostEqual(16.0346843f, (float)vector.Norm(1), 7);
AssertHelpers.AlmostEqual(16.0346843f, vector.Norm(1).Real, 7);
}
[Test]
public void CanComputeSquareNorm()
{
var vector = CreateVector(Data);
AssertHelpers.AlmostEqual(60f, (float)vector.Norm(2) * (float)vector.Norm(2), 6);
AssertHelpers.AlmostEqual(60f, vector.Norm(2).Real * vector.Norm(2).Real, 6);
}
[Test]
@ -65,14 +65,14 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32
public void CanComputeNormP(int p, float expected)
{
var vector = CreateVector(Data);
AssertHelpers.AlmostEqual(expected, (float)vector.Norm(p), 7);
AssertHelpers.AlmostEqual(expected, vector.Norm(p).Real, 7);
}
[Test]
public void CanComputeNormInfinity()
{
var vector = CreateVector(Data);
AssertHelpers.AlmostEqual(5.0990195, (float)vector.Norm(Double.PositiveInfinity), 7);
AssertHelpers.AlmostEqual(5.0990195, vector.Norm(Double.PositiveInfinity).Real, 7);
}
[Test]

132
src/UnitTests/LinearAlgebraTests/Double/UserDefinedVectorTests.cs

@ -3,9 +3,7 @@
// http://numerics.mathdotnet.com
// http://github.com/mathnet/mathnet-numerics
// http://mathnetnumerics.codeplex.com
//
// Copyright (c) 2009-2010 Math.NET
//
// Permission is hereby granted, free of charge, to any person
// obtaining a copy of this software and associated documentation
// files (the "Software"), to deal in the Software without
@ -14,10 +12,8 @@
// copies of the Software, and to permit persons to whom the
// Software is furnished to do so, subject to the following
// conditions:
//
// The above copyright notice and this permission notice shall be
// included in all copies or substantial portions of the Software.
//
// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES
// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
@ -30,14 +26,11 @@
namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Double
{
using System;
using System.Collections.Generic;
using Distributions;
using LinearAlgebra.Double;
using LinearAlgebra.Generic;
using Properties;
using Threading;
internal class UserDefinedVector : Vector<double>
internal class UserDefinedVector : Vector
{
private readonly double[] _data;
@ -75,125 +68,6 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Double
{
return new UserDefinedVector(size);
}
public override Vector<double> Negate()
{
var result = new UserDefinedVector(Count);
CommonParallel.For(
0,
_data.Length,
index => result[index] = -_data[index]);
return result;
}
public override Vector<double> Random(int length, IContinuousDistribution randomDistribution)
{
if (length < 1)
{
throw new ArgumentException(Resources.ArgumentMustBePositive, "length");
}
var v = (UserDefinedVector)CreateVector(length);
for (var index = 0; index < v._data.Length; index++)
{
v._data[index] = randomDistribution.Sample();
}
return v;
}
public override Vector<double> Random(int length, IDiscreteDistribution randomDistribution)
{
if (length < 1)
{
throw new ArgumentException(Resources.ArgumentMustBePositive, "length");
}
var v = (UserDefinedVector)CreateVector(length);
for (var index = 0; index < v._data.Length; index++)
{
v._data[index] = randomDistribution.Sample();
}
return v;
}
public override int MinimumIndex()
{
var index = 0;
var min = _data[0];
for (var i = 1; i < Count; i++)
{
if (min > _data[i])
{
index = i;
min = _data[i];
}
}
return index;
}
public override int MaximumIndex()
{
var index = 0;
var max = _data[0];
for (var i = 1; i < Count; i++)
{
if (max < _data[i])
{
index = i;
max = _data[i];
}
}
return index;
}
public override Vector<double> Normalize(double p)
{
if (p < 0.0)
{
throw new ArgumentOutOfRangeException("p");
}
var norm = Norm(p);
var clone = Clone();
if (norm == 0.0)
{
return clone;
}
clone.Multiply(1.0 / norm, clone);
return clone;
}
protected sealed override double AddT(double val1, double val2)
{
return val1 + val2;
}
protected sealed override double SubtractT(double val1, double val2)
{
return val1 - val2;
}
protected sealed override double MultiplyT(double val1, double val2)
{
return val1 * val2;
}
protected sealed override double DivideT(double val1, double val2)
{
return val1 / val2;
}
protected sealed override double AbsoluteT(double val1)
{
return Math.Abs(val1);
}
}
public class UserDefinedVectorTests : VectorTests
@ -214,4 +88,4 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Double
return vector;
}
}
}
}

34
src/UnitTests/LinearAlgebraTests/Single/Solvers/StopCriterium/ResidualStopCriteriumTest.cs

@ -12,7 +12,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single.Solvers.StopCrite
[ExpectedArgumentOutOfRangeException]
public void CreateWithNegativeMaximum()
{
new ResidualStopCriterium(-0.1);
new ResidualStopCriterium(-0.1f);
Assert.Fail();
}
@ -28,10 +28,10 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single.Solvers.StopCrite
[MultipleAsserts]
public void Create()
{
var criterium = new ResidualStopCriterium(1e-8, 50);
var criterium = new ResidualStopCriterium(1e-6f, 50);
Assert.IsNotNull(criterium, "There should be a criterium");
Assert.AreEqual(1e-8, criterium.Maximum, "Incorrect maximum");
Assert.AreEqual(1e-6f, criterium.Maximum, "Incorrect maximum");
Assert.AreEqual(50, criterium.MinimumIterationsBelowMaximum, "Incorrect iteration count");
}
@ -39,7 +39,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single.Solvers.StopCrite
[MultipleAsserts]
public void ResetMaximum()
{
var criterium = new ResidualStopCriterium(1e-8, 50);
var criterium = new ResidualStopCriterium(1e-6f, 50);
Assert.IsNotNull(criterium, "There should be a criterium");
criterium.ResetMaximumResidualToDefault();
@ -50,7 +50,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single.Solvers.StopCrite
[MultipleAsserts]
public void ResetMinimumIterationsBelowMaximum()
{
var criterium = new ResidualStopCriterium(1e-8, 50);
var criterium = new ResidualStopCriterium(1e-6f, 50);
Assert.IsNotNull(criterium, "There should be a criterium");
criterium.ResetMinimumIterationsBelowMaximumToDefault();
@ -61,7 +61,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single.Solvers.StopCrite
[ExpectedArgumentOutOfRangeException]
public void DetermineStatusWithIllegalIterationNumber()
{
var criterium = new ResidualStopCriterium(1e-8, 50);
var criterium = new ResidualStopCriterium(1e-6f, 50);
Assert.IsNotNull(criterium, "There should be a criterium");
criterium.DetermineStatus(-1,
@ -75,7 +75,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single.Solvers.StopCrite
[ExpectedArgumentNullException]
public void DetermineStatusWithNullSolutionVector()
{
var criterium = new ResidualStopCriterium(1e-8, 50);
var criterium = new ResidualStopCriterium(1e-6f, 50);
Assert.IsNotNull(criterium, "There should be a criterium");
criterium.DetermineStatus(1,
@ -89,7 +89,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single.Solvers.StopCrite
[ExpectedArgumentNullException]
public void DetermineStatusWithNullSourceVector()
{
var criterium = new ResidualStopCriterium(1e-8, 50);
var criterium = new ResidualStopCriterium(1e-6f, 50);
Assert.IsNotNull(criterium, "There should be a criterium");
criterium.DetermineStatus(1,
@ -103,7 +103,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single.Solvers.StopCrite
[ExpectedArgumentNullException]
public void DetermineStatusWithNullResidualVector()
{
var criterium = new ResidualStopCriterium(1e-8, 50);
var criterium = new ResidualStopCriterium(1e-6f, 50);
Assert.IsNotNull(criterium, "There should be a criterium");
criterium.DetermineStatus(1,
@ -116,7 +116,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single.Solvers.StopCrite
[ExpectedArgumentException]
public void DetermineStatusWithNonMatchingSolutionVector()
{
var criterium = new ResidualStopCriterium(1e-8, 50);
var criterium = new ResidualStopCriterium(1e-6f, 50);
Assert.IsNotNull(criterium, "There should be a criterium");
criterium.DetermineStatus(1,
@ -129,7 +129,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single.Solvers.StopCrite
[ExpectedArgumentException]
public void DetermineStatusWithNonMatchingSourceVector()
{
var criterium = new ResidualStopCriterium(1e-8, 50);
var criterium = new ResidualStopCriterium(1e-6f, 50);
Assert.IsNotNull(criterium, "There should be a criterium");
criterium.DetermineStatus(1,
@ -142,7 +142,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single.Solvers.StopCrite
[ExpectedArgumentException]
public void DetermineStatusWithNonMatchingResidualVector()
{
var criterium = new ResidualStopCriterium(1e-8, 50);
var criterium = new ResidualStopCriterium(1e-6f, 50);
Assert.IsNotNull(criterium, "There should be a criterium");
criterium.DetermineStatus(1,
@ -155,7 +155,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single.Solvers.StopCrite
[MultipleAsserts]
public void DetermineStatusWithSourceNaN()
{
var criterium = new ResidualStopCriterium(1e-3, 10);
var criterium = new ResidualStopCriterium(1e-3f, 10);
Assert.IsNotNull(criterium, "There should be a criterium");
var solution = new DenseVector(new[] { 1.0f, 1.0f, 2.0f });
@ -170,7 +170,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single.Solvers.StopCrite
[MultipleAsserts]
public void DetermineStatusWithResidualNaN()
{
var criterium = new ResidualStopCriterium(1e-3, 10);
var criterium = new ResidualStopCriterium(1e-3f, 10);
Assert.IsNotNull(criterium, "There should be a criterium");
var solution = new DenseVector(new[] { 1.0f, 1.0f, 2.0f });
@ -203,7 +203,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single.Solvers.StopCrite
[MultipleAsserts]
public void DetermineStatus()
{
var criterium = new ResidualStopCriterium(1e-3, 10);
var criterium = new ResidualStopCriterium(1e-3f, 10);
Assert.IsNotNull(criterium, "There should be a criterium");
// Note that the solution vector isn't actually being used so ...
@ -226,7 +226,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single.Solvers.StopCrite
[MultipleAsserts]
public void ResetCalculationState()
{
var criterium = new ResidualStopCriterium(1e-3, 10);
var criterium = new ResidualStopCriterium(1e-3f, 10);
Assert.IsNotNull(criterium, "There should be a criterium");
var solution = new DenseVector(new[] { 0.001f, 0.001f, 0.002f });
@ -244,7 +244,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single.Solvers.StopCrite
[MultipleAsserts]
public void Clone()
{
var criterium = new ResidualStopCriterium(1e-3, 10);
var criterium = new ResidualStopCriterium(1e-3f, 10);
Assert.IsNotNull(criterium, "There should be a criterium");
var clone = criterium.Clone();

159
src/UnitTests/LinearAlgebraTests/Single/UserDefinedVectorTests.cs

@ -35,10 +35,11 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single
using System.Linq;
using Distributions;
using LinearAlgebra.Generic;
using LinearAlgebra.Single;
using Properties;
using Threading;
internal class UserDefinedVector : Vector<float>
internal class UserDefinedVector : Vector
{
private readonly float[] _data;
@ -76,162 +77,6 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single
{
return new UserDefinedVector(size);
}
public override Vector<float> Negate()
{
var result = new UserDefinedVector(Count);
CommonParallel.For(
0,
_data.Length,
index => result[index] = -_data[index]);
return result;
}
public override Vector<float> Random(int length, IContinuousDistribution randomDistribution)
{
if (length < 1)
{
throw new ArgumentException(Resources.ArgumentMustBePositive, "length");
}
var v = (UserDefinedVector)CreateVector(length);
for (var index = 0; index < v._data.Length; index++)
{
v._data[index] = (float)randomDistribution.Sample();
}
return v;
}
public override Vector<float> Random(int length, IDiscreteDistribution randomDistribution)
{
if (length < 1)
{
throw new ArgumentException(Resources.ArgumentMustBePositive, "length");
}
var v = (UserDefinedVector)CreateVector(length);
for (var index = 0; index < v._data.Length; index++)
{
v._data[index] = randomDistribution.Sample();
}
return v;
}
public override int MinimumIndex()
{
var index = 0;
var min = _data[0];
for (var i = 1; i < Count; i++)
{
if (min > _data[i])
{
index = i;
min = _data[i];
}
}
return index;
}
public override int MaximumIndex()
{
var index = 0;
var max = _data[0];
for (var i = 1; i < Count; i++)
{
if (max < _data[i])
{
index = i;
max = _data[i];
}
}
return index;
}
public override double Norm(double p)
{
if (p < 0.0)
{
throw new ArgumentOutOfRangeException("p");
}
if (1.0 == p)
{
return CommonParallel.Aggregate(
0,
Count,
index => Math.Abs(this[index]));
}
if (2.0 == p)
{
return _data.Aggregate(0.0f, SpecialFunctions.Hypotenuse);
}
if (Double.IsPositiveInfinity(p))
{
return CommonParallel.Select(
0,
Count,
(index, localData) => Math.Max(localData, Math.Abs(_data[index])),
Math.Max);
}
var sum = CommonParallel.Aggregate(
0,
Count,
index => Math.Pow(Math.Abs(_data[index]), p));
return Math.Pow(sum, 1.0 / p);
}
public override Vector<float> Normalize(double p)
{
if (p < 0.0)
{
throw new ArgumentOutOfRangeException("p");
}
var norm = Norm(p);
var clone = Clone();
if (norm == 0.0)
{
return clone;
}
clone.Multiply(1.0f / (float)norm, clone);
return clone;
}
protected sealed override float AddT(float val1, float val2)
{
return val1 + val2;
}
protected sealed override float SubtractT(float val1, float val2)
{
return val1 - val2;
}
protected sealed override float MultiplyT(float val1, float val2)
{
return val1 * val2;
}
protected sealed override float DivideT(float val1, float val2)
{
return val1 / val2;
}
protected sealed override double AbsoluteT(float val1)
{
return Math.Abs(val1);
}
}
public class UserDefinedVectorTests : VectorTests

10
src/UnitTests/LinearAlgebraTests/Single/VectorTests.Norm.cs

@ -80,11 +80,11 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single
{
var vector = CreateVector(Data);
var result = vector.Normalize(2);
AssertHelpers.AlmostEqual(0.134839967f, result[0], 7);
AssertHelpers.AlmostEqual(0.269679934f, result[1], 7);
AssertHelpers.AlmostEqual(0.404519916f, result[2], 7);
AssertHelpers.AlmostEqual(0.539359868f, result[3], 7);
AssertHelpers.AlmostEqual(0.6741998f, result[4], 7);
AssertHelpers.AlmostEqual(0.134839967f, result[0], 6);
AssertHelpers.AlmostEqual(0.269679934f, result[1], 6);
AssertHelpers.AlmostEqual(0.404519916f, result[2], 6);
AssertHelpers.AlmostEqual(0.539359868f, result[3], 6);
AssertHelpers.AlmostEqual(0.6741998f, result[4], 6);
}
}
}
Loading…
Cancel
Save