integer_adm: adm_cm_avx2/avx512 not bit-exact with the C reference for scale-0 coefficients >= 15360; checkasm's ±8000 inputs never reach it

[lusoris]

Summary

adm_cm_avx2 (scale 0) is not bit-exact with the C reference adm_cm for large scale-0 coefficients. adm_cm_avx512 contains the first of the two constructs below (its cube shift is already arithmetic, _mm512_srai_epi64); on an AVX-512 host the widened checkasm run below fails for adm_cm_*_avx512 as well. checkasm does not catch it because check_adm.c fills the bands with values in ±8000 and the divergence starts at |coefficient| ≥ 15360. With the same test widened to ±16000, adm_cm_*_avx2 fails on 8 of 10 seeds with the stock 1e-4 tolerance (re-run on 6ec23e8: 10 of 10 seeds, and 60 of 60 for seeds 1 to 60), so it is the input range, not the tolerance, that hides it.

End-to-end effect on a synthetic 576x324 clip with hard 0/255 stripes (--feature adm, AVX2 vs --cpumask 4294967295), printed at full precision:

metric	C	AVX2
integer_adm_scale0 (frame 2)	0.93857384037524272	0.93857357475808256
integer_adm2 (frame 2)	0.93674829731911791	0.93674826962705859
integer_aim (frame 0)	0.13268025833079128	0.13267989648046769
integer_adm3 (frame 0)	0.90374926045970283	0.90374944138486468

Scales 1–3 (i4_adm_cm) are identical; the int64 accumulators of adm_cm differ. The deltas round away at the CLI's default 6 decimals, which is presumably why nobody noticed.

Root cause: two independent divergences in x86/adm_avx2.c

1. The threshold's centre term. The C reference (integer_adm.c, ADM_CM_THRESH_S_I_J and the other eight scale-0 threshold macros) computes

sum += (int16_t)(((ONE_BY_15 * abs((int32_t) src_ptr[j])) + 2048) >> 12);

With ONE_BY_15 = 8738 the value reaches 32768 at |src| = 15360 and the (int16_t) cast wraps it negative. ADM_CM_THRESH_S_I_J_avx256 (lines 242–244, and the _0_J / _H_M_1_J variants) computes the same term in int32 and never truncates:

src01 = _mm256_srai_epi32(_mm256_add_epi32(_mm256_mullo_epi32(_mm256_abs_epi32(src01), one_by_15), const_2048_32b), 12);

So for |src| ≥ 15360 the C threshold is smaller (or negative) and the C x is larger than the AVX2 x. The scalar tail macros copied into adm_avx2.c do carry the cast, which is why only the 8-wide interior differs.

2. The cube's shift. ADM_CM_ACCUM_ROUND in C truncates x_sq to int32_t and then applies >> (arithmetic) to the possibly negative product:

x_sq = (int32_t)((((int64_t)x * x) + add_shift_xsq) >> shift_xsq);
val  = (((int64_t)x_sq * x) + add_shift_xcub) >> shift_xcub;

ADM_CM_ACCUM_ROUND_avx256 (lines 661–662) gets the same truncated product via _mm256_mul_epi32 but shifts it with _mm256_srli_epi64, a logical shift. They differ whenever x ≥ 2^30, which is reachable once divergence 1 makes the threshold negative. (The file already has an sra_epi64 helper.)

Reproducer

Any of these, on master 6ec23e8 (line numbers in this report refer to it; the first measurements were on 86da14d):

1. libvmaf/test/checkasm/check_adm.c, fill_band: change % 16001) - 8000 to % 32001) - 16000. Then checkasm --test=adm --function='adm_cm*' --repeat=10.

band range	tolerance	adm_cm_*_avx2
±8000 (stock)	stock 1e-4 or 0	pass
±15000	0	pass
±16000	stock 1e-4	fail, 8/10 seeds
±16000	0	fail
±20000 / ±29000 / ±32000	either	fail

2. Instrument both adm_cm returns with the raw accum_h/v/d and run any clip with hard edges; the int64 accumulators differ for the scale-0 calls only.

Candidate fix (AVX2, 5 lines)

Replicating the C semantics makes everything bit-exact again: widened checkasm at ±32000 with tol = 0 passes 20/20 seeds, the stock suite passes 132/132, and the clip above is identical to the last digit in every integer_adm* / integer_aim value.

@@ ADM_CM_THRESH_S_I_J_avx256 (and the _0_J / _H_M_1_J variants, 3 sites)
     src01 = _mm256_srai_epi32(_mm256_add_epi32(_mm256_mullo_epi32(_mm256_abs_epi32(src01), one_by_15), const_2048_32b), 12); \
+    src01 = _mm256_srai_epi32(_mm256_slli_epi32(src01, 16), 16); /* replicate the (int16_t) cast in the C reference */ \
     (same for src11, src21)
