// <copyright file="SparseCompressedRowMatrixStorage.cs" company="Math.NET">
|
// Math.NET Numerics, part of the Math.NET Project
|
// http://numerics.mathdotnet.com
|
// http://github.com/mathnet/mathnet-numerics
|
//
|
// 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
|
// files (the "Software"), to deal in the Software without
|
// restriction, including without limitation the rights to use,
|
// copy, modify, merge, publish, distribute, sublicense, and/or sell
|
// copies of the Software, and to permit persons to whom the
|
// Software is furnished to do so, subject to the following
|
// conditions:
|
//
|
// The above copyright notice and this permission notice shall be
|
// included in all copies or substantial portions of the Software.
|
//
|
// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
|
// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES
|
// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
|
// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT
|
// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY,
|
// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
|
// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
|
// OTHER DEALINGS IN THE SOFTWARE.
|
// </copyright>
|
|
using System;
|
using System.Collections.Generic;
|
using System.Linq;
|
using System.Runtime.Serialization;
|
|
namespace IStation.Numerics.LinearAlgebra.Storage
|
{
|
[Serializable]
|
[DataContract(Namespace = "urn:IStation/Numerics/LinearAlgebra")]
|
public class SparseCompressedRowMatrixStorage<T> : MatrixStorage<T>
|
where T : struct, IEquatable<T>, IFormattable
|
{
|
// [ruegg] public fields are OK here
|
|
/// <summary>
|
/// The array containing the row indices of the existing rows. Element "i" of the array gives the index of the
|
/// element in the <see cref="Values"/> 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.
|
/// </summary>
|
[DataMember(Order = 1)]
|
public readonly int[] RowPointers;
|
|
/// <summary>
|
/// 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 <see cref="Values"/> array.
|
/// </summary>
|
[DataMember(Order = 2)]
|
public int[] ColumnIndices;
|
|
/// <summary>
|
/// Array that contains the non-zero elements of matrix. Values of the non-zero elements of matrix are mapped into the values
|
/// array using the row-major storage mapping described in a compressed sparse row (CSR) format.
|
/// </summary>
|
[DataMember(Order = 3)]
|
public T[] Values;
|
|
/// <summary>
|
/// Gets the number of non zero elements in the matrix.
|
/// </summary>
|
/// <value>The number of non zero elements.</value>
|
public int ValueCount => RowPointers[RowCount];
|
|
internal SparseCompressedRowMatrixStorage(int rows, int columns)
|
: base(rows, columns)
|
{
|
RowPointers = new int[rows + 1];
|
ColumnIndices = new int[0];
|
Values = new T[0];
|
}
|
|
internal SparseCompressedRowMatrixStorage(int rows, int columns, int[] rowPointers, int[] columnIndices, T[] values)
|
: base(rows, columns)
|
{
|
RowPointers = rowPointers;
|
ColumnIndices = columnIndices;
|
Values = values;
|
|
// Explicit zeros are not intentionally removed.
|
// Sort ColumnIndices.
|
NormalizeOrdering();
|
NormalizeDuplicates();
|
}
|
|
/// <summary>
|
/// True if the matrix storage format is dense.
|
/// </summary>
|
public override bool IsDense => false;
|
|
/// <summary>
|
/// True if all fields of this matrix can be set to any value.
|
/// False if some fields are fixed, like on a diagonal matrix.
|
/// </summary>
|
public override bool IsFullyMutable => true;
|
|
/// <summary>
|
/// True if the specified field can be set to any value.
|
/// False if the field is fixed, like an off-diagonal field on a diagonal matrix.
|
/// </summary>
|
public override bool IsMutableAt(int row, int column)
|
{
|
return true;
|
}
|
|
/// <summary>
|
/// Retrieves the requested element without range checking.
|
/// </summary>
|
/// <param name="row">
|
/// The row of the element.
|
/// </param>
|
/// <param name="column">
|
/// The column of the element.
|
/// </param>
|
/// <returns>
|
/// The requested element.
|
/// </returns>
|
/// <remarks>Not range-checked.</remarks>
|
public override T At(int row, int column)
|
{
|
var index = FindItem(row, column);
|
return index >= 0 ? Values[index] : Zero;
|
}
|
|
/// <summary>
|
/// Sets the element without range checking.
|
/// </summary>
|
/// <param name="row"> The row of the element. </param>
|
/// <param name="column"> The column of the element. </param>
|
/// <param name="value"> The value to set the element to. </param>
|
/// <remarks>WARNING: This method is not thread safe. Use "lock" with it and be sure to avoid deadlocks.</remarks>
|
public override void At(int row, int column, T value)
|
{
|
var index = FindItem(row, column);
|
if (index >= 0)
|
{
|
// Non-zero item found in matrix
|
if (Zero.Equals(value))
|
{
|
// Delete existing item
|
RemoveAtIndexUnchecked(index, row);
|
}
|
else
|
{
|
// Update item
|
Values[index] = value;
|
}
|
}
|
else
|
{
|
// Item not found. Add new value
|
if (Zero.Equals(value))
|
{
|
return;
|
}
|
|
index = ~index;
|
var valueCount = RowPointers[RowPointers.Length - 1];
|
|
// Check if the storage needs to be increased
|
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);
|
if (size > int.MaxValue)
|
{
|
throw new NotSupportedException("We only support sparse matrix with less than int.MaxValue elements.");
|
}
|
|
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);
|
|
// Add the value and the column index
|
Values[index] = value;
|
ColumnIndices[index] = column;
|
|
// 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++)
|
{
|
RowPointers[i] += 1;
|
}
|
}
|
}
|
|
/// <summary>
|
/// Delete value from internal storage
|
/// </summary>
|
/// <param name="itemIndex">Index of value in nonZeroValues array</param>
|
/// <param name="row">Row number of matrix</param>
|
/// <remarks>WARNING: This method is not thread safe. Use "lock" with it and be sure to avoid deadlocks</remarks>
|
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);
|
|
// Decrease value in Row
|
for (var i = row + 1; i < RowPointers.Length; i++)
|
{
|
RowPointers[i] -= 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))
|
{
|
Array.Resize(ref Values, valueCount);
|
Array.Resize(ref ColumnIndices, valueCount);
|
}
|
}
|
|
/// <summary>
|
/// Find item Index in nonZeroValues array
|
/// </summary>
|
/// <param name="row">Matrix row index</param>
|
/// <param name="column">Matrix column index</param>
|
/// <returns>Item index</returns>
|
/// <remarks>WARNING: This method is not thread safe. Use "lock" with it and be sure to avoid deadlocks</remarks>
|
public int FindItem(int row, int column)
|
{
|
// Determine bounds in columnIndices array where this item should be searched (using rowIndex)
|
return Array.BinarySearch(ColumnIndices, RowPointers[row], RowPointers[row + 1] - RowPointers[row], column);
|
}
|
|
/// <summary>
|
/// Calculates the amount with which to grow the storage array's if they need to be
|
/// increased in size.
|
/// </summary>
|
/// <returns>The amount grown.</returns>
|
int GrowthSize()
|
{
|
int delta;
|
if (Values.Length > 1024)
|
{
|
delta = Values.Length/4;
|
}
|
else
|
{
|
if (Values.Length > 256)
|
{
|
delta = 512;
|
}
|
else
|
{
|
delta = Values.Length > 64 ? 128 : 32;
|
}
|
}
|
|
return delta;
|
}
|
|
public void Normalize()
|
{
|
NormalizeOrdering();
|
NormalizeZeros();
|
}
|
|
public void NormalizeOrdering()
|
{
|
for (int i = 0; i < RowCount; i++)
|
{
|
int index = RowPointers[i];
|
int count = RowPointers[i + 1] - index;
|
if (count > 1)
|
{
|
Sorting.Sort(ColumnIndices, Values, index, count);
|
}
|
}
|
}
|
|
public void NormalizeZeros()
|
{
|
MapInplace(x => x, Zeros.AllowSkip);
|
}
|
|
/// <summary>
|
/// Eliminate duplicate entries by adding them together.
|
/// </summary>
|
public void NormalizeDuplicates()
|
{
|
var builder = BuilderInstance<T>.Matrix;
|
|
int valueCount = 0;
|
for (int i = 0; i < RowCount; i++)
|
{
|
int index = RowPointers[i];
|
int last = RowPointers[i + 1];
|
while (index < last)
|
{
|
var col = ColumnIndices[index];
|
var val = Values[index];
|
index++;
|
while (index < last)
|
{
|
if (ColumnIndices[index] == col)
|
{
|
val = builder.Add(val, Values[index]);
|
index++;
|
}
|
else
|
{
|
break;
|
}
|
}
|
ColumnIndices[valueCount] = col;
|
Values[valueCount] = val;
|
valueCount++;
|
}
|
RowPointers[i + 1] = valueCount;
|
}
|
|
// Remove extra space from arrays.
|
Array.Resize(ref Values, valueCount);
|
Array.Resize(ref ColumnIndices, valueCount);
|
}
|
|
/// <summary>
|
/// Fill zeros explicitly on the diagonal entries as required by the Intel MKL direct sparse solver.
|
/// </summary>
|
public void PopulateExplicitZerosOnDiagonal()
|
{
|
var delta = 0; // number of missing diagonal entries
|
|
for (int row = 0; row < RowCount; row++)
|
{
|
var found = false;
|
for (int j = RowPointers[row]; j < RowPointers[row + 1]; j++)
|
{
|
if (ColumnIndices[j] == row)
|
{
|
found = true;
|
break;
|
}
|
}
|
if (!found) delta++;
|
}
|
|
if (delta > 0)
|
{
|
var size = Values.Length + delta;
|
if (size > int.MaxValue)
|
{
|
throw new NotSupportedException("We only support sparse matrix with less than int.MaxValue elements.");
|
}
|
|
var newRowPointers = new int[RowCount + 1];
|
var newColumnIndices = new int[size];
|
var newValues = new T[size];
|
|
delta = 0;
|
for (int row = 0; row < RowCount; row++)
|
{
|
var found = false;
|
for (int j = RowPointers[row]; j < RowPointers[row + 1]; j++)
|
{
|
newColumnIndices[j + delta] = ColumnIndices[j];
|
newValues[j + delta] = Values[j];
|
if (ColumnIndices[j] == row)
|
{
|
found = true;
|
}
|
}
|
if (!found)
|
{
|
var start = RowPointers[row] + delta;
|
var end = RowPointers[row + 1] + delta;
|
var count = end - start + 1;
|
|
newColumnIndices[end] = row;
|
newValues[end] = Zero;
|
|
// Ordering may be not necessary
|
Sorting.Sort(newColumnIndices, newValues, start, count);
|
|
delta++;
|
}
|
newRowPointers[row + 1] = RowPointers[row + 1] + delta;
|
}
|
|
Array.Copy(newRowPointers, RowPointers, RowCount + 1);
|
ColumnIndices = newColumnIndices;
|
Values = newValues;
|
}
|
}
|
|
/// <summary>
|
/// Returns a hash code for this instance.
|
/// </summary>
|
/// <returns>
|
/// A hash code for this instance, suitable for use in hashing algorithms and data structures like a hash table.
|
/// </returns>
|
public override int GetHashCode()
|
{
|
var values = Values;
|
var hashNum = Math.Min(ValueCount, 25);
|
int hash = 17;
|
unchecked
|
{
|
for (var i = 0; i < hashNum; i++)
|
{
|
hash = hash*31 + values[i].GetHashCode();
|
}
|
}
|
return hash;
|
}
|
|
// CLEARING
|
|
public override void Clear()
|
{
|
Array.Clear(RowPointers, 0, RowPointers.Length);
|
}
|
|
internal override void ClearUnchecked(int rowIndex, int rowCount, int columnIndex, int columnCount)
|
{
|
if (rowIndex == 0 && columnIndex == 0 && rowCount == RowCount && columnCount == ColumnCount)
|
{
|
Clear();
|
return;
|
}
|
|
var valueCount = RowPointers[RowPointers.Length - 1];
|
|
for (int row = rowIndex + rowCount - 1; row >= rowIndex; row--)
|
{
|
var startIndex = RowPointers[row];
|
var endIndex = RowPointers[row + 1];
|
|
// empty row
|
if (startIndex == endIndex)
|
{
|
continue;
|
}
|
|
// multiple entries in row
|
var first = Array.BinarySearch(ColumnIndices, startIndex, endIndex - startIndex, columnIndex);
|
var last = Array.BinarySearch(ColumnIndices, startIndex, endIndex - startIndex, columnIndex + columnCount - 1);
|
if (first < 0) first = ~first;
|
if (last < 0) last = ~last - 1;
|
int count = last - first + 1;
|
|
if (count > 0)
|
{
|
// 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);
|
|
// Decrease value in Row
|
for (var k = row + 1; k < RowPointers.Length; k++)
|
{
|
RowPointers[k] -= 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))
|
{
|
Array.Resize(ref Values, valueCount);
|
Array.Resize(ref ColumnIndices, valueCount);
|
}
|
}
|
|
internal override void ClearRowsUnchecked(int[] rowIndices)
|
{
|
var rows = new bool[RowCount];
|
for (int i = 0; i < rowIndices.Length; i++)
|
{
|
rows[rowIndices[i]] = true;
|
}
|
MapIndexedInplace((i, j, x) => rows[i] ? Zero : x, Zeros.AllowSkip);
|
}
|
|
internal override void ClearColumnsUnchecked(int[] columnIndices)
|
{
|
var columns = new bool[ColumnCount];
|
for (int i = 0; i < columnIndices.Length; i++)
|
{
|
columns[columnIndices[i]] = true;
|
}
|
MapIndexedInplace((i, j, x) => columns[j] ? Zero : x, Zeros.AllowSkip);
|
}
|
|
// INITIALIZATION
|
|
public static SparseCompressedRowMatrixStorage<T> OfMatrix(MatrixStorage<T> matrix)
|
{
|
var storage = new SparseCompressedRowMatrixStorage<T>(matrix.RowCount, matrix.ColumnCount);
|
matrix.CopyToUnchecked(storage, ExistingData.AssumeZeros);
|
return storage;
|
}
|
|
public static SparseCompressedRowMatrixStorage<T> OfValue(int rows, int columns, T value)
|
{
|
if (Zero.Equals(value))
|
{
|
return new SparseCompressedRowMatrixStorage<T>(rows, columns);
|
}
|
|
var storage = new SparseCompressedRowMatrixStorage<T>(rows, columns);
|
|
var values = new T[rows * columns];
|
for (int i = 0; i < values.Length; i++)
|
{
|
values[i] = value;
|
}
|
|
var rowPointers = storage.RowPointers;
|
for (int i = 0; i <= rows; i++)
|
{
|
rowPointers[i] = i*columns;
|
}
|
|
var columnIndices = new int[values.Length];
|
for (int row = 0; row < rows; row++)
|
{
|
int offset = row*columns;
|
for (int col = 0; col < columns; col++)
|
{
|
columnIndices[offset + col] = col;
|
}
|
}
|
|
rowPointers[rows] = values.Length;
|
storage.ColumnIndices = columnIndices;
|
storage.Values = values;
|
return storage;
|
}
|
|
|
public static SparseCompressedRowMatrixStorage<T> OfInit(int rows, int columns, Func<int, int, T> init)
|
{
|
var storage = new SparseCompressedRowMatrixStorage<T>(rows, columns);
|
var rowPointers = storage.RowPointers;
|
var columnIndices = new List<int>();
|
var values = new List<T>();
|
|
for (int row = 0; row < rows; row++)
|
{
|
rowPointers[row] = values.Count;
|
for (int col = 0; col < columns; col++)
|
{
|
var x = init(row, col);
|
if (!Zero.Equals(x))
|
{
|
values.Add(x);
|
columnIndices.Add(col);
|
}
|
}
|
}
|
|
rowPointers[rows] = values.Count;
|
storage.ColumnIndices = columnIndices.ToArray();
|
storage.Values = values.ToArray();
|
return storage;
|
}
|
|
public static SparseCompressedRowMatrixStorage<T> OfDiagonalInit(int rows, int columns, Func<int, T> init)
|
{
|
int diagonalLength = Math.Min(rows, columns);
|
var storage = new SparseCompressedRowMatrixStorage<T>(rows, columns);
|
var rowPointers = storage.RowPointers;
|
var columnIndices = new List<int>(diagonalLength);
|
var values = new List<T>(diagonalLength);
|
|
for (int i = 0; i < diagonalLength; i++)
|
{
|
rowPointers[i] = values.Count;
|
var x = init(i);
|
if (!Zero.Equals(x))
|
{
|
values.Add(x);
|
columnIndices.Add(i);
|
}
|
}
|
|
rowPointers[rows] = values.Count;
|
storage.ColumnIndices = columnIndices.ToArray();
|
storage.Values = values.ToArray();
|
return storage;
|
}
|
|
public static SparseCompressedRowMatrixStorage<T> OfArray(T[,] array)
|
{
|
var storage = new SparseCompressedRowMatrixStorage<T>(array.GetLength(0), array.GetLength(1));
|
var rowPointers = storage.RowPointers;
|
var columnIndices = new List<int>();
|
var values = new List<T>();
|
|
for (int row = 0; row < storage.RowCount; row++)
|
{
|
rowPointers[row] = values.Count;
|
for (int col = 0; col < storage.ColumnCount; col++)
|
{
|
if (!Zero.Equals(array[row, col]))
|
{
|
values.Add(array[row, col]);
|
columnIndices.Add(col);
|
}
|
}
|
}
|
|
rowPointers[storage.RowCount] = values.Count;
|
storage.ColumnIndices = columnIndices.ToArray();
|
storage.Values = values.ToArray();
|
return storage;
|
}
|
|
public static SparseCompressedRowMatrixStorage<T> OfRowArrays(T[][] data)
|
{
|
if (data.Length <= 0)
|
{
|
throw new ArgumentOutOfRangeException(nameof(data), "Matrices can not be empty and must have at least one row and column.");
|
}
|
|
var storage = new SparseCompressedRowMatrixStorage<T>(data.Length, data[0].Length);
|
var rowPointers = storage.RowPointers;
|
var columnIndices = new List<int>();
|
var values = new List<T>();
|
|
for (int row = 0; row < storage.RowCount; row++)
|
{
|
rowPointers[row] = values.Count;
|
for (int col = 0; col < storage.ColumnCount; col++)
|
{
|
T x = data[row][col];
|
if (!Zero.Equals(x))
|
{
|
values.Add(x);
|
columnIndices.Add(col);
|
}
|
}
|
}
|
|
rowPointers[storage.RowCount] = values.Count;
|
storage.ColumnIndices = columnIndices.ToArray();
|
storage.Values = values.ToArray();
|
return storage;
|
}
|
|
public static SparseCompressedRowMatrixStorage<T> OfColumnArrays(T[][] data)
|
{
|
if (data.Length <= 0)
|
{
|
throw new ArgumentOutOfRangeException(nameof(data), "Matrices can not be empty and must have at least one row and column.");
|
}
|
|
var storage = new SparseCompressedRowMatrixStorage<T>(data[0].Length, data.Length);
|
var rowPointers = storage.RowPointers;
|
var columnIndices = new List<int>();
|
var values = new List<T>();
|
|
for (int row = 0; row < storage.RowCount; row++)
|
{
|
rowPointers[row] = values.Count;
|
for (int col = 0; col < storage.ColumnCount; col++)
|
{
|
T x = data[col][row];
|
if (!Zero.Equals(x))
|
{
|
values.Add(x);
|
columnIndices.Add(col);
|
}
|
}
|
}
|
|
rowPointers[storage.RowCount] = values.Count;
|
storage.ColumnIndices = columnIndices.ToArray();
|
storage.Values = values.ToArray();
|
return storage;
|
}
|
|
public static SparseCompressedRowMatrixStorage<T> OfRowVectors(VectorStorage<T>[] data)
|
{
|
if (data.Length <= 0)
|
{
|
throw new ArgumentOutOfRangeException(nameof(data), "Matrices can not be empty and must have at least one row and column.");
|
}
|
|
var storage = new SparseCompressedRowMatrixStorage<T>(data.Length, data[0].Length);
|
var rowPointers = storage.RowPointers;
|
var columnIndices = new List<int>();
|
var values = new List<T>();
|
|
// TODO PERF: Optimize for sparse and dense cases
|
for (int row = 0; row < storage.RowCount; row++)
|
{
|
var vector = data[row];
|
rowPointers[row] = values.Count;
|
for (int col = 0; col < storage.ColumnCount; col++)
|
{
|
var x = vector.At(col);
|
if (!Zero.Equals(x))
|
{
|
values.Add(x);
|
columnIndices.Add(col);
|
}
|
}
|
}
|
|
rowPointers[storage.RowCount] = values.Count;
|
storage.ColumnIndices = columnIndices.ToArray();
|
storage.Values = values.ToArray();
|
return storage;
|
}
|
|
public static SparseCompressedRowMatrixStorage<T> OfColumnVectors(VectorStorage<T>[] data)
|
{
|
if (data.Length <= 0)
|
{
|
throw new ArgumentOutOfRangeException(nameof(data), "Matrices can not be empty and must have at least one row and column.");
|
}
|
|
var storage = new SparseCompressedRowMatrixStorage<T>(data[0].Length, data.Length);
|
var rowPointers = storage.RowPointers;
|
var columnIndices = new List<int>();
|
var values = new List<T>();
|
|
// TODO PERF: Optimize for sparse and dense cases
|
for (int row = 0; row < storage.RowCount; row++)
|
{
|
rowPointers[row] = values.Count;
|
for (int col = 0; col < storage.ColumnCount; col++)
|
{
|
var x = data[col].At(row);
|
if (!Zero.Equals(x))
|
{
|
values.Add(x);
|
columnIndices.Add(col);
|
}
|
}
|
}
|
|
rowPointers[storage.RowCount] = values.Count;
|
storage.ColumnIndices = columnIndices.ToArray();
|
storage.Values = values.ToArray();
|
return storage;
|
}
|
|
public static SparseCompressedRowMatrixStorage<T> OfIndexedEnumerable(int rows, int columns, IEnumerable<Tuple<int, int, T>> data)
|
{
|
var trows = new List<Tuple<int, T>>[rows];
|
foreach (var item in data)
|
{
|
if (!Zero.Equals(item.Item3))
|
{
|
var row = trows[item.Item1] ?? (trows[item.Item1] = new List<Tuple<int, T>>());
|
row.Add(new Tuple<int, T>(item.Item2, item.Item3));
|
}
|
}
|
|
var storage = new SparseCompressedRowMatrixStorage<T>(rows, columns);
|
var rowPointers = storage.RowPointers;
|
var columnIndices = new List<int>();
|
var values = new List<T>();
|
|
int index = 0;
|
for (int row = 0; row < rows; row++)
|
{
|
rowPointers[row] = index;
|
var trow = trows[row];
|
if (trow != null)
|
{
|
trow.Sort();
|
foreach (var item in trow)
|
{
|
values.Add(item.Item2);
|
columnIndices.Add(item.Item1);
|
index++;
|
}
|
}
|
}
|
|
rowPointers[rows] = values.Count;
|
storage.ColumnIndices = columnIndices.ToArray();
|
storage.Values = values.ToArray();
|
return storage;
|
}
|
|
public static SparseCompressedRowMatrixStorage<T> OfRowEnumerables(int rows, int columns, IEnumerable<IEnumerable<T>> data)
|
{
|
var storage = new SparseCompressedRowMatrixStorage<T>(rows, columns);
|
var rowPointers = storage.RowPointers;
|
var columnIndices = new List<int>();
|
var values = new List<T>();
|
|
using (var rowIterator = data.GetEnumerator())
|
{
|
for (int row = 0; row < rows; row++)
|
{
|
if (!rowIterator.MoveNext()) throw new ArgumentOutOfRangeException(nameof(data), $"The given array has the wrong length. Should be {rows}.");
|
rowPointers[row] = values.Count;
|
using (var columnIterator = rowIterator.Current.GetEnumerator())
|
{
|
for (int col = 0; col < columns; col++)
|
{
|
if (!columnIterator.MoveNext()) throw new ArgumentOutOfRangeException(nameof(data), $"The given array has the wrong length. Should be {columns}.");
|
if (!Zero.Equals(columnIterator.Current))
|
{
|
values.Add(columnIterator.Current);
|
columnIndices.Add(col);
|
}
|
}
|
if (columnIterator.MoveNext()) throw new ArgumentOutOfRangeException(nameof(data), $"The given array has the wrong length. Should be {columns}.");
|
}
|
}
|
if (rowIterator.MoveNext()) throw new ArgumentOutOfRangeException(nameof(data), $"The given array has the wrong length. Should be {rows}.");
|
}
|
|
rowPointers[rows] = values.Count;
|
storage.ColumnIndices = columnIndices.ToArray();
|
storage.Values = values.ToArray();
|
return storage;
|
}
|
|
public static SparseCompressedRowMatrixStorage<T> OfColumnEnumerables(int rows, int columns, IEnumerable<IEnumerable<T>> data)
|
{
|
var trows = new List<Tuple<int, T>>[rows];
|
using (var columnIterator = data.GetEnumerator())
|
{
|
for (int column = 0; column < columns; column++)
|
{
|
if (!columnIterator.MoveNext()) throw new ArgumentOutOfRangeException(nameof(data), $"The given array has the wrong length. Should be {columns}.");
|
using (var rowIterator = columnIterator.Current.GetEnumerator())
|
{
|
for (int row = 0; row < rows; row++)
|
{
|
if (!rowIterator.MoveNext()) throw new ArgumentOutOfRangeException(nameof(data), $"The given array has the wrong length. Should be {rows}.");
|
if (!Zero.Equals(rowIterator.Current))
|
{
|
var trow = trows[row] ?? (trows[row] = new List<Tuple<int, T>>());
|
trow.Add(new Tuple<int, T>(column, rowIterator.Current));
|
}
|
}
|
}
|
}
|
}
|
|
var storage = new SparseCompressedRowMatrixStorage<T>(rows, columns);
|
var rowPointers = storage.RowPointers;
|
var columnIndices = new List<int>();
|
var values = new List<T>();
|
|
int index = 0;
|
for (int row = 0; row < rows; row++)
|
{
|
rowPointers[row] = index;
|
var trow = trows[row];
|
if (trow != null)
|
{
|
trow.Sort();
|
foreach (var item in trow)
|
{
|
values.Add(item.Item2);
|
columnIndices.Add(item.Item1);
|
index++;
|
}
|
}
|
}
|
|
rowPointers[rows] = values.Count;
|
storage.ColumnIndices = columnIndices.ToArray();
|
storage.Values = values.ToArray();
|
return storage;
|
}
|
|
public static SparseCompressedRowMatrixStorage<T> OfRowMajorEnumerable(int rows, int columns, IEnumerable<T> data)
|
{
|
var storage = new SparseCompressedRowMatrixStorage<T>(rows, columns);
|
var rowPointers = storage.RowPointers;
|
var columnIndices = new List<int>();
|
var values = new List<T>();
|
|
using (var iterator = data.GetEnumerator())
|
{
|
for (int row = 0; row < rows; row++)
|
{
|
rowPointers[row] = values.Count;
|
for (int col = 0; col < columns; col++)
|
{
|
iterator.MoveNext();
|
if (!Zero.Equals(iterator.Current))
|
{
|
values.Add(iterator.Current);
|
columnIndices.Add(col);
|
}
|
}
|
}
|
}
|
|
rowPointers[rows] = values.Count;
|
storage.ColumnIndices = columnIndices.ToArray();
|
storage.Values = values.ToArray();
|
return storage;
|
}
|
|
public static SparseCompressedRowMatrixStorage<T> OfColumnMajorList(int rows, int columns, IList<T> data)
|
{
|
if (rows*columns != data.Count)
|
{
|
throw new ArgumentException("Matrix dimensions must agree.");
|
}
|
|
var storage = new SparseCompressedRowMatrixStorage<T>(rows, columns);
|
var rowPointers = storage.RowPointers;
|
var columnIndices = new List<int>();
|
var values = new List<T>();
|
|
for (int row = 0; row < rows; row++)
|
{
|
rowPointers[row] = values.Count;
|
for (int col = 0; col < columns; col++)
|
{
|
var item = data[row + (col*rows)];
|
if (!Zero.Equals(item))
|
{
|
values.Add(item);
|
columnIndices.Add(col);
|
}
|
}
|
}
|
|
rowPointers[rows] = values.Count;
|
storage.ColumnIndices = columnIndices.ToArray();
|
storage.Values = values.ToArray();
|
return storage;
|
}
|
|
/// <summary>
|
/// Create a new sparse storage from a compressed sparse row (CSR) format.
|
/// This new storage will be independent from the given arrays.
|
/// A new memory block will be allocated for storing the matrix.
|
/// </summary>
|
/// <param name="rows">The number of rows.</param>
|
/// <param name="columns">The number of columns.</param>
|
/// <param name="valueCount">The number of stored values including explicit zeros.</param>
|
/// <param name="rowPointers">The row pointer array of the compressed sparse row format.</param>
|
/// <param name="columnIndices">The column index array of the compressed sparse row format.</param>
|
/// <param name="values">The data array of the compressed sparse row format.</param>
|
/// <returns>The sparse storage from the compressed sparse row format.</returns>
|
/// <remarks>Duplicate entries will be summed together and explicit zeros will be not intentionally removed.</remarks>
|
public static SparseCompressedRowMatrixStorage<T> OfCompressedSparseRowFormat(int rows, int columns, int valueCount, int[] rowPointers, int[] columnIndices, T[] values)
|
{
|
if (values == null) throw new NullReferenceException(nameof(values));
|
if (columnIndices == null) throw new NullReferenceException(nameof(columnIndices));
|
if (rowPointers == null) throw new NullReferenceException(nameof(rowPointers));
|
if (rowPointers.Length < rows) throw new Exception($"The given array has the wrong length. Should be {rows + 1}.");
|
if (valueCount != rowPointers[rows]) throw new Exception($"{nameof(valueCount)} should be same to {rowPointers[rows]}");
|
|
// copy arrays to new memory block.
|
var storage = new SparseCompressedRowMatrixStorage<T>(rows, columns);
|
var csrValues = new T[valueCount];
|
Array.Copy(values, csrValues, valueCount);
|
var csrColumnIndices = new int[valueCount];
|
Array.Copy(columnIndices, csrColumnIndices, valueCount);
|
Array.Copy(rowPointers, storage.RowPointers, rows + 1);
|
|
storage.ColumnIndices = csrColumnIndices;
|
storage.Values = csrValues;
|
storage.NormalizeOrdering();
|
storage.NormalizeDuplicates();
|
return storage;
|
}
|
|
/// <summary>
|
/// Create a new sparse matrix storage from a compressed sparse column (CSC) format.
|
/// This new storage will be independent from the given arrays.
|
/// A new memory block will be allocated.
|
/// </summary>
|
/// <param name="rows">The number of rows.</param>
|
/// <param name="columns">The number of columns.</param>
|
/// <param name="valueCount">The number of stored values including explicit zeros.</param>
|
/// <param name="rowIndices">The row index array of the compressed sparse column format.</param>
|
/// <param name="columnPointers">The column pointer array of the compressed sparse column format.</param>
|
/// <param name="values">The data array of the compressed sparse column format.</param>
|
/// <returns>The sparse storage from the compressed sparse column format.</returns>
|
/// <remarks>Duplicate entries will be summed together and explicit zeros will be not intentionally removed.</remarks>
|
public static SparseCompressedRowMatrixStorage<T> OfCompressedSparseColumnFormat(int rows, int columns, int valueCount, int[] rowIndices, int[] columnPointers, T[] values)
|
{
|
if (values == null) throw new NullReferenceException(nameof(values));
|
if (rowIndices == null) throw new NullReferenceException(nameof(rowIndices));
|
if (columnPointers == null) throw new NullReferenceException(nameof(columnPointers));
|
if (columnPointers.Length < columns) throw new Exception($"The given array has the wrong length. Should be {columns + 1}.");
|
if (valueCount != columnPointers[columns]) throw new Exception($"{nameof(valueCount)} should be same to {columnPointers[columns]}");
|
|
// convert from CSC to CSR
|
var storage = new SparseCompressedRowMatrixStorage<T>(rows, columns);
|
var csrValues = new T[valueCount];
|
var csrRowPointers = storage.RowPointers;
|
var csrColumnIndices = new int[valueCount];
|
|
for (int i = 0; i < columns; i++)
|
{
|
for (int j = columnPointers[i]; j < columnPointers[i + 1]; j++)
|
{
|
csrRowPointers[rowIndices[j] + 1]++;
|
}
|
}
|
|
for (int i = 1; i < rows + 1; i++)
|
{
|
csrRowPointers[i] += csrRowPointers[i - 1];
|
}
|
|
var curr = new int[rows];
|
for (int i = 0; i < columns; i++)
|
{
|
for (int j = columnPointers[i]; j < columnPointers[i + 1]; j++)
|
{
|
var loc = csrRowPointers[rowIndices[j]] + curr[rowIndices[j]];
|
curr[rowIndices[j]]++;
|
csrColumnIndices[loc] = i;
|
csrValues[loc] = values[j];
|
}
|
}
|
|
storage.ColumnIndices = csrColumnIndices;
|
storage.Values = csrValues;
|
storage.NormalizeOrdering();
|
storage.NormalizeDuplicates();
|
return storage;
|
}
|
|
/// <summary>
|
/// Create a new sparse storage from a coordinate (COO) format.
|
/// This new storage will be independent from the given arrays.
|
/// A new memory block will be allocated for storing the matrix.
|
/// </summary>
|
/// <param name="rows">The number of rows.</param>
|
/// <param name="columns">The number of columns.</param>
|
/// <param name="valueCount">The number of stored values including explicit zeros.</param>
|
/// <param name="rowIndices">The row index array of the coordinate format.</param>
|
/// <param name="columnIndices">The column index array of the coordinate format.</param>
|
/// <param name="values">The data array of the coordinate format.</param>
|
/// <returns>The sparse storage from the coordinate format.</returns>
|
/// <remarks>Duplicate entries will be summed together and
|
/// explicit zeros will be not intentionally removed.</remarks>
|
public static SparseCompressedRowMatrixStorage<T> OfCoordinateFormat(int rows, int columns, int valueCount, int[] rowIndices, int[] columnIndices, T[] values)
|
{
|
if (values == null) throw new NullReferenceException(nameof(values));
|
if (rowIndices == null) throw new NullReferenceException(nameof(rowIndices));
|
if (columnIndices == null) throw new NullReferenceException(nameof(columnIndices));
|
if (rowIndices.Length < valueCount || columnIndices.Length < valueCount || values.Length < valueCount)
|
{
|
throw new Exception($"The given array has the wrong length. Should be {valueCount}.");
|
}
|
|
// convert from COO to CSR
|
var storage = new SparseCompressedRowMatrixStorage<T>(rows, columns);
|
var csrRowPointers = storage.RowPointers;
|
var csrColumnIndices = new int[valueCount];
|
var csrValues = new T[valueCount];
|
|
for (int i = 0; i < valueCount; i++)
|
{
|
csrRowPointers[rowIndices[i]]++;
|
}
|
|
for (int i = 0, cumsum = 0; i < rows; i++)
|
{
|
var temp = csrRowPointers[i];
|
csrRowPointers[i] = cumsum;
|
cumsum += temp;
|
}
|
|
csrRowPointers[rows] = valueCount;
|
|
for (int i = 0; i < valueCount; i++)
|
{
|
var row = rowIndices[i];
|
var loc = csrRowPointers[row];
|
|
csrColumnIndices[loc] = columnIndices[i];
|
csrValues[loc] = values[i];
|
|
csrRowPointers[row]++;
|
}
|
|
for (int i = 0, last = 0; i <= rows; i++)
|
{
|
var temp = csrRowPointers[i];
|
csrRowPointers[i] = last;
|
last = temp;
|
}
|
|
storage.ColumnIndices = csrColumnIndices;
|
storage.Values = csrValues;
|
storage.NormalizeOrdering();
|
storage.NormalizeDuplicates();
|
return storage;
|
}
|
|
// MATRIX COPY
|
|
internal override void CopyToUnchecked(MatrixStorage<T> target, ExistingData existingData)
|
{
|
if (target is SparseCompressedRowMatrixStorage<T> sparseTarget)
|
{
|
CopyToUnchecked(sparseTarget);
|
return;
|
}
|
|
if (target is DenseColumnMajorMatrixStorage<T> denseTarget)
|
{
|
CopyToUnchecked(denseTarget, existingData);
|
return;
|
}
|
|
// FALL BACK
|
|
if (existingData == ExistingData.Clear)
|
{
|
target.Clear();
|
}
|
|
if (ValueCount != 0)
|
{
|
for (int row = 0; row < RowCount; row++)
|
{
|
var startIndex = RowPointers[row];
|
var endIndex = RowPointers[row + 1];
|
for (var j = startIndex; j < endIndex; j++)
|
{
|
target.At(row, ColumnIndices[j], Values[j]);
|
}
|
}
|
}
|
}
|
|
void CopyToUnchecked(SparseCompressedRowMatrixStorage<T> target)
|
{
|
target.Values = new T[ValueCount];
|
target.ColumnIndices = new int[ValueCount];
|
|
if (ValueCount != 0)
|
{
|
Array.Copy(Values, 0, target.Values, 0, ValueCount);
|
Buffer.BlockCopy(ColumnIndices, 0, target.ColumnIndices, 0, ValueCount*Constants.SizeOfInt);
|
Buffer.BlockCopy(RowPointers, 0, target.RowPointers, 0, (RowCount + 1)*Constants.SizeOfInt);
|
}
|
}
|
|
void CopyToUnchecked(DenseColumnMajorMatrixStorage<T> target, ExistingData existingData)
|
{
|
if (existingData == ExistingData.Clear)
|
{
|
target.Clear();
|
}
|
|
// TODO: proper implementation
|
|
if (ValueCount != 0)
|
{
|
for (int row = 0; row < RowCount; row++)
|
{
|
var startIndex = RowPointers[row];
|
var endIndex = RowPointers[row + 1];
|
for (var j = startIndex; j < endIndex; j++)
|
{
|
target.At(row, ColumnIndices[j], Values[j]);
|
}
|
}
|
}
|
}
|
|
internal override void CopySubMatrixToUnchecked(MatrixStorage<T> target,
|
int sourceRowIndex, int targetRowIndex, int rowCount,
|
int sourceColumnIndex, int targetColumnIndex, int columnCount,
|
ExistingData existingData)
|
{
|
if (target == null)
|
{
|
throw new ArgumentNullException(nameof(target));
|
}
|
|
if (target is SparseCompressedRowMatrixStorage<T> sparseTarget)
|
{
|
CopySubMatrixToUnchecked(sparseTarget,
|
sourceRowIndex, targetRowIndex, rowCount,
|
sourceColumnIndex, targetColumnIndex, columnCount,
|
existingData);
|
return;
|
}
|
|
// FALL BACK
|
|
if (existingData == ExistingData.Clear)
|
{
|
target.ClearUnchecked(targetRowIndex, rowCount, targetColumnIndex, columnCount);
|
}
|
|
for (int i = sourceRowIndex, row = 0; i < sourceRowIndex + rowCount; i++, row++)
|
{
|
var startIndex = RowPointers[i];
|
var endIndex = RowPointers[i + 1];
|
|
for (int j = startIndex; j < endIndex; j++)
|
{
|
// check if the column index is in the range
|
if ((ColumnIndices[j] >= sourceColumnIndex) && (ColumnIndices[j] < sourceColumnIndex + columnCount))
|
{
|
var column = ColumnIndices[j] - sourceColumnIndex;
|
target.At(targetRowIndex + row, targetColumnIndex + column, Values[j]);
|
}
|
}
|
}
|
}
|
|
void CopySubMatrixToUnchecked(SparseCompressedRowMatrixStorage<T> target,
|
int sourceRowIndex, int targetRowIndex, int rowCount,
|
int sourceColumnIndex, int targetColumnIndex, int columnCount,
|
ExistingData existingData)
|
{
|
var rowOffset = targetRowIndex - sourceRowIndex;
|
var columnOffset = targetColumnIndex - sourceColumnIndex;
|
|
// special case for empty target - much faster
|
if (target.ValueCount == 0)
|
{
|
// note: ValueCount is maximum resulting ValueCount (just using max to avoid internal copying)
|
// resulting arrays will likely be smaller - unless all values fit in the chosen range.
|
var values = new List<T>(ValueCount);
|
var columnIndices = new List<int>(ValueCount);
|
var rowPointers = target.RowPointers;
|
|
for (int i = sourceRowIndex; i < sourceRowIndex + rowCount; i++)
|
{
|
rowPointers[i + rowOffset] = values.Count;
|
|
var startIndex = RowPointers[i];
|
var endIndex = RowPointers[i + 1];
|
|
// note: we might be able to replace this loop with Array.Copy (perf)
|
for (int k = startIndex; k < endIndex; k++)
|
{
|
// check if the column index is in the range
|
if ((ColumnIndices[k] >= sourceColumnIndex) && (ColumnIndices[k] < sourceColumnIndex + columnCount))
|
{
|
values.Add(Values[k]);
|
columnIndices.Add(ColumnIndices[k] + columnOffset);
|
}
|
}
|
}
|
|
for (int i = targetRowIndex + rowCount; i < rowPointers.Length; i++)
|
{
|
rowPointers[i] = values.Count;
|
}
|
|
target.RowPointers[target.RowCount] = values.Count;
|
target.Values = values.ToArray();
|
target.ColumnIndices = columnIndices.ToArray();
|
|
return;
|
}
|
|
if (existingData == ExistingData.Clear)
|
{
|
target.ClearUnchecked(targetRowIndex, rowCount, targetColumnIndex, columnCount);
|
}
|
|
// NOTE: potential for more efficient implementation
|
for (int i = sourceRowIndex, row = 0; row < rowCount; i++, row++)
|
{
|
var startIndex = RowPointers[i];
|
var endIndex = RowPointers[i + 1];
|
|
for (int j = startIndex; j < endIndex; j++)
|
{
|
// check if the column index is in the range
|
if ((ColumnIndices[j] >= sourceColumnIndex) && (ColumnIndices[j] < sourceColumnIndex + columnCount))
|
{
|
var column = ColumnIndices[j] - sourceColumnIndex;
|
target.At(targetRowIndex + row, targetColumnIndex + column, Values[j]);
|
}
|
}
|
}
|
}
|
|
// ROW COPY
|
|
internal override void CopySubRowToUnchecked(VectorStorage<T> target, int rowIndex,
|
int sourceColumnIndex, int targetColumnIndex, int columnCount, ExistingData existingData)
|
{
|
// Determine bounds in columnIndices array where this item should be searched (using rowIndex)
|
var startIndexOfRow = RowPointers[rowIndex];
|
var endIndexOfRow = RowPointers[rowIndex + 1];
|
|
if (startIndexOfRow == endIndexOfRow)
|
{
|
if (existingData == ExistingData.Clear)
|
{
|
target.Clear(targetColumnIndex, columnCount);
|
}
|
return;
|
}
|
|
if (target is SparseVectorStorage<T> targetSparse)
|
{
|
if ((sourceColumnIndex == 0) && (targetColumnIndex == 0) && (columnCount == ColumnCount) && (ColumnCount == targetSparse.Length))
|
{
|
// rebuild of the values, indices, no clean necessary
|
targetSparse.ValueCount = endIndexOfRow - startIndexOfRow;
|
targetSparse.Values = new T[targetSparse.ValueCount];
|
targetSparse.Indices = new int[targetSparse.ValueCount];
|
Array.Copy(ColumnIndices, startIndexOfRow, targetSparse.Indices, 0, targetSparse.ValueCount);
|
Array.Copy(Values, startIndexOfRow, targetSparse.Values, 0, targetSparse.ValueCount);
|
}
|
else
|
{
|
int sourceStartPos = Array.BinarySearch(ColumnIndices, startIndexOfRow, endIndexOfRow - startIndexOfRow, sourceColumnIndex);
|
if (sourceStartPos < 0)
|
{
|
sourceStartPos = ~sourceStartPos;
|
}
|
int sourceEndPos = Array.BinarySearch(ColumnIndices, startIndexOfRow, endIndexOfRow - startIndexOfRow, sourceColumnIndex + columnCount);
|
if (sourceEndPos < 0)
|
{
|
sourceEndPos = ~sourceEndPos;
|
}
|
int positionsToCopy = sourceEndPos - sourceStartPos;
|
if (positionsToCopy > 0)
|
{
|
// rebuild the target (no clean necessary)
|
int targetStartPos = Array.BinarySearch(targetSparse.Indices,0, targetSparse.ValueCount, targetColumnIndex);
|
if (targetStartPos < 0)
|
{
|
targetStartPos = ~targetStartPos;
|
}
|
int targetEndPos = Array.BinarySearch(targetSparse.Indices,0,targetSparse.ValueCount, targetColumnIndex + columnCount);
|
if (targetEndPos < 0)
|
{
|
targetEndPos = Math.Max(~targetEndPos, targetStartPos);
|
}
|
int newValueCount = targetSparse.ValueCount - (targetEndPos - targetStartPos) + positionsToCopy;
|
T[] newValues = new T[newValueCount];
|
int[] newIndices = new int[newValueCount];
|
// copy before
|
Array.Copy(targetSparse.Indices, 0, newIndices, 0, targetStartPos);
|
Array.Copy(targetSparse.Values, 0, newValues, 0, targetStartPos);
|
// copy values themselves, with new positions
|
int shiftRight = targetColumnIndex - sourceColumnIndex;
|
for (int i = 0; i < positionsToCopy;++i)
|
{
|
newIndices[targetStartPos + i] = ColumnIndices[sourceStartPos + i] + shiftRight;
|
}
|
Array.Copy(Values, sourceStartPos, newValues, targetStartPos, positionsToCopy);
|
// copy after
|
Array.Copy(targetSparse.Indices, targetEndPos, newIndices, positionsToCopy + targetStartPos, targetSparse.ValueCount - targetEndPos);
|
Array.Copy(targetSparse.Values, targetEndPos, newValues, positionsToCopy + targetStartPos, targetSparse.ValueCount - targetEndPos);
|
targetSparse.Values = newValues;
|
targetSparse.Indices = newIndices;
|
targetSparse.ValueCount = newValueCount;
|
}
|
else
|
{
|
// although there are no values to copy, we still need to clean the existing values (if necessary)
|
if (existingData == ExistingData.Clear)
|
{
|
target.Clear(targetColumnIndex, columnCount);
|
}
|
}
|
}
|
return;
|
}
|
// FALLBACK
|
if (existingData == ExistingData.Clear)
|
{
|
target.Clear(targetColumnIndex, columnCount);
|
}
|
// If there are non-zero elements use base class implementation
|
for (int i = sourceColumnIndex, j = 0; i < sourceColumnIndex + columnCount; i++, j++)
|
{
|
var index = FindItem(rowIndex, i);
|
target.At(j, index >= 0 ? Values[index] : Zero);
|
}
|
}
|
|
// TRANSPOSE
|
|
internal override void TransposeToUnchecked(MatrixStorage<T> target, ExistingData existingData)
|
{
|
if (target is SparseCompressedRowMatrixStorage<T> sparseTarget)
|
{
|
TransposeToUnchecked(sparseTarget);
|
return;
|
}
|
|
if (target is DenseColumnMajorMatrixStorage<T> denseTarget)
|
{
|
TransposeToUnchecked(denseTarget, existingData);
|
return;
|
}
|
|
// FALL BACK
|
|
if (existingData == ExistingData.Clear)
|
{
|
target.Clear();
|
}
|
|
if (ValueCount != 0)
|
{
|
for (int row = 0; row < RowCount; row++)
|
{
|
var startIndex = RowPointers[row];
|
var endIndex = RowPointers[row + 1];
|
for (var j = startIndex; j < endIndex; j++)
|
{
|
target.At(ColumnIndices[j], row, Values[j]);
|
}
|
}
|
}
|
}
|
|
void TransposeToUnchecked(SparseCompressedRowMatrixStorage<T> target)
|
{
|
target.Values = new T[ValueCount];
|
target.ColumnIndices = new int[ValueCount];
|
var cx = target.Values;
|
var cp = target.RowPointers;
|
var ci = target.ColumnIndices;
|
|
// Column counts
|
int[] w = new int[ColumnCount];
|
for (int p = 0; p < RowPointers[RowCount]; p++)
|
{
|
w[ColumnIndices[p]]++;
|
}
|
|
// Column pointers
|
int nz = 0;
|
for (int i = 0; i < ColumnCount; i++)
|
{
|
cp[i] = nz;
|
nz += w[i];
|
w[i] = cp[i];
|
}
|
cp[ColumnCount] = nz;
|
|
for (int i = 0; i < RowCount; i++)
|
{
|
for (int p = RowPointers[i]; p < RowPointers[i + 1]; p++)
|
{
|
int j = w[ColumnIndices[p]]++;
|
|
// Place A(i,j) as entry C(j,i)
|
ci[j] = i;
|
cx[j] = Values[p];
|
}
|
}
|
}
|
|
void TransposeToUnchecked(DenseColumnMajorMatrixStorage<T> target, ExistingData existingData)
|
{
|
if (existingData == ExistingData.Clear)
|
{
|
target.Clear();
|
}
|
|
if (ValueCount != 0)
|
{
|
for (int row = 0; row < RowCount; row++)
|
{
|
var targetIndex = row * ColumnCount;
|
var startIndex = RowPointers[row];
|
var endIndex = RowPointers[row + 1];
|
for (var j = startIndex; j < endIndex; j++)
|
{
|
target.Data[targetIndex + ColumnIndices[j]] = Values[j];
|
}
|
}
|
}
|
}
|
|
internal override void TransposeSquareInplaceUnchecked()
|
{
|
var cx = new T[ValueCount]; //target.Values;
|
var cp = new int[RowCount + 1];
|
var ci = new int[ValueCount]; //target.ColumnIndices;
|
|
// Column counts
|
int[] w = new int[ColumnCount];
|
for (int p = 0; p < RowPointers[RowCount]; p++)
|
{
|
w[ColumnIndices[p]]++;
|
}
|
|
// Column pointers
|
int nz = 0;
|
for (int i = 0; i < ColumnCount; i++)
|
{
|
cp[i] = nz;
|
nz += w[i];
|
w[i] = cp[i];
|
}
|
cp[ColumnCount] = nz;
|
|
for (int i = 0; i < RowCount; i++)
|
{
|
for (int p = RowPointers[i]; p < RowPointers[i + 1]; p++)
|
{
|
int j = w[ColumnIndices[p]]++;
|
|
// Place A(i,j) as entry C(j,i)
|
ci[j] = i;
|
cx[j] = Values[p];
|
}
|
}
|
|
Array.Copy(cx, 0, Values, 0, ValueCount);
|
Buffer.BlockCopy(ci, 0, ColumnIndices, 0, ValueCount * Constants.SizeOfInt);
|
Buffer.BlockCopy(cp, 0, RowPointers, 0, (RowCount + 1) * Constants.SizeOfInt);
|
}
|
|
// EXTRACT
|
|
public override T[] ToRowMajorArray()
|
{
|
var ret = new T[RowCount*ColumnCount];
|
if (ValueCount != 0)
|
{
|
for (int row = 0; row < RowCount; row++)
|
{
|
var offset = row*ColumnCount;
|
var startIndex = RowPointers[row];
|
var endIndex = RowPointers[row + 1];
|
for (var j = startIndex; j < endIndex; j++)
|
{
|
ret[offset + ColumnIndices[j]] = Values[j];
|
}
|
}
|
}
|
return ret;
|
}
|
|
public override T[] ToColumnMajorArray()
|
{
|
var ret = new T[RowCount*ColumnCount];
|
if (ValueCount != 0)
|
{
|
for (int row = 0; row < RowCount; row++)
|
{
|
var startIndex = RowPointers[row];
|
var endIndex = RowPointers[row + 1];
|
for (var j = startIndex; j < endIndex; j++)
|
{
|
ret[(ColumnIndices[j])*RowCount + row] = Values[j];
|
}
|
}
|
}
|
return ret;
|
}
|
|
public override T[][] ToRowArrays()
|
{
|
var ret = new T[RowCount][];
|
if (ValueCount != 0)
|
{
|
for (int row = 0; row < RowCount; row++)
|
{
|
var array = new T[ColumnCount];
|
var startIndex = RowPointers[row];
|
var endIndex = RowPointers[row + 1];
|
for (var j = startIndex; j < endIndex; j++)
|
{
|
array[ColumnIndices[j]] = Values[j];
|
}
|
ret[row] = array;
|
}
|
}
|
return ret;
|
}
|
|
public override T[][] ToColumnArrays()
|
{
|
var ret = new T[ColumnCount][];
|
for (int j = 0; j < ColumnCount; j++)
|
{
|
ret[j] = new T[RowCount];
|
}
|
if (ValueCount != 0)
|
{
|
for (int row = 0; row < RowCount; row++)
|
{
|
var startIndex = RowPointers[row];
|
var endIndex = RowPointers[row + 1];
|
for (var j = startIndex; j < endIndex; j++)
|
{
|
ret[ColumnIndices[j]][row] = Values[j];
|
}
|
}
|
}
|
return ret;
|
}
|
|
public override T[,] ToArray()
|
{
|
var ret = new T[RowCount, ColumnCount];
|
if (ValueCount != 0)
|
{
|
for (int row = 0; row < RowCount; row++)
|
{
|
var startIndex = RowPointers[row];
|
var endIndex = RowPointers[row + 1];
|
for (var j = startIndex; j < endIndex; j++)
|
{
|
ret[row, ColumnIndices[j]] = Values[j];
|
}
|
}
|
}
|
return ret;
|
}
|
|
// ENUMERATION
|
|
public override IEnumerable<T> Enumerate()
|
{
|
int k = 0;
|
for (int row = 0; row < RowCount; row++)
|
{
|
for (int col = 0; col < ColumnCount; col++)
|
{
|
yield return k < RowPointers[row + 1] && ColumnIndices[k] == col
|
? Values[k++]
|
: Zero;
|
}
|
}
|
}
|
|
public override IEnumerable<Tuple<int, int, T>> EnumerateIndexed()
|
{
|
int k = 0;
|
for (int row = 0; row < RowCount; row++)
|
{
|
for (int col = 0; col < ColumnCount; col++)
|
{
|
yield return k < RowPointers[row + 1] && ColumnIndices[k] == col
|
? new Tuple<int, int, T>(row, col, Values[k++])
|
: new Tuple<int, int, T>(row, col, Zero);
|
}
|
}
|
}
|
|
public override IEnumerable<T> EnumerateNonZero()
|
{
|
return Values.Take(ValueCount).Where(x => !Zero.Equals(x));
|
}
|
|
public override IEnumerable<Tuple<int, int, T>> EnumerateNonZeroIndexed()
|
{
|
for (int row = 0; row < RowCount; row++)
|
{
|
var startIndex = RowPointers[row];
|
var endIndex = RowPointers[row + 1];
|
for (var j = startIndex; j < endIndex; j++)
|
{
|
if (!Zero.Equals(Values[j]))
|
{
|
yield return new Tuple<int, int, T>(row, ColumnIndices[j], Values[j]);
|
}
|
}
|
}
|
}
|
|
// FIND
|
|
public override Tuple<int, int, T> Find(Func<T, bool> predicate, Zeros zeros)
|
{
|
for (int row = 0; row < RowCount; row++)
|
{
|
var startIndex = RowPointers[row];
|
var endIndex = RowPointers[row + 1];
|
for (var j = startIndex; j < endIndex; j++)
|
{
|
if (predicate(Values[j]))
|
{
|
return new Tuple<int, int, T>(row, ColumnIndices[j], Values[j]);
|
}
|
}
|
}
|
if (zeros == Zeros.Include && ValueCount < (RowCount * ColumnCount))
|
{
|
if (predicate(Zero))
|
{
|
int k = 0;
|
for (int row = 0; row < RowCount; row++)
|
{
|
for (int col = 0; col < ColumnCount; col++)
|
{
|
if (k < RowPointers[row + 1] && ColumnIndices[k] == col)
|
{
|
k++;
|
}
|
else
|
{
|
return new Tuple<int, int, T>(row, col, Zero);
|
}
|
}
|
}
|
}
|
}
|
return null;
|
}
|
|
internal override Tuple<int, int, T, TOther> Find2Unchecked<TOther>(MatrixStorage<TOther> other, Func<T, TOther, bool> predicate, Zeros zeros)
|
{
|
if (other is DenseColumnMajorMatrixStorage<TOther> denseOther)
|
{
|
TOther[] otherData = denseOther.Data;
|
int k = 0;
|
for (int row = 0; row < RowCount; row++)
|
{
|
for (int col = 0; col < ColumnCount; col++)
|
{
|
bool available = k < RowPointers[row + 1] && ColumnIndices[k] == col;
|
if (predicate(available ? Values[k++] : Zero, otherData[col*RowCount + row]))
|
{
|
return new Tuple<int, int, T, TOther>(row, col, available ? Values[k - 1] : Zero, otherData[col*RowCount + row]);
|
}
|
}
|
}
|
return null;
|
}
|
|
if (other is DiagonalMatrixStorage<TOther> diagonalOther)
|
{
|
TOther[] otherData = diagonalOther.Data;
|
TOther otherZero = BuilderInstance<TOther>.Matrix.Zero;
|
|
// Full Scan
|
if (zeros == Zeros.Include && predicate(Zero, otherZero))
|
{
|
int k = 0;
|
for (int row = 0; row < RowCount; row++)
|
{
|
for (int col = 0; col < ColumnCount; col++)
|
{
|
bool available = k < RowPointers[row + 1] && ColumnIndices[k] == col;
|
if (predicate(available ? Values[k++] : Zero, row == col ? otherData[row] : otherZero))
|
{
|
return new Tuple<int, int, T, TOther>(row, col, available ? Values[k - 1] : Zero, row == col ? otherData[row] : otherZero);
|
}
|
}
|
}
|
return null;
|
}
|
|
// Sparse Scan
|
for (int row = 0; row < RowCount; row++)
|
{
|
bool diagonal = false;
|
var startIndex = RowPointers[row];
|
var endIndex = RowPointers[row + 1];
|
for (var j = startIndex; j < endIndex; j++)
|
{
|
if (ColumnIndices[j] == row)
|
{
|
diagonal = true;
|
if (predicate(Values[j], otherData[row]))
|
{
|
return new Tuple<int, int, T, TOther>(row, row, Values[j], otherData[row]);
|
}
|
}
|
else
|
{
|
if (predicate(Values[j], otherZero))
|
{
|
return new Tuple<int, int, T, TOther>(row, ColumnIndices[j], Values[j], otherZero);
|
}
|
}
|
}
|
if (!diagonal && row < ColumnCount)
|
{
|
if (predicate(Zero, otherData[row]))
|
{
|
return new Tuple<int, int, T, TOther>(row, row, Zero, otherData[row]);
|
}
|
}
|
}
|
return null;
|
}
|
|
if (other is SparseCompressedRowMatrixStorage<TOther> sparseOther)
|
{
|
int[] otherRowPointers = sparseOther.RowPointers;
|
int[] otherColumnIndices = sparseOther.ColumnIndices;
|
TOther[] otherValues = sparseOther.Values;
|
TOther otherZero = BuilderInstance<TOther>.Matrix.Zero;
|
|
if (zeros == Zeros.Include)
|
{
|
int k = 0, otherk = 0;
|
for (int row = 0; row < RowCount; row++)
|
{
|
for (int col = 0; col < ColumnCount; col++)
|
{
|
bool available = k < RowPointers[row + 1] && ColumnIndices[k] == col;
|
bool otherAvailable = otherk < otherRowPointers[row + 1] && otherColumnIndices[otherk] == col;
|
if (predicate(available ? Values[k++] : Zero, otherAvailable ? otherValues[otherk++] : otherZero))
|
{
|
return new Tuple<int, int, T, TOther>(row, col, available ? Values[k - 1] : Zero, otherAvailable ? otherValues[otherk - 1] : otherZero);
|
}
|
}
|
}
|
return null;
|
}
|
|
for (int row = 0; row < RowCount; row++)
|
{
|
var endIndex = RowPointers[row + 1];
|
var otherEndIndex = otherRowPointers[row + 1];
|
var k = RowPointers[row];
|
var otherk = otherRowPointers[row];
|
while (k < endIndex || otherk < otherEndIndex)
|
{
|
if (k == endIndex || otherk < otherEndIndex && ColumnIndices[k] > otherColumnIndices[otherk])
|
{
|
if (predicate(Zero, otherValues[otherk++]))
|
{
|
return new Tuple<int, int, T, TOther>(row, otherColumnIndices[otherk - 1], Zero, otherValues[otherk - 1]);
|
}
|
}
|
else if (otherk == otherEndIndex || ColumnIndices[k] < otherColumnIndices[otherk])
|
{
|
if (predicate(Values[k++], otherZero))
|
{
|
return new Tuple<int, int, T, TOther>(row, ColumnIndices[k - 1], Values[k - 1], otherZero);
|
}
|
}
|
else
|
{
|
if (predicate(Values[k++], otherValues[otherk++]))
|
{
|
return new Tuple<int, int, T, TOther>(row, ColumnIndices[k - 1], Values[k - 1], otherValues[otherk - 1]);
|
}
|
}
|
}
|
}
|
return null;
|
}
|
|
// FALL BACK
|
|
return base.Find2Unchecked(other, predicate, zeros);
|
}
|
|
// FUNCTIONAL COMBINATORS: MAP
|
|
public override void MapInplace(Func<T, T> f, Zeros zeros)
|
{
|
if (zeros == Zeros.Include || !Zero.Equals(f(Zero)))
|
{
|
var newRowPointers = RowPointers;
|
var newColumnIndices = new List<int>(ColumnIndices.Length);
|
var newValues = new List<T>(Values.Length);
|
|
int k = 0;
|
for (int row = 0; row < RowCount; row++)
|
{
|
newRowPointers[row] = newValues.Count;
|
for (int col = 0; col < ColumnCount; col++)
|
{
|
var item = k < RowPointers[row + 1] && ColumnIndices[k] == col ? f(Values[k++]) : f(Zero);
|
if (!Zero.Equals(item))
|
{
|
newValues.Add(item);
|
newColumnIndices.Add(col);
|
}
|
}
|
}
|
|
ColumnIndices = newColumnIndices.ToArray();
|
Values = newValues.ToArray();
|
newRowPointers[RowCount] = newValues.Count;
|
}
|
else
|
{
|
// we can safely do this in-place:
|
int nonZero = 0;
|
for (int row = 0; row < RowCount; row++)
|
{
|
var startIndex = RowPointers[row];
|
var endIndex = RowPointers[row + 1];
|
RowPointers[row] = nonZero;
|
for (var j = startIndex; j < endIndex; j++)
|
{
|
var item = f(Values[j]);
|
if (!Zero.Equals(item))
|
{
|
Values[nonZero] = item;
|
ColumnIndices[nonZero] = ColumnIndices[j];
|
nonZero++;
|
}
|
}
|
}
|
Array.Resize(ref ColumnIndices, nonZero);
|
Array.Resize(ref Values, nonZero);
|
RowPointers[RowCount] = nonZero;
|
}
|
}
|
|
public override void MapIndexedInplace(Func<int, int, T, T> f, Zeros zeros)
|
{
|
if (zeros == Zeros.Include || !Zero.Equals(f(0, 1, Zero)))
|
{
|
var newRowPointers = RowPointers;
|
var newColumnIndices = new List<int>(ColumnIndices.Length);
|
var newValues = new List<T>(Values.Length);
|
|
int k = 0;
|
for (int row = 0; row < RowCount; row++)
|
{
|
newRowPointers[row] = newValues.Count;
|
for (int col = 0; col < ColumnCount; col++)
|
{
|
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);
|
newColumnIndices.Add(col);
|
}
|
}
|
}
|
|
ColumnIndices = newColumnIndices.ToArray();
|
Values = newValues.ToArray();
|
newRowPointers[RowCount] = newValues.Count;
|
}
|
else
|
{
|
// we can safely do this in-place:
|
int nonZero = 0;
|
for (int row = 0; row < RowCount; row++)
|
{
|
var startIndex = RowPointers[row];
|
var endIndex = RowPointers[row + 1];
|
RowPointers[row] = nonZero;
|
for (var j = startIndex; j < endIndex; j++)
|
{
|
var item = f(row, ColumnIndices[j], Values[j]);
|
if (!Zero.Equals(item))
|
{
|
Values[nonZero] = item;
|
ColumnIndices[nonZero] = ColumnIndices[j];
|
nonZero++;
|
}
|
}
|
}
|
Array.Resize(ref ColumnIndices, nonZero);
|
Array.Resize(ref Values, nonZero);
|
RowPointers[RowCount] = nonZero;
|
}
|
}
|
|
internal override void MapToUnchecked<TU>(MatrixStorage<TU> target, Func<T, TU> f, Zeros zeros, ExistingData existingData)
|
{
|
var processZeros = zeros == Zeros.Include || !Zero.Equals(f(Zero));
|
|
if (target is SparseCompressedRowMatrixStorage<TU> sparseTarget)
|
{
|
var newRowPointers = sparseTarget.RowPointers;
|
var newColumnIndices = new List<int>(ColumnIndices.Length);
|
var newValues = new List<TU>(Values.Length);
|
|
if (processZeros)
|
{
|
int k = 0;
|
for (int row = 0; row < RowCount; row++)
|
{
|
newRowPointers[row] = newValues.Count;
|
for (int col = 0; col < ColumnCount; col++)
|
{
|
var item = k < RowPointers[row + 1] && ColumnIndices[k] == col ? f(Values[k++]) : f(Zero);
|
if (!Zero.Equals(item))
|
{
|
newValues.Add(item);
|
newColumnIndices.Add(col);
|
}
|
}
|
}
|
}
|
else
|
{
|
for (int row = 0; row < RowCount; row++)
|
{
|
newRowPointers[row] = newValues.Count;
|
var startIndex = RowPointers[row];
|
var endIndex = RowPointers[row + 1];
|
for (var j = startIndex; j < endIndex; j++)
|
{
|
var item = f(Values[j]);
|
if (!Zero.Equals(item))
|
{
|
newValues.Add(item);
|
newColumnIndices.Add(ColumnIndices[j]);
|
}
|
}
|
}
|
}
|
|
sparseTarget.ColumnIndices = newColumnIndices.ToArray();
|
sparseTarget.Values = newValues.ToArray();
|
newRowPointers[RowCount] = newValues.Count;
|
return;
|
}
|
|
// FALL BACK
|
|
if (existingData == ExistingData.Clear && !processZeros)
|
{
|
target.Clear();
|
}
|
|
if (processZeros)
|
{
|
for (int row = 0; row < RowCount; row++)
|
{
|
var index = RowPointers[row];
|
var endIndex = RowPointers[row + 1];
|
for (int j = 0; j < ColumnCount; j++)
|
{
|
if (index < endIndex && j == ColumnIndices[index])
|
{
|
target.At(row, j, f(Values[index]));
|
index = Math.Min(index + 1, endIndex);
|
}
|
else
|
{
|
target.At(row, j, f(Zero));
|
}
|
}
|
}
|
}
|
else
|
{
|
for (int row = 0; row < RowCount; row++)
|
{
|
var startIndex = RowPointers[row];
|
var endIndex = RowPointers[row + 1];
|
for (var j = startIndex; j < endIndex; j++)
|
{
|
target.At(row, ColumnIndices[j], f(Values[j]));
|
}
|
}
|
}
|
}
|
|
internal override void MapIndexedToUnchecked<TU>(MatrixStorage<TU> target, Func<int, int, T, TU> f, Zeros zeros, ExistingData existingData)
|
{
|
var processZeros = zeros == Zeros.Include || !Zero.Equals(f(0, 1, Zero));
|
|
if (target is SparseCompressedRowMatrixStorage<TU> sparseTarget)
|
{
|
var newRowPointers = sparseTarget.RowPointers;
|
var newColumnIndices = new List<int>(ColumnIndices.Length);
|
var newValues = new List<TU>(Values.Length);
|
|
if (processZeros)
|
{
|
int k = 0;
|
for (int row = 0; row < RowCount; row++)
|
{
|
newRowPointers[row] = newValues.Count;
|
for (int col = 0; col < ColumnCount; col++)
|
{
|
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);
|
newColumnIndices.Add(col);
|
}
|
}
|
}
|
}
|
else
|
{
|
for (int row = 0; row < RowCount; row++)
|
{
|
newRowPointers[row] = newValues.Count;
|
var startIndex = RowPointers[row];
|
var endIndex = RowPointers[row + 1];
|
for (var j = startIndex; j < endIndex; j++)
|
{
|
var item = f(row, ColumnIndices[j], Values[j]);
|
if (!Zero.Equals(item))
|
{
|
newValues.Add(item);
|
newColumnIndices.Add(ColumnIndices[j]);
|
}
|
}
|
}
|
}
|
|
sparseTarget.ColumnIndices = newColumnIndices.ToArray();
|
sparseTarget.Values = newValues.ToArray();
|
newRowPointers[RowCount] = newValues.Count;
|
return;
|
}
|
|
// FALL BACK
|
|
if (existingData == ExistingData.Clear && !processZeros)
|
{
|
target.Clear();
|
}
|
|
if (processZeros)
|
{
|
for (int row = 0; row < RowCount; row++)
|
{
|
var index = RowPointers[row];
|
var endIndex = RowPointers[row + 1];
|
for (int j = 0; j < ColumnCount; j++)
|
{
|
if (index < endIndex && j == ColumnIndices[index])
|
{
|
target.At(row, j, f(row, j, Values[index]));
|
index = Math.Min(index + 1, endIndex);
|
}
|
else
|
{
|
target.At(row, j, f(row, j, Zero));
|
}
|
}
|
}
|
}
|
else
|
{
|
for (int row = 0; row < RowCount; row++)
|
{
|
var startIndex = RowPointers[row];
|
var endIndex = RowPointers[row + 1];
|
for (var j = startIndex; j < endIndex; j++)
|
{
|
target.At(row, ColumnIndices[j], f(row, ColumnIndices[j], Values[j]));
|
}
|
}
|
}
|
}
|
|
internal override void MapSubMatrixIndexedToUnchecked<TU>(MatrixStorage<TU> target, Func<int, int, T, TU> f,
|
int sourceRowIndex, int targetRowIndex, int rowCount,
|
int sourceColumnIndex, int targetColumnIndex, int columnCount,
|
Zeros zeros, ExistingData existingData)
|
{
|
if (target is SparseCompressedRowMatrixStorage<TU> sparseTarget)
|
{
|
MapSubMatrixIndexedToUnchecked(sparseTarget, f, sourceRowIndex, targetRowIndex, rowCount, sourceColumnIndex, targetColumnIndex, columnCount, zeros, existingData);
|
return;
|
}
|
|
// FALL BACK
|
|
var processZeros = zeros == Zeros.Include || !Zero.Equals(f(0, 1, Zero));
|
if (existingData == ExistingData.Clear && !processZeros)
|
{
|
target.ClearUnchecked(targetRowIndex, rowCount, targetColumnIndex, columnCount);
|
}
|
|
if (processZeros)
|
{
|
for (int sr = sourceRowIndex, tr = targetRowIndex; sr < sourceRowIndex + rowCount; sr++, tr++)
|
{
|
var index = RowPointers[sr];
|
var endIndex = RowPointers[sr + 1];
|
|
// move forward to our sub-range
|
for (; ColumnIndices[index] < sourceColumnIndex && index < endIndex; index++)
|
{
|
}
|
for (int sc = sourceColumnIndex, tc = targetColumnIndex; sc < sourceColumnIndex + columnCount; sc++, tc++)
|
{
|
if (index < endIndex && sc == ColumnIndices[index])
|
{
|
target.At(tr, tc, f(tr, tc, Values[index]));
|
index = Math.Min(index + 1, endIndex);
|
}
|
else
|
{
|
target.At(tr, tc, f(tr, tc, Zero));
|
}
|
}
|
}
|
}
|
else
|
{
|
int columnOffset = targetColumnIndex - sourceColumnIndex;
|
for (int sr = sourceRowIndex, tr = targetRowIndex; sr < sourceRowIndex + rowCount; sr++, tr++)
|
{
|
var startIndex = RowPointers[sr];
|
var endIndex = RowPointers[sr + 1];
|
for (int k = startIndex; k < endIndex; k++)
|
{
|
// check if the column index is in the range
|
if ((ColumnIndices[k] >= sourceColumnIndex) && (ColumnIndices[k] < sourceColumnIndex + columnCount))
|
{
|
int tc = ColumnIndices[k] + columnOffset;
|
target.At(tr, tc, f(tr, tc, Values[k]));
|
}
|
}
|
}
|
}
|
}
|
|
void MapSubMatrixIndexedToUnchecked<TU>(SparseCompressedRowMatrixStorage<TU> target, Func<int, int, T, TU> f,
|
int sourceRowIndex, int targetRowIndex, int rowCount,
|
int sourceColumnIndex, int targetColumnIndex, int columnCount,
|
Zeros zeros, ExistingData existingData)
|
where TU : struct, IEquatable<TU>, IFormattable
|
{
|
var processZeros = zeros == Zeros.Include || !Zero.Equals(f(0, 1, Zero));
|
if (existingData == ExistingData.Clear && !processZeros)
|
{
|
target.ClearUnchecked(targetRowIndex, rowCount, targetColumnIndex, columnCount);
|
}
|
|
var rowOffset = targetRowIndex - sourceRowIndex;
|
var columnOffset = targetColumnIndex - sourceColumnIndex;
|
var zero = Matrix<TU>.Zero;
|
|
// special case for empty target - much faster
|
if (target.ValueCount == 0)
|
{
|
var values = new List<TU>(ValueCount);
|
var columnIndices = new List<int>(ValueCount);
|
var rowPointers = target.RowPointers;
|
|
if (processZeros)
|
{
|
for (int sr = sourceRowIndex; sr < sourceRowIndex + rowCount; sr++)
|
{
|
int tr = sr + rowOffset;
|
rowPointers[tr] = values.Count;
|
|
var index = RowPointers[sr];
|
var endIndex = RowPointers[sr + 1];
|
|
// move forward to our sub-range
|
for (; ColumnIndices[index] < sourceColumnIndex && index < endIndex; index++)
|
{
|
}
|
for (int sc = sourceColumnIndex, tc = targetColumnIndex; sc < sourceColumnIndex + columnCount; sc++, tc++)
|
{
|
if (index < endIndex && sc == ColumnIndices[index])
|
{
|
TU item = f(tr, tc, Values[index]);
|
if (!zero.Equals(item))
|
{
|
values.Add(item);
|
columnIndices.Add(tc);
|
}
|
index = Math.Min(index + 1, endIndex);
|
}
|
else
|
{
|
TU item = f(tr, tc, Zero);
|
if (!zero.Equals(item))
|
{
|
values.Add(item);
|
columnIndices.Add(tc);
|
}
|
}
|
}
|
}
|
}
|
else
|
{
|
for (int sr = sourceRowIndex; sr < sourceRowIndex + rowCount; sr++)
|
{
|
int tr = sr + rowOffset;
|
rowPointers[tr] = values.Count;
|
|
var startIndex = RowPointers[sr];
|
var endIndex = RowPointers[sr + 1];
|
|
for (int k = startIndex; k < endIndex; k++)
|
{
|
// check if the column index is in the range
|
if ((ColumnIndices[k] >= sourceColumnIndex) && (ColumnIndices[k] < sourceColumnIndex + columnCount))
|
{
|
int tc = ColumnIndices[k] + columnOffset;
|
TU item = f(tr, tc, Values[k]);
|
if (!zero.Equals(item))
|
{
|
values.Add(item);
|
columnIndices.Add(tc);
|
}
|
}
|
}
|
}
|
}
|
|
for (int i = targetRowIndex + rowCount; i < rowPointers.Length; i++)
|
{
|
rowPointers[i] = values.Count;
|
}
|
|
target.RowPointers[target.RowCount] = values.Count;
|
target.Values = values.ToArray();
|
target.ColumnIndices = columnIndices.ToArray();
|
return;
|
}
|
|
// TODO: proper general sparse case - the following is essentially a fall back, not leveraging the target data structure
|
|
if (processZeros)
|
{
|
for (int sr = sourceRowIndex, tr = targetRowIndex; sr < sourceRowIndex + rowCount; sr++, tr++)
|
{
|
var index = RowPointers[sr];
|
var endIndex = RowPointers[sr + 1];
|
|
// move forward to our sub-range
|
for (; ColumnIndices[index] < sourceColumnIndex && index < endIndex; index++)
|
{
|
}
|
for (int sc = sourceColumnIndex, tc = targetColumnIndex; sc < sourceColumnIndex + columnCount; sc++, tc++)
|
{
|
if (index < endIndex && sc == ColumnIndices[index])
|
{
|
target.At(tr, tc, f(tr, tc, Values[index]));
|
index = Math.Min(index + 1, endIndex);
|
}
|
else
|
{
|
target.At(tr, tc, f(tr, tc, Zero));
|
}
|
}
|
}
|
}
|
else
|
{
|
for (int sr = sourceRowIndex, tr = targetRowIndex; sr < sourceRowIndex + rowCount; sr++, tr++)
|
{
|
var startIndex = RowPointers[sr];
|
var endIndex = RowPointers[sr + 1];
|
for (int k = startIndex; k < endIndex; k++)
|
{
|
// check if the column index is in the range
|
if ((ColumnIndices[k] >= sourceColumnIndex) && (ColumnIndices[k] < sourceColumnIndex + columnCount))
|
{
|
int tc = ColumnIndices[k] + columnOffset;
|
target.At(tr, tc, f(tr, tc, Values[k]));
|
}
|
}
|
}
|
}
|
}
|
|
// FUNCTIONAL COMBINATORS: FOLD
|
|
internal override void FoldByRowUnchecked<TU>(TU[] target, Func<TU, T, TU> f, Func<TU, int, TU> finalize, TU[] state, Zeros zeros)
|
{
|
if (zeros == Zeros.AllowSkip)
|
{
|
for (int row = 0; row < RowCount; row++)
|
{
|
var startIndex = RowPointers[row];
|
var endIndex = RowPointers[row + 1];
|
TU s = state[row];
|
for (var j = startIndex; j < endIndex; j++)
|
{
|
s = f(s, Values[j]);
|
}
|
target[row] = finalize(s, endIndex - startIndex);
|
}
|
}
|
else
|
{
|
for (int row = 0; row < RowCount; row++)
|
{
|
var index = RowPointers[row];
|
var endIndex = RowPointers[row + 1];
|
TU s = state[row];
|
for (int j = 0; j < ColumnCount; j++)
|
{
|
if (index < endIndex && j == ColumnIndices[index])
|
{
|
s = f(s, Values[index]);
|
index = Math.Min(index + 1, endIndex);
|
}
|
else
|
{
|
s = f(s, Zero);
|
}
|
}
|
target[row] = finalize(s, ColumnCount);
|
}
|
}
|
}
|
|
internal override void FoldByColumnUnchecked<TU>(TU[] target, Func<TU, T, TU> f, Func<TU, int, TU> finalize, TU[] state, Zeros zeros)
|
{
|
if (!ReferenceEquals(state, target))
|
{
|
Array.Copy(state, 0, target, 0, state.Length);
|
}
|
if (zeros == Zeros.AllowSkip)
|
{
|
int[] count = new int[ColumnCount];
|
for (int row = 0; row < RowCount; row++)
|
{
|
var startIndex = RowPointers[row];
|
var endIndex = RowPointers[row + 1];
|
for (var j = startIndex; j < endIndex; j++)
|
{
|
var column = ColumnIndices[j];
|
target[column] = f(target[column], Values[j]);
|
count[column]++;
|
}
|
}
|
for (int j = 0; j < ColumnCount; j++)
|
{
|
target[j] = finalize(target[j], count[j]);
|
}
|
}
|
else
|
{
|
for (int row = 0; row < RowCount; row++)
|
{
|
var index = RowPointers[row];
|
var endIndex = RowPointers[row + 1];
|
for (int j = 0; j < ColumnCount; j++)
|
{
|
if (index < endIndex && j == ColumnIndices[index])
|
{
|
target[j] = f(target[j], Values[index]);
|
index = Math.Min(index + 1, endIndex);
|
}
|
else
|
{
|
target[j] = f(target[j], Zero);
|
}
|
}
|
}
|
for (int j = 0; j < ColumnCount; j++)
|
{
|
target[j] = finalize(target[j], RowCount);
|
}
|
}
|
}
|
|
internal override TState Fold2Unchecked<TOther, TState>(MatrixStorage<TOther> other, Func<TState, T, TOther, TState> f, TState state, Zeros zeros)
|
{
|
if (other is DenseColumnMajorMatrixStorage<TOther> denseOther)
|
{
|
TOther[] otherData = denseOther.Data;
|
int k = 0;
|
for (int row = 0; row < RowCount; row++)
|
{
|
for (int col = 0; col < ColumnCount; col++)
|
{
|
bool available = k < RowPointers[row + 1] && ColumnIndices[k] == col;
|
state = f(state, available ? Values[k++] : Zero, otherData[col*RowCount + row]);
|
}
|
}
|
return state;
|
}
|
|
if (other is DiagonalMatrixStorage<TOther> diagonalOther)
|
{
|
TOther[] otherData = diagonalOther.Data;
|
TOther otherZero = BuilderInstance<TOther>.Matrix.Zero;
|
|
if (zeros == Zeros.Include)
|
{
|
int k = 0;
|
for (int row = 0; row < RowCount; row++)
|
{
|
for (int col = 0; col < ColumnCount; col++)
|
{
|
bool available = k < RowPointers[row + 1] && ColumnIndices[k] == col;
|
state = f(state, available ? Values[k++] : Zero, row == col ? otherData[row] : otherZero);
|
}
|
}
|
return state;
|
}
|
|
for (int row = 0; row < RowCount; row++)
|
{
|
bool diagonal = false;
|
|
var startIndex = RowPointers[row];
|
var endIndex = RowPointers[row + 1];
|
for (var j = startIndex; j < endIndex; j++)
|
{
|
if (ColumnIndices[j] == row)
|
{
|
diagonal = true;
|
state = f(state, Values[j], otherData[row]);
|
}
|
else
|
{
|
state = f(state, Values[j], otherZero);
|
}
|
}
|
|
if (!diagonal && row < ColumnCount)
|
{
|
state = f(state, Zero, otherData[row]);
|
}
|
}
|
return state;
|
}
|
|
if (other is SparseCompressedRowMatrixStorage<TOther> sparseOther)
|
{
|
int[] otherRowPointers = sparseOther.RowPointers;
|
int[] otherColumnIndices = sparseOther.ColumnIndices;
|
TOther[] otherValues = sparseOther.Values;
|
TOther otherZero = BuilderInstance<TOther>.Matrix.Zero;
|
|
if (zeros == Zeros.Include)
|
{
|
int k = 0, otherk = 0;
|
for (int row = 0; row < RowCount; row++)
|
{
|
for (int col = 0; col < ColumnCount; col++)
|
{
|
bool available = k < RowPointers[row + 1] && ColumnIndices[k] == col;
|
bool otherAvailable = otherk < otherRowPointers[row + 1] && otherColumnIndices[otherk] == col;
|
state = f(state, available ? Values[k++] : Zero, otherAvailable ? otherValues[otherk++] : otherZero);
|
}
|
}
|
return state;
|
}
|
|
for (int row = 0; row < RowCount; row++)
|
{
|
var startIndex = RowPointers[row];
|
var endIndex = RowPointers[row + 1];
|
var otherStartIndex = otherRowPointers[row];
|
var otherEndIndex = otherRowPointers[row + 1];
|
|
var j1 = startIndex;
|
var j2 = otherStartIndex;
|
|
while (j1 < endIndex || j2 < otherEndIndex)
|
{
|
if (j1 == endIndex || j2 < otherEndIndex && ColumnIndices[j1] > otherColumnIndices[j2])
|
{
|
state = f(state, Zero, otherValues[j2++]);
|
}
|
else if (j2 == otherEndIndex || ColumnIndices[j1] < otherColumnIndices[j2])
|
{
|
state = f(state, Values[j1++], otherZero);
|
}
|
else
|
{
|
state = f(state, Values[j1++], otherValues[j2++]);
|
}
|
}
|
}
|
return state;
|
}
|
|
// FALL BACK
|
|
return base.Fold2Unchecked(other, f, state, zeros);
|
}
|
}
|
}
|
