perf: lightning-fast vectorized PDF calculators

.github/scripts/check-kernel-isolation.py

@@ -0,0 +1,132 @@
+"""Check that the SIMD kernel objects share no code with the library.
+
+GCC and Clang objects must define no shared functions. MSVC objects may
+define them with baseline instructions only.
+"""
+
+import re
+import subprocess
+import sys
+from pathlib import Path
+
+OBJECT_NAME = re.compile(
+    r"CalculationKernel(SSE2|AVX2|AVX512|NEON)(?:\.cpp)?\.(?:o|obj)$"
+)
+COFF_FUNCTION = re.compile(
+    r"SECT[0-9A-F]+\s+notype \(\)\s+External\s+\| (\S+)"
+)
+INSTRUCTION = re.compile(
+    r"^\s+[0-9A-F]+:\s+((?:[0-9A-F]{2}\s+)+)([a-z?][a-z0-9?]*)(?:\s|$)"
+)
+# Instructions above the baseline that have no VEX prefix.
+EXTRA_ISA = {"tzcnt", "lzcnt", "popcnt", "crc32"}
+
+
+def run(*command):
+    return subprocess.run(
+        command, check=True, capture_output=True, text=True, errors="replace"
+    ).stdout
+
+
+def private_symbol(name, isa):
+    namespace = f"xsimd_{isa.lower()}"
+    if f"{namespace}::" in name or f"@{namespace}@@" in name:
+        return True
+    return "(anonymous namespace)" in name or "?A0x" in name
+
+
+def kernel_entry(name, isa):
+    entry = f"{isa.lower()}CalculationKernel"
+    return name == f"diffpy::srreal::kernels::{entry}()" or name.startswith(
+        f"?{entry}@kernels@srreal@diffpy@@"
+    )
+
+
+def check_nm(path, isa):
+    output = run(
+        "nm",
+        "--defined-only",
+        "--extern-only",
+        "--demangle",
+        "--format=posix",
+        str(path),
+    )
+    functions = []
+    for line in output.splitlines():
+        match = re.fullmatch(r"(.+) ([A-Za-z?]) [0-9a-fA-F]+(?: .*)?", line)
+        if not match:
+            raise ValueError(f"Unrecognized nm output: {line}")
+        name, kind = match.groups()
+        if kind in {"T", "W", "i", "I"}:
+            functions.append(name)
+    if not any(kernel_entry(name, isa) for name in functions):
+        raise ValueError(f"Missing {isa} entry point in {path}")
+    return [
+        f"shared code symbol: {name}"
+        for name in functions
+        if not private_symbol(name, isa) and not kernel_entry(name, isa)
+    ]
+
+
+def check_coff(path, isa):
+    symbols = run("dumpbin", "/nologo", "/symbols", str(path))
+    functions = COFF_FUNCTION.findall(symbols)
+    if not any(kernel_entry(name, isa) for name in functions):
+        raise ValueError(f"Missing {isa} entry point in {path} (LTO object?)")
+    shared = {
+        name
+        for name in functions
+        if not private_symbol(name, isa) and not kernel_entry(name, isa)
+    }
+    output = run("dumpbin", "/nologo", "/disasm", str(path))
+    bodies = {}
+    current = None
+    for line in output.splitlines():
+        # A function starts with its decorated name.
+        if line and not line[0].isspace() and line.endswith(":"):
+            current = line.split(" ", 1)[0].rstrip(":")
+            bodies[current] = []
+        elif current is not None:
+            match = INSTRUCTION.match(line)
+            if match:
+                bodies[current].append(match.groups())
+    errors = []
+    for name in sorted(shared):
+        if not bodies.get(name):
+            raise ValueError(f"No disassembly for shared function {name}")
+        for raw, mnemonic in bodies[name]:
+            # C4 and C5 start VEX instructions, 62 starts EVEX ones.
+            if (
+                raw.split()[0] in {"C4", "C5", "62"}
+                or mnemonic.startswith("v")
+                or mnemonic in EXTRA_ISA
+                or "?" in mnemonic
+            ):
+                errors.append(f"non-baseline instruction {mnemonic}: {name}")
+                break
+    return errors
+
+
+def main():
+    build_dir = Path(sys.argv[1])
+    failures = []
+    found = False
+    for path in sorted(build_dir.rglob("CalculationKernel*")):
+        match = OBJECT_NAME.fullmatch(path.name)
+        if not match:
+            continue
+        found = True
+        check = check_coff if path.suffix == ".obj" else check_nm
+        errors = check(path, match[1])
+        if errors:
+            failures.extend(f"{path}: {error}" for error in errors)
+        else:
+            print(f"PASS {path}")
+    if not found:
+        failures.append(f"No SIMD kernel objects in {build_dir}")
+    if failures:
+        raise SystemExit("\n".join(failures))
+
+
+if __name__ == "__main__":
+    main()

