|
|
|
@ -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); |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>Calculates the induced L1 norm of this matrix.</summary>
|
|
|
|
/// <returns>The maximum absolute column sum of the matrix.</returns>
|
|
|
|
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; |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>Calculates the induced infinity norm of this matrix.</summary>
|
|
|
|
/// <returns>The maximum absolute row sum of the matrix.</returns>
|
|
|
|
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; |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>Calculates the entry-wise Frobenius norm of this matrix.</summary>
|
|
|
|
/// <returns>The square root of the sum of the squared values.</returns>
|
|
|
|
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); |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
/// 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)
|
|
|
|
/// </summary>
|
|
|
|
public override Vector<double> 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<double>.Build.Dense(ret); |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
/// 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)
|
|
|
|
/// </summary>
|
|
|
|
public override Vector<double> 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<double>.Build.Dense(ret); |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
/// 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.
|
|
|
|
/// </summary>
|
|
|
|
public override sealed Matrix<double> NormalizeRows(double norm) |
|
|
|
/// <returns>The conjugate transpose of this matrix.</returns>
|
|
|
|
public override sealed Matrix<double> ConjugateTranspose() |
|
|
|
{ |
|
|
|
var norminv = ((DenseVectorStorage<double>)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(); |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
/// 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.
|
|
|
|
/// </summary>
|
|
|
|
public override sealed Matrix<double> NormalizeColumns(double norm) |
|
|
|
/// <param name="result">The result of the conjugation.</param>
|
|
|
|
protected override sealed void DoConjugate(Matrix<double> result) |
|
|
|
{ |
|
|
|
var norminv = ((DenseVectorStorage<double>)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; |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
/// Calculates the value sum of each row vector.
|
|
|
|
/// </summary>
|
|
|
|
public override Vector<double> RowSums() |
|
|
|
{ |
|
|
|
var ret = new double[RowCount]; |
|
|
|
Storage.FoldByRowUnchecked(ret, (s, x) => s + x, (x, c) => x, ret, Zeros.AllowSkip); |
|
|
|
return Vector<double>.Build.Dense(ret); |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
/// Calculates the absolute value sum of each row vector.
|
|
|
|
/// </summary>
|
|
|
|
public override Vector<double> RowAbsoluteSums() |
|
|
|
{ |
|
|
|
var ret = new double[RowCount]; |
|
|
|
Storage.FoldByRowUnchecked(ret, (s, x) => s + Math.Abs(x), (x, c) => x, ret, Zeros.AllowSkip); |
|
|
|
return Vector<double>.Build.Dense(ret); |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
/// Calculates the value sum of each column vector.
|
|
|
|
/// </summary>
|
|
|
|
public override Vector<double> ColumnSums() |
|
|
|
{ |
|
|
|
var ret = new double[ColumnCount]; |
|
|
|
Storage.FoldByColumnUnchecked(ret, (s, x) => s + x, (x, c) => x, ret, Zeros.AllowSkip); |
|
|
|
return Vector<double>.Build.Dense(ret); |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
/// Calculates the absolute value sum of each column vector.
|
|
|
|
/// </summary>
|
|
|
|
public override Vector<double> ColumnAbsoluteSums() |
|
|
|
{ |
|
|
|
var ret = new double[ColumnCount]; |
|
|
|
Storage.FoldByColumnUnchecked(ret, (s, x) => s + Math.Abs(x), (x, c) => x, ret, Zeros.AllowSkip); |
|
|
|
return Vector<double>.Build.Dense(ret); |
|
|
|
CopyTo(result); |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
/// Returns the conjugate transpose of this matrix.
|
|
|
|
/// Negate each element of this matrix and place the results into the result matrix.
|
|
|
|
/// </summary>
|
|
|
|
/// <returns>The conjugate transpose of this matrix.</returns>
|
|
|
|
public override sealed Matrix<double> ConjugateTranspose() |
|
|
|
/// <param name="result">The result of the negation.</param>
|
|
|
|
protected override void DoNegate(Matrix<double> result) |
|
|
|
{ |
|
|
|
return Transpose(); |
|
|
|
Map(x => -x, result, Zeros.AllowSkip); |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
@ -260,13 +97,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double |
|
|
|
/// <param name="result">The matrix to store the result of the addition.</param>
|
|
|
|
protected override void DoAdd(double scalar, Matrix<double> 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); |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
@ -278,13 +109,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double |
|
|
|
/// <exception cref="ArgumentOutOfRangeException">If the two matrices don't have the same dimensions.</exception>
|
|
|
|
protected override void DoAdd(Matrix<double> other, Matrix<double> 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); |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
@ -294,13 +119,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double |
|
|
|
/// <param name="result">The matrix to store the result of the subtraction.</param>
|
|
|
|
protected override void DoSubtract(double scalar, Matrix<double> 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); |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
@ -312,13 +131,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double |
|
|
|
/// <exception cref="ArgumentOutOfRangeException">If the two matrices don't have the same dimensions.</exception>
|
|
|
|
protected override void DoSubtract(Matrix<double> other, Matrix<double> 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); |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
@ -328,13 +141,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double |
|
|
|
/// <param name="result">The matrix to store the result of the multiplication.</param>
|
|
|
|
protected override void DoMultiply(double scalar, Matrix<double> 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); |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
@ -362,23 +169,17 @@ namespace MathNet.Numerics.LinearAlgebra.Double |
|
|
|
/// <param name="result">The matrix to store the result of the division.</param>
|
|
|
|
protected override void DoDivide(double divisor, Matrix<double> result) |
|
|
|
{ |
|
|
|
DoMultiply(1.0/divisor, result); |
|
|
|
Map(x => x/divisor, result, divisor == 0.0 ? Zeros.Include : Zeros.AllowSkip); |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
/// Divides a scalar by each element of the matrix and stores the result in the result matrix.
|
|
|
|
/// </summary>
|
|
|
|
/// <param name="dividend">The scalar to add.</param>
|
|
|
|
/// <param name="dividend">The scalar to divide by each element of the matrix.</param>
|
|
|
|
/// <param name="result">The matrix to store the result of the division.</param>
|
|
|
|
protected override void DoDivideByThis(double dividend, Matrix<double> 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); |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
@ -492,35 +293,6 @@ namespace MathNet.Numerics.LinearAlgebra.Double |
|
|
|
DoTransposeThisAndMultiply(rightSide, result); |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
/// Negate each element of this matrix and place the results into the result matrix.
|
|
|
|
/// </summary>
|
|
|
|
/// <param name="result">The result of the negation.</param>
|
|
|
|
protected override void DoNegate(Matrix<double> result) |
|
|
|
{ |
|
|
|
for (var i = 0; i < RowCount; i++) |
|
|
|
{ |
|
|
|
for (var j = 0; j < ColumnCount; j++) |
|
|
|
{ |
|
|
|
result.At(i, j, -At(i, j)); |
|
|
|
} |
|
|
|
} |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
/// Complex conjugates each element of this matrix and place the results into the result matrix.
|
|
|
|
/// </summary>
|
|
|
|
/// <param name="result">The result of the conjugation.</param>
|
|
|
|
protected override sealed void DoConjugate(Matrix<double> result) |
|
|
|
{ |
|
|
|
if (ReferenceEquals(this, result)) |
|
|
|
{ |
|
|
|
return; |
|
|
|
} |
|
|
|
|
|
|
|
CopyTo(result); |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
/// 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 |
|
|
|
/// <param name="result">Matrix to store the results in.</param>
|
|
|
|
protected override void DoModulus(double divisor, Matrix<double> 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); |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
@ -546,13 +312,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double |
|
|
|
/// <param name="result">A vector to store the results in.</param>
|
|
|
|
protected override void DoModulusByThis(double dividend, Matrix<double> 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); |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
@ -563,13 +323,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double |
|
|
|
/// <param name="result">Matrix to store the results in.</param>
|
|
|
|
protected override void DoRemainder(double divisor, Matrix<double> 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); |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
@ -580,13 +334,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double |
|
|
|
/// <param name="result">A vector to store the results in.</param>
|
|
|
|
protected override void DoRemainderByThis(double dividend, Matrix<double> 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); |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
@ -596,13 +344,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double |
|
|
|
/// <param name="result">The matrix to store the result of the pointwise multiplication.</param>
|
|
|
|
protected override void DoPointwiseMultiply(Matrix<double> other, Matrix<double> 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); |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
@ -612,13 +354,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double |
|
|
|
/// <param name="result">The matrix to store the result of the pointwise division.</param>
|
|
|
|
protected override void DoPointwiseDivide(Matrix<double> divisor, Matrix<double> 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); |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
@ -628,7 +364,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double |
|
|
|
/// <param name="result">The vector to store the result of the pointwise power.</param>
|
|
|
|
protected override void DoPointwisePower(double exponent, Matrix<double> result) |
|
|
|
{ |
|
|
|
Map(x => Math.Pow(x, exponent), result, Zeros.AllowSkip); |
|
|
|
Map(x => Math.Pow(x, exponent), result, exponent > 0.0 ? Zeros.AllowSkip : Zeros.Include); |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
@ -639,13 +375,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double |
|
|
|
/// <param name="result">The result of the modulus.</param>
|
|
|
|
protected override void DoPointwiseModulus(Matrix<double> divisor, Matrix<double> 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); |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
@ -656,13 +386,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double |
|
|
|
/// <param name="result">The result of the modulus.</param>
|
|
|
|
protected override void DoPointwiseRemainder(Matrix<double> divisor, Matrix<double> 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); |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
@ -683,6 +407,192 @@ namespace MathNet.Numerics.LinearAlgebra.Double |
|
|
|
Map(Math.Log, result, Zeros.Include); |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>Calculates the induced L1 norm of this matrix.</summary>
|
|
|
|
/// <returns>The maximum absolute column sum of the matrix.</returns>
|
|
|
|
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; |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>Calculates the induced infinity norm of this matrix.</summary>
|
|
|
|
/// <returns>The maximum absolute row sum of the matrix.</returns>
|
|
|
|
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; |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>Calculates the entry-wise Frobenius norm of this matrix.</summary>
|
|
|
|
/// <returns>The square root of the sum of the squared values.</returns>
|
|
|
|
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); |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
/// 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)
|
|
|
|
/// </summary>
|
|
|
|
public override Vector<double> 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<double>.Build.Dense(ret); |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
/// 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)
|
|
|
|
/// </summary>
|
|
|
|
public override Vector<double> 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<double>.Build.Dense(ret); |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
/// 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)
|
|
|
|
/// </summary>
|
|
|
|
public override sealed Matrix<double> NormalizeRows(double norm) |
|
|
|
{ |
|
|
|
var norminv = ((DenseVectorStorage<double>)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; |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
/// 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)
|
|
|
|
/// </summary>
|
|
|
|
public override sealed Matrix<double> NormalizeColumns(double norm) |
|
|
|
{ |
|
|
|
var norminv = ((DenseVectorStorage<double>)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; |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
/// Calculates the value sum of each row vector.
|
|
|
|
/// </summary>
|
|
|
|
public override Vector<double> RowSums() |
|
|
|
{ |
|
|
|
var ret = new double[RowCount]; |
|
|
|
Storage.FoldByRowUnchecked(ret, (s, x) => s + x, (x, c) => x, ret, Zeros.AllowSkip); |
|
|
|
return Vector<double>.Build.Dense(ret); |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
/// Calculates the absolute value sum of each row vector.
|
|
|
|
/// </summary>
|
|
|
|
public override Vector<double> RowAbsoluteSums() |
|
|
|
{ |
|
|
|
var ret = new double[RowCount]; |
|
|
|
Storage.FoldByRowUnchecked(ret, (s, x) => s + Math.Abs(x), (x, c) => x, ret, Zeros.AllowSkip); |
|
|
|
return Vector<double>.Build.Dense(ret); |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
/// Calculates the value sum of each column vector.
|
|
|
|
/// </summary>
|
|
|
|
public override Vector<double> ColumnSums() |
|
|
|
{ |
|
|
|
var ret = new double[ColumnCount]; |
|
|
|
Storage.FoldByColumnUnchecked(ret, (s, x) => s + x, (x, c) => x, ret, Zeros.AllowSkip); |
|
|
|
return Vector<double>.Build.Dense(ret); |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
/// Calculates the absolute value sum of each column vector.
|
|
|
|
/// </summary>
|
|
|
|
public override Vector<double> ColumnAbsoluteSums() |
|
|
|
{ |
|
|
|
var ret = new double[ColumnCount]; |
|
|
|
Storage.FoldByColumnUnchecked(ret, (s, x) => s + Math.Abs(x), (x, c) => x, ret, Zeros.AllowSkip); |
|
|
|
return Vector<double>.Build.Dense(ret); |
|
|
|
} |
|
|
|
|
|
|
|
/// <summary>
|
|
|
|
/// Computes the trace of this matrix.
|
|
|
|
/// </summary>
|
|
|
|
|
Christoph Ruegg
12 years ago