codigo - matriz dispersa c++

¿Cuál es la mejor manera de crear una matriz dispersa en C++?

Question

Estoy trabajando en un proyecto que requiere la manipulación de matrices enormes, específicamente suma piramidal para un cálculo de cópula.

En resumen, necesito hacer un seguimiento de un número relativamente pequeño de valores (generalmente un valor de 1, y en casos excepcionales, más de 1) en un mar de ceros en la matriz (matriz multidimensional).

Una matriz dispersa permite al usuario almacenar una pequeña cantidad de valores y asumir que todos los registros indefinidos son un valor preestablecido. Como físicamente no es posible almacenar todos los valores en la memoria, necesito almacenar solo los pocos elementos que no sean cero. Esto podría ser varios millones de entradas.

La velocidad es una gran prioridad, y también me gustaría elegir dinámicamente el número de variables en la clase en tiempo de ejecución.

Actualmente trabajo en un sistema que usa un árbol de búsqueda binario (b-tree) para almacenar entradas. ¿Alguien sabe de un mejor sistema?

Answer 1

Aquí hay una implementación relativamente simple que debe proporcionar una búsqueda rápida razonable (usando una tabla hash) así como una iteración rápida sobre elementos distintos de cero en una fila / columna.

// Copyright 2014 Leo Osvald
//
// Licensed under the Apache License, Version 2.0 (the "License");
// you may not use this file except in compliance with the License.
// You may obtain a copy of the License at
//
// http://www.apache.org/licenses/LICENSE-2.0
//
// Unless required by applicable law or agreed to in writing, software
// distributed under the License is distributed on an "AS IS" BASIS,
// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
// See the License for the specific language governing permissions and
// limitations under the License.
#ifndef UTIL_IMMUTABLE_SPARSE_MATRIX_HPP_
#define UTIL_IMMUTABLE_SPARSE_MATRIX_HPP_
#include <algorithm>
#include <limits>
#include <map>
#include <type_traits>
#include <unordered_map>
#include <utility>
#include <vector>
// A simple time-efficient implementation of an immutable sparse matrix
// Provides efficient iteration of non-zero elements by rows/cols,
// e.g. to iterate over a range [row_from, row_to) x [col_from, col_to):
// for (int row = row_from; row < row_to; ++row) {
// for (auto col_range = sm.nonzero_col_range(row, col_from, col_to);
// col_range.first != col_range.second; ++col_range.first) {
// int col = *col_range.first;
// // use sm(row, col)
// ...
// }
template<typename T = double, class Coord = int>
class SparseMatrix {
struct PointHasher;
typedef std::map< Coord, std::vector<Coord> > NonZeroList;
typedef std::pair<Coord, Coord> Point;
public:
typedef T ValueType;
typedef Coord CoordType;
typedef typename NonZeroList::mapped_type::const_iterator CoordIter;
typedef std::pair<CoordIter, CoordIter> CoordIterRange;
SparseMatrix() = default;
// Reads a matrix stored in MatrixMarket-like format, i.e.:
// <num_rows> <num_cols> <num_entries>
// <row_1> <col_1> <val_1>
// ...
// Note: the header (lines starting with ''%'' are ignored).
template<class InputStream, size_t max_line_length = 1024>
void Init(InputStream& is) {
rows_.clear(), cols_.clear();
values_.clear();
// skip the header (lines beginning with ''%'', if any)
decltype(is.tellg()) offset = 0;
for (char buf[max_line_length + 1];
is.getline(buf, sizeof(buf)) && buf[0] == ''%''; )
offset = is.tellg();
is.seekg(offset);
size_t n;
is >> row_count_ >> col_count_ >> n;
values_.reserve(n);
while (n--) {
Coord row, col;
typename std::remove_cv<T>::type val;
is >> row >> col >> val;
values_[Point(--row, --col)] = val;
rows_[col].push_back(row);
cols_[row].push_back(col);
}
SortAndShrink(rows_);
SortAndShrink(cols_);
}
const T& operator()(const Coord& row, const Coord& col) const {
static const T kZero = T();
auto it = values_.find(Point(row, col));
if (it != values_.end())
return it->second;
return kZero;
}
CoordIterRange
nonzero_col_range(Coord row, Coord col_from, Coord col_to) const {
CoordIterRange r;
GetRange(cols_, row, col_from, col_to, &r);
return r;
}
CoordIterRange
nonzero_row_range(Coord col, Coord row_from, Coord row_to) const {
CoordIterRange r;
GetRange(rows_, col, row_from, row_to, &r);
return r;
}
Coord row_count() const { return row_count_; }
Coord col_count() const { return col_count_; }
size_t nonzero_count() const { return values_.size(); }
size_t element_count() const { return size_t(row_count_) * col_count_; }
private:
typedef std::unordered_map<Point,
typename std::remove_cv<T>::type,
PointHasher> ValueMap;
struct PointHasher {
size_t operator()(const Point& p) const {
return p.first << (std::numeric_limits<Coord>::digits >> 1) ^ p.second;
}
};
static void SortAndShrink(NonZeroList& list) {
for (auto& it : list) {
auto& indices = it.second;
indices.shrink_to_fit();
std::sort(indices.begin(), indices.end());
}
// insert a sentinel vector to handle the case of all zeroes
if (list.empty())
list.emplace(Coord(), std::vector<Coord>(Coord()));
}
static void GetRange(const NonZeroList& list, Coord i, Coord from, Coord to,
CoordIterRange* r) {
auto lr = list.equal_range(i);
if (lr.first == lr.second) {
r->first = r->second = list.begin()->second.end();
return;
}
auto begin = lr.first->second.begin(), end = lr.first->second.end();
r->first = lower_bound(begin, end, from);
r->second = lower_bound(r->first, end, to);
}
ValueMap values_;
NonZeroList rows_, cols_;
Coord row_count_, col_count_;
};
#endif /* UTIL_IMMUTABLE_SPARSE_MATRIX_HPP_ */