.github/workflows/matrix-and-codecov.yml

@@ -75,15 +75,18 @@ jobs:
           cmake --build --preset release --parallel ${{ matrix.jobs }} || exit /b 1
           ctest --test-dir build\release --output-on-failure --no-tests=error
 
+      - name: Check SIMD kernel isolation
+        run: python .github/scripts/check-kernel-isolation.py build
+
       - name: Generate coverage report
-        if: runner.os == 'Linux'
+        if: matrix.coverage == 'ON'
         run: >-
           gcovr --gcov-executable "$CONDA_PREFIX/bin/gcov"
           --root . --filter 'src/diffpy/'
           --xml coverage.xml
 
       - name: Save coverage report
-        if: runner.os == 'Linux'
+        if: matrix.coverage == 'ON'
         uses: actions/upload-artifact@v4
         with:
           name: coverage-report
@@ -92,7 +95,7 @@ jobs:
 
       - name: Upload coverage to Codecov
         # Fork push builds have no upstream Codecov token; keep their artifact.
-        if: runner.os == 'Linux' && github.repository == 'diffpy/libdiffpy'
+        if: matrix.coverage == 'ON' && github.repository == 'diffpy/libdiffpy'
         uses: codecov/codecov-action@v5
         with:
           files: coverage.xml

CMakeLists.txt

@@ -35,6 +35,8 @@ endif()
 option(LIBDIFFPY_BUILD_OBJECT_LIB "Build an OBJECT library for embedding" OFF)
 option(LIBDIFFPY_INSTALL "Install libdiffpy and its CMake package" ON)
 option(LIBDIFFPY_ENABLE_OBJCRYST "Build with the ObjCryst structure adapter" ON)
+option(LIBDIFFPY_LTO
+  "Link the optimized shared library with link-time optimization" ON)
 option(LIBDIFFPY_BUILD_TESTS "Build the Catch2 unit tests" OFF)
 option(LIBDIFFPY_TEST_INSTALLED
   "Build only the unit tests and link them with the installed libdiffpy" OFF)
@@ -294,6 +296,27 @@ endif()
 if(MSVC)
   # Required for the Boost.Serialization code, including embedded builds.
   target_compile_options(diffpy PRIVATE /bigobj)
