diff --git a/src/Numerics/LinearAlgebra/Complex/Matrix.cs b/src/Numerics/LinearAlgebra/Complex/Matrix.cs
index 8a8594ea..0b111ca4 100644
--- a/src/Numerics/LinearAlgebra/Complex/Matrix.cs
+++ b/src/Numerics/LinearAlgebra/Complex/Matrix.cs
@@ -4,7 +4,7 @@
// http://github.com/mathnet/mathnet-numerics
// http://mathnetnumerics.codeplex.com
//
-// Copyright (c) 2009-2013 Math.NET
+// Copyright (c) 2009-2015 Math.NET
//
// Permission is hereby granted, free of charge, to any person
// obtaining a copy of this software and associated documentation
@@ -65,201 +65,33 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
MapInplace(x => x.Magnitude < threshold ? Complex.Zero : x, Zeros.AllowSkip);
}
- /// Calculates the induced L1 norm of this matrix.
- /// The maximum absolute column sum of the matrix.
- public override double L1Norm()
- {
- var norm = 0d;
- for (var j = 0; j < ColumnCount; j++)
- {
- var s = 0d;
- for (var i = 0; i < RowCount; i++)
- {
- s += At(i, j).Magnitude;
- }
- norm = Math.Max(norm, s);
- }
- return norm;
- }
-
- /// Calculates the induced infinity norm of this matrix.
- /// The maximum absolute row sum of the matrix.
- public override double InfinityNorm()
- {
- var norm = 0d;
- for (var i = 0; i < RowCount; i++)
- {
- var s = 0d;
- for (var j = 0; j < ColumnCount; j++)
- {
- s += At(i, j).Magnitude;
- }
- norm = Math.Max(norm, s);
- }
- return norm;
- }
-
- /// Calculates the entry-wise Frobenius norm of this matrix.
- /// The square root of the sum of the squared values.
- public override double FrobeniusNorm()
- {
- var transpose = ConjugateTranspose();
- var aat = this*transpose;
- var norm = 0d;
- for (var i = 0; i < RowCount; i++)
- {
- norm += aat.At(i, i).Magnitude;
- }
- return Math.Sqrt(norm);
- }
-
///
- /// Calculates the p-norms of all row vectors.
- /// Typical values for p are 1.0 (L1, Manhattan norm), 2.0 (L2, Euclidean norm) and positive infinity (infinity norm)
- ///
- public override Vector RowNorms(double norm)
- {
- if (norm <= 0.0)
- {
- throw new ArgumentOutOfRangeException("norm", Resources.ArgumentMustBePositive);
- }
-
- var ret = new double[RowCount];
- if (norm == 2.0)
- {
- Storage.FoldByRowUnchecked(ret, (s, x) => s + x.MagnitudeSquared(), (x, c) => Math.Sqrt(x), ret, Zeros.AllowSkip);
- }
- else if (norm == 1.0)
- {
- Storage.FoldByRowUnchecked(ret, (s, x) => s + x.Magnitude, (x, c) => x, ret, Zeros.AllowSkip);
- }
- else if (double.IsPositiveInfinity(norm))
- {
- Storage.FoldByRowUnchecked(ret, (s, x) => Math.Max(s, x.Magnitude), (x, c) => x, ret, Zeros.AllowSkip);
- }
- else
- {
- double invnorm = 1.0/norm;
- Storage.FoldByRowUnchecked(ret, (s, x) => s + Math.Pow(x.Magnitude, norm), (x, c) => Math.Pow(x, invnorm), ret, Zeros.AllowSkip);
- }
- return Vector.Build.Dense(ret);
- }
-
- ///
- /// Calculates the p-norms of all column vectors.
- /// Typical values for p are 1.0 (L1, Manhattan norm), 2.0 (L2, Euclidean norm) and positive infinity (infinity norm)
- ///
- public override Vector ColumnNorms(double norm)
- {
- if (norm <= 0.0)
- {
- throw new ArgumentOutOfRangeException("norm", Resources.ArgumentMustBePositive);
- }
-
- var ret = new double[ColumnCount];
- if (norm == 2.0)
- {
- Storage.FoldByColumnUnchecked(ret, (s, x) => s + x.MagnitudeSquared(), (x, c) => Math.Sqrt(x), ret, Zeros.AllowSkip);
- }
- else if (norm == 1.0)
- {
- Storage.FoldByColumnUnchecked(ret, (s, x) => s + x.Magnitude, (x, c) => x, ret, Zeros.AllowSkip);
- }
- else if (double.IsPositiveInfinity(norm))
- {
- Storage.FoldByColumnUnchecked(ret, (s, x) => Math.Max(s, x.Magnitude), (x, c) => x, ret, Zeros.AllowSkip);
- }
- else
- {
- double invnorm = 1.0/norm;
- Storage.FoldByColumnUnchecked(ret, (s, x) => s + Math.Pow(x.Magnitude, norm), (x, c) => Math.Pow(x, invnorm), ret, Zeros.AllowSkip);
- }
- return Vector.Build.Dense(ret);
- }
-
- ///
- /// Normalizes all row vectors to a unit p-norm.
- /// Typical values for p are 1.0 (L1, Manhattan norm), 2.0 (L2, Euclidean norm) and positive infinity (infinity norm)
- ///
- public override sealed Matrix NormalizeRows(double norm)
- {
- var norminv = ((DenseVectorStorage)RowNorms(norm).Storage).Data;
- for (int i = 0; i < norminv.Length; i++)
- {
- norminv[i] = norminv[i] == 0d ? 1d : 1d/norminv[i];
- }
-
- var result = Build.SameAs(this, RowCount, ColumnCount);
- Storage.MapIndexedTo(result.Storage, (i, j, x) => norminv[i]*x, Zeros.AllowSkip, ExistingData.AssumeZeros);
- return result;
- }
-
- ///
- /// Normalizes all column vectors to a unit p-norm.
- /// Typical values for p are 1.0 (L1, Manhattan norm), 2.0 (L2, Euclidean norm) and positive infinity (infinity norm)
- ///
- public override sealed Matrix NormalizeColumns(double norm)
- {
- var norminv = ((DenseVectorStorage)ColumnNorms(norm).Storage).Data;
- for (int i = 0; i < norminv.Length; i++)
- {
- norminv[i] = norminv[i] == 0d ? 1d : 1d/norminv[i];
- }
-
- var result = Build.SameAs(this, RowCount, ColumnCount);
- Storage.MapIndexedTo(result.Storage, (i, j, x) => norminv[j]*x, Zeros.AllowSkip, ExistingData.AssumeZeros);
- return result;
- }
-
- ///
- /// Calculates the value sum of each row vector.
- ///
- public override Vector RowSums()
- {
- var ret = new Complex[RowCount];
- Storage.FoldByRowUnchecked(ret, (s, x) => s + x, (x, c) => x, ret, Zeros.AllowSkip);
- return Vector.Build.Dense(ret);
- }
-
- ///
- /// Calculates the absolute value sum of each row vector.
- ///
- public override Vector RowAbsoluteSums()
- {
- var ret = new Complex[RowCount];
- Storage.FoldByRowUnchecked(ret, (s, x) => s + x.Magnitude, (x, c) => x, ret, Zeros.AllowSkip);
- return Vector.Build.Dense(ret);
- }
-
- ///
- /// Calculates the value sum of each column vector.
+ /// Returns the conjugate transpose of this matrix.
///
- public override Vector ColumnSums()
+ /// The conjugate transpose of this matrix.
+ public override sealed Matrix ConjugateTranspose()
{
- var ret = new Complex[ColumnCount];
- Storage.FoldByColumnUnchecked(ret, (s, x) => s + x, (x, c) => x, ret, Zeros.AllowSkip);
- return Vector.Build.Dense(ret);
+ var ret = Transpose();
+ ret.MapInplace(c => c.Conjugate(), Zeros.AllowSkip);
+ return ret;
}
///
- /// Calculates the absolute value sum of each column vector.
+ /// Complex conjugates each element of this matrix and place the results into the result matrix.
///
- public override Vector ColumnAbsoluteSums()
+ /// The result of the conjugation.
+ protected override void DoConjugate(Matrix result)
{
- var ret = new Complex[ColumnCount];
- Storage.FoldByColumnUnchecked(ret, (s, x) => s + x.Magnitude, (x, c) => x, ret, Zeros.AllowSkip);
- return Vector.Build.Dense(ret);
+ Map(Complex.Conjugate, result, Zeros.AllowSkip);
}
///
- /// Returns the conjugate transpose of this matrix.
+ /// Negate each element of this matrix and place the results into the result matrix.
///
- /// The conjugate transpose of this matrix.
- public override sealed Matrix ConjugateTranspose()
+ /// The result of the negation.
+ protected override void DoNegate(Matrix result)
{
- var ret = Transpose();
- ret.MapInplace(c => c.Conjugate(), Zeros.AllowSkip);
- return ret;
+ Map(Complex.Negate, result, Zeros.AllowSkip);
}
///
@@ -269,13 +101,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
/// The matrix to store the result of the addition.
protected override void DoAdd(Complex scalar, Matrix result)
{
- for (var i = 0; i < RowCount; i++)
- {
- for (var j = 0; j < ColumnCount; j++)
- {
- result.At(i, j, At(i, j) + scalar);
- }
- }
+ Map(x => x + scalar, result, Zeros.Include);
}
///
@@ -287,13 +113,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
/// If the two matrices don't have the same dimensions.
protected override void DoAdd(Matrix other, Matrix result)
{
- for (var i = 0; i < RowCount; i++)
- {
- for (var j = 0; j < ColumnCount; j++)
- {
- result.At(i, j, At(i, j) + other.At(i, j));
- }
- }
+ Map2(Complex.Add, other, result, Zeros.AllowSkip);
}
///
@@ -303,13 +123,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
/// The matrix to store the result of the subtraction.
protected override void DoSubtract(Complex scalar, Matrix result)
{
- for (var i = 0; i < RowCount; i++)
- {
- for (var j = 0; j < ColumnCount; j++)
- {
- result.At(i, j, At(i, j) - scalar);
- }
- }
+ Map(x => x - scalar, result, Zeros.Include);
}
///
@@ -321,13 +135,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
/// If the two matrices don't have the same dimensions.
protected override void DoSubtract(Matrix other, Matrix result)
{
- for (var i = 0; i < RowCount; i++)
- {
- for (var j = 0; j < ColumnCount; j++)
- {
- result.At(i, j, At(i, j) - other.At(i, j));
- }
- }
+ Map2(Complex.Subtract, other, result, Zeros.AllowSkip);
}
///
@@ -337,13 +145,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
/// The matrix to store the result of the multiplication.
protected override void DoMultiply(Complex scalar, Matrix result)
{
- for (var i = 0; i < RowCount; i++)
- {
- for (var j = 0; j < ColumnCount; j++)
- {
- result.At(i, j, At(i, j)*scalar);
- }
- }
+ Map(x => x*scalar, result, Zeros.AllowSkip);
}
///
@@ -392,23 +194,17 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
/// The matrix to store the result of the division.
protected override void DoDivide(Complex divisor, Matrix result)
{
- DoMultiply(1.0/divisor, result);
+ Map(x => x/divisor, result, divisor.IsZero() ? Zeros.Include : Zeros.AllowSkip);
}
///
/// Divides a scalar by each element of the matrix and stores the result in the result matrix.
///
- /// The scalar to add.
+ /// The scalar to divide by each element of the matrix.
/// The matrix to store the result of the division.
protected override void DoDivideByThis(Complex dividend, Matrix result)
{
- for (var i = 0; i < RowCount; i++)
- {
- for (var j = 0; j < ColumnCount; j++)
- {
- result.At(i, j, dividend/At(i, j));
- }
- }
+ Map(x => dividend/x, result, Zeros.Include);
}
///
@@ -531,36 +327,6 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
}
}
- ///
- /// Negate each element of this matrix and place the results into the result matrix.
- ///
- /// The result of the negation.
- protected override void DoNegate(Matrix result)
- {
- for (var i = 0; i < RowCount; i++)
- {
- for (var j = 0; j < ColumnCount; j++)
- {
- result.At(i, j, -At(i, j));
- }
- }
- }
-
- ///
- /// Complex conjugates each element of this matrix and place the results into the result matrix.
- ///
- /// The result of the conjugation.
- protected override void DoConjugate(Matrix result)
- {
- for (var i = 0; i < RowCount; i++)
- {
- for (var j = 0; j < ColumnCount; j++)
- {
- result.At(i, j, At(i, j).Conjugate());
- }
- }
- }
-
///
/// Pointwise multiplies this matrix with another matrix and stores the result into the result matrix.
///
@@ -568,13 +334,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
/// The matrix to store the result of the pointwise multiplication.
protected override void DoPointwiseMultiply(Matrix other, Matrix result)
{
- for (var j = 0; j < ColumnCount; j++)
- {
- for (var i = 0; i < RowCount; i++)
- {
- result.At(i, j, At(i, j)*other.At(i, j));
- }
- }
+ Map2(Complex.Multiply, other, result, Zeros.AllowSkip);
}
///
@@ -584,13 +344,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
/// The matrix to store the result of the pointwise division.
protected override void DoPointwiseDivide(Matrix divisor, Matrix result)
{
- for (var j = 0; j < ColumnCount; j++)
- {
- for (var i = 0; i < RowCount; i++)
- {
- result.At(i, j, At(i, j)/divisor.At(i, j));
- }
- }
+ Map2(Complex.Divide, divisor, result, Zeros.Include);
}
///
@@ -600,7 +354,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
/// The vector to store the result of the pointwise power.
protected override void DoPointwisePower(Complex exponent, Matrix result)
{
- Map(x => x.Power(exponent), result, Zeros.AllowSkip);
+ Map(x => x.Power(exponent), result, Zeros.Include);
}
///
@@ -687,6 +441,194 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
Map(Complex.Log, result, Zeros.Include);
}
+
+
+ /// Calculates the induced L1 norm of this matrix.
+ /// The maximum absolute column sum of the matrix.
+ public override double L1Norm()
+ {
+ var norm = 0d;
+ for (var j = 0; j < ColumnCount; j++)
+ {
+ var s = 0d;
+ for (var i = 0; i < RowCount; i++)
+ {
+ s += At(i, j).Magnitude;
+ }
+ norm = Math.Max(norm, s);
+ }
+ return norm;
+ }
+
+ /// Calculates the induced infinity norm of this matrix.
+ /// The maximum absolute row sum of the matrix.
+ public override double InfinityNorm()
+ {
+ var norm = 0d;
+ for (var i = 0; i < RowCount; i++)
+ {
+ var s = 0d;
+ for (var j = 0; j < ColumnCount; j++)
+ {
+ s += At(i, j).Magnitude;
+ }
+ norm = Math.Max(norm, s);
+ }
+ return norm;
+ }
+
+ /// Calculates the entry-wise Frobenius norm of this matrix.
+ /// The square root of the sum of the squared values.
+ public override double FrobeniusNorm()
+ {
+ var transpose = ConjugateTranspose();
+ var aat = this*transpose;
+ var norm = 0d;
+ for (var i = 0; i < RowCount; i++)
+ {
+ norm += aat.At(i, i).Magnitude;
+ }
+ return Math.Sqrt(norm);
+ }
+
+ ///
+ /// Calculates the p-norms of all row vectors.
+ /// Typical values for p are 1.0 (L1, Manhattan norm), 2.0 (L2, Euclidean norm) and positive infinity (infinity norm)
+ ///
+ public override Vector RowNorms(double norm)
+ {
+ if (norm <= 0.0)
+ {
+ throw new ArgumentOutOfRangeException("norm", Resources.ArgumentMustBePositive);
+ }
+
+ var ret = new double[RowCount];
+ if (norm == 2.0)
+ {
+ Storage.FoldByRowUnchecked(ret, (s, x) => s + x.MagnitudeSquared(), (x, c) => Math.Sqrt(x), ret, Zeros.AllowSkip);
+ }
+ else if (norm == 1.0)
+ {
+ Storage.FoldByRowUnchecked(ret, (s, x) => s + x.Magnitude, (x, c) => x, ret, Zeros.AllowSkip);
+ }
+ else if (double.IsPositiveInfinity(norm))
+ {
+ Storage.FoldByRowUnchecked(ret, (s, x) => Math.Max(s, x.Magnitude), (x, c) => x, ret, Zeros.AllowSkip);
+ }
+ else
+ {
+ double invnorm = 1.0/norm;
+ Storage.FoldByRowUnchecked(ret, (s, x) => s + Math.Pow(x.Magnitude, norm), (x, c) => Math.Pow(x, invnorm), ret, Zeros.AllowSkip);
+ }
+ return Vector.Build.Dense(ret);
+ }
+
+ ///
+ /// Calculates the p-norms of all column vectors.
+ /// Typical values for p are 1.0 (L1, Manhattan norm), 2.0 (L2, Euclidean norm) and positive infinity (infinity norm)
+ ///
+ public override Vector ColumnNorms(double norm)
+ {
+ if (norm <= 0.0)
+ {
+ throw new ArgumentOutOfRangeException("norm", Resources.ArgumentMustBePositive);
+ }
+
+ var ret = new double[ColumnCount];
+ if (norm == 2.0)
+ {
+ Storage.FoldByColumnUnchecked(ret, (s, x) => s + x.MagnitudeSquared(), (x, c) => Math.Sqrt(x), ret, Zeros.AllowSkip);
+ }
+ else if (norm == 1.0)
+ {
+ Storage.FoldByColumnUnchecked(ret, (s, x) => s + x.Magnitude, (x, c) => x, ret, Zeros.AllowSkip);
+ }
+ else if (double.IsPositiveInfinity(norm))
+ {
+ Storage.FoldByColumnUnchecked(ret, (s, x) => Math.Max(s, x.Magnitude), (x, c) => x, ret, Zeros.AllowSkip);
+ }
+ else
+ {
+ double invnorm = 1.0/norm;
+ Storage.FoldByColumnUnchecked(ret, (s, x) => s + Math.Pow(x.Magnitude, norm), (x, c) => Math.Pow(x, invnorm), ret, Zeros.AllowSkip);
+ }
+ return Vector.Build.Dense(ret);
+ }
+
+ ///
+ /// Normalizes all row vectors to a unit p-norm.
+ /// Typical values for p are 1.0 (L1, Manhattan norm), 2.0 (L2, Euclidean norm) and positive infinity (infinity norm)
+ ///
+ public override sealed Matrix NormalizeRows(double norm)
+ {
+ var norminv = ((DenseVectorStorage)RowNorms(norm).Storage).Data;
+ for (int i = 0; i < norminv.Length; i++)
+ {
+ norminv[i] = norminv[i] == 0d ? 1d : 1d/norminv[i];
+ }
+
+ var result = Build.SameAs(this, RowCount, ColumnCount);
+ Storage.MapIndexedTo(result.Storage, (i, j, x) => norminv[i]*x, Zeros.AllowSkip, ExistingData.AssumeZeros);
+ return result;
+ }
+
+ ///
+ /// Normalizes all column vectors to a unit p-norm.
+ /// Typical values for p are 1.0 (L1, Manhattan norm), 2.0 (L2, Euclidean norm) and positive infinity (infinity norm)
+ ///
+ public override sealed Matrix NormalizeColumns(double norm)
+ {
+ var norminv = ((DenseVectorStorage)ColumnNorms(norm).Storage).Data;
+ for (int i = 0; i < norminv.Length; i++)
+ {
+ norminv[i] = norminv[i] == 0d ? 1d : 1d/norminv[i];
+ }
+
+ var result = Build.SameAs(this, RowCount, ColumnCount);
+ Storage.MapIndexedTo(result.Storage, (i, j, x) => norminv[j]*x, Zeros.AllowSkip, ExistingData.AssumeZeros);
+ return result;
+ }
+
+ ///
+ /// Calculates the value sum of each row vector.
+ ///
+ public override Vector RowSums()
+ {
+ var ret = new Complex[RowCount];
+ Storage.FoldByRowUnchecked(ret, (s, x) => s + x, (x, c) => x, ret, Zeros.AllowSkip);
+ return Vector.Build.Dense(ret);
+ }
+
+ ///
+ /// Calculates the absolute value sum of each row vector.
+ ///
+ public override Vector RowAbsoluteSums()
+ {
+ var ret = new Complex[RowCount];
+ Storage.FoldByRowUnchecked(ret, (s, x) => s + x.Magnitude, (x, c) => x, ret, Zeros.AllowSkip);
+ return Vector.Build.Dense(ret);
+ }
+
+ ///
+ /// Calculates the value sum of each column vector.
+ ///
+ public override Vector ColumnSums()
+ {
+ var ret = new Complex[ColumnCount];
+ Storage.FoldByColumnUnchecked(ret, (s, x) => s + x, (x, c) => x, ret, Zeros.AllowSkip);
+ return Vector.Build.Dense(ret);
+ }
+
+ ///
+ /// Calculates the absolute value sum of each column vector.
+ ///
+ public override Vector ColumnAbsoluteSums()
+ {
+ var ret = new Complex[ColumnCount];
+ Storage.FoldByColumnUnchecked(ret, (s, x) => s + x.Magnitude, (x, c) => x, ret, Zeros.AllowSkip);
+ return Vector.Build.Dense(ret);
+ }
+
///
/// Computes the trace of this matrix.
///
diff --git a/src/Numerics/LinearAlgebra/Complex/Vector.cs b/src/Numerics/LinearAlgebra/Complex/Vector.cs
index 6cda7da6..a057e15f 100644
--- a/src/Numerics/LinearAlgebra/Complex/Vector.cs
+++ b/src/Numerics/LinearAlgebra/Complex/Vector.cs
@@ -4,7 +4,7 @@
// http://github.com/mathnet/mathnet-numerics
// http://mathnetnumerics.codeplex.com
//
-// Copyright (c) 2009-2013 Math.NET
+// Copyright (c) 2009-2015 Math.NET
//
// Permission is hereby granted, free of charge, to any person
// obtaining a copy of this software and associated documentation
@@ -63,6 +63,24 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
MapInplace(x => x.Magnitude < threshold ? Complex.Zero : x, Zeros.AllowSkip);
}
+ ///
+ /// Conjugates vector and save result to
+ ///
+ /// Target vector
+ protected override void DoConjugate(Vector result)
+ {
+ Map(Complex.Conjugate, result, Zeros.AllowSkip);
+ }
+
+ ///
+ /// Negates vector and saves result to
+ ///
+ /// Target vector
+ protected override void DoNegate(Vector result)
+ {
+ Map(Complex.Negate, result, Zeros.AllowSkip);
+ }
+
///
/// Adds a scalar to each element of the vector and stores the result in the result vector.
///
@@ -74,10 +92,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
///
protected override void DoAdd(Complex scalar, Vector result)
{
- for (var index = 0; index < Count; index++)
- {
- result.At(index, At(index) + scalar);
- }
+ Map(x => x + scalar, result, Zeros.Include);
}
///
@@ -91,10 +106,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
///
protected override void DoAdd(Vector other, Vector result)
{
- for (var index = 0; index < Count; index++)
- {
- result.At(index, At(index) + other.At(index));
- }
+ Map2(Complex.Add, other, result, Zeros.AllowSkip);
}
///
@@ -108,7 +120,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
///
protected override void DoSubtract(Complex scalar, Vector result)
{
- DoAdd(-scalar, result);
+ Map(x => x - scalar, result, Zeros.Include);
}
///
@@ -122,10 +134,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
///
protected override void DoSubtract(Vector other, Vector result)
{
- for (var index = 0; index < Count; index++)
- {
- result.At(index, At(index) - other.At(index));
- }
+ Map2(Complex.Subtract, other, result, Zeros.AllowSkip);
}
///
@@ -139,10 +148,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
///
protected override void DoMultiply(Complex scalar, Vector result)
{
- for (var index = 0; index < Count; index++)
- {
- result.At(index, At(index) * scalar);
- }
+ Map(x => x*scalar, result, Zeros.AllowSkip);
}
///
@@ -156,7 +162,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
///
protected override void DoDivide(Complex divisor, Vector result)
{
- DoMultiply(1 / divisor, result);
+ Map(x => x/divisor, result, divisor.IsZero() ? Zeros.Include : Zeros.AllowSkip);
}
///
@@ -166,10 +172,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
/// The vector to store the result of the division.
protected override void DoDivideByThis(Complex dividend, Vector result)
{
- for (var index = 0; index < Count; index++)
- {
- result.At(index, dividend / At(index));
- }
+ Map(x => dividend/x, result, Zeros.Include);
}
///
@@ -179,10 +182,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
/// The vector to store the result of the pointwise multiplication.
protected override void DoPointwiseMultiply(Vector other, Vector result)
{
- for (var index = 0; index < Count; index++)
- {
- result.At(index, At(index) * other.At(index));
- }
+ Map2(Complex.Multiply, other, result, Zeros.AllowSkip);
}
///
@@ -192,10 +192,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
/// The vector to store the result of the pointwise division.
protected override void DoPointwiseDivide(Vector divisor, Vector result)
{
- for (var index = 0; index < Count; index++)
- {
- result.At(index, At(index) / divisor.At(index));
- }
+ Map2(Complex.Divide, divisor, result, Zeros.Include);
}
///
@@ -205,7 +202,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
/// The vector to store the result of the pointwise power.
protected override void DoPointwisePower(Complex exponent, Vector result)
{
- Map(x => x.Power(exponent), result, Zeros.AllowSkip);
+ Map(x => x.Power(exponent), result, Zeros.Include);
}
///
@@ -453,30 +450,6 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
return Math.Pow(sum, 1.0/p);
}
- ///
- /// Conjugates vector and save result to
- ///
- /// Target vector
- protected override void DoConjugate(Vector result)
- {
- for (var index = 0; index < Count; index++)
- {
- result.At(index, At(index).Conjugate());
- }
- }
-
- ///
- /// Negates vector and saves result to
- ///
- /// Target vector
- protected override void DoNegate(Vector result)
- {
- for (var index = 0; index < Count; index++)
- {
- result.At(index, -At(index));
- }
- }
-
///
/// Returns the index of the absolute maximum element.
///
diff --git a/src/Numerics/LinearAlgebra/Complex32/Matrix.cs b/src/Numerics/LinearAlgebra/Complex32/Matrix.cs
index cc9d5c95..a4501582 100644
--- a/src/Numerics/LinearAlgebra/Complex32/Matrix.cs
+++ b/src/Numerics/LinearAlgebra/Complex32/Matrix.cs
@@ -4,7 +4,7 @@
// http://github.com/mathnet/mathnet-numerics
// http://mathnetnumerics.codeplex.com
//
-// Copyright (c) 2009-2013 Math.NET
+// Copyright (c) 2009-2015 Math.NET
//
// Permission is hereby granted, free of charge, to any person
// obtaining a copy of this software and associated documentation
@@ -60,201 +60,33 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
MapInplace(x => x.Magnitude < threshold ? Complex32.Zero : x, Zeros.AllowSkip);
}
- /// Calculates the induced L1 norm of this matrix.
- /// The maximum absolute column sum of the matrix.
- public override double L1Norm()
- {
- var norm = 0d;
- for (var j = 0; j < ColumnCount; j++)
- {
- var s = 0d;
- for (var i = 0; i < RowCount; i++)
- {
- s += At(i, j).Magnitude;
- }
- norm = Math.Max(norm, s);
- }
- return norm;
- }
-
- /// Calculates the induced infinity norm of this matrix.
- /// The maximum absolute row sum of the matrix.
- public override double InfinityNorm()
- {
- var norm = 0d;
- for (var i = 0; i < RowCount; i++)
- {
- var s = 0d;
- for (var j = 0; j < ColumnCount; j++)
- {
- s += At(i, j).Magnitude;
- }
- norm = Math.Max(norm, s);
- }
- return norm;
- }
-
- /// Calculates the entry-wise Frobenius norm of this matrix.
- /// The square root of the sum of the squared values.
- public override double FrobeniusNorm()
- {
- var transpose = ConjugateTranspose();
- var aat = this*transpose;
- var norm = 0d;
- for (var i = 0; i < RowCount; i++)
- {
- norm += aat.At(i, i).Magnitude;
- }
- return Math.Sqrt(norm);
- }
-
- ///
- /// Calculates the p-norms of all row vectors.
- /// Typical values for p are 1.0 (L1, Manhattan norm), 2.0 (L2, Euclidean norm) and positive infinity (infinity norm)
- ///
- public override Vector RowNorms(double norm)
- {
- if (norm <= 0.0)
- {
- throw new ArgumentOutOfRangeException("norm", Resources.ArgumentMustBePositive);
- }
-
- var ret = new double[RowCount];
- if (norm == 2.0)
- {
- Storage.FoldByRowUnchecked(ret, (s, x) => s + x.MagnitudeSquared, (x, c) => Math.Sqrt(x), ret, Zeros.AllowSkip);
- }
- else if (norm == 1.0)
- {
- Storage.FoldByRowUnchecked(ret, (s, x) => s + x.Magnitude, (x, c) => x, ret, Zeros.AllowSkip);
- }
- else if (double.IsPositiveInfinity(norm))
- {
- Storage.FoldByRowUnchecked(ret, (s, x) => Math.Max(s, x.Magnitude), (x, c) => x, ret, Zeros.AllowSkip);
- }
- else
- {
- double invnorm = 1.0/norm;
- Storage.FoldByRowUnchecked(ret, (s, x) => s + Math.Pow(x.Magnitude, norm), (x, c) => Math.Pow(x, invnorm), ret, Zeros.AllowSkip);
- }
- return Vector.Build.Dense(ret);
- }
-
- ///
- /// Calculates the p-norms of all column vectors.
- /// Typical values for p are 1.0 (L1, Manhattan norm), 2.0 (L2, Euclidean norm) and positive infinity (infinity norm)
- ///
- public override Vector ColumnNorms(double norm)
- {
- if (norm <= 0.0)
- {
- throw new ArgumentOutOfRangeException("norm", Resources.ArgumentMustBePositive);
- }
-
- var ret = new double[ColumnCount];
- if (norm == 2.0)
- {
- Storage.FoldByColumnUnchecked(ret, (s, x) => s + x.MagnitudeSquared, (x, c) => Math.Sqrt(x), ret, Zeros.AllowSkip);
- }
- else if (norm == 1.0)
- {
- Storage.FoldByColumnUnchecked(ret, (s, x) => s + x.Magnitude, (x, c) => x, ret, Zeros.AllowSkip);
- }
- else if (double.IsPositiveInfinity(norm))
- {
- Storage.FoldByColumnUnchecked(ret, (s, x) => Math.Max(s, x.Magnitude), (x, c) => x, ret, Zeros.AllowSkip);
- }
- else
- {
- double invnorm = 1.0/norm;
- Storage.FoldByColumnUnchecked(ret, (s, x) => s + Math.Pow(x.Magnitude, norm), (x, c) => Math.Pow(x, invnorm), ret, Zeros.AllowSkip);
- }
- return Vector.Build.Dense(ret);
- }
-
- ///
- /// Normalizes all row vectors to a unit p-norm.
- /// Typical values for p are 1.0 (L1, Manhattan norm), 2.0 (L2, Euclidean norm) and positive infinity (infinity norm)
- ///
- public override sealed Matrix NormalizeRows(double norm)
- {
- var norminv = ((DenseVectorStorage)RowNorms(norm).Storage).Data;
- for (int i = 0; i < norminv.Length; i++)
- {
- norminv[i] = norminv[i] == 0d ? 1d : 1d/norminv[i];
- }
-
- var result = Build.SameAs(this, RowCount, ColumnCount);
- Storage.MapIndexedTo(result.Storage, (i, j, x) => ((float)norminv[i])*x, Zeros.AllowSkip, ExistingData.AssumeZeros);
- return result;
- }
-
///
- /// Normalizes all column vectors to a unit p-norm.
- /// Typical values for p are 1.0 (L1, Manhattan norm), 2.0 (L2, Euclidean norm) and positive infinity (infinity norm)
- ///
- public override sealed Matrix NormalizeColumns(double norm)
- {
- var norminv = ((DenseVectorStorage)ColumnNorms(norm).Storage).Data;
- for (int i = 0; i < norminv.Length; i++)
- {
- norminv[i] = norminv[i] == 0d ? 1d : 1d/norminv[i];
- }
-
- var result = Build.SameAs(this, RowCount, ColumnCount);
- Storage.MapIndexedTo(result.Storage, (i, j, x) => ((float)norminv[j])*x, Zeros.AllowSkip, ExistingData.AssumeZeros);
- return result;
- }
-
- ///
- /// Calculates the value sum of each row vector.
- ///
- public override Vector RowSums()
- {
- var ret = new Complex32[RowCount];
- Storage.FoldByRowUnchecked(ret, (s, x) => s + x, (x, c) => x, ret, Zeros.AllowSkip);
- return Vector.Build.Dense(ret);
- }
-
- ///
- /// Calculates the absolute value sum of each row vector.
- ///
- public override Vector RowAbsoluteSums()
- {
- var ret = new Complex32[RowCount];
- Storage.FoldByRowUnchecked(ret, (s, x) => s + x.Magnitude, (x, c) => x, ret, Zeros.AllowSkip);
- return Vector.Build.Dense(ret);
- }
-
- ///
- /// Calculates the value sum of each column vector.
+ /// Returns the conjugate transpose of this matrix.
///
- public override Vector ColumnSums()
+ /// The conjugate transpose of this matrix.
+ public override sealed Matrix ConjugateTranspose()
{
- var ret = new Complex32[ColumnCount];
- Storage.FoldByColumnUnchecked(ret, (s, x) => s + x, (x, c) => x, ret, Zeros.AllowSkip);
- return Vector.Build.Dense(ret);
+ var ret = Transpose();
+ ret.MapInplace(c => c.Conjugate(), Zeros.AllowSkip);
+ return ret;
}
///
- /// Calculates the absolute value sum of each column vector.
+ /// Complex conjugates each element of this matrix and place the results into the result matrix.
///
- public override Vector ColumnAbsoluteSums()
+ /// The result of the conjugation.
+ protected override void DoConjugate(Matrix result)
{
- var ret = new Complex32[ColumnCount];
- Storage.FoldByColumnUnchecked(ret, (s, x) => s + x.Magnitude, (x, c) => x, ret, Zeros.AllowSkip);
- return Vector.Build.Dense(ret);
+ Map(Complex32.Conjugate, result, Zeros.AllowSkip);
}
///
- /// Returns the conjugate transpose of this matrix.
+ /// Negate each element of this matrix and place the results into the result matrix.
///
- /// The conjugate transpose of this matrix.
- public override sealed Matrix ConjugateTranspose()
+ /// The result of the negation.
+ protected override void DoNegate(Matrix result)
{
- var ret = Transpose();
- ret.MapInplace(c => c.Conjugate(), Zeros.AllowSkip);
- return ret;
+ Map(Complex32.Negate, result, Zeros.AllowSkip);
}
///
@@ -264,13 +96,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
/// The matrix to store the result of the addition.
protected override void DoAdd(Complex32 scalar, Matrix result)
{
- for (var i = 0; i < RowCount; i++)
- {
- for (var j = 0; j < ColumnCount; j++)
- {
- result.At(i, j, At(i, j) + scalar);
- }
- }
+ Map(x => x + scalar, result, Zeros.Include);
}
///
@@ -282,13 +108,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
/// If the two matrices don't have the same dimensions.
protected override void DoAdd(Matrix other, Matrix result)
{
- for (var i = 0; i < RowCount; i++)
- {
- for (var j = 0; j < ColumnCount; j++)
- {
- result.At(i, j, At(i, j) + other.At(i, j));
- }
- }
+ Map2(Complex32.Add, other, result, Zeros.AllowSkip);
}
///
@@ -298,13 +118,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
/// The matrix to store the result of the subtraction.
protected override void DoSubtract(Complex32 scalar, Matrix result)
{
- for (var i = 0; i < RowCount; i++)
- {
- for (var j = 0; j < ColumnCount; j++)
- {
- result.At(i, j, At(i, j) - scalar);
- }
- }
+ Map(x => x - scalar, result, Zeros.Include);
}
///
@@ -316,13 +130,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
/// If the two matrices don't have the same dimensions.
protected override void DoSubtract(Matrix other, Matrix result)
{
- for (var i = 0; i < RowCount; i++)
- {
- for (var j = 0; j < ColumnCount; j++)
- {
- result.At(i, j, At(i, j) - other.At(i, j));
- }
- }
+ Map2(Complex32.Subtract, other, result, Zeros.AllowSkip);
}
///
@@ -332,13 +140,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
/// The matrix to store the result of the multiplication.
protected override void DoMultiply(Complex32 scalar, Matrix result)
{
- for (var i = 0; i < RowCount; i++)
- {
- for (var j = 0; j < ColumnCount; j++)
- {
- result.At(i, j, At(i, j)*scalar);
- }
- }
+ Map(x => x*scalar, result, Zeros.AllowSkip);
}
///
@@ -366,23 +168,17 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
/// The matrix to store the result of the division.
protected override void DoDivide(Complex32 divisor, Matrix result)
{
- DoMultiply(1.0f/divisor, result);
+ Map(x => x/divisor, result, divisor.IsZero() ? Zeros.Include : Zeros.AllowSkip);
}
///
/// Divides a scalar by each element of the matrix and stores the result in the result matrix.
///
- /// The scalar to add.
+ /// The scalar to divide by each element of the matrix.
/// The matrix to store the result of the division.
protected override void DoDivideByThis(Complex32 dividend, Matrix result)
{
- for (var i = 0; i < RowCount; i++)
- {
- for (var j = 0; j < ColumnCount; j++)
- {
- result.At(i, j, dividend/At(i, j));
- }
- }
+ Map(x => dividend/x, result, Zeros.Include);
}
///
@@ -526,36 +322,6 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
}
}
- ///
- /// Negate each element of this matrix and place the results into the result matrix.
- ///
- /// The result of the negation.
- protected override void DoNegate(Matrix result)
- {
- for (var i = 0; i < RowCount; i++)
- {
- for (var j = 0; j < ColumnCount; j++)
- {
- result.At(i, j, -At(i, j));
- }
- }
- }
-
- ///
- /// Complex conjugates each element of this matrix and place the results into the result matrix.
- ///
- /// The result of the conjugation.
- protected override void DoConjugate(Matrix result)
- {
- for (var i = 0; i < RowCount; i++)
- {
- for (var j = 0; j < ColumnCount; j++)
- {
- result.At(i, j, At(i, j).Conjugate());
- }
- }
- }
-
///
/// Pointwise multiplies this matrix with another matrix and stores the result into the result matrix.
///
@@ -563,13 +329,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
/// The matrix to store the result of the pointwise multiplication.
protected override void DoPointwiseMultiply(Matrix other, Matrix result)
{
- for (var j = 0; j < ColumnCount; j++)
- {
- for (var i = 0; i < RowCount; i++)
- {
- result.At(i, j, At(i, j)*other.At(i, j));
- }
- }
+ Map2(Complex32.Multiply, other, result, Zeros.AllowSkip);
}
///
@@ -579,13 +339,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
/// The matrix to store the result of the pointwise division.
protected override void DoPointwiseDivide(Matrix divisor, Matrix result)
{
- for (var j = 0; j < ColumnCount; j++)
- {
- for (var i = 0; i < RowCount; i++)
- {
- result.At(i, j, At(i, j)/divisor.At(i, j));
- }
- }
+ Map2(Complex32.Divide, divisor, result, Zeros.Include);
}
///
@@ -595,7 +349,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
/// The vector to store the result of the pointwise power.
protected override void DoPointwisePower(Complex32 exponent, Matrix result)
{
- Map(x => x.Power(exponent), result, Zeros.AllowSkip);
+ Map(x => x.Power(exponent), result, Zeros.Include);
}
///
@@ -682,6 +436,192 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
Map(Complex32.Log, result, Zeros.Include);
}
+ /// Calculates the induced L1 norm of this matrix.
+ /// The maximum absolute column sum of the matrix.
+ public override double L1Norm()
+ {
+ var norm = 0d;
+ for (var j = 0; j < ColumnCount; j++)
+ {
+ var s = 0d;
+ for (var i = 0; i < RowCount; i++)
+ {
+ s += At(i, j).Magnitude;
+ }
+ norm = Math.Max(norm, s);
+ }
+ return norm;
+ }
+
+ /// Calculates the induced infinity norm of this matrix.
+ /// The maximum absolute row sum of the matrix.
+ public override double InfinityNorm()
+ {
+ var norm = 0d;
+ for (var i = 0; i < RowCount; i++)
+ {
+ var s = 0d;
+ for (var j = 0; j < ColumnCount; j++)
+ {
+ s += At(i, j).Magnitude;
+ }
+ norm = Math.Max(norm, s);
+ }
+ return norm;
+ }
+
+ /// Calculates the entry-wise Frobenius norm of this matrix.
+ /// The square root of the sum of the squared values.
+ public override double FrobeniusNorm()
+ {
+ var transpose = ConjugateTranspose();
+ var aat = this*transpose;
+ var norm = 0d;
+ for (var i = 0; i < RowCount; i++)
+ {
+ norm += aat.At(i, i).Magnitude;
+ }
+ return Math.Sqrt(norm);
+ }
+
+ ///
+ /// Calculates the p-norms of all row vectors.
+ /// Typical values for p are 1.0 (L1, Manhattan norm), 2.0 (L2, Euclidean norm) and positive infinity (infinity norm)
+ ///
+ public override Vector RowNorms(double norm)
+ {
+ if (norm <= 0.0)
+ {
+ throw new ArgumentOutOfRangeException("norm", Resources.ArgumentMustBePositive);
+ }
+
+ var ret = new double[RowCount];
+ if (norm == 2.0)
+ {
+ Storage.FoldByRowUnchecked(ret, (s, x) => s + x.MagnitudeSquared, (x, c) => Math.Sqrt(x), ret, Zeros.AllowSkip);
+ }
+ else if (norm == 1.0)
+ {
+ Storage.FoldByRowUnchecked(ret, (s, x) => s + x.Magnitude, (x, c) => x, ret, Zeros.AllowSkip);
+ }
+ else if (double.IsPositiveInfinity(norm))
+ {
+ Storage.FoldByRowUnchecked(ret, (s, x) => Math.Max(s, x.Magnitude), (x, c) => x, ret, Zeros.AllowSkip);
+ }
+ else
+ {
+ double invnorm = 1.0/norm;
+ Storage.FoldByRowUnchecked(ret, (s, x) => s + Math.Pow(x.Magnitude, norm), (x, c) => Math.Pow(x, invnorm), ret, Zeros.AllowSkip);
+ }
+ return Vector.Build.Dense(ret);
+ }
+
+ ///
+ /// Calculates the p-norms of all column vectors.
+ /// Typical values for p are 1.0 (L1, Manhattan norm), 2.0 (L2, Euclidean norm) and positive infinity (infinity norm)
+ ///
+ public override Vector ColumnNorms(double norm)
+ {
+ if (norm <= 0.0)
+ {
+ throw new ArgumentOutOfRangeException("norm", Resources.ArgumentMustBePositive);
+ }
+
+ var ret = new double[ColumnCount];
+ if (norm == 2.0)
+ {
+ Storage.FoldByColumnUnchecked(ret, (s, x) => s + x.MagnitudeSquared, (x, c) => Math.Sqrt(x), ret, Zeros.AllowSkip);
+ }
+ else if (norm == 1.0)
+ {
+ Storage.FoldByColumnUnchecked(ret, (s, x) => s + x.Magnitude, (x, c) => x, ret, Zeros.AllowSkip);
+ }
+ else if (double.IsPositiveInfinity(norm))
+ {
+ Storage.FoldByColumnUnchecked(ret, (s, x) => Math.Max(s, x.Magnitude), (x, c) => x, ret, Zeros.AllowSkip);
+ }
+ else
+ {
+ double invnorm = 1.0/norm;
+ Storage.FoldByColumnUnchecked(ret, (s, x) => s + Math.Pow(x.Magnitude, norm), (x, c) => Math.Pow(x, invnorm), ret, Zeros.AllowSkip);
+ }
+ return Vector.Build.Dense(ret);
+ }
+
+ ///
+ /// Normalizes all row vectors to a unit p-norm.
+ /// Typical values for p are 1.0 (L1, Manhattan norm), 2.0 (L2, Euclidean norm) and positive infinity (infinity norm)
+ ///
+ public override sealed Matrix NormalizeRows(double norm)
+ {
+ var norminv = ((DenseVectorStorage)RowNorms(norm).Storage).Data;
+ for (int i = 0; i < norminv.Length; i++)
+ {
+ norminv[i] = norminv[i] == 0d ? 1d : 1d/norminv[i];
+ }
+
+ var result = Build.SameAs(this, RowCount, ColumnCount);
+ Storage.MapIndexedTo(result.Storage, (i, j, x) => ((float)norminv[i])*x, Zeros.AllowSkip, ExistingData.AssumeZeros);
+ return result;
+ }
+
+ ///
+ /// Normalizes all column vectors to a unit p-norm.
+ /// Typical values for p are 1.0 (L1, Manhattan norm), 2.0 (L2, Euclidean norm) and positive infinity (infinity norm)
+ ///
+ public override sealed Matrix NormalizeColumns(double norm)
+ {
+ var norminv = ((DenseVectorStorage)ColumnNorms(norm).Storage).Data;
+ for (int i = 0; i < norminv.Length; i++)
+ {
+ norminv[i] = norminv[i] == 0d ? 1d : 1d/norminv[i];
+ }
+
+ var result = Build.SameAs(this, RowCount, ColumnCount);
+ Storage.MapIndexedTo(result.Storage, (i, j, x) => ((float)norminv[j])*x, Zeros.AllowSkip, ExistingData.AssumeZeros);
+ return result;
+ }
+
+ ///
+ /// Calculates the value sum of each row vector.
+ ///
+ public override Vector RowSums()
+ {
+ var ret = new Complex32[RowCount];
+ Storage.FoldByRowUnchecked(ret, (s, x) => s + x, (x, c) => x, ret, Zeros.AllowSkip);
+ return Vector.Build.Dense(ret);
+ }
+
+ ///
+ /// Calculates the absolute value sum of each row vector.
+ ///
+ public override Vector RowAbsoluteSums()
+ {
+ var ret = new Complex32[RowCount];
+ Storage.FoldByRowUnchecked(ret, (s, x) => s + x.Magnitude, (x, c) => x, ret, Zeros.AllowSkip);
+ return Vector.Build.Dense(ret);
+ }
+
+ ///
+ /// Calculates the value sum of each column vector.
+ ///
+ public override Vector ColumnSums()
+ {
+ var ret = new Complex32[ColumnCount];
+ Storage.FoldByColumnUnchecked(ret, (s, x) => s + x, (x, c) => x, ret, Zeros.AllowSkip);
+ return Vector.Build.Dense(ret);
+ }
+
+ ///
+ /// Calculates the absolute value sum of each column vector.
+ ///
+ public override Vector ColumnAbsoluteSums()
+ {
+ var ret = new Complex32[ColumnCount];
+ Storage.FoldByColumnUnchecked(ret, (s, x) => s + x.Magnitude, (x, c) => x, ret, Zeros.AllowSkip);
+ return Vector.Build.Dense(ret);
+ }
+
///
/// Computes the trace of this matrix.
///
diff --git a/src/Numerics/LinearAlgebra/Complex32/Vector.cs b/src/Numerics/LinearAlgebra/Complex32/Vector.cs
index 865296d5..41d93489 100644
--- a/src/Numerics/LinearAlgebra/Complex32/Vector.cs
+++ b/src/Numerics/LinearAlgebra/Complex32/Vector.cs
@@ -4,7 +4,7 @@
// http://github.com/mathnet/mathnet-numerics
// http://mathnetnumerics.codeplex.com
//
-// Copyright (c) 2009-2013 Math.NET
+// Copyright (c) 2009-2015 Math.NET
//
// Permission is hereby granted, free of charge, to any person
// obtaining a copy of this software and associated documentation
@@ -58,6 +58,24 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
MapInplace(x => x.Magnitude < threshold ? Complex32.Zero : x, Zeros.AllowSkip);
}
+ ///
+ /// Conjugates vector and save result to
+ ///
+ /// Target vector
+ protected override void DoConjugate(Vector result)
+ {
+ Map(Complex32.Conjugate, result, Zeros.AllowSkip);
+ }
+
+ ///
+ /// Negates vector and saves result to
+ ///
+ /// Target vector
+ protected override void DoNegate(Vector result)
+ {
+ Map(Complex32.Negate, result, Zeros.AllowSkip);
+ }
+
///
/// Adds a scalar to each element of the vector and stores the result in the result vector.
///
@@ -69,10 +87,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
///
protected override void DoAdd(Complex32 scalar, Vector result)
{
- for (var index = 0; index < Count; index++)
- {
- result.At(index, At(index) + scalar);
- }
+ Map(x => x + scalar, result, Zeros.Include);
}
///
@@ -86,10 +101,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
///
protected override void DoAdd(Vector other, Vector result)
{
- for (var index = 0; index < Count; index++)
- {
- result.At(index, At(index) + other.At(index));
- }
+ Map2(Complex32.Add, other, result, Zeros.AllowSkip);
}
///
@@ -103,7 +115,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
///
protected override void DoSubtract(Complex32 scalar, Vector result)
{
- DoAdd(-scalar, result);
+ Map(x => x - scalar, result, Zeros.Include);
}
///
@@ -117,10 +129,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
///
protected override void DoSubtract(Vector other, Vector result)
{
- for (var index = 0; index < Count; index++)
- {
- result.At(index, At(index) - other.At(index));
- }
+ Map2(Complex32.Subtract, other, result, Zeros.AllowSkip);
}
///
@@ -134,10 +143,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
///
protected override void DoMultiply(Complex32 scalar, Vector result)
{
- for (var index = 0; index < Count; index++)
- {
- result.At(index, At(index) * scalar);
- }
+ Map(x => x*scalar, result, Zeros.AllowSkip);
}
///
@@ -151,7 +157,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
///
protected override void DoDivide(Complex32 divisor, Vector result)
{
- DoMultiply(1 / divisor, result);
+ Map(x => x/divisor, result, divisor.IsZero() ? Zeros.Include : Zeros.AllowSkip);
}
///
@@ -161,10 +167,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
/// The vector to store the result of the division.
protected override void DoDivideByThis(Complex32 dividend, Vector result)
{
- for (var index = 0; index < Count; index++)
- {
- result.At(index, dividend / At(index));
- }
+ Map(x => dividend/x, result, Zeros.Include);
}
///
@@ -174,10 +177,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
/// The vector to store the result of the pointwise multiplication.
protected override void DoPointwiseMultiply(Vector other, Vector result)
{
- for (var index = 0; index < Count; index++)
- {
- result.At(index, At(index) * other.At(index));
- }
+ Map2(Complex32.Multiply, other, result, Zeros.AllowSkip);
}
///
@@ -187,10 +187,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
/// The vector to store the result of the pointwise division.
protected override void DoPointwiseDivide(Vector divisor, Vector result)
{
- for (var index = 0; index < Count; index++)
- {
- result.At(index, At(index) / divisor.At(index));
- }
+ Map2(Complex32.Divide, divisor, result, Zeros.Include);
}
///
@@ -200,7 +197,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
/// The vector to store the result of the pointwise power.
protected override void DoPointwisePower(Complex32 exponent, Vector result)
{
- Map(x => x.Power(exponent), result, Zeros.AllowSkip);
+ Map(x => x.Power(exponent), result, Zeros.Include);
}
///
@@ -448,30 +445,6 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
return Math.Pow(sum, 1.0/p);
}
- ///
- /// Conjugates vector and save result to
- ///
- /// Target vector
- protected override void DoConjugate(Vector result)
- {
- for (var index = 0; index < Count; index++)
- {
- result.At(index, At(index).Conjugate());
- }
- }
-
- ///
- /// Negates vector and saves result to
- ///
- /// Target vector
- protected override void DoNegate(Vector result)
- {
- for (var index = 0; index < Count; index++)
- {
- result.At(index, -At(index));
- }
- }
-
///
/// Returns the index of the absolute maximum element.
///
diff --git a/src/Numerics/LinearAlgebra/Double/Matrix.cs b/src/Numerics/LinearAlgebra/Double/Matrix.cs
index 75dd70a0..edf8c603 100644
--- a/src/Numerics/LinearAlgebra/Double/Matrix.cs
+++ b/src/Numerics/LinearAlgebra/Double/Matrix.cs
@@ -4,7 +4,7 @@
// http://github.com/mathnet/mathnet-numerics
// http://mathnetnumerics.codeplex.com
//
-// Copyright (c) 2009-2013 Math.NET
+// Copyright (c) 2009-2015 Math.NET
//
// Permission is hereby granted, free of charge, to any person
// obtaining a copy of this software and associated documentation
@@ -58,199 +58,36 @@ namespace MathNet.Numerics.LinearAlgebra.Double
MapInplace(x => Math.Abs(x) < threshold ? 0d : x, Zeros.AllowSkip);
}
- /// Calculates the induced L1 norm of this matrix.
- /// The maximum absolute column sum of the matrix.
- public override double L1Norm()
- {
- var norm = 0d;
- for (var j = 0; j < ColumnCount; j++)
- {
- var s = 0d;
- for (var i = 0; i < RowCount; i++)
- {
- s += Math.Abs(At(i, j));
- }
- norm = Math.Max(norm, s);
- }
- return norm;
- }
-
- /// Calculates the induced infinity norm of this matrix.
- /// The maximum absolute row sum of the matrix.
- public override double InfinityNorm()
- {
- var norm = 0d;
- for (var i = 0; i < RowCount; i++)
- {
- var s = 0d;
- for (var j = 0; j < ColumnCount; j++)
- {
- s += Math.Abs(At(i, j));
- }
- norm = Math.Max(norm, s);
- }
- return norm;
- }
-
- /// Calculates the entry-wise Frobenius norm of this matrix.
- /// The square root of the sum of the squared values.
- public override double FrobeniusNorm()
- {
- var transpose = Transpose();
- var aat = this*transpose;
- var norm = 0d;
- for (var i = 0; i < RowCount; i++)
- {
- norm += aat.At(i, i);
- }
- return Math.Sqrt(norm);
- }
-
- ///
- /// Calculates the p-norms of all row vectors.
- /// Typical values for p are 1.0 (L1, Manhattan norm), 2.0 (L2, Euclidean norm) and positive infinity (infinity norm)
- ///
- public override Vector RowNorms(double norm)
- {
- if (norm <= 0.0)
- {
- throw new ArgumentOutOfRangeException("norm", Resources.ArgumentMustBePositive);
- }
-
- var ret = new double[RowCount];
- if (norm == 2.0)
- {
- Storage.FoldByRowUnchecked(ret, (s, x) => s + x*x, (x, c) => Math.Sqrt(x), ret, Zeros.AllowSkip);
- }
- else if (norm == 1.0)
- {
- Storage.FoldByRowUnchecked(ret, (s, x) => s + Math.Abs(x), (x, c) => x, ret, Zeros.AllowSkip);
- }
- else if (double.IsPositiveInfinity(norm))
- {
- Storage.FoldByRowUnchecked(ret, (s, x) => Math.Max(s, Math.Abs(x)), (x, c) => x, ret, Zeros.AllowSkip);
- }
- else
- {
- double invnorm = 1.0/norm;
- Storage.FoldByRowUnchecked(ret, (s, x) => s + Math.Pow(Math.Abs(x), norm), (x, c) => Math.Pow(x, invnorm), ret, Zeros.AllowSkip);
- }
- return Vector.Build.Dense(ret);
- }
-
- ///
- /// Calculates the p-norms of all column vectors.
- /// Typical values for p are 1.0 (L1, Manhattan norm), 2.0 (L2, Euclidean norm) and positive infinity (infinity norm)
- ///
- public override Vector ColumnNorms(double norm)
- {
- if (norm <= 0.0)
- {
- throw new ArgumentOutOfRangeException("norm", Resources.ArgumentMustBePositive);
- }
-
- var ret = new double[ColumnCount];
- if (norm == 2.0)
- {
- Storage.FoldByColumnUnchecked(ret, (s, x) => s + x*x, (x, c) => Math.Sqrt(x), ret, Zeros.AllowSkip);
- }
- else if (norm == 1.0)
- {
- Storage.FoldByColumnUnchecked(ret, (s, x) => s + Math.Abs(x), (x, c) => x, ret, Zeros.AllowSkip);
- }
- else if (double.IsPositiveInfinity(norm))
- {
- Storage.FoldByColumnUnchecked(ret, (s, x) => Math.Max(s, Math.Abs(x)), (x, c) => x, ret, Zeros.AllowSkip);
- }
- else
- {
- double invnorm = 1.0/norm;
- Storage.FoldByColumnUnchecked(ret, (s, x) => s + Math.Pow(Math.Abs(x), norm), (x, c) => Math.Pow(x, invnorm), ret, Zeros.AllowSkip);
- }
- return Vector.Build.Dense(ret);
- }
-
///
- /// Normalizes all row vectors to a unit p-norm.
- /// Typical values for p are 1.0 (L1, Manhattan norm), 2.0 (L2, Euclidean norm) and positive infinity (infinity norm)
+ /// Returns the conjugate transpose of this matrix.
///
- public override sealed Matrix NormalizeRows(double norm)
+ /// The conjugate transpose of this matrix.
+ public override sealed Matrix ConjugateTranspose()
{
- var norminv = ((DenseVectorStorage)RowNorms(norm).Storage).Data;
- for (int i = 0; i < norminv.Length; i++)
- {
- norminv[i] = norminv[i] == 0d ? 1d : 1d/norminv[i];
- }
-
- var result = Build.SameAs(this, RowCount, ColumnCount);
- Storage.MapIndexedTo(result.Storage, (i, j, x) => norminv[i]*x, Zeros.AllowSkip, ExistingData.AssumeZeros);
- return result;
+ return Transpose();
}
///
- /// Normalizes all column vectors to a unit p-norm.
- /// Typical values for p are 1.0 (L1, Manhattan norm), 2.0 (L2, Euclidean norm) and positive infinity (infinity norm)
+ /// Complex conjugates each element of this matrix and place the results into the result matrix.
///
- public override sealed Matrix NormalizeColumns(double norm)
+ /// The result of the conjugation.
+ protected override sealed void DoConjugate(Matrix result)
{
- var norminv = ((DenseVectorStorage)ColumnNorms(norm).Storage).Data;
- for (int i = 0; i < norminv.Length; i++)
+ if (ReferenceEquals(this, result))
{
- norminv[i] = norminv[i] == 0d ? 1d : 1d/norminv[i];
+ return;
}
- var result = Build.SameAs(this, RowCount, ColumnCount);
- Storage.MapIndexedTo(result.Storage, (i, j, x) => norminv[j]*x, Zeros.AllowSkip, ExistingData.AssumeZeros);
- return result;
- }
-
- ///
- /// Calculates the value sum of each row vector.
- ///
- public override Vector RowSums()
- {
- var ret = new double[RowCount];
- Storage.FoldByRowUnchecked(ret, (s, x) => s + x, (x, c) => x, ret, Zeros.AllowSkip);
- return Vector.Build.Dense(ret);
- }
-
- ///
- /// Calculates the absolute value sum of each row vector.
- ///
- public override Vector RowAbsoluteSums()
- {
- var ret = new double[RowCount];
- Storage.FoldByRowUnchecked(ret, (s, x) => s + Math.Abs(x), (x, c) => x, ret, Zeros.AllowSkip);
- return Vector.Build.Dense(ret);
- }
-
- ///
- /// Calculates the value sum of each column vector.
- ///
- public override Vector ColumnSums()
- {
- var ret = new double[ColumnCount];
- Storage.FoldByColumnUnchecked(ret, (s, x) => s + x, (x, c) => x, ret, Zeros.AllowSkip);
- return Vector.Build.Dense(ret);
- }
-
- ///
- /// Calculates the absolute value sum of each column vector.
- ///
- public override Vector ColumnAbsoluteSums()
- {
- var ret = new double[ColumnCount];
- Storage.FoldByColumnUnchecked(ret, (s, x) => s + Math.Abs(x), (x, c) => x, ret, Zeros.AllowSkip);
- return Vector.Build.Dense(ret);
+ CopyTo(result);
}
///
- /// Returns the conjugate transpose of this matrix.
+ /// Negate each element of this matrix and place the results into the result matrix.
///
- /// The conjugate transpose of this matrix.
- public override sealed Matrix ConjugateTranspose()
+ /// The result of the negation.
+ protected override void DoNegate(Matrix result)
{
- return Transpose();
+ Map(x => -x, result, Zeros.AllowSkip);
}
///
@@ -260,13 +97,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double
/// The matrix to store the result of the addition.
protected override void DoAdd(double scalar, Matrix result)
{
- for (var i = 0; i < RowCount; i++)
- {
- for (var j = 0; j < ColumnCount; j++)
- {
- result.At(i, j, At(i, j) + scalar);
- }
- }
+ Map(x => x + scalar, result, Zeros.Include);
}
///
@@ -278,13 +109,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double
/// If the two matrices don't have the same dimensions.
protected override void DoAdd(Matrix other, Matrix result)
{
- for (var i = 0; i < RowCount; i++)
- {
- for (var j = 0; j < ColumnCount; j++)
- {
- result.At(i, j, At(i, j) + other.At(i, j));
- }
- }
+ Map2((x, y) => x + y, other, result, Zeros.AllowSkip);
}
///
@@ -294,13 +119,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double
/// The matrix to store the result of the subtraction.
protected override void DoSubtract(double scalar, Matrix result)
{
- for (var i = 0; i < RowCount; i++)
- {
- for (var j = 0; j < ColumnCount; j++)
- {
- result.At(i, j, At(i, j) - scalar);
- }
- }
+ Map(x => x - scalar, result, Zeros.Include);
}
///
@@ -312,13 +131,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double
/// If the two matrices don't have the same dimensions.
protected override void DoSubtract(Matrix other, Matrix result)
{
- for (var i = 0; i < RowCount; i++)
- {
- for (var j = 0; j < ColumnCount; j++)
- {
- result.At(i, j, At(i, j) - other.At(i, j));
- }
- }
+ Map2((x, y) => x - y, other, result, Zeros.AllowSkip);
}
///
@@ -328,13 +141,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double
/// The matrix to store the result of the multiplication.
protected override void DoMultiply(double scalar, Matrix result)
{
- for (var i = 0; i < RowCount; i++)
- {
- for (var j = 0; j < ColumnCount; j++)
- {
- result.At(i, j, At(i, j)*scalar);
- }
- }
+ Map(x => x*scalar, result, Zeros.AllowSkip);
}
///
@@ -362,23 +169,17 @@ namespace MathNet.Numerics.LinearAlgebra.Double
/// The matrix to store the result of the division.
protected override void DoDivide(double divisor, Matrix result)
{
- DoMultiply(1.0/divisor, result);
+ Map(x => x/divisor, result, divisor == 0.0 ? Zeros.Include : Zeros.AllowSkip);
}
///
/// Divides a scalar by each element of the matrix and stores the result in the result matrix.
///
- /// The scalar to add.
+ /// The scalar to divide by each element of the matrix.
/// The matrix to store the result of the division.
protected override void DoDivideByThis(double dividend, Matrix result)
{
- for (var i = 0; i < RowCount; i++)
- {
- for (var j = 0; j < ColumnCount; j++)
- {
- result.At(i, j, dividend/At(i, j));
- }
- }
+ Map(x => dividend/x, result, Zeros.Include);
}
///
@@ -492,35 +293,6 @@ namespace MathNet.Numerics.LinearAlgebra.Double
DoTransposeThisAndMultiply(rightSide, result);
}
- ///
- /// Negate each element of this matrix and place the results into the result matrix.
- ///
- /// The result of the negation.
- protected override void DoNegate(Matrix result)
- {
- for (var i = 0; i < RowCount; i++)
- {
- for (var j = 0; j < ColumnCount; j++)
- {
- result.At(i, j, -At(i, j));
- }
- }
- }
-
- ///
- /// Complex conjugates each element of this matrix and place the results into the result matrix.
- ///
- /// The result of the conjugation.
- protected override sealed void DoConjugate(Matrix result)
- {
- if (ReferenceEquals(this, result))
- {
- return;
- }
-
- CopyTo(result);
- }
-
///
/// Computes the canonical modulus, where the result has the sign of the divisor,
/// for the given divisor each element of the matrix.
@@ -529,13 +301,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double
/// Matrix to store the results in.
protected override void DoModulus(double divisor, Matrix result)
{
- for (var row = 0; row < RowCount; row++)
- {
- for (var column = 0; column < ColumnCount; column++)
- {
- result.At(row, column, Euclid.Modulus(At(row, column), divisor));
- }
- }
+ Map(x => Euclid.Modulus(x, divisor), result, Zeros.Include);
}
///
@@ -546,13 +312,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double
/// A vector to store the results in.
protected override void DoModulusByThis(double dividend, Matrix result)
{
- for (var row = 0; row < RowCount; row++)
- {
- for (var column = 0; column < ColumnCount; column++)
- {
- result.At(row, column, Euclid.Modulus(dividend, At(row, column)));
- }
- }
+ Map(x => Euclid.Modulus(dividend, x), result, Zeros.Include);
}
///
@@ -563,13 +323,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double
/// Matrix to store the results in.
protected override void DoRemainder(double divisor, Matrix result)
{
- for (var row = 0; row < RowCount; row++)
- {
- for (var column = 0; column < ColumnCount; column++)
- {
- result.At(row, column, At(row, column)%divisor);
- }
- }
+ Map(x => Euclid.Remainder(x, divisor), result, Zeros.Include);
}
///
@@ -580,13 +334,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double
/// A vector to store the results in.
protected override void DoRemainderByThis(double dividend, Matrix result)
{
- for (var row = 0; row < RowCount; row++)
- {
- for (var column = 0; column < ColumnCount; column++)
- {
- result.At(row, column, dividend%At(row, column));
- }
- }
+ Map(x => Euclid.Remainder(dividend, x), result, Zeros.Include);
}
///
@@ -596,13 +344,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double
/// The matrix to store the result of the pointwise multiplication.
protected override void DoPointwiseMultiply(Matrix other, Matrix result)
{
- for (var j = 0; j < ColumnCount; j++)
- {
- for (var i = 0; i < RowCount; i++)
- {
- result.At(i, j, At(i, j)*other.At(i, j));
- }
- }
+ Map2((x, y) => x*y, other, result, Zeros.AllowSkip);
}
///
@@ -612,13 +354,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double
/// The matrix to store the result of the pointwise division.
protected override void DoPointwiseDivide(Matrix divisor, Matrix result)
{
- for (var j = 0; j < ColumnCount; j++)
- {
- for (var i = 0; i < RowCount; i++)
- {
- result.At(i, j, At(i, j)/divisor.At(i, j));
- }
- }
+ Map2((x, y) => x/y, divisor, result, Zeros.Include);
}
///
@@ -628,7 +364,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double
/// The vector to store the result of the pointwise power.
protected override void DoPointwisePower(double exponent, Matrix result)
{
- Map(x => Math.Pow(x, exponent), result, Zeros.AllowSkip);
+ Map(x => Math.Pow(x, exponent), result, exponent > 0.0 ? Zeros.AllowSkip : Zeros.Include);
}
///
@@ -639,13 +375,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double
/// The result of the modulus.
protected override void DoPointwiseModulus(Matrix divisor, Matrix result)
{
- for (var j = 0; j < ColumnCount; j++)
- {
- for (var i = 0; i < RowCount; i++)
- {
- result.At(i, j, Euclid.Modulus(At(i, j), divisor.At(i, j)));
- }
- }
+ Map2(Euclid.Modulus, divisor, result, Zeros.Include);
}
///
@@ -656,13 +386,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double
/// The result of the modulus.
protected override void DoPointwiseRemainder(Matrix divisor, Matrix result)
{
- for (var j = 0; j < ColumnCount; j++)
- {
- for (var i = 0; i < RowCount; i++)
- {
- result.At(i, j, At(i, j)%divisor.At(i, j));
- }
- }
+ Map2(Euclid.Remainder, divisor, result, Zeros.Include);
}
///
@@ -683,6 +407,192 @@ namespace MathNet.Numerics.LinearAlgebra.Double
Map(Math.Log, result, Zeros.Include);
}
+ /// Calculates the induced L1 norm of this matrix.
+ /// The maximum absolute column sum of the matrix.
+ public override double L1Norm()
+ {
+ var norm = 0d;
+ for (var j = 0; j < ColumnCount; j++)
+ {
+ var s = 0d;
+ for (var i = 0; i < RowCount; i++)
+ {
+ s += Math.Abs(At(i, j));
+ }
+ norm = Math.Max(norm, s);
+ }
+ return norm;
+ }
+
+ /// Calculates the induced infinity norm of this matrix.
+ /// The maximum absolute row sum of the matrix.
+ public override double InfinityNorm()
+ {
+ var norm = 0d;
+ for (var i = 0; i < RowCount; i++)
+ {
+ var s = 0d;
+ for (var j = 0; j < ColumnCount; j++)
+ {
+ s += Math.Abs(At(i, j));
+ }
+ norm = Math.Max(norm, s);
+ }
+ return norm;
+ }
+
+ /// Calculates the entry-wise Frobenius norm of this matrix.
+ /// The square root of the sum of the squared values.
+ public override double FrobeniusNorm()
+ {
+ var transpose = Transpose();
+ var aat = this*transpose;
+ var norm = 0d;
+ for (var i = 0; i < RowCount; i++)
+ {
+ norm += aat.At(i, i);
+ }
+ return Math.Sqrt(norm);
+ }
+
+ ///
+ /// Calculates the p-norms of all row vectors.
+ /// Typical values for p are 1.0 (L1, Manhattan norm), 2.0 (L2, Euclidean norm) and positive infinity (infinity norm)
+ ///
+ public override Vector RowNorms(double norm)
+ {
+ if (norm <= 0.0)
+ {
+ throw new ArgumentOutOfRangeException("norm", Resources.ArgumentMustBePositive);
+ }
+
+ var ret = new double[RowCount];
+ if (norm == 2.0)
+ {
+ Storage.FoldByRowUnchecked(ret, (s, x) => s + x*x, (x, c) => Math.Sqrt(x), ret, Zeros.AllowSkip);
+ }
+ else if (norm == 1.0)
+ {
+ Storage.FoldByRowUnchecked(ret, (s, x) => s + Math.Abs(x), (x, c) => x, ret, Zeros.AllowSkip);
+ }
+ else if (double.IsPositiveInfinity(norm))
+ {
+ Storage.FoldByRowUnchecked(ret, (s, x) => Math.Max(s, Math.Abs(x)), (x, c) => x, ret, Zeros.AllowSkip);
+ }
+ else
+ {
+ double invnorm = 1.0/norm;
+ Storage.FoldByRowUnchecked(ret, (s, x) => s + Math.Pow(Math.Abs(x), norm), (x, c) => Math.Pow(x, invnorm), ret, Zeros.AllowSkip);
+ }
+ return Vector.Build.Dense(ret);
+ }
+
+ ///
+ /// Calculates the p-norms of all column vectors.
+ /// Typical values for p are 1.0 (L1, Manhattan norm), 2.0 (L2, Euclidean norm) and positive infinity (infinity norm)
+ ///
+ public override Vector ColumnNorms(double norm)
+ {
+ if (norm <= 0.0)
+ {
+ throw new ArgumentOutOfRangeException("norm", Resources.ArgumentMustBePositive);
+ }
+
+ var ret = new double[ColumnCount];
+ if (norm == 2.0)
+ {
+ Storage.FoldByColumnUnchecked(ret, (s, x) => s + x*x, (x, c) => Math.Sqrt(x), ret, Zeros.AllowSkip);
+ }
+ else if (norm == 1.0)
+ {
+ Storage.FoldByColumnUnchecked(ret, (s, x) => s + Math.Abs(x), (x, c) => x, ret, Zeros.AllowSkip);
+ }
+ else if (double.IsPositiveInfinity(norm))
+ {
+ Storage.FoldByColumnUnchecked(ret, (s, x) => Math.Max(s, Math.Abs(x)), (x, c) => x, ret, Zeros.AllowSkip);
+ }
+ else
+ {
+ double invnorm = 1.0/norm;
+ Storage.FoldByColumnUnchecked(ret, (s, x) => s + Math.Pow(Math.Abs(x), norm), (x, c) => Math.Pow(x, invnorm), ret, Zeros.AllowSkip);
+ }
+ return Vector.Build.Dense(ret);
+ }
+
+ ///
+ /// Normalizes all row vectors to a unit p-norm.
+ /// Typical values for p are 1.0 (L1, Manhattan norm), 2.0 (L2, Euclidean norm) and positive infinity (infinity norm)
+ ///
+ public override sealed Matrix NormalizeRows(double norm)
+ {
+ var norminv = ((DenseVectorStorage)RowNorms(norm).Storage).Data;
+ for (int i = 0; i < norminv.Length; i++)
+ {
+ norminv[i] = norminv[i] == 0d ? 1d : 1d/norminv[i];
+ }
+
+ var result = Build.SameAs(this, RowCount, ColumnCount);
+ Storage.MapIndexedTo(result.Storage, (i, j, x) => norminv[i]*x, Zeros.AllowSkip, ExistingData.AssumeZeros);
+ return result;
+ }
+
+ ///
+ /// Normalizes all column vectors to a unit p-norm.
+ /// Typical values for p are 1.0 (L1, Manhattan norm), 2.0 (L2, Euclidean norm) and positive infinity (infinity norm)
+ ///
+ public override sealed Matrix NormalizeColumns(double norm)
+ {
+ var norminv = ((DenseVectorStorage)ColumnNorms(norm).Storage).Data;
+ for (int i = 0; i < norminv.Length; i++)
+ {
+ norminv[i] = norminv[i] == 0d ? 1d : 1d/norminv[i];
+ }
+
+ var result = Build.SameAs(this, RowCount, ColumnCount);
+ Storage.MapIndexedTo(result.Storage, (i, j, x) => norminv[j]*x, Zeros.AllowSkip, ExistingData.AssumeZeros);
+ return result;
+ }
+
+ ///
+ /// Calculates the value sum of each row vector.
+ ///
+ public override Vector RowSums()
+ {
+ var ret = new double[RowCount];
+ Storage.FoldByRowUnchecked(ret, (s, x) => s + x, (x, c) => x, ret, Zeros.AllowSkip);
+ return Vector.Build.Dense(ret);
+ }
+
+ ///
+ /// Calculates the absolute value sum of each row vector.
+ ///
+ public override Vector RowAbsoluteSums()
+ {
+ var ret = new double[RowCount];
+ Storage.FoldByRowUnchecked(ret, (s, x) => s + Math.Abs(x), (x, c) => x, ret, Zeros.AllowSkip);
+ return Vector.Build.Dense(ret);
+ }
+
+ ///
+ /// Calculates the value sum of each column vector.
+ ///
+ public override Vector ColumnSums()
+ {
+ var ret = new double[ColumnCount];
+ Storage.FoldByColumnUnchecked(ret, (s, x) => s + x, (x, c) => x, ret, Zeros.AllowSkip);
+ return Vector