diff --git a/c/glm.c b/c/glm.c index 861dcf0..b3451c0 100644 --- a/c/glm.c +++ b/c/glm.c @@ -614,18 +614,30 @@ static inline int32_t dot_i8i8(const int8_t *w, const int8_t *x, int I){ #elif defined(__ARM_NEON) /* ARM: SDOT nativo se disponibile (Apple Silicon: sempre); altrimenti vmull/vpadal. * Stesso bound anti-overflow del trucco AVX2: coppie <= 128*127*2 = 32512 < 32767. */ +#if defined(__ARM_FEATURE_DOTPROD) + /* 4 accumulatori indipendenti: SDOT ha latenza ~3-4 cicli, con un solo acc la + * catena seriale strozza il core a ~26 GB/s di pesi; con 4 lane indipendenti il + * dot diventa memory-bound (misurato su M4: 26 -> 63 GB/s per core, 2.4x). */ + int32x4_t a0=vdupq_n_s32(0),a1=vdupq_n_s32(0),a2=vdupq_n_s32(0),a3=vdupq_n_s32(0); + for(;i+64<=I;i+=64){ + a0=vdotq_s32(a0,vld1q_s8(w+i), vld1q_s8(x+i)); + a1=vdotq_s32(a1,vld1q_s8(w+i+16),vld1q_s8(x+i+16)); + a2=vdotq_s32(a2,vld1q_s8(w+i+32),vld1q_s8(x+i+32)); + a3=vdotq_s32(a3,vld1q_s8(w+i+48),vld1q_s8(x+i+48)); + } + int32x4_t acc=vaddq_s32(vaddq_s32(a0,a1),vaddq_s32(a2,a3)); + for(;i+16<=I;i+=16) acc=vdotq_s32(acc,vld1q_s8(w+i),vld1q_s8(x+i)); + sum=vaddvq_s32(acc); +#else int32x4_t acc=vdupq_n_s32(0); for(;i+16<=I;i+=16){ int8x16_t wv=vld1q_s8(w+i), xv=vld1q_s8(x+i); -#if defined(__ARM_FEATURE_DOTPROD) - acc=vdotq_s32(acc,wv,xv); -#else int16x8_t p=vmull_s8(vget_low_s8(wv),vget_low_s8(xv)); p=vmlal_s8(p,vget_high_s8(wv),vget_high_s8(xv)); acc=vpadalq_s16(acc,p); -#endif } sum=vaddvq_s32(acc); +#endif #elif defined(__VSX__) /* POWER8: vec_msum (s8 x u8 -> s32) somma i prodotti byte DIRETTAMENTE in lane * s32, 16 byte/iter: il bound anti-saturazione a 16 bit di maddubs qui non serve. @@ -705,6 +717,28 @@ static inline int32_t dot_i4i8(const uint8_t *w4, const int8_t *x, int I){ sum=hsum256_i32(acc); #elif defined(__ARM_NEON) const uint8x16_t m4q=vdupq_n_u8(0x0F); const int8x16_t b8q=vdupq_n_s8(8); +#if defined(__ARM_FEATURE_DOTPROD) + /* 4 accumulatori indipendenti (vedi dot_i8i8): spezza la catena seriale su acc. + * Misurato su M4: 12.4 -> 29.9 GB/s di pesi per core (2.4x). */ + int32x4_t a0=vdupq_n_s32(0),a1=vdupq_n_s32(0),a2=vdupq_n_s32(0),a3=vdupq_n_s32(0); + for(;i+64<=I;i+=64){ + uint8x16_t byA=vld1q_u8(w4+(i>>1)), byB=vld1q_u8(w4+(i>>1)+16); + uint8x16x2_t zA=vzipq_u8(vandq_u8(byA,m4q), vshrq_n_u8(byA,4)); /* nibble in ordine */ + uint8x16x2_t zB=vzipq_u8(vandq_u8(byB,m4q), vshrq_n_u8(byB,4)); + a0=vdotq_s32(a0,vsubq_s8(vreinterpretq_s8_u8(zA.val[0]),b8q),vld1q_s8(x+i)); + a1=vdotq_s32(a1,vsubq_s8(vreinterpretq_s8_u8(zA.val[1]),b8q),vld1q_s8(x+i+16)); + a2=vdotq_s32(a2,vsubq_s8(vreinterpretq_s8_u8(zB.val[0]),b8q),vld1q_s8(x+i+32)); + a3=vdotq_s32(a3,vsubq_s8(vreinterpretq_s8_u8(zB.val[1]),b8q),vld1q_s8(x+i+48)); + } + int32x4_t acc=vaddq_s32(vaddq_s32(a0,a1),vaddq_s32(a2,a3)); + for(;i+32<=I;i+=32){ + uint8x16_t by=vld1q_u8(w4+(i>>1)); /* 16 byte = 32 nibble */ + uint8x16x2_t z=vzipq_u8(vandq_u8(by,m4q), vshrq_n_u8(by,4)); /* nibble in ordine */ + acc=vdotq_s32(acc,vsubq_s8(vreinterpretq_s8_u8(z.val[0]),b8q),vld1q_s8(x+i)); + acc=vdotq_s32(acc,vsubq_s8(vreinterpretq_s8_u8(z.val[1]),b8q),vld1q_s8(x+i+16)); + } + sum=vaddvq_s32(acc); +#else int32x4_t acc=vdupq_n_s32(0); for(;i+32<=I;i+=32){ uint8x16_t by=vld1q_u8(w4+(i>>1)); /* 16 byte = 32 nibble */ @@ -712,18 +746,15 @@ static inline int32_t dot_i4i8(const uint8_t *w4, const int8_t *x, int I){ int8x16_t w0=vsubq_s8(vreinterpretq_s8_u8(z.val[0]),b8q); int8x16_t w1=vsubq_s8(vreinterpretq_s8_u8(z.val[1]),b8q); int8x16_t x0=vld1q_s8(x+i), x1=vld1q_s8(x+i+16); -#if defined(__ARM_FEATURE_DOTPROD) - acc=vdotq_s32(acc,w0,x0); acc=vdotq_s32(acc,w1,x1); -#else int16x8_t p=vmull_s8(vget_low_s8(w0),vget_low_s8(x0)); /* |w|<=8: nessun overflow */ p=vmlal_s8(p,vget_high_s8(w0),vget_high_s8(x0)); acc=vpadalq_s16(acc,p); p=vmull_s8(vget_low_s8(w1),vget_low_s8(x1)); p=vmlal_s8(p,vget_high_s8(w1),vget_high_s8(x1)); acc=vpadalq_s16(acc,p); -#endif } sum=vaddvq_s32(acc); +#endif #elif defined(__VSX__) /* 16 byte = 32 nibble. vec_mergeh/vec_mergel su ppc64le (GCC) interallacciano come * unpacklo/unpackhi x86 (verificato empiricamente su POWER8): i nibble escono in