C中的二元向量和矩阵操作

Question

我试图在C中实现一个数据结构,这将允许我有效地操作**二进制**矩阵(仅包含1或0).我将解释我必须对此矩阵应用哪些操作,并想知道使用哪种最佳数据结构？

操作在字段F_2中完成(这意味着1 + 1 = 0,其他操作保持不变).我有一个k*n矩阵(k< n)调用H.最多k= 2325和n= 3009.

我必须对此矩阵执行的操作是:

我将仅使用行交换和行添加来部分对角化它.一旦完成,我将不再使用行操作,并将在此矩阵上运行大量(!)列添加(我的意思是"很多"是关于((nk)/ 2)³列添加)

我正在考虑矩阵的数据结构:

对于矩阵系数,我考虑在一个单个unsigned int中一次存储多个位的序列.例如,我可以将序列存储(11001011)到uint8_t 203(从二进制转换为十进制)

• 这是个好主意吗 ？

如果我这样做,我有两个选择:

我可以使用uint16_t或uint64_t系数在许多4*4或8*8子矩阵中分割我的矩阵H.

• 这是一个很好的选择(在时间效率方面),如果是,是否更好地使用uint16_t或uint64_t？

另外我想存储每一行中的多个uint32_t或uint64_t,然后操作我的部分对角化.接下来切换到将矩阵编码为n列向量以处理剩余操作的结构.

• 你认为这更有效吗？

无论我使用什么方法,我都必须有效地访问nunsigned int(uint16,32或64)的第一位.我怎么做 ？

Answer 1

为获得最佳性能,请使用行指针数组进行行交换和行添加.使用<stdint.h>和最小支持字长的快速无符号整数类型 - 我建议使用uint_fast32_t,除非您打算在16位或8位处理器上运行它.

完成所有行交换和行添加后,转置数组.虽然此操作"慢",但以下列操作将如此快以抵消转置成本.

考虑以下:

#include <stdint.h>
#include <limits.h>

typedef uint_fast32_t  word;
#define WORD_BITS      (CHAR_BIT * sizeof (word))

typedef struct {
    int    rows;  /* Number of row vectors */
    int    cols;  /* Number of defined bits in each row */
    int    words; /* Number of words per row vector */
    word **row;   /* Array of pointers */
} row_matrix;

typedef struct {
    int    rows;  /* Number of defined bits in each column */
    int    cols;  /* Number of col vectors */
    int    words; /* Number of words per col vector */
    word **col;
} col_matrix;

虽然您可以使用单一类型来描述两种矩阵形式,但使用单独的类型可以使代码和函数更易于维护.你最终会得到一些重复的代码,但与清晰直观的类型相比,这是一个小问题.

在32位系统上,uint_fast32_t通常是32位类型.在64位系统上,它通常是64位.该WORD_BITS宏展开为一个比特的数量word-它并不总是32!

对位进行编号的最简单方法是将矩阵中最左边的位指定为位0,并将这些位存储在每个字的最低有效位中.如果你有row_matrix *rm,那么行row,列的位col是

!!(rm->row[row][col / WORD_BITS] & ((word)1U << (col % WORD_BITS)))

这!!是not-not运算符:如果参数非零,则产生1,否则它产生0.因为我们从字中屏蔽了一个位,所以"bit is set"值将是2的幂(1, 2,4,8,16,32,64等).

要设置该位,请使用

rm->row[row][col / WORD_BITS] |= (word)1U << (col % WORD_BITS);

要清除一点,您需要使用除目标位1之外的所有掩码进行二进制AND运算.使用not运算符很容易实现~:

rm->row[row][col / WORD_BITS] &= ~((word)1U << (col % WORD_BITS));

相应的操作col_matrix *cm是

!!(cm->col[col][row / WORD_BITS] & ((word)1U << (row % WORD_BITS)))
cm->col[col][row / WORD_BITS] |= (word)1U << (row % WORD_BITS);
cm->col[col][row / WORD_BITS] &= ~((word)1U << (row % WORD_BITS));

尽管除法/和模数(或余数)%通常是慢速的(与加法,减法甚至乘法相比),但WORD_BITS在所有广泛使用的架构中,这将是两个编译时常量的幂.我所知道的所有编译器都会将上述内容转换为快速位移和二进制AND运算符.

要添加行srcrow到行dstrow,您只需对所有单词执行二进制异或:

{
    const word *const src = rm->row[srcrow];
    word *const       dst = rm->row[dstrow];
    int               w = rm->words;
    while (w-->0)
        dst[w] ^= src[w];
}

类似地,对于列矩阵,

{
    const word *const src = cm->col[srccol];
    word *const       dst = cm->col[dstcol];
    int               w = cm->words;
    while (w-->0)
        dst[w] ^= src[w];
}

请注意,如果组合两行以上,则可以非常有效地执行此操作; 它会比连续添加更快.Intel和AMD CPU非常擅长预测上述模式,因此您可以使用多个源行/列.此外,目的地不必参与结果,但如果我猜对了你正在实施的算法,我猜你想要它.

如果您知道目标体系结构具有SSE2或更高版本,甚至AVX,则可以分别使用emmintrin.h或者immintrin.h头文件用于编译器内置类型和运算符,这些类型和运算符允许您一次分别对XOR 128位和256位进行异或; 有时会给你很大的提升.

由于向量类型需要C标准称之为"过度对齐",因此您还需要包含mm_malloc.h,使用_mm_malloc()和_mm_free()分配数据的行/列向量word- 显然是words向上舍入以便您可以访问行/列作为合适的整数字类型(__m128i对于SSE*,__m256i对于AVX).

就个人而言,我总是首先实现非版本化版本,然后是单元测试的一些"讨厌"测试用例,然后才能看到它的矢量化.这样做的好处是,您可以将非版本化版本作为初步版本提供给那些将使用它的人,并且您可以比较矢量化和非矢量化案例之间的测试用例结果,以查看其中一个或另一个是否有错误.

转置操作非常简单,虽然我建议使用三重循环:最内层循环一个字中的位.此外,您可能想要检查哪个顺序 - 行或列主要 - 最适合外循环; 根据矩阵大小,您可能会看到巨大的差异.(这是由于缓存行为:您希望CPU能够预测访问模式,而不必重新加载相同的缓存行.在最好的情况下,在最近几年的AMD和Intel x86- 64个处理器,如果两个矩阵都适合缓存,则可以接近缓存速度.)

所有上述内容都可以在单个头文件中实现 - 如果目标体系结构支持SSE2/AVX,甚至包括矢量化版本 - 因此实现起来应该不会太难.

有问题吗？