性能分析表明，这个函数是我应用程序的瓶颈所在：
static inline int countEqualChars(const char* string1, const char* string2, int size) {
int r = 0;
for (int j = 0; j < size; ++j) {
if (string1[j] == string2[j]) {
++r;
}
}
return r;
}
即使使用-O3和-march=native，G++ 4.7.2仍然无法将此函数向量化(我检查了汇编器输出)。我不是SSE和其它相关技术的专家，但我认为一次比较多个字符应该更快。有任何加速的想法吗？目标架构是x86-64。
pcmpeqb 比较两个 16 字节向量，并生成一个向量，在它们不同的地方为零，在它们匹配的地方为 -1。使用此方法一次比较 16 个字节，并将结果添加到累加器向量中（确保最多积累 255 个向量比较的结果以避免溢出）。完成后，累加器中有 16 个结果。对它们求和并取反即可获得相等元素的数量。向量化的编译器标志：
-ftree-vectorize
-ftree-vectorize -march=<your_architecture>（使用计算机上所有可用的指令集扩展，而不仅仅是基线，如x86-64的SSE2）。使用-march=native来优化编译器正在运行的机器。）-march=<foo>还设置了-mtune=<foo>，这也是好的。
使用SSEx内部函数：
使用16字节（根据您实际要使用的向量大小）填充和对齐缓冲区
使用_mm_set1_epi8(0)创建一个累加器countU8
对于所有n/16个输入（子）向量，执行以下操作：
使用_mm_load_si128或_mm_loadu_si128（对于非对齐加载）从两个字符串中加载16个字符
_mm_cmpeq_epi8并行比较八位字节。每个匹配都产生0xFF（-1），否则为0x00。
使用_mm_sub_epi8将上述结果向量从countU8减去（减去-1->+1）
在255个周期后，必须将16个8位计数器提取为较大的整数类型，以防止溢出。请参见此良好答案中的展开和水平添加方式：https://dev59.com/FGTWa4cB1Zd3GeqPF7Sn#10930706
代码:
#include <iostream>
#include <vector>
#include <cassert>
#include <cstdint>
#include <climits>
#include <cstring>
#include <emmintrin.h>
#ifdef __SSE2__
#if !defined(UINTPTR_MAX) || !defined(UINT64_MAX) || !defined(UINT32_MAX)
# error "Limit macros are not defined"
#endif
#if UINTPTR_MAX == UINT64_MAX
#define PTR_64
#elif UINTPTR_MAX == UINT32_MAX
#define PTR_32
#else
# error "Current UINTPTR_MAX is not supported"
#endif
template<typename T>
void print_vector(std::ostream& out,const __m128i& vec)
{
static_assert(sizeof(vec) % sizeof(T) == 0,"Invalid element size");
std::cout << '{';
const T* const end = reinterpret_cast<const T*>(&vec)-1;
const T* const upper = end+(sizeof(vec)/sizeof(T));
for(const T* elem = upper;
elem != end;
--elem
)
{
if(elem != upper)
std::cout << ',';
std::cout << +(*elem);
}
std::cout << '}' << std::endl;
}
#define PRINT_VECTOR(_TYPE,_VEC) do{ std::cout << #_VEC << " : "; print_vector<_TYPE>(std::cout,_VEC); } while(0)
///@note SSE2 required (macro: __SSE2__)
///@warning Not tested!
size_t counteq_epi8(const __m128i* a_in,const __m128i* b_in,size_t count)
{
assert(a_in != nullptr && (uintptr_t(a_in) % 16) == 0);
assert(b_in != nullptr && (uintptr_t(b_in) % 16) == 0);
//assert(count > 0);
/*
//maybe not so good with all that branching and additional loop variables
__m128i accumulatorU8 = _mm_set1_epi8(0);
__m128i sum2xU64 = _mm_set1_epi8(0);
for(size_t i = 0;i < count;++i)
{
//this operation could also be unrolled, where multiple result registers would be accumulated
accumulatorU8 = _mm_sub_epi8(accumulatorU8,_mm_cmpeq_epi8(*a_in++,*b_in++));
if(i % 255 == 0)
{
//before overflow of uint8, the counter will be extracted
__m128i sum2xU16 = _mm_sad_epu8(accumulatorU8,_mm_set1_epi8(0));
sum2xU64 = _mm_add_epi64(sum2xU64,sum2xU16);
//reset accumulatorU8
accumulatorU8 = _mm_set1_epi8(0);
}
}
//blindly accumulate remaining values
__m128i sum2xU16 = _mm_sad_epu8(accumulatorU8,_mm_set1_epi8(0));
sum2xU64 = _mm_add_epi64(sum2xU64,sum2xU16);
//do a horizontal addition of the two counter values
sum2xU64 = _mm_add_epi64(sum2xU64,_mm_srli_si128(sum2xU64,64/8));
#if defined PTR_64
return _mm_cvtsi128_si64(sum2xU64);
#elif defined PTR_32
return _mm_cvtsi128_si32(sum2xU64);
#else
# error "macro PTR_(32|64) is not set"
#endif
*/
__m128i sum2xU64 = _mm_set1_epi32(0);
while(count--)
{
__m128i matches = _mm_sub_epi8(_mm_set1_epi32(0),_mm_cmpeq_epi8(*a_in++,*b_in++));
__m128i sum2xU16 = _mm_sad_epu8(matches,_mm_set1_epi32(0));
sum2xU64 = _mm_add_epi64(sum2xU64,sum2xU16);
#ifndef NDEBUG
PRINT_VECTOR(uint16_t,sum2xU64);
#endif
}
//do a horizontal addition of the two counter values
sum2xU64 = _mm_add_epi64(sum2xU64,_mm_srli_si128(sum2xU64,64/8));
#ifndef NDEBUG
std::cout << "----------------------------------------" << std::endl;
PRINT_VECTOR(uint16_t,sum2xU64);
#endif
#if !defined(UINTPTR_MAX) || !defined(UINT64_MAX) || !defined(UINT32_MAX)
# error "Limit macros are not defined"
#endif
#if defined PTR_64
return _mm_cvtsi128_si64(sum2xU64);
#elif defined PTR_32
return _mm_cvtsi128_si32(sum2xU64);
#else
# error "macro PTR_(32|64) is not set"
#endif
}
#endif
int main(int argc, char* argv[])
{
std::vector<__m128i> a(64); // * 16 bytes
std::vector<__m128i> b(a.size());
const size_t nBytes = a.size() * sizeof(std::vector<__m128i>::value_type);
char* const a_out = reinterpret_cast<char*>(a.data());
char* const b_out = reinterpret_cast<char*>(b.data());
memset(a_out,0,nBytes);
memset(b_out,0,nBytes);
a_out[1023] = 1;
b_out[1023] = 1;
size_t equalBytes = counteq_epi8(a.data(),b.data(),a.size());
std::cout << "equalBytes = " << equalBytes << std::endl;
return 0;
}
我最快的针对大型和小型数组的SSE实现：
size_t counteq_epi8(const __m128i* a_in,const __m128i* b_in,size_t count)
{
assert((count > 0 ? a_in != nullptr : true) && (uintptr_t(a_in) % sizeof(__m128i)) == 0);
assert((count > 0 ? b_in != nullptr : true) && (uintptr_t(b_in) % sizeof(__m128i)) == 0);
//assert(count > 0);
const size_t maxInnerLoops = 255;
const size_t nNestedLoops = count / maxInnerLoops;
const size_t nRemainderLoops = count % maxInnerLoops;
const __m128i zero = _mm_setzero_si128();
__m128i sum16xU8 = zero;
__m128i sum2xU64 = zero;
for(size_t i = 0;i < nNestedLoops;++i)
{
for(size_t j = 0;j < maxInnerLoops;++j)
{
sum16xU8 = _mm_sub_epi8(sum16xU8,_mm_cmpeq_epi8(*a_in++,*b_in++));
}
sum2xU64 = _mm_add_epi64(sum2xU64,_mm_sad_epu8(sum16xU8,zero));
sum16xU8 = zero;
}
for(size_t j = 0;j < nRemainderLoops;++j)
{
sum16xU8 = _mm_sub_epi8(sum16xU8,_mm_cmpeq_epi8(*a_in++,*b_in++));
}
sum2xU64 = _mm_add_epi64(sum2xU64,_mm_sad_epu8(sum16xU8,zero));
sum2xU64 = _mm_add_epi64(sum2xU64,_mm_srli_si128(sum2xU64,64/8));
#if UINTPTR_MAX == UINT64_MAX
return _mm_cvtsi128_si64(sum2xU64);
#elif UINTPTR_MAX == UINT32_MAX
return _mm_cvtsi128_si32(sum2xU64);
#else
# error "macro PTR_(32|64) is not set"
#endif
}
psadbw 进行更高效的水平加法，接着移动高到低的64位再进行加法运算。 - Stephen Canonassert（equalBytes == nBytes）; 与您的实现失败（似乎缺少剩余的“size％something”循环）。 关于内存加载：我们可以尝试对某些缓存行进行_mm_prefetch() ，但我认为这只对较旧的Pentium（II / III）处理器产生相当大的影响。 - Sam当前gcc的自动向量化是帮助编译器理解代码易于向量化的问题。对于您的情况：如果您去除条件并以更加命令式的方式重写代码，它将理解向量化请求：
static inline int count(const char* string1, const char* string2, int size) {
int r = 0;
bool b;
for (int j = 0; j < size; ++j) {
b = (string1[j] == string2[j]);
r += b;
}
return r;
}
movdqa 16(%rsp), %xmm1
movl $.LC2, %esi
pxor %xmm2, %xmm2
movzbl 416(%rsp), %edx
movdqa .LC1(%rip), %xmm3
pcmpeqb 224(%rsp), %xmm1
cmpb %dl, 208(%rsp)
movzbl 417(%rsp), %eax
movl $1, %edi
pand %xmm3, %xmm1
movdqa %xmm1, %xmm5
sete %dl
movdqa %xmm1, %xmm4
movzbl %dl, %edx
punpcklbw %xmm2, %xmm5
punpckhbw %xmm2, %xmm4
pxor %xmm1, %xmm1
movdqa %xmm5, %xmm6
movdqa %xmm5, %xmm0
movdqa %xmm4, %xmm5
punpcklwd %xmm1, %xmm6
(etc.)
#include <omp.h>
#include <stdalign.h>
#include <stddef.h>
#include <stdint.h>
#define ALIGNMENT 16U
size_t countEqualBytes(const size_t n, const uint8_t a[n], const uint8_t b[n]) {
size_t sum = 0;
size_t i = 0;
if (n >= 32U) {
const size_t sentinel = n - 31U;
// #pragma omp parallel for reduction(+:sum) schedule(static)
for (i = 0; i < sentinel; i += 32U) {
sum += (size_t)((a[i] == b[i]) +
(a[i + 1] == b[i + 1]) +
(a[i + 2] == b[i + 2]) +
(a[i + 3] == b[i + 3]) +
(a[i + 4] == b[i + 4]) +
(a[i + 5] == b[i + 5]) +
(a[i + 6] == b[i + 6]) +
(a[i + 7] == b[i + 7]) +
(a[i + 8] == b[i + 8]) +
(a[i + 9] == b[i + 9]) +
(a[i + 10] == b[i + 10]) +
(a[i + 11] == b[i + 11]) +
(a[i + 12] == b[i + 12]) +
(a[i + 13] == b[i + 13]) +
(a[i + 14] == b[i + 14]) +
(a[i + 15] == b[i + 15]) +
(a[i + 16] == b[i + 16]) +
(a[i + 17] == b[i + 17]) +
(a[i + 18] == b[i + 18]) +
(a[i + 19] == b[i + 19]) +
(a[i + 20] == b[i + 20]) +
(a[i + 21] == b[i + 21]) +
(a[i + 22] == b[i + 22]) +
(a[i + 23] == b[i + 23]) +
(a[i + 24] == b[i + 24]) +
(a[i + 25] == b[i + 25]) +
(a[i + 26] == b[i + 26]) +
(a[i + 27] == b[i + 27]) +
(a[i + 28] == b[i + 28]) +
(a[i + 29] == b[i + 29]) +
(a[i + 30] == b[i + 30]) +
(a[i + 31] == b[i + 31]));
}
}
for (; i<n; i++) {
sum += (a[i] != b[i]);
}
return sum;
}
-std=c17 -O3 -march=x86-64-v4选项能够将此关键循环编译为：.LBB0_5: # =>This Inner Loop Header: Depth=1
vmovdqu ymm0, ymmword ptr [rsi + r10]
vmovdqu ymm1, ymmword ptr [rsi + r10 + 32]
vmovdqu ymm2, ymmword ptr [rsi + r10 + 64]
vmovdqu ymm3, ymmword ptr [rsi + r10 + 96]
vpcmpeqb k0, ymm0, ymmword ptr [rdx + r10]
kmovd ebx, k0
popcnt ebx, ebx
add rbx, rax
vpcmpeqb k0, ymm1, ymmword ptr [rdx + r10 + 32]
kmovd eax, k0
popcnt eax, eax
add rax, rbx
vpcmpeqb k0, ymm2, ymmword ptr [rdx + r10 + 64]
kmovd ebx, k0
popcnt ebx, ebx
add rbx, rax
vpcmpeqb k0, ymm3, ymmword ptr [rdx + r10 + 96]
kmovd eax, k0
popcnt eax, eax
add rax, rbx
sub r10, -128
add r9, -4
jne .LBB0_5
.LBB0_8: # =>This Inner Loop Header: Depth=1
vmovdqu ymm0, ymmword ptr [r10 + rbx]
vpcmpeqb k0, ymm0, ymmword ptr [r9 + rbx]
kmovd ecx, k0
popcnt ecx, ecx
add rax, rcx
add rbx, 32
cmp r8, rbx
jne .LBB0_8
-march=x86-64-v4编译来生成AVX-512代码。另外，如果你想避免缓存行分割，你需要将数据按照寄存器宽度对齐。如果你使用YMM或ZMM寄存器，那么对齐宽度就是32或64；除非你使用XMM寄存器（例如SSE或128位AVX指令），否则16字节对齐没有什么特殊之处。 - Peter Cordesvpsubb），而不是5个（加载/加载+比较/移位掩码/人口计数/加法）或者没有微融合的6个。要让编译器生成这样的良好汇编代码，仍然需要使用内部函数（请参见其他答案）。 - Peter Cordesshuffles)，最终得到vpsadbw指令，所以比这个方法还要慢。你的回答仅展示了AVX-512，并忽略了一个Godbolt链接(https://godbolt.org/z/3deoqMPbj)，以便人们可以更快地查看其他选项。是的，在AVX2上它使用了`vpmovmskb`和`popcnt`，但没有展开你的循环。 - Peter Cordes
-march=native，这意味着他的CPU支持的任何-mfoo标志。 - user784668