@@ ADM_CM_ACCUM_ROUND_avx256
-    x_sq_lo = _mm256_srli_epi64(_mm256_add_epi64(_mm256_mul_epi32(x_sq_lo, x), _mm256_set1_epi64x(add_shift_xcub)), shift_xcub); \
-    x_sq_hi = _mm256_srli_epi64(_mm256_add_epi64(_mm256_mul_epi32(x_sq_hi, _mm256_srli_epi64(x, 32)), _mm256_set1_epi64x(add_shift_xcub)), shift_xcub); \
+    x_sq_lo = sra_epi64(_mm256_add_epi64(_mm256_mul_epi32(x_sq_lo, x), _mm256_set1_epi64x(add_shift_xcub)), _mm256_set1_epi64x(shift_xcub)); \
+    x_sq_hi = sra_epi64(_mm256_add_epi64(_mm256_mul_epi32(x_sq_hi, _mm256_srli_epi64(x, 32)), _mm256_set1_epi64x(add_shift_xcub)), _mm256_set1_epi64x(shift_xcub)); \

adm_avx512.c lines 372–374 need the same treatment; its accumulate macro already shifts arithmetically.

Whether the right fix is to replicate the C wrap (above; keeps every published score) or to remove the (int16_t) cast from the reference (arguably the intended arithmetic, but it changes scores on extreme-contrast content and touches AVX2/AVX-512 together) is your call. Either way the SIMD paths and the reference should agree.

Also observed, not root-caused

adm_csf_*_avx2 fails checkasm's exact 2-D compare (check_adm.c:297) for every size at ±32000 and passes at ±15000. Likely the same family of int16 truncation differences. Widening fill_band to the full int16 domain would surface both, and tol = 0 for adm_cm / i4_adm_cm seems appropriate for an integer kernel.

Environment

i9-12900K (AVX2, no AVX-512), Ubuntu 24.04 container, gcc 13.3.0, meson 1.3.2, -Dbuildtype=release -Denable_float=false -Denable_checkasm=true, master @ 86da14d. Re-checked on 2026-10-01 against 6ec23e8 on a Ryzen 9 9950X3D (AVX2 and AVX-512), gcc 16.2.1, meson 1.12.1, -Dc_args=-Wno-error=incompatible-pointer-types (#1630): the widened checkasm table below reproduces for adm_cm_*_avx2 (10 of 10 seeds for every range from ±16000 up, none for ±8000 and ±15000, with the stock tolerance and with tol = 0), the candidate fix passes 20 of 20 seeds at ±32000 with tol = 0, and adm_csf_*_avx2 fails at ±32000 and passes at ±15000. The exact stripe clip of the end-to-end table was not regenerated; on random 0/255 clips (576x324, three frames, 2%, 20% and 50% of the samples inverted in the distorted clip) the AVX2 path (--cpumask 16) and --cpumask -1 differ in integer_adm_scale0, integer_adm2, integer_aim and integer_adm3 by up to 1.1e-3, 1.6e-3 and 5.9e-3 respectively, and are identical to the last digit with the candidate fix.

[lusoris]

#1602 fixes this; I should have linked it when I opened the issue. It takes the other route: the centre tap stays in 32 bits in the scalar, x86 edge and CUDA code, as anematode suggested in review there, and a0c54d6 makes the AVX2 cube shift arithmetic. Three corrections to the description, checked against 86da14d:

• AVX-512 has only the first divergence. ADM_CM_ACCUM_ROUND_avx512 already shifts the cube with _mm512_srai_epi64 (adm_avx512.c:764-765); only its centre term (adm_avx512.c:361-363) differs from C.
• The extractor does not reach the second divergence, even with the wrap. The diagonal centre term cannot wrap (csf_a for d is at most 22916 * 49417 >> 17, about 8640), so the threshold stays at or above -65536. A coefficient at the dwt2 bound of 22916 then gives x of at most about 9.0e8 for h and v and 1.40e9 for d. Those are below the 2^30 and about 1.52e9 at which x_sq leaves int32_t. The end-to-end differences in the table come from the centre term. Only checkasm reaches the cube shift, because it fills csf_f with signed values. With check_adm_cm's bands at ±22916 and the stock tolerance, adm_cm_*_avx2 fails on 20 of 20 seeds on 86da14d, on 13 of 20 with only the centre-tap change, and on none with both.
• Neither route keeps every published score. Where the wrap is reached, the default AVX2 dispatch and --cpumask -1 already disagree, so one of them has to change. With feature/adm: keep contrast masking thresholds in range #1602, the three reference pairs and the 16-bit src01 pair are identical to 86da14d at 17 significant digits under both.

I will keep this open until #1602 is merged or declined.