diff --git a/src/Numerics/LinearAlgebra/Complex/DiagonalMatrix.cs b/src/Numerics/LinearAlgebra/Complex/DiagonalMatrix.cs index cfe8a807..35ec6452 100644 --- a/src/Numerics/LinearAlgebra/Complex/DiagonalMatrix.cs +++ b/src/Numerics/LinearAlgebra/Complex/DiagonalMatrix.cs @@ -1174,68 +1174,11 @@ namespace MathNet.Numerics.LinearAlgebra.Complex /// is not positive. public override Matrix SubMatrix(int rowIndex, int rowCount, int columnIndex, int columnCount) { - if (rowIndex >= RowCount || rowIndex < 0) - { - throw new ArgumentOutOfRangeException("rowIndex"); - } - - if (columnIndex >= ColumnCount || columnIndex < 0) - { - throw new ArgumentOutOfRangeException("columnIndex"); - } - - if (rowCount < 1) - { - throw new ArgumentException(Resources.ArgumentMustBePositive, "rowCount"); - } - - if (columnCount < 1) - { - throw new ArgumentException(Resources.ArgumentMustBePositive, "columnCount"); - } - - var colMax = columnIndex + columnCount; - var rowMax = rowIndex + rowCount; - - if (rowMax > RowCount) - { - throw new ArgumentOutOfRangeException("rowCount"); - } + // TODO: if rowIndex == columnIndex, use a diagonal matrix instead of a sparse one - if (colMax > ColumnCount) - { - throw new ArgumentOutOfRangeException("columnCount"); - } - - var result = new SparseMatrix(rowCount, columnCount); - - if (rowIndex > columnIndex && columnIndex + columnCount > rowIndex) - { - int columnInit = rowIndex - columnIndex; - int end = Math.Min(columnCount, rowCount + columnInit); - for (var i = 0; columnInit + i < end; i++) - { - result[i, columnInit + i] = _data[rowIndex + i]; - } - } - else if (rowIndex < columnIndex && rowIndex + rowCount > columnIndex) - { - int rowInit = columnIndex - rowIndex; - int end = Math.Min(columnCount + rowInit, rowCount); - for (var i = 0; rowInit + i < end; i++) - { - result[rowInit + i, i] = _data[columnIndex + i]; - } - } - else - { - for (var i = 0; i < Math.Min(columnCount, rowCount); i++) - { - result[i, i] = _data[rowIndex + i]; - } - } - - return result; + var storage = new SparseCompressedRowMatrixStorage(rowCount, columnCount, Complex.Zero); + _storage.CopySubMatrixTo(storage, rowIndex, 0, rowCount, columnIndex, 0, columnCount, true); + return new SparseMatrix(storage); } /// diff --git a/src/Numerics/LinearAlgebra/Complex/SparseMatrix.cs b/src/Numerics/LinearAlgebra/Complex/SparseMatrix.cs index 415891c5..b2eae281 100644 --- a/src/Numerics/LinearAlgebra/Complex/SparseMatrix.cs +++ b/src/Numerics/LinearAlgebra/Complex/SparseMatrix.cs @@ -56,6 +56,12 @@ namespace MathNet.Numerics.LinearAlgebra.Complex get { return _storage.ValueCount; } } + internal SparseMatrix(SparseCompressedRowMatrixStorage storage) + : base(storage.RowCount, storage.ColumnCount) + { + _storage = storage; + } + /// /// Initializes a new instance of the class. /// diff --git a/src/Numerics/LinearAlgebra/Complex32/DiagonalMatrix.cs b/src/Numerics/LinearAlgebra/Complex32/DiagonalMatrix.cs index 121c807f..586eb060 100644 --- a/src/Numerics/LinearAlgebra/Complex32/DiagonalMatrix.cs +++ b/src/Numerics/LinearAlgebra/Complex32/DiagonalMatrix.cs @@ -1179,68 +1179,11 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 /// is not positive. public override Matrix SubMatrix(int rowIndex, int rowCount, int columnIndex, int columnCount) { - if (rowIndex >= RowCount || rowIndex < 0) - { - throw new ArgumentOutOfRangeException("rowIndex"); - } - - if (columnIndex >= ColumnCount || columnIndex < 0) - { - throw new ArgumentOutOfRangeException("columnIndex"); - } - - if (rowCount < 1) - { - throw new ArgumentException(Resources.ArgumentMustBePositive, "rowCount"); - } - - if (columnCount < 1) - { - throw new ArgumentException(Resources.ArgumentMustBePositive, "columnCount"); - } - - var colMax = columnIndex + columnCount; - var rowMax = rowIndex + rowCount; - - if (rowMax > RowCount) - { - throw new ArgumentOutOfRangeException("rowCount"); - } + // TODO: if rowIndex == columnIndex, use a diagonal matrix instead of a sparse one - if (colMax > ColumnCount) - { - throw new ArgumentOutOfRangeException("columnCount"); - } - - var result = new SparseMatrix(rowCount, columnCount); - - if (rowIndex > columnIndex && columnIndex + columnCount > rowIndex) - { - int columnInit = rowIndex - columnIndex; - int end = Math.Min(columnCount, rowCount + columnInit); - for (var i = 0; columnInit + i < end; i++) - { - result[i, columnInit + i] = _data[rowIndex + i]; - } - } - else if (rowIndex < columnIndex && rowIndex + rowCount > columnIndex) - { - int rowInit = columnIndex - rowIndex; - int end = Math.Min(columnCount + rowInit, rowCount); - for (var i = 0; rowInit + i < end; i++) - { - result[rowInit + i, i] = _data[columnIndex + i]; - } - } - else - { - for (var i = 0; i < Math.Min(columnCount, rowCount); i++) - { - result[i, i] = _data[rowIndex + i]; - } - } - - return result; + var storage = new SparseCompressedRowMatrixStorage(rowCount, columnCount, Complex32.Zero); + _storage.CopySubMatrixTo(storage, rowIndex, 0, rowCount, columnIndex, 0, columnCount, true); + return new SparseMatrix(storage); } /// diff --git a/src/Numerics/LinearAlgebra/Complex32/SparseMatrix.cs b/src/Numerics/LinearAlgebra/Complex32/SparseMatrix.cs index 06d5ab8f..86b2e01a 100644 --- a/src/Numerics/LinearAlgebra/Complex32/SparseMatrix.cs +++ b/src/Numerics/LinearAlgebra/Complex32/SparseMatrix.cs @@ -56,6 +56,12 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 get { return _storage.ValueCount; } } + internal SparseMatrix(SparseCompressedRowMatrixStorage storage) + : base(storage.RowCount, storage.ColumnCount) + { + _storage = storage; + } + /// /// Initializes a new instance of the class. /// diff --git a/src/Numerics/LinearAlgebra/Double/DiagonalMatrix.cs b/src/Numerics/LinearAlgebra/Double/DiagonalMatrix.cs index 92f938c3..11af89cf 100644 --- a/src/Numerics/LinearAlgebra/Double/DiagonalMatrix.cs +++ b/src/Numerics/LinearAlgebra/Double/DiagonalMatrix.cs @@ -1168,68 +1168,11 @@ namespace MathNet.Numerics.LinearAlgebra.Double /// is not positive. public override Matrix SubMatrix(int rowIndex, int rowCount, int columnIndex, int columnCount) { - if (rowIndex >= RowCount || rowIndex < 0) - { - throw new ArgumentOutOfRangeException("rowIndex"); - } - - if (columnIndex >= ColumnCount || columnIndex < 0) - { - throw new ArgumentOutOfRangeException("columnIndex"); - } - - if (rowCount < 1) - { - throw new ArgumentException(Resources.ArgumentMustBePositive, "rowCount"); - } - - if (columnCount < 1) - { - throw new ArgumentException(Resources.ArgumentMustBePositive, "columnCount"); - } - - var colMax = columnIndex + columnCount; - var rowMax = rowIndex + rowCount; - - if (rowMax > RowCount) - { - throw new ArgumentOutOfRangeException("rowCount"); - } + // TODO: if rowIndex == columnIndex, use a diagonal matrix instead of a sparse one - if (colMax > ColumnCount) - { - throw new ArgumentOutOfRangeException("columnCount"); - } - - var result = new SparseMatrix(rowCount, columnCount); - - if (rowIndex > columnIndex && columnIndex + columnCount > rowIndex) - { - int columnInit = rowIndex - columnIndex; - int end = Math.Min(columnCount, rowCount + columnInit); - for (var i = 0; columnInit + i < end; i++) - { - result[i, columnInit + i] = _data[rowIndex + i]; - } - } - else if (rowIndex < columnIndex && rowIndex + rowCount > columnIndex) - { - int rowInit = columnIndex - rowIndex; - int end = Math.Min(columnCount + rowInit, rowCount); - for (var i = 0; rowInit + i < end; i++) - { - result[rowInit + i, i] = _data[columnIndex + i]; - } - } - else - { - for (var i = 0; i < Math.Min(columnCount, rowCount); i++) - { - result[i, i] = _data[rowIndex + i]; - } - } - - return result; + var storage = new SparseCompressedRowMatrixStorage(rowCount, columnCount, 0d); + _storage.CopySubMatrixTo(storage, rowIndex, 0, rowCount, columnIndex, 0, columnCount, true); + return new SparseMatrix(storage); } /// diff --git a/src/Numerics/LinearAlgebra/Double/SparseMatrix.cs b/src/Numerics/LinearAlgebra/Double/SparseMatrix.cs index 1cffbb80..05f1c8d4 100644 --- a/src/Numerics/LinearAlgebra/Double/SparseMatrix.cs +++ b/src/Numerics/LinearAlgebra/Double/SparseMatrix.cs @@ -55,6 +55,12 @@ namespace MathNet.Numerics.LinearAlgebra.Double get { return _storage.ValueCount; } } + internal SparseMatrix(SparseCompressedRowMatrixStorage storage) + : base(storage.RowCount, storage.ColumnCount) + { + _storage = storage; + } + /// /// Initializes a new instance of the class. /// diff --git a/src/Numerics/LinearAlgebra/Single/DiagonalMatrix.cs b/src/Numerics/LinearAlgebra/Single/DiagonalMatrix.cs index 9ecd9b52..c602f74b 100644 --- a/src/Numerics/LinearAlgebra/Single/DiagonalMatrix.cs +++ b/src/Numerics/LinearAlgebra/Single/DiagonalMatrix.cs @@ -1173,68 +1173,11 @@ namespace MathNet.Numerics.LinearAlgebra.Single /// is not positive. public override Matrix SubMatrix(int rowIndex, int rowCount, int columnIndex, int columnCount) { - if (rowIndex >= RowCount || rowIndex < 0) - { - throw new ArgumentOutOfRangeException("rowIndex"); - } - - if (columnIndex >= ColumnCount || columnIndex < 0) - { - throw new ArgumentOutOfRangeException("columnIndex"); - } - - if (rowCount < 1) - { - throw new ArgumentException(Resources.ArgumentMustBePositive, "rowCount"); - } - - if (columnCount < 1) - { - throw new ArgumentException(Resources.ArgumentMustBePositive, "columnCount"); - } - - var colMax = columnIndex + columnCount; - var rowMax = rowIndex + rowCount; - - if (rowMax > RowCount) - { - throw new ArgumentOutOfRangeException("rowCount"); - } + // TODO: if rowIndex == columnIndex, use a diagonal matrix instead of a sparse one - if (colMax > ColumnCount) - { - throw new ArgumentOutOfRangeException("columnCount"); - } - - var result = new SparseMatrix(rowCount, columnCount); - - if (rowIndex > columnIndex && columnIndex + columnCount > rowIndex) - { - int columnInit = rowIndex - columnIndex; - int end = Math.Min(columnCount, rowCount + columnInit); - for (var i = 0; columnInit + i < end; i++) - { - result[i, columnInit + i] = _data[rowIndex + i]; - } - } - else if (rowIndex < columnIndex && rowIndex + rowCount > columnIndex) - { - int rowInit = columnIndex - rowIndex; - int end = Math.Min(columnCount + rowInit, rowCount); - for (var i = 0; rowInit + i < end; i++) - { - result[rowInit + i, i] = _data[columnIndex + i]; - } - } - else - { - for (var i = 0; i < Math.Min(columnCount, rowCount); i++) - { - result[i, i] = _data[rowIndex + i]; - } - } - - return result; + var storage = new SparseCompressedRowMatrixStorage(rowCount, columnCount, 0f); + _storage.CopySubMatrixTo(storage, rowIndex, 0, rowCount, columnIndex, 0, columnCount, true); + return new SparseMatrix(storage); } /// diff --git a/src/Numerics/LinearAlgebra/Single/SparseMatrix.cs b/src/Numerics/LinearAlgebra/Single/SparseMatrix.cs index 99668cd7..bacd08f7 100644 --- a/src/Numerics/LinearAlgebra/Single/SparseMatrix.cs +++ b/src/Numerics/LinearAlgebra/Single/SparseMatrix.cs @@ -55,6 +55,12 @@ namespace MathNet.Numerics.LinearAlgebra.Single get { return _storage.ValueCount; } } + internal SparseMatrix(SparseCompressedRowMatrixStorage storage) + : base(storage.RowCount, storage.ColumnCount) + { + _storage = storage; + } + /// /// Initializes a new instance of the class. /// diff --git a/src/Numerics/LinearAlgebra/Storage/ArgumentValidation.cs b/src/Numerics/LinearAlgebra/Storage/ArgumentValidation.cs new file mode 100644 index 00000000..58810168 --- /dev/null +++ b/src/Numerics/LinearAlgebra/Storage/ArgumentValidation.cs @@ -0,0 +1,77 @@ +using System; +using MathNet.Numerics.Properties; + +namespace MathNet.Numerics.LinearAlgebra.Storage +{ + // ReSharper disable UnusedParameter.Global + internal static class ArgumentValidation + { + public static void CopySubMatrixTo( + int sourceRowCount, int sourceColumnCount, + int targetRowCount, int targetColumnCount, + int sourceRowIndex, int targetRowIndex, int rowCount, + int sourceColumnIndex, int targetColumnIndex, int columnCount) + { + if (rowCount < 1) + { + throw new ArgumentException(Resources.ArgumentMustBePositive, "rowCount"); + } + + if (columnCount < 1) + { + throw new ArgumentException(Resources.ArgumentMustBePositive, "columnCount"); + } + + // Verify Source + + if (sourceRowIndex >= sourceRowCount || sourceRowIndex < 0) + { + throw new ArgumentOutOfRangeException("sourceRowIndex"); + } + + if (sourceColumnIndex >= sourceColumnCount || sourceColumnIndex < 0) + { + throw new ArgumentOutOfRangeException("sourceColumnIndex"); + } + + var sourceRowMax = sourceRowIndex + rowCount; + var sourceColumnMax = sourceColumnIndex + columnCount; + + if (sourceRowMax > sourceRowCount) + { + throw new ArgumentOutOfRangeException("rowCount"); + } + + if (sourceColumnMax > sourceColumnCount) + { + throw new ArgumentOutOfRangeException("columnCount"); + } + + // Verify Target + + if (targetRowIndex >= targetRowCount || targetRowIndex < 0) + { + throw new ArgumentOutOfRangeException("targetRowIndex"); + } + + if (targetColumnIndex >= targetColumnCount || targetColumnIndex < 0) + { + throw new ArgumentOutOfRangeException("targetColumnIndex"); + } + + var targetRowMax = targetRowIndex + rowCount; + var targetColumnMax = targetColumnIndex + columnCount; + + if (targetRowMax > targetRowCount) + { + throw new ArgumentOutOfRangeException("rowCount"); + } + + if (targetColumnMax > targetColumnCount) + { + throw new ArgumentOutOfRangeException("columnCount"); + } + } + } + // ReSharper restore UnusedParameter.Global +} diff --git a/src/Numerics/LinearAlgebra/Storage/DenseColumnMajorMatrixStorage.cs b/src/Numerics/LinearAlgebra/Storage/DenseColumnMajorMatrixStorage.cs index aba3e5a3..3d7d65cb 100644 --- a/src/Numerics/LinearAlgebra/Storage/DenseColumnMajorMatrixStorage.cs +++ b/src/Numerics/LinearAlgebra/Storage/DenseColumnMajorMatrixStorage.cs @@ -121,43 +121,6 @@ namespace MathNet.Numerics.LinearAlgebra.Storage int sourceRowIndex, int targetRowIndex, int rowCount, int sourceColumnIndex, int targetColumnIndex, int columnCount) { - if (rowCount < 1) - { - throw new ArgumentException(Resources.ArgumentMustBePositive, "rowCount"); - } - - if (columnCount < 1) - { - throw new ArgumentException(Resources.ArgumentMustBePositive, "columnCount"); - } - - // Verify Source - - if (sourceRowIndex >= RowCount || sourceRowIndex < 0) - { - throw new ArgumentOutOfRangeException("sourceRowIndex"); - } - - if (sourceColumnIndex >= ColumnCount || sourceColumnIndex < 0) - { - throw new ArgumentOutOfRangeException("sourceColumnIndex"); - } - - var sourceRowMax = sourceRowIndex + rowCount; - var sourceColumnMax = sourceColumnIndex + columnCount; - - if (sourceRowMax > RowCount) - { - throw new ArgumentOutOfRangeException("rowCount"); - } - - if (sourceColumnMax > ColumnCount) - { - throw new ArgumentOutOfRangeException("columnCount"); - } - - // Verify Target - if (target == null) { throw new ArgumentNullException("target"); @@ -168,32 +131,12 @@ namespace MathNet.Numerics.LinearAlgebra.Storage throw new NotSupportedException(); } - if (targetRowIndex >= target.RowCount || targetRowIndex < 0) - { - throw new ArgumentOutOfRangeException("targetRowIndex"); - } - - if (targetColumnIndex >= target.ColumnCount || targetColumnIndex < 0) - { - throw new ArgumentOutOfRangeException("targetColumnIndex"); - } - - var targetRowMax = targetRowIndex + rowCount; - var targetColumnMax = targetColumnIndex + columnCount; - - if (targetRowMax > target.RowCount) - { - throw new ArgumentOutOfRangeException("rowCount"); - } - - if (targetColumnMax > target.ColumnCount) - { - throw new ArgumentOutOfRangeException("columnCount"); - } - - // Copy + ArgumentValidation.CopySubMatrixTo(RowCount, ColumnCount, + target.RowCount, target.ColumnCount, + sourceRowIndex, targetRowIndex, rowCount, + sourceColumnIndex, targetColumnIndex, columnCount); - for (int j = sourceColumnIndex, jj = targetColumnIndex; j < sourceColumnMax; j++, jj++) + for (int j = sourceColumnIndex, jj = targetColumnIndex; j < sourceColumnIndex + columnCount; j++, jj++) { //Buffer.BlockCopy(Data, j*RowCount + sourceRowIndex, target.Data, jj*target.RowCount + targetRowIndex, rowCount * System.Runtime.InteropServices.Marshal.SizeOf(typeof(T))); Array.Copy(Data, j*RowCount + sourceRowIndex, target.Data, jj*target.RowCount + targetRowIndex, rowCount); diff --git a/src/Numerics/LinearAlgebra/Storage/SparseDiagonalMatrixStorage.cs b/src/Numerics/LinearAlgebra/Storage/SparseDiagonalMatrixStorage.cs index 20e3d722..e23cf8a8 100644 --- a/src/Numerics/LinearAlgebra/Storage/SparseDiagonalMatrixStorage.cs +++ b/src/Numerics/LinearAlgebra/Storage/SparseDiagonalMatrixStorage.cs @@ -174,5 +174,113 @@ namespace MathNet.Numerics.LinearAlgebra.Storage target.Data[i*(target.RowCount + 1)] = Data[i]; } } + + public void CopySubMatrixTo(DenseColumnMajorMatrixStorage target, + int sourceRowIndex, int targetRowIndex, int rowCount, + int sourceColumnIndex, int targetColumnIndex, int columnCount, + bool skipClearing = false) + { + if (target == null) + { + throw new ArgumentNullException("target"); + } + + ArgumentValidation.CopySubMatrixTo(RowCount, ColumnCount, + target.RowCount, target.ColumnCount, + sourceRowIndex, targetRowIndex, rowCount, + sourceColumnIndex, targetColumnIndex, columnCount); + + if (!skipClearing) + { + target.Clear(); + } + + if (sourceRowIndex > sourceColumnIndex && sourceColumnIndex + columnCount > sourceRowIndex) + { + // column by column, but skip resulting zero columns at the beginning + + int columnInit = sourceRowIndex - sourceColumnIndex; + int offset = (columnInit + targetColumnIndex) * target.RowCount + targetRowIndex; + int step = target.RowCount + 1; + int end = Math.Min(columnCount - columnInit, rowCount) + sourceRowIndex; + + for (int i = sourceRowIndex, j = offset; i < end; i++, j += step) + { + target.Data[j] = Data[i]; + } + } + else if (sourceRowIndex < sourceColumnIndex && sourceRowIndex + rowCount > sourceColumnIndex) + { + // row by row, but skip resulting zero rows at the beginning + + int rowInit = sourceColumnIndex - sourceRowIndex; + int offset = targetColumnIndex*target.RowCount + rowInit + targetRowIndex; + int step = target.RowCount + 1; + int end = Math.Min(columnCount, rowCount - rowInit) + sourceColumnIndex; + + for (int i = sourceColumnIndex, j = offset; i < end; i++, j += step) + { + target.Data[j] = Data[i]; + } + } + else + { + int offset = targetColumnIndex*target.RowCount + targetRowIndex; + int step = target.RowCount + 1; + var end = Math.Min(columnCount, rowCount) + sourceRowIndex; + + for (int i = sourceRowIndex, j = offset; i < end; i++, j += step) + { + target.Data[j] = Data[i]; + } + } + } + + public void CopySubMatrixTo(SparseCompressedRowMatrixStorage target, + int sourceRowIndex, int targetRowIndex, int rowCount, + int sourceColumnIndex, int targetColumnIndex, int columnCount, + bool skipClearing = false) + { + if (target == null) + { + throw new ArgumentNullException("target"); + } + + ArgumentValidation.CopySubMatrixTo(RowCount, ColumnCount, + target.RowCount, target.ColumnCount, + sourceRowIndex, targetRowIndex, rowCount, + sourceColumnIndex, targetColumnIndex, columnCount); + + if (!skipClearing) + { + target.Clear(); + } + + if (sourceRowIndex > sourceColumnIndex && sourceColumnIndex + columnCount > sourceRowIndex) + { + // column by column, but skip resulting zero columns at the beginning + int columnInit = sourceRowIndex - sourceColumnIndex; + for (var i = 0; i < Math.Min(columnCount - columnInit, rowCount); i++) + { + target.At(i + targetRowIndex, columnInit + i + targetColumnIndex, Data[sourceRowIndex + i]); + } + } + else if (sourceRowIndex < sourceColumnIndex && sourceRowIndex + rowCount > sourceColumnIndex) + { + // row by row, but skip resulting zero rows at the beginning + int rowInit = sourceColumnIndex - sourceRowIndex; + for (var i = 0; i < Math.Min(columnCount, rowCount - rowInit); i++) + { + target.At(rowInit + i + targetRowIndex, i + targetColumnIndex, Data[sourceColumnIndex + i]); + } + } + else + { + for (var i = 0; i < Math.Min(columnCount, rowCount); i++) + { + target.At(i + targetRowIndex, i + targetColumnIndex, Data[sourceRowIndex + i]); + } + } + } } } diff --git a/src/Numerics/Numerics.csproj b/src/Numerics/Numerics.csproj index 24e2c91b..d0f376b1 100644 --- a/src/Numerics/Numerics.csproj +++ b/src/Numerics/Numerics.csproj @@ -122,6 +122,7 @@ + diff --git a/src/Portable/Portable.csproj b/src/Portable/Portable.csproj index da3a7848..b6ab9e48 100644 --- a/src/Portable/Portable.csproj +++ b/src/Portable/Portable.csproj @@ -888,6 +888,9 @@ LinearAlgebra\Single\Vector.cs + + LinearAlgebra\Storage\ArgumentValidation.cs + LinearAlgebra\Storage\DenseColumnMajorMatrixStorage.cs