CPU: even-odd (symmetric) factorization for the AVX basis-apply — benchmark + RFC
- Dominant language
- C
- Stars
- 265
- Forks
- 78
- Avg merge
- 2d 8h
- Merged PRs (30d)
- 11
Description
Following a suggestion from @jedbrown to open an issue describing thoughts and plans — here's the RFC for adding an even-odd (symmetric) factorization path to the AVX tensor contraction kernel.
## Proposal
Add an even-odd factorization code path to `CeedTensorContractApply_Avx` in `backends/avx/ceed-avx-tensor.c`. For bases whose 1D interpolation/gradient matrix has the point symmetry of GLL/GL nodes (`B[i,j] = ±B[n-1-i, n-1-j]`), this precomputes symmetric and antisymmetric halves and uses them to halve the FMA count in the inner contraction loop.
## Evidence
I've benchmarked this in [mf-kernels](https://github.com/mohitt31/mf-kernels), a standalone set of five 1D tensor-contraction kernel variants (naive, pitfall, explicit AVX2, register-blocked AVX2, and even-odd AVX2). The isolated effect of even-odd vs register-blocked AVX2 on the bare 1D contraction:
| p | blocked (GFLOP/s) | even-odd (GFLOP/s) | speedup |
|---|---:|---:|---:|
| 5 | 20.82 | 21.73 | 1.04× |
| 7 | 21.58 | 25.92 | 1.20× |
| 9 | 21.09 | 28.36 | 1.34× |
(AMD EPYC Zen 3, single core, AVX2+FMA, gcc 13.3, `-O3 -march=native`. GFLOP/s uses the standard algorithm's FLOP count so the algorithmic saving shows as a wall-clock gain, not a smaller denominator.)
The gain grows with order because the ratio of halved FMAs to fixed overhead improves. At p=9, even-odd reaches ~28 GFLOP/s on the bare contraction.
The same repo also has a full 3D sum-factorized operator benchmark (BP1 mass + BP3 Poisson, p=1..8, structured hex mesh, ~0.5M DOFs) against MFEM native partial assembly and libCEED's `/cpu/self/avx/blocked` and `/cpu/self/xsmm/blocked` backends. This establishes the contraction kernels in a realistic operator context and confirms correctness (worst-case max-abs difference vs MFEM's own output: 3.8e-15 across all variants and orders).
## Where this does NOT help — honest caveats
1. **Even-odd does not beat register-blocking in the full operator.** Once gather/scatter, element restriction, and the pointwise quadrature-data apply are included, the operator becomes partly memory-bound and the contraction's FMA savings wash out. In my full-operator benchmark, the register-blocked variant matches or leads even-odd at every order on both BP1 and BP3 — they only converge at p=8 BP3 (33.9 vs 33.7 M DOF/s). The CeedTensorContract is only one component of the full CeedOperator pipeline.
2. **LIBXSMM leads at high order.** libCEED with LIBXSMM outperforms both the mf-kernels variants and the plain AVX backend at higher orders (BP1 p=8: LIBXSMM 124.9 vs mf-kernels blocked 96.2 M DOF/s). This proposal improves the plain AVX backend's contraction kernel only, not the LIBXSMM path.
3. **Scope is `CeedTensorContract` only.** Specifically `CeedTensorContract_Avx_Blocked` and the dispatch in `CeedTensorContractApply_Avx`. Not the full CeedOperator, not restriction, not the LIBXSMM backend.
4. **Portable SIMD intrinsics** (AVX2/FMA via `immintrin.h`, matching the existing backend style), not vendor-specific. This is unrelated to the SVE backend work in #1984.
## Proposed plan
1. **This RFC** — confirm the scoping makes sense.
2. **WIP PR** adding the even-odd path to `backends/avx/ceed-avx-tensor.c`:
- Precompute the symmetric/antisymmetric halves of the basis matrix at `CeedBasis` setup time.
- Gate activation on (a) detected symmetry in the 1D basis matrix and (b) a minimum order threshold so it only activates where it actually wins.
- Existing contraction behavior is untouched for non-symmetric bases or low orders.
3. **Reproducible benchmark**: BP1 + BP3, p=1..8, single-core, comparing the even-odd path against the current AVX backend (`/cpu/self/avx/blocked`).
4. **Correctness**: match current backend output to ~1e-15.
5. `make format` + `make prove-all` clean.
## Question
Does this scoping look right? Is there a particular order range or operator type that's most important for the AVX backend's downstream users?
Contributor guide
Research direction
Start in backends/avx/ceed-avx-tensor.c by reading CeedTensorContract_Avx_Blocked and the dispatch in CeedTensorContractApply_Avx. Review the existing AVX backend style and the linked mf-kernels benchmark before evaluating the proposed symmetry detection and order threshold. Done means the scoped path is correct, benchmarked against the current AVX backend, and make format plus make prove-all pass.
Written by the indexing model from the issue text.
Assessment
- Tech stack
- c
- Domain
- hpc, performance
- Issue type
- Feature
- Difficulty
- 5/5
- Estimated time
- Over a week
- Activity status
- Active
- Clarity
- Clearly specified
- Newbie friendliness
- 38/100