在C ++中创建稀疏数组的最佳方式是什么?

我正在研究一个需要处理巨大矩阵的项目,特别是用于copula计算的金字塔加法。

总之,我需要在矩阵(多维数组)中的零海中追踪相对较少的值(通常值为1,在极少数情况下超过1)。

稀疏数组允许用户存储少量值,并假定所有未定义的记录都是预设值。由于在物理上不可能将所有值存储在内存中,所以我只需要存储少数几个非零元素。这可能是数百万条目。

速度是一个重要的优先事项,我也想在运行时动态选择类中的变量数量。

我目前正在使用二叉搜索树(b-tree)存储条目的系统。有谁知道更好的系统?

0

11 答案

哈希表具有快速插入和查找功能。你可以编写一个简单的哈希函数,因为你知道你只需要处理整数对作为键。

0
额外

对于C ++,地图运行良好。数百万个物体不会成为问题。我的电脑上有1000万件物品需要4.4秒左右,大约57兆。

我的测试应用程序如下所示:

#include 
#include 
#include  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; } 

现在要动态选择变量的数量,最简单的解决方案是将索引表示为字符串,然后将字符串用作映射关键字。例如,位于[23] [55]的项目可以通过“23,55”字符串表示。我们也可以扩展这个解决方案以获得更高的维例如对于三维,任意索引将看起来像“34,45,56”。这种技术的简单实现如下:

std::map data 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;
0
额外
由于它没有被固定很长时间,我冒昧地以正确的实施取而代之。
额外 作者 celtschk,
从这个获取元素范围的性能或检查范围是否完全在数组中?
额外 作者 Ivan G.,
运算符<�的实现不正确。考虑Triple {1,2,3}和Triple {3,2,1},两者都不会少于另一个。一个正确的实现将检查x然后y然后z顺序而不是一次全部。
额外 作者 Whanhee,
只有一百万个元素的秒?这似乎很糟糕。你应该考虑使用哈希函数和 unordered_map 。查看 github.com/victorprad/InfiniTAM 他们使用散列((x * 73856093u) ^(y * 19349669u)^(z * 83492791u)),并且可以将数百万个样本以良好的帧率整合到稀疏的3D网格中。
额外 作者 masterxilo,

实现稀疏矩阵的最好方法是不要实现它们 - 至少不是你自己的。我会建议BLAS(我认为它是LAPACK的一部分),它可以处理真正巨大的矩阵。

0
额外
LAPACK是一个密集矩阵库。标准BLAS也适用于密集矩阵。有一个稀疏的BLAS包(通过NIST),但这与标准BLAS不同。
额外 作者 codehippo,

这是一个相对简单的实现,应该提供一个合理的快速查找(使用散列表)以及对行/列中的非零元素进行快速迭代。

// 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 
#include 
#include  #include  #include  #include  #include  // 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 class SparseMatrix { struct PointHasher; typedef std::map< Coord, std::vector > NonZeroList; typedef std::pair Point; public: typedef T ValueType; typedef Coord CoordType; typedef typename NonZeroList::mapped_type::const_iterator CoordIter; typedef std::pair CoordIterRange; SparseMatrix() = default; // Reads a matrix stored in MatrixMarket-like format, i.e.: //    //    // ... // Note: the header (lines starting with '%' are ignored). template 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::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::type, PointHasher> ValueMap; struct PointHasher { size_t operator()(const Point& p) const { return p.first << (std::numeric_limits::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())); } 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_ */ 

为简单起见,它是 immutable ,但您可以使其变为可变的;如果你想要一个合理高效的“插入”(将零改为非零),务必将 std :: vector 更改为 std :: set

0
额外

完整的解决方案列表可以在维基百科中找到。为方便起见,我引用了以下相关部分。

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

密钥字典(DOK)

