diff --git a/src/Numerics/LinearAlgebra/Complex/SparseMatrix.cs b/src/Numerics/LinearAlgebra/Complex/SparseMatrix.cs index fc366b6a..5c8338d4 100644 --- a/src/Numerics/LinearAlgebra/Complex/SparseMatrix.cs +++ b/src/Numerics/LinearAlgebra/Complex/SparseMatrix.cs @@ -379,13 +379,11 @@ namespace MathNet.Numerics.LinearAlgebra.Complex var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var row = 0; row < result.RowCount; row++) { - var startIndex = rowPointers[row]; - var endIndex = row < rowPointers.Length - 1 ? rowPointers[row + 1] : valueCount; - for (var j = startIndex; j < endIndex; j++) + var endIndex = rowPointers[row + 1]; + for (var j = rowPointers[row]; j < endIndex; j++) { if (row >= columnIndices[j]) { @@ -446,13 +444,11 @@ namespace MathNet.Numerics.LinearAlgebra.Complex var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var row = 0; row < result.RowCount; row++) { - var startIndex = rowPointers[row]; - var endIndex = row < rowPointers.Length - 1 ? rowPointers[row + 1] : valueCount; - for (var j = startIndex; j < endIndex; j++) + var endIndex = rowPointers[row + 1]; + for (var j = rowPointers[row]; j < endIndex; j++) { if (row <= columnIndices[j]) { @@ -514,13 +510,11 @@ namespace MathNet.Numerics.LinearAlgebra.Complex var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var row = 0; row < result.RowCount; row++) { - var startIndex = rowPointers[row]; - var endIndex = row < rowPointers.Length - 1 ? rowPointers[row + 1] : valueCount; - for (var j = startIndex; j < endIndex; j++) + var endIndex = rowPointers[row + 1]; + for (var j = rowPointers[row]; j < endIndex; j++) { if (row > columnIndices[j]) { @@ -582,13 +576,11 @@ namespace MathNet.Numerics.LinearAlgebra.Complex var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var row = 0; row < result.RowCount; row++) { - var startIndex = rowPointers[row]; - var endIndex = row < rowPointers.Length - 1 ? rowPointers[row + 1] : valueCount; - for (var j = startIndex; j < endIndex; j++) + var endIndex = rowPointers[row + 1]; + for (var j = rowPointers[row]; j < endIndex; j++) { if (row < columnIndices[j]) { @@ -607,25 +599,21 @@ namespace MathNet.Numerics.LinearAlgebra.Complex var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; var ret = new SparseCompressedRowMatrixStorage(ColumnCount, RowCount) { - ColumnIndices = new int[valueCount], - Values = new Complex[valueCount] + ColumnIndices = new int[_storage.ValueCount], + Values = new Complex[_storage.ValueCount] }; // Do an 'inverse' CopyTo iterate over the rows - for (var i = 0; i < rowPointers.Length; i++) + for (var i = 0; i < RowCount; i++) { - // Get the begin / end index for the current row var startIndex = rowPointers[i]; - var endIndex = i < rowPointers.Length - 1 ? rowPointers[i + 1] : NonZerosCount; + var endIndex = rowPointers[i + 1]; - // Get the values for the current row if (startIndex == endIndex) { - // Begin and end are equal. There are no values in the row, Move to the next row continue; } @@ -645,16 +633,13 @@ namespace MathNet.Numerics.LinearAlgebra.Complex var aat = (SparseCompressedRowMatrixStorage) (this*ConjugateTranspose()).Storage; var norm = 0d; - for (var i = 0; i < aat.RowPointers.Length; i++) + for (var i = 0; i < aat.RowCount; i++) { - // Get the begin / end index for the current row var startIndex = aat.RowPointers[i]; - var endIndex = i < aat.RowPointers.Length - 1 ? aat.RowPointers[i + 1] : aat.ValueCount; + var endIndex = aat.RowPointers[i + 1]; - // Get the values for the current row if (startIndex == endIndex) { - // Begin and end are equal. There are no values in the row, Move to the next row continue; } @@ -677,19 +662,15 @@ namespace MathNet.Numerics.LinearAlgebra.Complex { var rowPointers = _storage.RowPointers; var values = _storage.Values; - var valueCount = _storage.ValueCount; var norm = 0d; - for (var i = 0; i < rowPointers.Length; i++) + for (var i = 0; i < RowCount; i++) { - // Get the begin / end index for the current row var startIndex = rowPointers[i]; - var endIndex = i < rowPointers.Length - 1 ? rowPointers[i + 1] : valueCount; + var endIndex = rowPointers[i + 1]; - // Get the values for the current row if (startIndex == endIndex) { - // Begin and end are equal. There are no values in the row, Move to the next row continue; } @@ -716,21 +697,22 @@ namespace MathNet.Numerics.LinearAlgebra.Complex /// public static SparseMatrix Identity(int order) { - var m = new SparseCompressedRowMatrixStorage(order, order) + var storage = new SparseCompressedRowMatrixStorage(order, order) { - ValueCount = order, Values = new Complex[order], ColumnIndices = new int[order] }; for (var i = 0; i < order; i++) { - m.Values[i] = 1d; - m.ColumnIndices[i] = i; - m.RowPointers[i] = i; + storage.Values[i] = 1d; + storage.ColumnIndices[i] = i; + storage.RowPointers[i] = i; } - return new SparseMatrix(m); + storage.RowPointers[order] = order; + + return new SparseMatrix(storage); } #endregion @@ -781,11 +763,8 @@ namespace MathNet.Numerics.LinearAlgebra.Complex var leftStorage = left._storage; for (var i = 0; i < leftStorage.RowCount; i++) { - // Get the begin / end index for the current row - var startIndex = leftStorage.RowPointers[i]; - var endIndex = i < leftStorage.RowPointers.Length - 1 ? leftStorage.RowPointers[i + 1] : leftStorage.ValueCount; - - for (var j = startIndex; j < endIndex; j++) + var endIndex = leftStorage.RowPointers[i + 1]; + for (var j = leftStorage.RowPointers[i]; j < endIndex; j++) { var columnIndex = leftStorage.ColumnIndices[j]; var resVal = leftStorage.Values[j] + result.At(i, columnIndex); @@ -824,11 +803,8 @@ namespace MathNet.Numerics.LinearAlgebra.Complex { for (var i = 0; i < otherStorage.RowCount; i++) { - // Get the begin / end index for the current row - var startIndex = otherStorage.RowPointers[i]; - var endIndex = i < otherStorage.RowPointers.Length - 1 ? otherStorage.RowPointers[i + 1] : otherStorage.ValueCount; - - for (var j = startIndex; j < endIndex; j++) + var endIndex = otherStorage.RowPointers[i + 1]; + for (var j = otherStorage.RowPointers[i]; j < endIndex; j++) { var columnIndex = otherStorage.ColumnIndices[j]; var resVal = sparseResult.At(i, columnIndex) - otherStorage.Values[j]; @@ -848,15 +824,11 @@ namespace MathNet.Numerics.LinearAlgebra.Complex var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var i = 0; i < RowCount; i++) { - // Get the begin / end index for the current row - var startIndex = rowPointers[i]; - var endIndex = i < rowPointers.Length - 1 ? rowPointers[i + 1] : valueCount; - - for (var j = startIndex; j < endIndex; j++) + var endIndex = rowPointers[i + 1]; + for (var j = rowPointers[i]; j < endIndex; j++) { var columnIndex = columnIndices[j]; var resVal = sparseResult.At(i, columnIndex) + values[j]; @@ -935,13 +907,12 @@ namespace MathNet.Numerics.LinearAlgebra.Complex var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var row = 0; row < RowCount; row++) { - // Get the begin / end index for the current row var startIndex = rowPointers[row]; - var endIndex = row < rowPointers.Length - 1 ? rowPointers[row + 1] : valueCount; + var endIndex = rowPointers[row + 1]; + if (startIndex == endIndex) { continue; @@ -973,13 +944,12 @@ namespace MathNet.Numerics.LinearAlgebra.Complex var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var row = 0; row < RowCount; row++) { - // Get the begin / end index for the current row var startIndex = rowPointers[row]; - var endIndex = row < rowPointers.Length - 1 ? rowPointers[row + 1] : valueCount; + var endIndex = rowPointers[row + 1]; + if (startIndex == endIndex) { continue; @@ -1015,15 +985,14 @@ namespace MathNet.Numerics.LinearAlgebra.Complex var rowPointers = _storage.RowPointers; var values = _storage.Values; - var valueCount = _storage.ValueCount; var otherStorage = otherSparse._storage; for (var j = 0; j < RowCount; j++) { - // Get the begin / end index for the row var startIndexOther = otherStorage.RowPointers[j]; - var endIndexOther = j < otherStorage.RowPointers.Length - 1 ? otherStorage.RowPointers[j + 1] : otherStorage.ValueCount; + var endIndexOther = otherStorage.RowPointers[j + 1]; + if (startIndexOther == endIndexOther) { continue; @@ -1032,9 +1001,10 @@ namespace MathNet.Numerics.LinearAlgebra.Complex for (var i = 0; i < RowCount; i++) { // Multiply row of matrix A on row of matrix B - // Get the begin / end index for the row + var startIndexThis = rowPointers[i]; - var endIndexThis = i < rowPointers.Length - 1 ? rowPointers[i + 1] : valueCount; + var endIndexThis = rowPointers[i + 1]; + if (startIndexThis == endIndexThis) { continue; @@ -1077,15 +1047,11 @@ namespace MathNet.Numerics.LinearAlgebra.Complex var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var i = 0; i < RowCount; i++) { - // Get the begin / end index for the current row - var startIndex = rowPointers[i]; - var endIndex = i < rowPointers.Length - 1 ? rowPointers[i + 1] : valueCount; - - for (var j = startIndex; j < endIndex; j++) + var endIndex = rowPointers[i + 1]; + for (var j = rowPointers[i]; j < endIndex; j++) { var resVal = values[j]*other.At(i, columnIndices[j]); if (!resVal.IsZero()) @@ -1108,15 +1074,11 @@ namespace MathNet.Numerics.LinearAlgebra.Complex var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var i = 0; i < RowCount; i++) { - // Get the begin / end index for the current row - var startIndex = rowPointers[i]; - var endIndex = i < rowPointers.Length - 1 ? rowPointers[i + 1] : valueCount; - - for (var j = startIndex; j < endIndex; j++) + var endIndex = rowPointers[i + 1]; + for (var j = rowPointers[i]; j < endIndex; j++) { if (!values[j].IsZero()) { @@ -1146,15 +1108,11 @@ namespace MathNet.Numerics.LinearAlgebra.Complex var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var i = 0; i < RowCount; i++) { - // Get the begin / end index for the current row - var startIndex = rowPointers[i]; - var endIndex = i < rowPointers.Length - 1 ? rowPointers[i + 1] : valueCount; - - for (var j = startIndex; j < endIndex; j++) + var endIndex = rowPointers[i + 1]; + for (var j = rowPointers[i]; j < endIndex; j++) { if (!values[j].IsZero()) { @@ -1178,7 +1136,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex // todo: we might be able to speed this up by caching one half of the matrix var rowPointers = _storage.RowPointers; - for (var row = 0; row < RowCount - 1; row++) + for (var row = 0; row < RowCount; row++) { var start = rowPointers[row]; var end = rowPointers[row + 1]; @@ -1194,16 +1152,6 @@ namespace MathNet.Numerics.LinearAlgebra.Complex } } - var lastRow = rowPointers.Length - 1; - - if (rowPointers[lastRow] < NonZerosCount) - { - if (!CheckIfOppositesAreEqual(rowPointers[lastRow], _storage.ValueCount, lastRow)) - { - return false; - } - } - return true; } } diff --git a/src/Numerics/LinearAlgebra/Complex32/SparseMatrix.cs b/src/Numerics/LinearAlgebra/Complex32/SparseMatrix.cs index db7951ac..c1efa898 100644 --- a/src/Numerics/LinearAlgebra/Complex32/SparseMatrix.cs +++ b/src/Numerics/LinearAlgebra/Complex32/SparseMatrix.cs @@ -374,13 +374,11 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var row = 0; row < result.RowCount; row++) { - var startIndex = rowPointers[row]; - var endIndex = row < rowPointers.Length - 1 ? rowPointers[row + 1] : valueCount; - for (var j = startIndex; j < endIndex; j++) + var endIndex = rowPointers[row + 1]; + for (var j = rowPointers[row]; j < endIndex; j++) { if (row >= columnIndices[j]) { @@ -441,13 +439,11 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var row = 0; row < result.RowCount; row++) { - var startIndex = rowPointers[row]; - var endIndex = row < rowPointers.Length - 1 ? rowPointers[row + 1] : valueCount; - for (var j = startIndex; j < endIndex; j++) + var endIndex = rowPointers[row + 1]; + for (var j = rowPointers[row]; j < endIndex; j++) { if (row <= columnIndices[j]) { @@ -509,13 +505,11 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var row = 0; row < result.RowCount; row++) { - var startIndex = rowPointers[row]; - var endIndex = row < rowPointers.Length - 1 ? rowPointers[row + 1] : valueCount; - for (var j = startIndex; j < endIndex; j++) + var endIndex = rowPointers[row + 1]; + for (var j = rowPointers[row]; j < endIndex; j++) { if (row > columnIndices[j]) { @@ -577,13 +571,11 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var row = 0; row < result.RowCount; row++) { - var startIndex = rowPointers[row]; - var endIndex = row < rowPointers.Length - 1 ? rowPointers[row + 1] : valueCount; - for (var j = startIndex; j < endIndex; j++) + var endIndex = rowPointers[row + 1]; + for (var j = rowPointers[row]; j < endIndex; j++) { if (row < columnIndices[j]) { @@ -602,25 +594,21 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; var ret = new SparseCompressedRowMatrixStorage(ColumnCount, RowCount) { - ColumnIndices = new int[valueCount], - Values = new Complex32[valueCount] + ColumnIndices = new int[_storage.ValueCount], + Values = new Complex32[_storage.ValueCount] }; // Do an 'inverse' CopyTo iterate over the rows - for (var i = 0; i < rowPointers.Length; i++) + for (var i = 0; i < RowCount; i++) { - // Get the begin / end index for the current row var startIndex = rowPointers[i]; - var endIndex = i < rowPointers.Length - 1 ? rowPointers[i + 1] : NonZerosCount; + var endIndex = rowPointers[i + 1]; - // Get the values for the current row if (startIndex == endIndex) { - // Begin and end are equal. There are no values in the row, Move to the next row continue; } @@ -640,16 +628,13 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 var aat = (SparseCompressedRowMatrixStorage) (this*ConjugateTranspose()).Storage; var norm = 0f; - for (var i = 0; i < aat.RowPointers.Length; i++) + for (var i = 0; i < aat.RowCount; i++) { - // Get the begin / end index for the current row var startIndex = aat.RowPointers[i]; - var endIndex = i < aat.RowPointers.Length - 1 ? aat.RowPointers[i + 1] : aat.ValueCount; + var endIndex = aat.RowPointers[i + 1]; - // Get the values for the current row if (startIndex == endIndex) { - // Begin and end are equal. There are no values in the row, Move to the next row continue; } @@ -672,19 +657,15 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 { var rowPointers = _storage.RowPointers; var values = _storage.Values; - var valueCount = _storage.ValueCount; var norm = 0f; - for (var i = 0; i < rowPointers.Length; i++) + for (var i = 0; i < RowCount; i++) { - // Get the begin / end index for the current row var startIndex = rowPointers[i]; - var endIndex = i < rowPointers.Length - 1 ? rowPointers[i + 1] : valueCount; + var endIndex = rowPointers[i + 1]; - // Get the values for the current row if (startIndex == endIndex) { - // Begin and end are equal. There are no values in the row, Move to the next row continue; } @@ -711,21 +692,22 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 /// public static SparseMatrix Identity(int order) { - var m = new SparseCompressedRowMatrixStorage(order, order) + var storage = new SparseCompressedRowMatrixStorage(order, order) { - ValueCount = order, Values = new Complex32[order], ColumnIndices = new int[order] }; for (var i = 0; i < order; i++) { - m.Values[i] = 1f; - m.ColumnIndices[i] = i; - m.RowPointers[i] = i; + storage.Values[i] = 1f; + storage.ColumnIndices[i] = i; + storage.RowPointers[i] = i; } - return new SparseMatrix(m); + storage.RowPointers[order] = order; + + return new SparseMatrix(storage); } #endregion @@ -776,11 +758,8 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 var leftStorage = left._storage; for (var i = 0; i < leftStorage.RowCount; i++) { - // Get the begin / end index for the current row - var startIndex = leftStorage.RowPointers[i]; - var endIndex = i < leftStorage.RowPointers.Length - 1 ? leftStorage.RowPointers[i + 1] : leftStorage.ValueCount; - - for (var j = startIndex; j < endIndex; j++) + var endIndex = leftStorage.RowPointers[i + 1]; + for (var j = leftStorage.RowPointers[i]; j < endIndex; j++) { var columnIndex = leftStorage.ColumnIndices[j]; var resVal = leftStorage.Values[j] + result.At(i, columnIndex); @@ -818,11 +797,8 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 { for (var i = 0; i < otherStorage.RowCount; i++) { - // Get the begin / end index for the current row - var startIndex = otherStorage.RowPointers[i]; - var endIndex = i < otherStorage.RowPointers.Length - 1 ? otherStorage.RowPointers[i + 1] : otherStorage.ValueCount; - - for (var j = startIndex; j < endIndex; j++) + var endIndex = otherStorage.RowPointers[i + 1]; + for (var j = otherStorage.RowPointers[i]; j < endIndex; j++) { var columnIndex = otherStorage.ColumnIndices[j]; var resVal = sparseResult.At(i, columnIndex) - otherStorage.Values[j]; @@ -842,15 +818,11 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var i = 0; i < RowCount; i++) { - // Get the begin / end index for the current row - var startIndex = rowPointers[i]; - var endIndex = i < rowPointers.Length - 1 ? rowPointers[i + 1] : valueCount; - - for (var j = startIndex; j < endIndex; j++) + var endIndex = rowPointers[i + 1]; + for (var j = rowPointers[i]; j < endIndex; j++) { var columnIndex = columnIndices[j]; var resVal = sparseResult.At(i, columnIndex) + values[j]; @@ -929,13 +901,12 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var row = 0; row < RowCount; row++) { - // Get the begin / end index for the current row var startIndex = rowPointers[row]; - var endIndex = row < rowPointers.Length - 1 ? rowPointers[row + 1] : valueCount; + var endIndex = rowPointers[row + 1]; + if (startIndex == endIndex) { continue; @@ -967,13 +938,12 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var row = 0; row < RowCount; row++) { - // Get the begin / end index for the current row var startIndex = rowPointers[row]; - var endIndex = row < rowPointers.Length - 1 ? rowPointers[row + 1] : valueCount; + var endIndex = rowPointers[row + 1]; + if (startIndex == endIndex) { continue; @@ -1009,15 +979,14 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 var rowPointers = _storage.RowPointers; var values = _storage.Values; - var valueCount = _storage.ValueCount; var otherStorage = otherSparse._storage; for (var j = 0; j < RowCount; j++) { - // Get the begin / end index for the row var startIndexOther = otherStorage.RowPointers[j]; - var endIndexOther = j < otherStorage.RowPointers.Length - 1 ? otherStorage.RowPointers[j + 1] : otherStorage.ValueCount; + var endIndexOther = otherStorage.RowPointers[j + 1]; + if (startIndexOther == endIndexOther) { continue; @@ -1026,9 +995,10 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 for (var i = 0; i < RowCount; i++) { // Multiply row of matrix A on row of matrix B - // Get the begin / end index for the row + var startIndexThis = rowPointers[i]; - var endIndexThis = i < rowPointers.Length - 1 ? rowPointers[i + 1] : valueCount; + var endIndexThis = rowPointers[i + 1]; + if (startIndexThis == endIndexThis) { continue; @@ -1071,15 +1041,11 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var i = 0; i < RowCount; i++) { - // Get the begin / end index for the current row - var startIndex = rowPointers[i]; - var endIndex = i < rowPointers.Length - 1 ? rowPointers[i + 1] : valueCount; - - for (var j = startIndex; j < endIndex; j++) + var endIndex = rowPointers[i + 1]; + for (var j = rowPointers[i]; j < endIndex; j++) { var resVal = values[j]*other.At(i, columnIndices[j]); if (!resVal.IsZero()) @@ -1102,15 +1068,11 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var i = 0; i < RowCount; i++) { - // Get the begin / end index for the current row - var startIndex = rowPointers[i]; - var endIndex = i < rowPointers.Length - 1 ? rowPointers[i + 1] : valueCount; - - for (var j = startIndex; j < endIndex; j++) + var endIndex = rowPointers[i + 1]; + for (var j = rowPointers[i]; j < endIndex; j++) { if (!values[j].IsZero()) { @@ -1140,15 +1102,11 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var i = 0; i < RowCount; i++) { - // Get the begin / end index for the current row - var startIndex = rowPointers[i]; - var endIndex = i < rowPointers.Length - 1 ? rowPointers[i + 1] : valueCount; - - for (var j = startIndex; j < endIndex; j++) + var endIndex = rowPointers[i + 1]; + for (var j = rowPointers[i]; j < endIndex; j++) { if (!values[j].IsZero()) { @@ -1172,7 +1130,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 // todo: we might be able to speed this up by caching one half of the matrix var rowPointers = _storage.RowPointers; - for (var row = 0; row < RowCount - 1; row++) + for (var row = 0; row < RowCount; row++) { var start = rowPointers[row]; var end = rowPointers[row + 1]; @@ -1188,16 +1146,6 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 } } - var lastRow = rowPointers.Length - 1; - - if (rowPointers[lastRow] < NonZerosCount) - { - if (!CheckIfOppositesAreEqual(rowPointers[lastRow], _storage.ValueCount, lastRow)) - { - return false; - } - } - return true; } } diff --git a/src/Numerics/LinearAlgebra/Double/Solvers/MILU0Preconditioner.cs b/src/Numerics/LinearAlgebra/Double/Solvers/MILU0Preconditioner.cs index e14043c3..89899fab 100644 --- a/src/Numerics/LinearAlgebra/Double/Solvers/MILU0Preconditioner.cs +++ b/src/Numerics/LinearAlgebra/Double/Solvers/MILU0Preconditioner.cs @@ -116,15 +116,10 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Solvers // Original matrix compressed sparse row storage. int[] ja = csr.ColumnIndices; double[] a = csr.Values; + int[] ia = csr.RowPointers; - // Make row pointers conform to standard CSR definition (ie. length n+1) - // TODO: We plan to migrate to this definition everywhere - var ia = new int[n + 1]; - Array.Copy(csr.RowPointers, ia, n); - ia[n] = csr.ValueCount; - - _alu = new double[csr.ValueCount + 1]; - _jlu = new int[csr.ValueCount + 1]; + _alu = new double[ia.Length]; + _jlu = new int[ia.Length]; _diag = new int[n]; int code = Compute(n, a, ja, ia, _alu, _jlu, _diag, UseModified); diff --git a/src/Numerics/LinearAlgebra/Double/SparseMatrix.cs b/src/Numerics/LinearAlgebra/Double/SparseMatrix.cs index 736a9cb5..cb34bc75 100644 --- a/src/Numerics/LinearAlgebra/Double/SparseMatrix.cs +++ b/src/Numerics/LinearAlgebra/Double/SparseMatrix.cs @@ -372,13 +372,11 @@ namespace MathNet.Numerics.LinearAlgebra.Double var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var row = 0; row < result.RowCount; row++) { - var startIndex = rowPointers[row]; - var endIndex = row < rowPointers.Length - 1 ? rowPointers[row + 1] : valueCount; - for (var j = startIndex; j < endIndex; j++) + var endIndex = rowPointers[row + 1]; + for (var j = rowPointers[row]; j < endIndex; j++) { if (row >= columnIndices[j]) { @@ -439,13 +437,11 @@ namespace MathNet.Numerics.LinearAlgebra.Double var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var row = 0; row < result.RowCount; row++) { - var startIndex = rowPointers[row]; - var endIndex = row < rowPointers.Length - 1 ? rowPointers[row + 1] : valueCount; - for (var j = startIndex; j < endIndex; j++) + var endIndex = rowPointers[row + 1]; + for (var j = rowPointers[row]; j < endIndex; j++) { if (row <= columnIndices[j]) { @@ -507,13 +503,11 @@ namespace MathNet.Numerics.LinearAlgebra.Double var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var row = 0; row < result.RowCount; row++) { - var startIndex = rowPointers[row]; - var endIndex = row < rowPointers.Length - 1 ? rowPointers[row + 1] : valueCount; - for (var j = startIndex; j < endIndex; j++) + var endIndex = rowPointers[row + 1]; + for (var j = rowPointers[row]; j < endIndex; j++) { if (row > columnIndices[j]) { @@ -575,13 +569,11 @@ namespace MathNet.Numerics.LinearAlgebra.Double var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var row = 0; row < result.RowCount; row++) { - var startIndex = rowPointers[row]; - var endIndex = row < rowPointers.Length - 1 ? rowPointers[row + 1] : valueCount; - for (var j = startIndex; j < endIndex; j++) + var endIndex = rowPointers[row + 1]; + for (var j = rowPointers[row]; j < endIndex; j++) { if (row < columnIndices[j]) { @@ -600,22 +592,19 @@ namespace MathNet.Numerics.LinearAlgebra.Double var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; var ret = new SparseCompressedRowMatrixStorage(ColumnCount, RowCount) { - ColumnIndices = new int[valueCount], - Values = new double[valueCount] + ColumnIndices = new int[_storage.ValueCount], + Values = new double[_storage.ValueCount] }; // Do an 'inverse' CopyTo iterate over the rows - for (var i = 0; i < rowPointers.Length; i++) + for (var i = 0; i < RowCount; i++) { - // Get the begin / end index for the current row var startIndex = rowPointers[i]; - var endIndex = i < rowPointers.Length - 1 ? rowPointers[i + 1] : valueCount; + var endIndex = rowPointers[i + 1]; - // Get the values for the current row if (startIndex == endIndex) { // Begin and end are equal. There are no values in the row, Move to the next row @@ -638,13 +627,11 @@ namespace MathNet.Numerics.LinearAlgebra.Double var aat = (SparseCompressedRowMatrixStorage) (this*Transpose()).Storage; var norm = 0d; - for (var i = 0; i < aat.RowPointers.Length; i++) + for (var i = 0; i < aat.RowCount; i++) { - // Get the begin / end index for the current row var startIndex = aat.RowPointers[i]; - var endIndex = i < aat.RowPointers.Length - 1 ? aat.RowPointers[i + 1] : aat.ValueCount; + var endIndex = aat.RowPointers[i + 1]; - // Get the values for the current row if (startIndex == endIndex) { // Begin and end are equal. There are no values in the row, Move to the next row @@ -669,16 +656,13 @@ namespace MathNet.Numerics.LinearAlgebra.Double { var rowPointers = _storage.RowPointers; var values = _storage.Values; - var valueCount = _storage.ValueCount; var norm = 0d; - for (var i = 0; i < rowPointers.Length; i++) + for (var i = 0; i < RowCount; i++) { - // Get the begin / end index for the current row var startIndex = rowPointers[i]; - var endIndex = i < rowPointers.Length - 1 ? rowPointers[i + 1] : valueCount; + var endIndex = rowPointers[i + 1]; - // Get the values for the current row if (startIndex == endIndex) { // Begin and end are equal. There are no values in the row, Move to the next row @@ -708,21 +692,22 @@ namespace MathNet.Numerics.LinearAlgebra.Double /// public static SparseMatrix Identity(int order) { - var m = new SparseCompressedRowMatrixStorage(order, order) + var storage = new SparseCompressedRowMatrixStorage(order, order) { - ValueCount = order, Values = new double[order], ColumnIndices = new int[order] }; for (var i = 0; i < order; i++) { - m.Values[i] = 1d; - m.ColumnIndices[i] = i; - m.RowPointers[i] = i; + storage.Values[i] = 1d; + storage.ColumnIndices[i] = i; + storage.RowPointers[i] = i; } - return new SparseMatrix(m); + storage.RowPointers[order] = order; + + return new SparseMatrix(storage); } #endregion @@ -773,11 +758,8 @@ namespace MathNet.Numerics.LinearAlgebra.Double var leftStorage = left._storage; for (var i = 0; i < leftStorage.RowCount; i++) { - // Get the begin / end index for the current row - var startIndex = leftStorage.RowPointers[i]; - var endIndex = i < leftStorage.RowPointers.Length - 1 ? leftStorage.RowPointers[i + 1] : leftStorage.ValueCount; - - for (var j = startIndex; j < endIndex; j++) + var endIndex = leftStorage.RowPointers[i + 1]; + for (var j = leftStorage.RowPointers[i]; j < endIndex; j++) { var columnIndex = leftStorage.ColumnIndices[j]; var resVal = leftStorage.Values[j] + result.At(i, columnIndex); @@ -816,11 +798,8 @@ namespace MathNet.Numerics.LinearAlgebra.Double { for (var i = 0; i < otherStorage.RowCount; i++) { - // Get the begin / end index for the current row - var startIndex = otherStorage.RowPointers[i]; - var endIndex = i < otherStorage.RowPointers.Length - 1 ? otherStorage.RowPointers[i + 1] : otherStorage.ValueCount; - - for (var j = startIndex; j < endIndex; j++) + var endIndex = otherStorage.RowPointers[i + 1]; + for (var j = otherStorage.RowPointers[i]; j < endIndex; j++) { var columnIndex = otherStorage.ColumnIndices[j]; var resVal = sparseResult.At(i, columnIndex) - otherStorage.Values[j]; @@ -840,15 +819,11 @@ namespace MathNet.Numerics.LinearAlgebra.Double var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var i = 0; i < RowCount; i++) { - // Get the begin / end index for the current row - var startIndex = rowPointers[i]; - var endIndex = i < rowPointers.Length - 1 ? rowPointers[i + 1] : valueCount; - - for (var j = startIndex; j < endIndex; j++) + var endIndex = rowPointers[i + 1]; + for (var j = rowPointers[i]; j < endIndex; j++) { var columnIndex = columnIndices[j]; var resVal = sparseResult.At(i, columnIndex) + values[j]; @@ -871,7 +846,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double return; } - if (scalar == 0.0 || _storage.ValueCount == 0) + if (scalar == 0.0 || NonZerosCount == 0) { result.Clear(); return; @@ -927,13 +902,12 @@ namespace MathNet.Numerics.LinearAlgebra.Double var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var row = 0; row < RowCount; row++) { - // Get the begin / end index for the current row var startIndex = rowPointers[row]; - var endIndex = row < rowPointers.Length - 1 ? rowPointers[row + 1] : valueCount; + var endIndex = rowPointers[row + 1]; + if (startIndex == endIndex) { continue; @@ -965,13 +939,12 @@ namespace MathNet.Numerics.LinearAlgebra.Double var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var row = 0; row < RowCount; row++) { - // Get the begin / end index for the current row var startIndex = rowPointers[row]; - var endIndex = row < rowPointers.Length - 1 ? rowPointers[row + 1] : valueCount; + var endIndex = rowPointers[row + 1]; + if (startIndex == endIndex) { continue; @@ -1007,15 +980,14 @@ namespace MathNet.Numerics.LinearAlgebra.Double var rowPointers = _storage.RowPointers; var values = _storage.Values; - var valueCount = _storage.ValueCount; var otherStorage = otherSparse._storage; for (var j = 0; j < RowCount; j++) { - // Get the begin / end index for the row var startIndexOther = otherStorage.RowPointers[j]; - var endIndexOther = j < otherStorage.RowPointers.Length - 1 ? otherStorage.RowPointers[j + 1] : otherStorage.ValueCount; + var endIndexOther = otherStorage.RowPointers[j + 1]; + if (startIndexOther == endIndexOther) { continue; @@ -1023,10 +995,9 @@ namespace MathNet.Numerics.LinearAlgebra.Double for (var i = 0; i < RowCount; i++) { - // Multiply row of matrix A on row of matrix B - // Get the begin / end index for the row var startIndexThis = rowPointers[i]; - var endIndexThis = i < rowPointers.Length - 1 ? rowPointers[i + 1] : valueCount; + var endIndexThis = rowPointers[i + 1]; + if (startIndexThis == endIndexThis) { continue; @@ -1069,15 +1040,11 @@ namespace MathNet.Numerics.LinearAlgebra.Double var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var i = 0; i < RowCount; i++) { - // Get the begin / end index for the current row - var startIndex = rowPointers[i]; - var endIndex = i < rowPointers.Length - 1 ? rowPointers[i + 1] : valueCount; - - for (var j = startIndex; j < endIndex; j++) + var endIndex = rowPointers[i + 1]; + for (var j = rowPointers[i]; j < endIndex; j++) { var resVal = values[j]*other.At(i, columnIndices[j]); if (resVal != 0d) @@ -1100,15 +1067,11 @@ namespace MathNet.Numerics.LinearAlgebra.Double var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var i = 0; i < RowCount; i++) { - // Get the begin / end index for the current row - var startIndex = rowPointers[i]; - var endIndex = i < rowPointers.Length - 1 ? rowPointers[i + 1] : valueCount; - - for (var j = startIndex; j < endIndex; j++) + var endIndex = rowPointers[i + 1]; + for (var j = rowPointers[i]; j < endIndex; j++) { if (values[j] != 0d) { @@ -1138,15 +1101,11 @@ namespace MathNet.Numerics.LinearAlgebra.Double var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var i = 0; i < RowCount; i++) { - // Get the begin / end index for the current row - var startIndex = rowPointers[i]; - var endIndex = i < rowPointers.Length - 1 ? rowPointers[i + 1] : valueCount; - - for (var j = startIndex; j < endIndex; j++) + var endIndex = rowPointers[i + 1]; + for (var j = rowPointers[i]; j < endIndex; j++) { if (values[j] != 0d) { @@ -1196,7 +1155,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double // todo: we might be able to speed this up by caching one half of the matrix var rowPointers = _storage.RowPointers; - for (var row = 0; row < RowCount - 1; row++) + for (var row = 0; row < RowCount; row++) { var start = rowPointers[row]; var end = rowPointers[row + 1]; @@ -1212,16 +1171,6 @@ namespace MathNet.Numerics.LinearAlgebra.Double } } - var lastRow = rowPointers.Length - 1; - - if (rowPointers[lastRow] < _storage.ValueCount) - { - if (!CheckIfOppositesAreEqual(rowPointers[lastRow], _storage.ValueCount, lastRow)) - { - return false; - } - } - return true; } } diff --git a/src/Numerics/LinearAlgebra/Single/SparseMatrix.cs b/src/Numerics/LinearAlgebra/Single/SparseMatrix.cs index 595b8b09..97f0e667 100644 --- a/src/Numerics/LinearAlgebra/Single/SparseMatrix.cs +++ b/src/Numerics/LinearAlgebra/Single/SparseMatrix.cs @@ -372,13 +372,11 @@ namespace MathNet.Numerics.LinearAlgebra.Single var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var row = 0; row < result.RowCount; row++) { - var startIndex = rowPointers[row]; - var endIndex = row < rowPointers.Length - 1 ? rowPointers[row + 1] : valueCount; - for (var j = startIndex; j < endIndex; j++) + var endIndex = rowPointers[row + 1]; + for (var j = rowPointers[row]; j < endIndex; j++) { if (row >= columnIndices[j]) { @@ -439,13 +437,11 @@ namespace MathNet.Numerics.LinearAlgebra.Single var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var row = 0; row < result.RowCount; row++) { - var startIndex = rowPointers[row]; - var endIndex = row < rowPointers.Length - 1 ? rowPointers[row + 1] : valueCount; - for (var j = startIndex; j < endIndex; j++) + var endIndex = rowPointers[row + 1]; + for (var j = rowPointers[row]; j < endIndex; j++) { if (row <= columnIndices[j]) { @@ -507,13 +503,11 @@ namespace MathNet.Numerics.LinearAlgebra.Single var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var row = 0; row < result.RowCount; row++) { - var startIndex = rowPointers[row]; - var endIndex = row < rowPointers.Length - 1 ? rowPointers[row + 1] : valueCount; - for (var j = startIndex; j < endIndex; j++) + var endIndex = rowPointers[row + 1]; + for (var j = rowPointers[row]; j < endIndex; j++) { if (row > columnIndices[j]) { @@ -575,13 +569,11 @@ namespace MathNet.Numerics.LinearAlgebra.Single var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var row = 0; row < result.RowCount; row++) { - var startIndex = rowPointers[row]; - var endIndex = row < rowPointers.Length - 1 ? rowPointers[row + 1] : valueCount; - for (var j = startIndex; j < endIndex; j++) + var endIndex = rowPointers[row + 1]; + for (var j = rowPointers[row]; j < endIndex; j++) { if (row < columnIndices[j]) { @@ -600,20 +592,19 @@ namespace MathNet.Numerics.LinearAlgebra.Single var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; var ret = new SparseCompressedRowMatrixStorage(ColumnCount, RowCount) { - ColumnIndices = new int[valueCount], - Values = new float[valueCount] + ColumnIndices = new int[_storage.ValueCount], + Values = new float[_storage.ValueCount] }; // Do an 'inverse' CopyTo iterate over the rows - for (var i = 0; i < rowPointers.Length; i++) + for (var i = 0; i < RowCount; i++) { // Get the begin / end index for the current row var startIndex = rowPointers[i]; - var endIndex = i < rowPointers.Length - 1 ? rowPointers[i + 1] : valueCount; + var endIndex = rowPointers[i + 1]; // Get the values for the current row if (startIndex == endIndex) @@ -638,11 +629,11 @@ namespace MathNet.Numerics.LinearAlgebra.Single var aat = (SparseCompressedRowMatrixStorage) (this*Transpose()).Storage; var norm = 0f; - for (var i = 0; i < aat.RowPointers.Length; i++) + for (var i = 0; i < aat.RowCount; i++) { // Get the begin / end index for the current row var startIndex = aat.RowPointers[i]; - var endIndex = i < aat.RowPointers.Length - 1 ? aat.RowPointers[i + 1] : aat.ValueCount; + var endIndex = aat.RowPointers[i + 1]; // Get the values for the current row if (startIndex == endIndex) @@ -669,14 +660,13 @@ namespace MathNet.Numerics.LinearAlgebra.Single { var rowPointers = _storage.RowPointers; var values = _storage.Values; - var valueCount = _storage.ValueCount; var norm = 0f; - for (var i = 0; i < rowPointers.Length; i++) + for (var i = 0; i < RowCount; i++) { // Get the begin / end index for the current row var startIndex = rowPointers[i]; - var endIndex = i < rowPointers.Length - 1 ? rowPointers[i + 1] : valueCount; + var endIndex = rowPointers[i + 1]; // Get the values for the current row if (startIndex == endIndex) @@ -708,21 +698,22 @@ namespace MathNet.Numerics.LinearAlgebra.Single /// public static SparseMatrix Identity(int order) { - var mStorage = new SparseCompressedRowMatrixStorage(order, order) + var storage = new SparseCompressedRowMatrixStorage(order, order) { - ValueCount = order, Values = new float[order], ColumnIndices = new int[order] }; for (var i = 0; i < order; i++) { - mStorage.Values[i] = 1f; - mStorage.ColumnIndices[i] = i; - mStorage.RowPointers[i] = i; + storage.Values[i] = 1f; + storage.ColumnIndices[i] = i; + storage.RowPointers[i] = i; } - return new SparseMatrix(mStorage); + storage.RowPointers[order] = order; + + return new SparseMatrix(storage); } #endregion @@ -773,11 +764,8 @@ namespace MathNet.Numerics.LinearAlgebra.Single var leftStorage = left._storage; for (var i = 0; i < leftStorage.RowCount; i++) { - // Get the begin / end index for the current row - var startIndex = leftStorage.RowPointers[i]; - var endIndex = i < leftStorage.RowPointers.Length - 1 ? leftStorage.RowPointers[i + 1] : leftStorage.ValueCount; - - for (var j = startIndex; j < endIndex; j++) + var endIndex = leftStorage.RowPointers[i + 1]; + for (var j = leftStorage.RowPointers[i]; j < endIndex; j++) { var columnIndex = leftStorage.ColumnIndices[j]; var resVal = leftStorage.Values[j] + result.At(i, columnIndex); @@ -815,11 +803,8 @@ namespace MathNet.Numerics.LinearAlgebra.Single { for (var i = 0; i < otherStorage.RowCount; i++) { - // Get the begin / end index for the current row - var startIndex = otherStorage.RowPointers[i]; - var endIndex = i < otherStorage.RowPointers.Length - 1 ? otherStorage.RowPointers[i + 1] : otherStorage.ValueCount; - - for (var j = startIndex; j < endIndex; j++) + var endIndex = otherStorage.RowPointers[i + 1]; + for (var j = otherStorage.RowPointers[i]; j < endIndex; j++) { var columnIndex = otherStorage.ColumnIndices[j]; var resVal = sparseResult.At(i, columnIndex) - otherStorage.Values[j]; @@ -839,15 +824,11 @@ namespace MathNet.Numerics.LinearAlgebra.Single var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var i = 0; i < RowCount; i++) { - // Get the begin / end index for the current row - var startIndex = rowPointers[i]; - var endIndex = i < rowPointers.Length - 1 ? rowPointers[i + 1] : valueCount; - - for (var j = startIndex; j < endIndex; j++) + var endIndex = rowPointers[i + 1]; + for (var j = rowPointers[i]; j < endIndex; j++) { var columnIndex = columnIndices[j]; var resVal = sparseResult.At(i, columnIndex) + values[j]; @@ -926,13 +907,12 @@ namespace MathNet.Numerics.LinearAlgebra.Single var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var row = 0; row < RowCount; row++) { - // Get the begin / end index for the current row var startIndex = rowPointers[row]; - var endIndex = row < rowPointers.Length - 1 ? rowPointers[row + 1] : valueCount; + var endIndex = rowPointers[row + 1]; + if (startIndex == endIndex) { continue; @@ -964,13 +944,12 @@ namespace MathNet.Numerics.LinearAlgebra.Single var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var row = 0; row < RowCount; row++) { - // Get the begin / end index for the current row var startIndex = rowPointers[row]; - var endIndex = row < rowPointers.Length - 1 ? rowPointers[row + 1] : valueCount; + var endIndex = rowPointers[row + 1]; + if (startIndex == endIndex) { continue; @@ -1006,15 +985,14 @@ namespace MathNet.Numerics.LinearAlgebra.Single var rowPointers = _storage.RowPointers; var values = _storage.Values; - var valueCount = _storage.ValueCount; var otherStorage = otherSparse._storage; for (var j = 0; j < RowCount; j++) { - // Get the begin / end index for the row var startIndexOther = otherStorage.RowPointers[j]; - var endIndexOther = j < otherStorage.RowPointers.Length - 1 ? otherStorage.RowPointers[j + 1] : otherStorage.ValueCount; + var endIndexOther = otherStorage.RowPointers[j + 1]; + if (startIndexOther == endIndexOther) { continue; @@ -1023,9 +1001,10 @@ namespace MathNet.Numerics.LinearAlgebra.Single for (var i = 0; i < RowCount; i++) { // Multiply row of matrix A on row of matrix B - // Get the begin / end index for the row + var startIndexThis = rowPointers[i]; - var endIndexThis = i < rowPointers.Length - 1 ? rowPointers[i + 1] : valueCount; + var endIndexThis = rowPointers[i + 1]; + if (startIndexThis == endIndexThis) { continue; @@ -1068,15 +1047,11 @@ namespace MathNet.Numerics.LinearAlgebra.Single var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var i = 0; i < RowCount; i++) { - // Get the begin / end index for the current row - var startIndex = rowPointers[i]; - var endIndex = i < rowPointers.Length - 1 ? rowPointers[i + 1] : valueCount; - - for (var j = startIndex; j < endIndex; j++) + var endIndex = rowPointers[i + 1]; + for (var j = rowPointers[i]; j < endIndex; j++) { var resVal = values[j]*other.At(i, columnIndices[j]); if (resVal != 0f) @@ -1099,15 +1074,11 @@ namespace MathNet.Numerics.LinearAlgebra.Single var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var i = 0; i < RowCount; i++) { - // Get the begin / end index for the current row - var startIndex = rowPointers[i]; - var endIndex = i < rowPointers.Length - 1 ? rowPointers[i + 1] : valueCount; - - for (var j = startIndex; j < endIndex; j++) + var endIndex = rowPointers[i + 1]; + for (var j = rowPointers[i]; j < endIndex; j++) { if (values[j] != 0f) { @@ -1137,15 +1108,11 @@ namespace MathNet.Numerics.LinearAlgebra.Single var rowPointers = _storage.RowPointers; var columnIndices = _storage.ColumnIndices; var values = _storage.Values; - var valueCount = _storage.ValueCount; for (var i = 0; i < RowCount; i++) { - // Get the begin / end index for the current row - var startIndex = rowPointers[i]; - var endIndex = i < rowPointers.Length - 1 ? rowPointers[i + 1] : valueCount; - - for (var j = startIndex; j < endIndex; j++) + var endIndex = rowPointers[i + 1]; + for (var j = rowPointers[i]; j < endIndex; j++) { if (values[j] != 0f) { @@ -1195,7 +1162,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single // todo: we might be able to speed this up by caching one half of the matrix var rowPointers = _storage.RowPointers; - for (var row = 0; row < RowCount - 1; row++) + for (var row = 0; row < RowCount; row++) { var start = rowPointers[row]; var end = rowPointers[row + 1]; @@ -1211,16 +1178,6 @@ namespace MathNet.Numerics.LinearAlgebra.Single } } - var lastRow = rowPointers.Length - 1; - - if (rowPointers[lastRow] < NonZerosCount) - { - if (!CheckIfOppositesAreEqual(rowPointers[lastRow], NonZerosCount, lastRow)) - { - return false; - } - } - return true; } } diff --git a/src/Numerics/LinearAlgebra/Storage/SparseCompressedRowMatrixStorage.cs b/src/Numerics/LinearAlgebra/Storage/SparseCompressedRowMatrixStorage.cs index 6e97cc20..dc109b9d 100644 --- a/src/Numerics/LinearAlgebra/Storage/SparseCompressedRowMatrixStorage.cs +++ b/src/Numerics/LinearAlgebra/Storage/SparseCompressedRowMatrixStorage.cs @@ -42,14 +42,16 @@ namespace MathNet.Numerics.LinearAlgebra.Storage // [ruegg] public fields are OK here /// - /// The array containing the row indices of the existing rows. Element "j" of the array gives the index of the - /// element in the array that is first non-zero element in a row "j" + /// The array containing the row indices of the existing rows. Element "i" of the array gives the index of the + /// element in the array that is first non-zero element in a row "i". + /// The last value is equal to ValueCount, so that the number of non-zero entries in row "i" is always + /// given by RowPointers[i+i] - RowPointers[i]. This array thus has length RowCount+1. /// public readonly int[] RowPointers; /// - /// An array containing the column indices of the non-zero values. Element "I" of the array - /// is the number of the column in matrix that contains the I-th value in the array. + /// An array containing the column indices of the non-zero values. Element "j" of the array + /// is the number of the column in matrix that contains the j-th value in the array. /// public int[] ColumnIndices; @@ -63,15 +65,17 @@ namespace MathNet.Numerics.LinearAlgebra.Storage /// Gets the number of non zero elements in the matrix. /// /// The number of non zero elements. - public int ValueCount; + public int ValueCount + { + get { return RowPointers[RowCount]; } + } internal SparseCompressedRowMatrixStorage(int rows, int columns) : base(rows, columns) { - RowPointers = new int[rows]; + RowPointers = new int[rows + 1]; ColumnIndices = new int[0]; Values = new T[0]; - ValueCount = 0; } /// @@ -152,34 +156,32 @@ namespace MathNet.Numerics.LinearAlgebra.Storage } index = ~index; + var valueCount = RowPointers[RowPointers.Length - 1]; // Check if the storage needs to be increased - if ((ValueCount == Values.Length) && (ValueCount < ((long)RowCount * ColumnCount))) + if ((valueCount == Values.Length) && (valueCount < ((long)RowCount * ColumnCount))) { // Value array is completely full so we increase the size // Determine the increase in size. We will not grow beyond the size of the matrix - var size = Math.Min(Values.Length + GrowthSize(), (long)RowCount * ColumnCount); + var size = Math.Min(Values.Length + GrowthSize(), (long) RowCount*ColumnCount); if (size > int.MaxValue) { throw new NotSupportedException(Resources.TooManyElements); } - Array.Resize(ref Values, (int)size); - Array.Resize(ref ColumnIndices, (int)size); + Array.Resize(ref Values, (int) size); + Array.Resize(ref ColumnIndices, (int) size); } // Move all values (with a position larger than index) in the value array to the next position // move all values (with a position larger than index) in the columIndices array to the next position - Array.Copy(Values, index, Values, index + 1, ValueCount - index); - Array.Copy(ColumnIndices, index, ColumnIndices, index + 1, ValueCount - index); + Array.Copy(Values, index, Values, index + 1, valueCount - index); + Array.Copy(ColumnIndices, index, ColumnIndices, index + 1, valueCount - index); // Add the value and the column index Values[index] = value; ColumnIndices[index] = column; - // increase the number of non-zero numbers by one - ValueCount += 1; - // add 1 to all the row indices for rows bigger than rowIndex // so that they point to the correct part of the value array again. for (var i = row + 1; i < RowPointers.Length; i++) @@ -197,10 +199,12 @@ namespace MathNet.Numerics.LinearAlgebra.Storage /// WARNING: This method is not thread safe. Use "lock" with it and be sure to avoid deadlocks void RemoveAtIndexUnchecked(int itemIndex, int row) { + var valueCount = RowPointers[RowPointers.Length - 1]; + // Move all values (with a position larger than index) in the value array to the previous position // move all values (with a position larger than index) in the columIndices array to the previous position - Array.Copy(Values, itemIndex + 1, Values, itemIndex, ValueCount - itemIndex - 1); - Array.Copy(ColumnIndices, itemIndex + 1, ColumnIndices, itemIndex, ValueCount - itemIndex - 1); + Array.Copy(Values, itemIndex + 1, Values, itemIndex, valueCount - itemIndex - 1); + Array.Copy(ColumnIndices, itemIndex + 1, ColumnIndices, itemIndex, valueCount - itemIndex - 1); // Decrease value in Row for (var i = row + 1; i < RowPointers.Length; i++) @@ -208,14 +212,14 @@ namespace MathNet.Numerics.LinearAlgebra.Storage RowPointers[i] -= 1; } - ValueCount -= 1; + valueCount -= 1; // Check whether we need to shrink the arrays. This is reasonable to do if // there are a lot of non-zero elements and storage is two times bigger - if ((ValueCount > 1024) && (ValueCount < Values.Length / 2)) + if ((valueCount > 1024) && (valueCount < Values.Length / 2)) { - Array.Resize(ref Values, ValueCount); - Array.Resize(ref ColumnIndices, ValueCount); + Array.Resize(ref Values, valueCount); + Array.Resize(ref ColumnIndices, valueCount); } } @@ -229,9 +233,7 @@ namespace MathNet.Numerics.LinearAlgebra.Storage public int FindItem(int row, int column) { // Determin bounds in columnIndices array where this item should be searched (using rowIndex) - var startIndex = RowPointers[row]; - var endIndex = row < RowPointers.Length - 1 ? RowPointers[row + 1] : ValueCount; - return Array.BinarySearch(ColumnIndices, startIndex, endIndex - startIndex, column); + return Array.BinarySearch(ColumnIndices, RowPointers[row], RowPointers[row + 1] - RowPointers[row], column); } /// @@ -244,7 +246,7 @@ namespace MathNet.Numerics.LinearAlgebra.Storage int delta; if (Values.Length > 1024) { - delta = Values.Length / 4; + delta = Values.Length/4; } else { @@ -263,7 +265,6 @@ namespace MathNet.Numerics.LinearAlgebra.Storage public override void Clear() { - ValueCount = 0; Array.Clear(RowPointers, 0, RowPointers.Length); } @@ -275,10 +276,12 @@ namespace MathNet.Numerics.LinearAlgebra.Storage return; } + var valueCount = RowPointers[RowPointers.Length - 1]; + for (int row = rowIndex + rowCount - 1; row >= rowIndex; row--) { var startIndex = RowPointers[row]; - var endIndex = row < RowPointers.Length - 1 ? RowPointers[row + 1] : ValueCount; + var endIndex = RowPointers[row + 1]; // empty row if (startIndex == endIndex) @@ -297,8 +300,8 @@ namespace MathNet.Numerics.LinearAlgebra.Storage { // Move all values (with a position larger than index) in the value array to the previous position // move all values (with a position larger than index) in the columIndices array to the previous position - Array.Copy(Values, first + count, Values, first, ValueCount - first - count); - Array.Copy(ColumnIndices, first + count, ColumnIndices, first, ValueCount - first - count); + Array.Copy(Values, first + count, Values, first, valueCount - first - count); + Array.Copy(ColumnIndices, first + count, ColumnIndices, first, valueCount - first - count); // Decrease value in Row for (var k = row + 1; k < RowPointers.Length; k++) @@ -306,16 +309,16 @@ namespace MathNet.Numerics.LinearAlgebra.Storage RowPointers[k] -= count; } - ValueCount -= count; + valueCount -= count; } } // Check whether we need to shrink the arrays. This is reasonable to do if // there are a lot of non-zero elements and storage is two times bigger - if ((ValueCount > 1024) && (ValueCount < Values.Length / 2)) + if ((valueCount > 1024) && (valueCount < Values.Length / 2)) { - Array.Resize(ref Values, ValueCount); - Array.Resize(ref ColumnIndices, ValueCount); + Array.Resize(ref Values, valueCount); + Array.Resize(ref ColumnIndices, valueCount); } } @@ -350,7 +353,7 @@ namespace MathNet.Numerics.LinearAlgebra.Storage if (ValueCount != sparse.ValueCount) { - // TODO: this is not always correct + // TODO: this is only correct if normalized return false; } @@ -382,7 +385,7 @@ namespace MathNet.Numerics.LinearAlgebra.Storage { for (var i = 0; i < hashNum; i++) { - hash = hash * 31 + values[i].GetHashCode(); + hash = hash*31 + values[i].GetHashCode(); } } return hash; @@ -418,9 +421,9 @@ namespace MathNet.Numerics.LinearAlgebra.Storage } } + rowPointers[rows] = values.Count; storage.ColumnIndices = columnIndices.ToArray(); storage.Values = values.ToArray(); - storage.ValueCount = values.Count; return storage; } @@ -442,9 +445,9 @@ namespace MathNet.Numerics.LinearAlgebra.Storage } } + rowPointers[rows] = values.Count; storage.ColumnIndices = columnIndices.ToArray(); storage.Values = values.ToArray(); - storage.ValueCount = values.Count; return storage; } @@ -460,7 +463,7 @@ namespace MathNet.Numerics.LinearAlgebra.Storage rowPointers[row] = values.Count; for (int col = 0; col < storage.ColumnCount; col++) { - if (!Zero.Equals(array[row,col])) + if (!Zero.Equals(array[row, col])) { values.Add(array[row, col]); columnIndices.Add(col); @@ -468,9 +471,9 @@ namespace MathNet.Numerics.LinearAlgebra.Storage } } + rowPointers[storage.RowCount] = values.Count; storage.ColumnIndices = columnIndices.ToArray(); storage.Values = values.ToArray(); - storage.ValueCount = values.Count; return storage; } @@ -495,9 +498,9 @@ namespace MathNet.Numerics.LinearAlgebra.Storage } } + rowPointers[storage.RowCount] = values.Count; storage.ColumnIndices = columnIndices.ToArray(); storage.Values = values.ToArray(); - storage.ValueCount = values.Count; return storage; } @@ -522,9 +525,9 @@ namespace MathNet.Numerics.LinearAlgebra.Storage } } + rowPointers[storage.RowCount] = values.Count; storage.ColumnIndices = columnIndices.ToArray(); storage.Values = values.ToArray(); - storage.ValueCount = values.Count; return storage; } @@ -551,9 +554,9 @@ namespace MathNet.Numerics.LinearAlgebra.Storage } } + rowPointers[storage.RowCount] = values.Count; storage.ColumnIndices = columnIndices.ToArray(); storage.Values = values.ToArray(); - storage.ValueCount = values.Count; return storage; } @@ -579,9 +582,9 @@ namespace MathNet.Numerics.LinearAlgebra.Storage } } + rowPointers[storage.RowCount] = values.Count; storage.ColumnIndices = columnIndices.ToArray(); storage.Values = values.ToArray(); - storage.ValueCount = values.Count; return storage; } @@ -619,9 +622,9 @@ namespace MathNet.Numerics.LinearAlgebra.Storage } } + rowPointers[rows] = values.Count; storage.ColumnIndices = columnIndices.ToArray(); storage.Values = values.ToArray(); - storage.ValueCount = values.Count; return storage; } @@ -654,9 +657,10 @@ namespace MathNet.Numerics.LinearAlgebra.Storage } if (rowIterator.MoveNext()) throw new ArgumentOutOfRangeException("data", string.Format(Resources.ArgumentArrayWrongLength, rows)); } + + rowPointers[rows] = values.Count; storage.ColumnIndices = columnIndices.ToArray(); storage.Values = values.ToArray(); - storage.ValueCount = values.Count; return storage; } @@ -705,9 +709,9 @@ namespace MathNet.Numerics.LinearAlgebra.Storage } } + rowPointers[rows] = values.Count; storage.ColumnIndices = columnIndices.ToArray(); storage.Values = values.ToArray(); - storage.ValueCount = values.Count; return storage; } @@ -735,15 +739,15 @@ namespace MathNet.Numerics.LinearAlgebra.Storage } } + rowPointers[rows] = values.Count; storage.ColumnIndices = columnIndices.ToArray(); storage.Values = values.ToArray(); - storage.ValueCount = values.Count; return storage; } public static SparseCompressedRowMatrixStorage OfColumnMajorList(int rows, int columns, IList data) { - if (rows * columns != data.Count) + if (rows*columns != data.Count) { throw new ArgumentOutOfRangeException(Resources.ArgumentMatrixDimensions); } @@ -767,9 +771,9 @@ namespace MathNet.Numerics.LinearAlgebra.Storage } } + rowPointers[rows] = values.Count; storage.ColumnIndices = columnIndices.ToArray(); storage.Values = values.ToArray(); - storage.ValueCount = values.Count; return storage; } @@ -803,7 +807,7 @@ namespace MathNet.Numerics.LinearAlgebra.Storage for (int row = 0; row < RowCount; row++) { var startIndex = RowPointers[row]; - var endIndex = row < RowPointers.Length - 1 ? RowPointers[row + 1] : ValueCount; + var endIndex = RowPointers[row + 1]; for (var j = startIndex; j < endIndex; j++) { target.At(row, ColumnIndices[j], Values[j]); @@ -814,15 +818,14 @@ namespace MathNet.Numerics.LinearAlgebra.Storage void CopyToUnchecked(SparseCompressedRowMatrixStorage target) { - target.ValueCount = ValueCount; target.Values = new T[ValueCount]; target.ColumnIndices = new int[ValueCount]; if (ValueCount != 0) { Array.Copy(Values, target.Values, ValueCount); - Buffer.BlockCopy(ColumnIndices, 0, target.ColumnIndices, 0, ValueCount * Constants.SizeOfInt); - Buffer.BlockCopy(RowPointers, 0, target.RowPointers, 0, RowCount * Constants.SizeOfInt); + Buffer.BlockCopy(ColumnIndices, 0, target.ColumnIndices, 0, ValueCount*Constants.SizeOfInt); + Buffer.BlockCopy(RowPointers, 0, target.RowPointers, 0, (RowCount + 1)*Constants.SizeOfInt); } } @@ -838,7 +841,7 @@ namespace MathNet.Numerics.LinearAlgebra.Storage for (int row = 0; row < RowCount; row++) { var startIndex = RowPointers[row]; - var endIndex = row < RowPointers.Length - 1 ? RowPointers[row + 1] : ValueCount; + var endIndex = RowPointers[row + 1]; for (var j = startIndex; j < endIndex; j++) { target.At(row, ColumnIndices[j], Values[j]); @@ -860,7 +863,7 @@ namespace MathNet.Numerics.LinearAlgebra.Storage var sparseTarget = target as SparseCompressedRowMatrixStorage; if (sparseTarget != null) { - CopySubMatrixToUnchecked(sparseTarget, + CopySubMatrixToUnchecked(sparseTarget, sourceRowIndex, targetRowIndex, rowCount, sourceColumnIndex, targetColumnIndex, columnCount, skipClearing); @@ -877,7 +880,7 @@ namespace MathNet.Numerics.LinearAlgebra.Storage for (int i = sourceRowIndex, row = 0; i < sourceRowIndex + rowCount; i++, row++) { var startIndex = RowPointers[i]; - var endIndex = i < RowPointers.Length - 1 ? RowPointers[i + 1] : ValueCount; + var endIndex = RowPointers[i + 1]; for (int j = startIndex; j < endIndex; j++) { @@ -913,7 +916,7 @@ namespace MathNet.Numerics.LinearAlgebra.Storage rowPointers[i + rowOffset] = values.Count; var startIndex = RowPointers[i]; - var endIndex = i < RowPointers.Length - 1 ? RowPointers[i + 1] : ValueCount; + var endIndex = RowPointers[i + 1]; // note: we might be able to replace this loop with Array.Copy (perf) for (int j = startIndex; j < endIndex; j++) @@ -927,12 +930,12 @@ namespace MathNet.Numerics.LinearAlgebra.Storage } } - for(int i=targetRowIndex + rowCount; i(row, col, Values[k++]) : new Tuple(row, col, Zero); } @@ -1087,7 +1090,7 @@ namespace MathNet.Numerics.LinearAlgebra.Storage for (int row = 0; row < RowCount; row++) { var startIndex = RowPointers[row]; - var endIndex = row < RowPointers.Length - 1 ? RowPointers[row + 1] : ValueCount; + var endIndex = RowPointers[row + 1]; for (var j = startIndex; j < endIndex; j++) { if (!Zero.Equals(Values[j])) @@ -1102,7 +1105,7 @@ namespace MathNet.Numerics.LinearAlgebra.Storage public override void MapInplace(Func f, bool forceMapZeros = false) { - var newRowPointers = new int[RowCount]; + var newRowPointers = new int[RowCount+1]; var newColumnIndices = new List(); var newValues = new List(); @@ -1114,9 +1117,7 @@ namespace MathNet.Numerics.LinearAlgebra.Storage newRowPointers[row] = newValues.Count; for (int col = 0; col < ColumnCount; col++) { - var item = k < (row < RowPointers.Length - 1 ? RowPointers[row + 1] : ValueCount) && (ColumnIndices[k]) == col - ? f(Values[k++]) - : f(Zero); + var item = k < RowPointers[row + 1] && ColumnIndices[k] == col ? f(Values[k++]) : f(Zero); if (!Zero.Equals(item)) { newValues.Add(item); @@ -1131,7 +1132,7 @@ namespace MathNet.Numerics.LinearAlgebra.Storage { newRowPointers[row] = newValues.Count; var startIndex = RowPointers[row]; - var endIndex = row < RowPointers.Length - 1 ? RowPointers[row + 1] : ValueCount; + var endIndex = RowPointers[row + 1]; for (var j = startIndex; j < endIndex; j++) { var item = f(Values[j]); @@ -1146,17 +1147,17 @@ namespace MathNet.Numerics.LinearAlgebra.Storage ColumnIndices = newColumnIndices.ToArray(); Values = newValues.ToArray(); - ValueCount = newValues.Count; - Array.Copy(newRowPointers, RowPointers, RowCount); + newRowPointers[RowCount] = newValues.Count; + Array.Copy(newRowPointers, RowPointers, newRowPointers.Length); } public override void MapIndexedInplace(Func f, bool forceMapZeros = false) { - var newRowPointers = new int[RowCount]; + var newRowPointers = new int[RowCount+1]; var newColumnIndices = new List(); var newValues = new List(); - if (forceMapZeros || !Zero.Equals(f(0,0,Zero))) + if (forceMapZeros || !Zero.Equals(f(0, 0, Zero))) { int k = 0; for (int row = 0; row < RowCount; row++) @@ -1164,9 +1165,7 @@ namespace MathNet.Numerics.LinearAlgebra.Storage newRowPointers[row] = newValues.Count; for (int col = 0; col < ColumnCount; col++) { - var item = k < (row < RowPointers.Length - 1 ? RowPointers[row + 1] : ValueCount) && (ColumnIndices[k]) == col - ? f(row, col, Values[k++]) - : f(row, col, Zero); + var item = k < RowPointers[row + 1] && ColumnIndices[k] == col ? f(row, col, Values[k++]) : f(row, col, Zero); if (!Zero.Equals(item)) { newValues.Add(item); @@ -1181,7 +1180,7 @@ namespace MathNet.Numerics.LinearAlgebra.Storage { newRowPointers[row] = newValues.Count; var startIndex = RowPointers[row]; - var endIndex = row < RowPointers.Length - 1 ? RowPointers[row + 1] : ValueCount; + var endIndex = RowPointers[row + 1]; for (var j = startIndex; j < endIndex; j++) { var item = f(row, ColumnIndices[j], Values[j]); @@ -1196,8 +1195,8 @@ namespace MathNet.Numerics.LinearAlgebra.Storage ColumnIndices = newColumnIndices.ToArray(); Values = newValues.ToArray(); - ValueCount = newValues.Count; - Array.Copy(newRowPointers, RowPointers, RowCount); + newRowPointers[RowCount] = newValues.Count; + Array.Copy(newRowPointers, RowPointers, newRowPointers.Length); } } }