Para simplificar, es immutable , pero puedes hacerlo mutable; asegúrese de cambiar std::vector a std::set si desea una "inserción" eficiente y razonable (cambiando un cero a un valor distinto de cero).

Answer 2

Boost tiene una implementación con plantilla de BLAS llamada uBLAS que contiene una matriz dispersa.

http://www.boost.org/doc/libs/1_36_0/libs/numeric/ublas/doc/index.htm

Answer 3

Como solo son importantes los valores con [a] [b] [c] ... [w] [x] [y] [z], solo almacenamos el índice, no el valor 1, que está en todas partes, siempre la misma + no forma de hash it. Al notar que la maldición de la dimensionalidad está presente, sugiera ir con alguna herramienta establecida NIST o Boost, al menos lea las fuentes para evitar el error innecesario.

Si el trabajo necesita capturar las distribuciones de dependencia temporal y las tendencias paramétricas de conjuntos de datos desconocidos, entonces un mapa o árbol B con raíz unificada probablemente no sea práctico. Podemos almacenar solo los índices, hash si el orden (sensibilidad para la presentación) puede subordinarse a la reducción del dominio del tiempo en el tiempo de ejecución, para todos los 1 valores. Como los valores distintos de cero distintos de uno son pocos, un candidato obvio para ellos es cualquier estructura de datos que pueda encontrar fácilmente y comprender. Si el conjunto de datos es verdaderamente vasto, del tamaño de un universo, sugiero algún tipo de ventana deslizante que administre el archivo / disco / persistente-yo usted mismo, moviendo partes de los datos en el alcance según sea necesario. (escribir un código que pueda entender) Si no está comprometido a proporcionar una solución real a un grupo de trabajo, no hacerlo lo deja a merced de los sistemas operativos de grado de consumidor que tienen el único objetivo de llevarse su almuerzo lejos de usted.

Answer 4

La lista completa de soluciones se puede encontrar en la wikipedia. Para mayor comodidad, he citado secciones relevantes de la siguiente manera.

https://en.wikipedia.org/wiki/Sparse_matrix#Dictionary_of_keys_.28DOK.29

Diccionario de llaves (DOK)

DOK consiste en un diccionario que correlaciona (fila, columna) -pairs con el valor de los elementos. Los elementos que faltan en el diccionario se toman como cero. El formato es bueno para construir incrementalmente una matriz dispersa en orden aleatorio, pero es deficiente para iterar sobre valores distintos de cero en orden lexicográfico. Uno normalmente construye una matriz en este formato y luego se convierte a otro formato más eficiente para el procesamiento. [1]

Lista de listas (LIL)

LIL almacena una lista por fila, con cada entrada que contiene el índice de la columna y el valor. Por lo general, estas entradas se mantienen ordenadas por índice de columna para una búsqueda más rápida. Este es otro formato bueno para la construcción de matriz incremental. [2]

Lista de coordenadas (COO)

COO almacena una lista de tuplas (filas, columnas, valores). Idealmente, las entradas se ordenan (por índice de fila, luego por índice de columna) para mejorar los tiempos de acceso aleatorio. Este es otro formato que es bueno para la construcción de matriz incremental. [3]

Fila escasa comprimida (formato CSR, CRS o Yale)