+else()
+  # Let GCC inline calls between the exported functions of the library.
+  include(CheckCXXCompilerFlag)
+  check_cxx_compiler_flag(-fno-semantic-interposition
+    LIBDIFFPY_HAS_NO_SEMANTIC_INTERPOSITION)
+  if(LIBDIFFPY_HAS_NO_SEMANTIC_INTERPOSITION)
+    target_compile_options(diffpy PRIVATE -fno-semantic-interposition)
+  endif()
+  # LTO inlines the short per-pair calls between source files.  Static and
+  # OBJECT libraries skip it, their objects would need the same compiler.
+  if(LIBDIFFPY_LTO AND _libdiffpy_library_type STREQUAL "SHARED_LIBRARY")
+    include(CheckIPOSupported)
+    check_ipo_supported(RESULT _ipo_supported OUTPUT _ipo_error LANGUAGES CXX)
+    if(_ipo_supported)
+      set_target_properties(diffpy PROPERTIES
+        INTERPROCEDURAL_OPTIMIZATION_RELEASE ON
+        INTERPROCEDURAL_OPTIMIZATION_RELWITHDEBINFO ON)
+    else()
+      message(STATUS "Link-time optimization is not supported, skipping it.")
+    endif()
+  endif()
 endif()
 target_compile_features(diffpy PUBLIC cxx_std_23)
 target_include_directories(diffpy PUBLIC
@@ -319,6 +342,7 @@ target_link_libraries(diffpy
 if(LIBDIFFPY_HAS_OBJCRYST)
   target_link_libraries(diffpy PUBLIC ObjCryst::ObjCryst)
 endif()
+include(CalculationKernels)
 
 # Installation ----------------------------------------------------------------
 

cmake/CalculationKernels.cmake

@@ -0,0 +1,113 @@
+# Calculation kernels.  The scalar kernel is always built.  With LIBDIFFPY_SIMD
+# and xsimd, all SIMD kernels of the target architecture are built too and
+# the library selects the best one the CPU supports at runtime, see
+# src/diffpy/srreal/kernels/CalculationKernelDispatch.cpp.
+option(LIBDIFFPY_SIMD
+  "Build the SIMD calculation kernels of the target architecture when xsimd \
+is found" ON)
+
+# Target architecture.  A universal macOS build matches no architecture
+# below and gets only the scalar kernel.
+set(_kernel_arch "${CMAKE_SYSTEM_PROCESSOR}")
+if(MSVC AND CMAKE_CXX_COMPILER_ARCHITECTURE_ID)
+  set(_kernel_arch "${CMAKE_CXX_COMPILER_ARCHITECTURE_ID}")
+endif()
+if(APPLE AND CMAKE_OSX_ARCHITECTURES)
+  set(_kernel_arch "${CMAKE_OSX_ARCHITECTURES}")
+endif()
+string(TOLOWER "${_kernel_arch}" _kernel_arch)
+if(_kernel_arch MATCHES "^(x86_64|amd64|x64)$")
+  set(_kernel_arch "x86_64")
+elseif(_kernel_arch MATCHES "^(aarch64|arm64)$")
+  set(_kernel_arch "arm64")
+else()
+  set(_kernel_arch "")
+endif()
+
+# The SIMD kernels need xsimd.
+set(_kernel_simd FALSE)
+if(LIBDIFFPY_SIMD AND _kernel_arch)
+  if(NOT TARGET xsimd)
+    find_package(xsimd CONFIG QUIET)
+  endif()
+  if(TARGET xsimd)
+    set(_kernel_simd TRUE)
+  else()
+    message(STATUS "xsimd not found, building only the scalar calculation "
+      "kernel.  Install xsimd for the SIMD kernels.")
+  endif()
+endif()
+
+# The main source glob is not recursive, add the kernels explicitly.
+set(_kernel_dir "${PROJECT_SOURCE_DIR}/src/diffpy/srreal/kernels")
+target_sources(diffpy PRIVATE
+  "${_kernel_dir}/CalculationKernelScalar.cpp"
+  "${_kernel_dir}/CalculationKernelDispatch.cpp")
+set(_kernel_names scalar)
+
+# All kernels need IEEE arithmetic to give consistent results.
+if(MSVC)
+  set(_kernel_math_options /fp:precise)
+else()
+  set(_kernel_math_options -fno-fast-math -ffp-contract=off)
+endif()
+set_source_files_properties("${_kernel_dir}/CalculationKernelScalar.cpp"
+  PROPERTIES COMPILE_OPTIONS "${_kernel_math_options}")
+
+# Inline functions compiled for the SIMD instruction sets must not be shared
+# with the rest of the library, see CalculationKernelSIMD.ipp.  MSVC compiles
+# the kernels for the SSE2 baseline, its intrinsics need no /arch:AVX*.  With
+# GCC and Clang -O2 inlines the functions also in unoptimized builds.  /GL-
+# and -fno-lto keep LTO from sharing them.
+set(_simd_isolation_options)
+if(_kernel_simd AND MSVC AND _kernel_arch STREQUAL "x86_64")
+  list(APPEND _simd_isolation_options /arch:SSE2 /GL-)
+elseif(_kernel_simd AND NOT MSVC)
+  include(CheckCXXCompilerFlag)
+  check_cxx_compiler_flag(-fno-lto LIBDIFFPY_HAS_NO_LTO)
+  if(LIBDIFFPY_HAS_NO_LTO)
+    list(APPEND _simd_isolation_options -fno-lto)
+  endif()
+  list(APPEND _simd_isolation_options
+    "$<$<NOT:$<CONFIG:Release,RelWithDebInfo>>:-O2>")
+endif()
+
+# Add the SIMD kernel in CalculationKernel<name>.cpp, compiled with the
+# instruction set flags in ARGN.  CalculationKernelDispatch.cpp uses the
+# kernel when LIBDIFFPY_KERNEL_<name> is defined.
+function(_libdiffpy_add_simd_kernel name)
+  set(_source "${_kernel_dir}/CalculationKernel${name}.cpp")
+  set(_options ${ARGN} ${_kernel_math_options} ${_simd_isolation_options})
+  # Each kernel sets up xsimd for itself.
+  set_source_files_properties("${_source}"
+    PROPERTIES COMPILE_OPTIONS "${_options}"
+      SKIP_PRECOMPILE_HEADERS ON
+      SKIP_UNITY_BUILD_INCLUSION ON)
+  target_sources(diffpy PRIVATE "${_source}")
+  target_compile_definitions(diffpy PRIVATE LIBDIFFPY_KERNEL_${name}=1)
+  string(TOLOWER "${name}" _lower)
+  set(_kernel_names ${_kernel_names} ${_lower} PARENT_SCOPE)
+endfunction()
+
+if(_kernel_simd AND _kernel_arch STREQUAL "x86_64")
+  _libdiffpy_add_simd_kernel(SSE2)
+  if(MSVC)
+    _libdiffpy_add_simd_kernel(AVX2)
+    _libdiffpy_add_simd_kernel(AVX512)
+  else()
+    _libdiffpy_add_simd_kernel(AVX2 -mavx2 -mfma)
+    _libdiffpy_add_simd_kernel(AVX512
+      -mavx512f -mavx512cd -mavx512dq -mfma)
+  endif()
+elseif(_kernel_simd AND _kernel_arch STREQUAL "arm64")
+  _libdiffpy_add_simd_kernel(NEON)
+endif()
+list(JOIN _kernel_names ", " _kernel_names)
+message(STATUS "Calculation kernels: ${_kernel_names}")
+
+if(_kernel_simd)
+  target_compile_definitions(diffpy PRIVATE LIBDIFFPY_SIMD_KERNELS=1)
+  # Include path only, keeps xsimd out of the static library's link interface.
+  target_include_directories(diffpy SYSTEM PRIVATE
+    $<TARGET_PROPERTY:xsimd,INTERFACE_INCLUDE_DIRECTORIES>)
+endif()

docs/source/development.rst

@@ -45,7 +45,21 @@ To use the same options for every build, set them in a ``CMakeUserPresets.json``
 file in the repository root. To use another generator, configure a separate
 directory with ``cmake -S . -B <build-dir> -G "<generator>"``. To disable the
 ObjCryst adapter, pass ``-DLIBDIFFPY_ENABLE_OBJCRYST=OFF``. For gprof
-instrumentation, use ``-DLIBDIFFPY_PROFILE=ON``.
+instrumentation, use ``-DLIBDIFFPY_PROFILE=ON``. To build the optimized
+shared library without link-time optimization, pass ``-DLIBDIFFPY_LTO=OFF``.
+
+Calculation kernels
+-------------------
+
+``PDFCalculator`` and ``DebyePDFCalculator`` run their inner loops in a
+calculation kernel. The scalar kernel is always built. With the default
+``-DLIBDIFFPY_SIMD=ON`` the SIMD kernels of the target architecture are built
+too when `xsimd <https://github.com/xtensor-stack/xsimd>`_ is found;
+``-DLIBDIFFPY_SIMD=OFF`` builds only the scalar kernel.
+
+The SIMD kernels are ``avx512``, ``avx2`` and ``sse2`` on x86-64 and ``neon``
+on arm64. On first use the library picks the best kernel that the CPU
+supports, for example ``avx2`` on a CPU without AVX-512.
 
 Embed libdiffpy in another CMake project
 ----------------------------------------