DOK由一个字典映射(行,列) - 对 -   元素的价值。字典中缺少的元素   被认为是零。这种格式对增量而言是很好的   以随机顺序构造稀疏矩阵,但迭代不好   按照字典顺序排列非零值。一个典型   以这种格式构造一个矩阵,然后再转换为另一个矩阵   高效的处理格式。[1]

列表清单(LIL)

LIL每行存储一个列表,每个条目包含列   指数和价值。通常情况下,这些条目将按照   列索引更快查找。这是另一种适合的格式   增量矩阵构造。[2]

协调清单(COO)

COO存储(行,列,值)元组列表。理想情况下,条目   (按行索引,然后是列索引)进行排序以改进随机访问   倍。这是另一种适用于增量矩阵的格式   结构。[3]

压缩稀疏行(CSR,CRS或耶鲁格式)

压缩稀疏行(CSR)或压缩行存储(CRS)格式   通过三个(一维)阵列表示矩阵M,即   分别包含非零值,行的范围和列   指数。这与COO类似,但是因此压缩了行索引   名字。这种格式允许快速行访问和矩阵向量   乘法(Mx)。
0
额外

Eigen is a C++ linear algebra library that has an implementation of a sparse matrix. It even supports matrix operations and solvers (LU factorization etc) that are optimized for sparse matrices.

0
额外

就像一个建议:使用字符串作为索引的方法实际上很慢。一个更有效率,但其他等效的解决方案将是使用矢量/数组。绝对不需要在字符串中编写索引。

typedef vector index_t;

struct index_cmp_t : binary_function {
    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 data;
index_t i(dims);
i[0] = 1;
i[1] = 2;
// ? etc.
data[i] = 42;

但是,使用 map 实际上并不是非常高效,因为使用平衡二叉搜索树进行实现。在这种情况下,更好的执行数据结构将是(随机)哈希表。

0
额外
@Andrew不完全。这个答案早于C ++ 11,因此早于 std :: unordered_map 。现在我推荐使用后者,但它不是我的答案的替代品,因为 std :: vector 似乎不能为 std :: hash提供合适的实现</代码>。也就是说,除非索引实际上是可变大小的(这不太可能),否则一个 std :: unordered_map <�?>> 实际上应该是不可用的(尽管接口绝对可以改进)。
额外 作者 Konrad Rudolph,
又名。 unordered_map
额外 作者 Andrew,

Boost有一个名为uBLAS的BLAS模板实现,它包含一个稀疏矩阵。

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

0
额外

因为只有[a] [b] [c] ... [w] [x] [y] [z]的值是重要的,所以我们只存储指数本身,而不是存储在任何地方的值1 - 总是相同+无法散列它。注意到维度的诅咒是存在的,建议使用一些既定的工具NIST或Boost,至少读取来源以避免不必要的失误。

如果工作需要捕获未知数据集的时间依赖性分布和参数趋势,那么具有单值根的Map或B-Tree可能不太实际。对于所有的1个值,我们只能自己存储指令,如果在运行时对订单(对于表示的敏感性)可以从属于时间域的缩减则进行散列。由于除了一个以外的非零值很少,因此显而易见的候选对象就是您可以轻松掌握和理解的任何数据结构。如果数据集真的是庞大的宇宙大小,我建议你自己管理文件/磁盘/永久性io的某种滑动窗口,根据需要将部分数据移动到范围内。 (编写可以理解的代码)如果您承诺为工作组提供实际的解决方案,如果不这样做,您将受到消费级操作系统的摆布,这些操作系统的唯一目标是让您的午餐远离您。

0
额外

索引比较中的小细节。你需要做一个词典比较,否则:

a= (1, 2, 1); b= (2, 1, 2);
(a

编辑:所以比较应该可能是:

return lhs.x
0
额外

我会建议做一些事情:

typedef std::tuple coord_t;
typedef boost::hash coord_hash_t;
typedef std::unordered_map 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;
    /* ... */
}

为了帮助保持数据稀疏,您可能需要编写 unorderd_map 的子类,其迭代器会自动跳过(并擦除)值为0的任何项目。

0
额外