El formato de fila dispersa comprimida (CSR) o almacenamiento de fila comprimida (CRS) representa una matriz M por tres matrices (unidimensionales), que respectivamente contienen valores distintos de cero, las extensiones de filas y los índices de columna. Es similar a COO, pero comprime los índices de filas, de ahí el nombre. Este formato permite un rápido acceso a filas y multiplicaciones de matriz-vector (Mx).

Answer 5

La mejor manera de implementar matrices dispersas es no implementarlas, al menos no por su cuenta. Sugeriría a BLAS (que creo que es parte de LAPACK) que puede manejar matrices realmente grandes.

Answer 6

Las tablas hash tienen una inserción rápida y miran hacia arriba. Podrías escribir una función hash simple ya que sabes que tratarías solo pares enteros como las teclas.

Answer 7

Para C ++, un mapa funciona bien. Varios millones de objetos no serán un problema. 10 millones de elementos tomaron aproximadamente 4.4 segundos y aproximadamente 57 meg en mi computadora.

Mi aplicación de prueba es la siguiente:

#include <stdio.h>
#include <stdlib.h>
#include <map>
class triple {
public:
int x;
int y;
int z;
bool operator<(const triple &other) const {
if (x < other.x) return true;
if (other.x < x) return false;
if (y < other.y) return true;
if (other.y < y) return false;
return z < other.z;
}
};
int main(int, char**)
{
std::map<triple,int> data;
triple point;
int i;
for (i = 0; i < 10000000; ++i) {
point.x = rand();
point.y = rand();
point.z = rand();
//printf("%d %d %d %d/n", i, point.x, point.y, point.z);
data[point] = i;
}
return 0;
}

Ahora, para elegir dinámicamente el número de variables, la solución más fácil es representar el índice como una cadena y luego usar una cadena como clave para el mapa. Por ejemplo, un artículo ubicado en [23] [55] se puede representar mediante una cadena de "23,55". También podemos extender esta solución para dimensiones más altas; por ejemplo, para tres dimensiones, un índice arbitrario se verá como "34,45,56". Una implementación simple de esta técnica es la siguiente:

std::map data<string,int> data;
char ix[100];
sprintf(ix, "%d,%d", x, y); // 2 vars
data[ix] = i;
sprintf(ix, "%d,%d,%d", x, y, z); // 3 vars
data[ix] = i;

Answer 8

Pequeño detalle en la comparación del índice. Necesitas hacer una comparación lexicográfica, de lo contrario:

a= (1, 2, 1); b= (2, 1, 2);
(a<b) == (b<a) is true, but b!=a

Editar: Entonces la comparación debería ser probablemente:

return lhs.x<rhs.x
? true
: lhs.x==rhs.x
? lhs.y<rhs.y
? true
: lhs.y==rhs.y
? lhs.z<rhs.z
: false
: false

Answer 9

Solo como consejo: el método que usa cadenas como índices es realmente muy lento. Una solución mucho más eficiente pero por lo demás equivalente sería usar vectores / matrices. No hay necesidad de escribir los índices en una cadena.

typedef vector<size_t> index_t;
struct index_cmp_t : binary_function<index_t, index_t, bool> {
bool operator ()(index_t const& a, index_t const& b) const {
for (index_t::size_type i = 0; i < a.size(); ++i)
if (a[i] != b[i])
return a[i] < b[i];
return false;
}
};
map<index_t, int, index_cmp_t> data;
index_t i(dims);
i[0] = 1;
i[1] = 2;
// … etc.
data[i] = 42;

Sin embargo, usar un map no es realmente muy eficiente debido a la implementación en términos de un árbol de búsqueda binario equilibrado. En este caso, las estructuras de datos de mucho mejor rendimiento serían una tabla hash (aleatoria).

Answer 10

Sugeriría hacer algo como:

typedef std::tuple<int, int, int> coord_t;
typedef boost::hash<coord_t> coord_hash_t;
typedef std::unordered_map<coord_hash_t, int, c_hash_t> sparse_array_t;
sparse_array_t the_data;
the_data[ { x, y, z } ] = 1; /* list-initialization is cool */
for( const auto& element : the_data ) {
int xx, yy, zz, val;
std::tie( std::tie( xx, yy, zz ), val ) = element;
/* ... */
}

Para ayudar a mantener sus datos dispersos, es posible que desee escribir una subclase de unorderd_map , cuyos iteradores salten automáticamente (y borren) cualquier elemento con un valor de 0.

Answer 11

Eigen es una biblioteca de álgebra lineal de C ++ que tiene una implementation de una matriz dispersa. Incluso admite operaciones de matriz y solucionadores (factorización de LU, etc.) que están optimizados para matrices dispersas.