Understanding which operations use vector units or many cores explains benchmarks and avoids oversubscription.
NumPy speeds up work at three levels, and they behave differently.
SIMD inside one core. Many ufunc and sorting loops are written with universal intrinsics (and Google Highway for some routines) that compile to SSE/AVX2/AVX-512 on x86, NEON/SVE on ARM, and VSX on POWER. Wheels contain a baseline plus extra variants, and the best one for the running CPU is chosen at import time (runtime dispatch). np.show_runtime() prints what was detected and enabled; the environment variable NPY_DISABLE_CPU_FEATURES can switch features off for debugging. SIMD paths apply mainly to contiguous, aligned, native-dtype data, another reason layout matters.
Multiple cores via BLAS. NumPy's own element-wise operations run on one thread. But matrix multiplication and most of np.linalg call a BLAS/LAPACK library (OpenBLAS in PyPI wheels, Accelerate on recent macOS wheels, MKL in some distributions), which is multithreaded. np.show_config() shows which one you have. Thread count is set with environment variables such as OPENBLAS_NUM_THREADS or OMP_NUM_THREADS before import, or at runtime with the threadpoolctl package. Running a process pool where each worker also starts a full BLAS thread pool oversubscribes the CPU; limit BLAS to one thread per worker.
Your own threads. NumPy releases the GIL during most array loops, so Python threads working on large arrays do run in parallel. With small arrays the per-call Python overhead dominates and threads don't help. NumPy 2.1+ supports the free-threaded CPython build (3.13t), removing the GIL from Python-level code as well. NumPy does not lock arrays: concurrent writes to the same array from several threads are a data race, regardless of the GIL.
Example
import numpy as np
import time
from concurrent.futures import ThreadPoolExecutor
np.show_runtime() # SIMD features found, BLAS info (prints a report)
a = np.random.default_rng(0).random((2000, 2000))
t = time.perf_counter(); a @ a; print("matmul (BLAS, threaded):", round(time.perf_counter() - t, 3))
# NumPy releases the GIL inside loops: threads can overlap big array work
chunks = np.array_split(np.random.default_rng(1).random(20_000_000), 4)
with ThreadPoolExecutor(max_workers=4) as ex:
total = sum(ex.map(lambda c: np.sqrt(c).sum(), chunks))
print(round(total))
# Limit BLAS threads at runtime (pip install threadpoolctl):
# from threadpoolctl import threadpool_limits
# with threadpool_limits(limits=1, user_api="blas"):
# a @ aRemember
- Ufunc loops use CPU-specific SIMD chosen at import;
np.show_runtime()shows it - Element-wise NumPy is single-threaded; matmul and linalg use multithreaded BLAS
- Set BLAS threads with env vars or threadpoolctl to avoid oversubscription
- NumPy releases the GIL in loops, but arrays aren't locked against concurrent writes
Watch out for
- Running 16 worker processes that each use a 16-thread BLAS pool and wondering why it is slower than one process.