Based on the presentation "What is ISPC?" by Dimo Yordanov, together with the accompanying ML kernel project.
Modern CPUs ship with wide SIMD (Single Instruction, Multiple Data) units - SSE, AVX on x86, NEON on ARM - that can process 4, 8, 16 or more numbers per instruction. Extracting that throughput from ordinary C++ is notoriously awkward. The solution for this was pioneered by Matt Pharr ISPC (the Intel Implicit SPMD Program Compiler) exists to close that gap: you write code that looks almost exactly like C, and the compiler maps it onto the machine's vector lanes for you auto-magicaly.
This article walks through what ISPC is, the programming model it introduces, and a concrete project that uses it to build the numerical core of a small neural network - activation functions, matrix multiplication, a loss function, and full backpropagation - all callable from a C++ training loop.
Before ISPC, a programmer who wanted to use SIMD had three unsatisfying options:
_mm256_add_ps and friends) are a middle ground, but they are verbose, tied to a specific instruction set, and force the programmer to think in terms of registers rather than algorithms.ISPC offers a fourth path: a high-level, C-like language where vectorization is part of the execution model rather than something you hope the optimizer notices.
ISPC is a compiler for a small language that is, deliberately, almost C. If you know C or C++, you can read ISPC immediately. Its defining idea is SPMD - Single Program, Multiple Data: you write one program as if it operated on a single scalar element, and the compiler runs many copies of that program called gangs in parallel across the SIMD lanes of the CPU.
This is similar to the mental model GPU programmers use in CUDA or OpenCL, brought to the CPU's vector units. The programmer reasons about "one element" while the hardware processes a whole vector's worth of elements at once.
ISPC originated at Intel as a research project and is now an open-source tool with broad instruction-set support.
The SPMD-on-CPU model is a good fit for any workload that is regular and data-parallel:
Getting the compiler is straightforward on the common platforms:
sudo apt install ispcbrew install ispcThe central runtime concept in ISPC is the gang - a group of program instances that execute together, in lockstep, mapped onto the lanes of a SIMD register. When you launch ISPC code, the work is distributed across the gang: each instance handles a different element, and they all step through the same instructions at the same time. The gang's size is bounded by the target's vector width, which is also where the model's limitations come from - divergent control flow, for example, must be handled by masking lanes on and off rather than truly branching independently.
The most important distinction in ISPC, and the one that determines whether your code is fast, is uniform vs. varying:
uniform value is scalar: it holds a single value shared by the whole gang and lives in an ordinary scalar register. Loop bounds, indices that are the same for every lane, and accumulators you intend to reduce are typically uniform.varying value is a vector: it holds a different value for each lane and lives in a SIMD register (XMM/YMM/ZMM on x86). A varying variable is the default for per-element data.Efficient ISPC comes from using uniform and varying correctly. Casts between them are allowed but only in one direction - you can widen a uniform into a varying, but turning a varying back into a uniform requires a reduction (for example summing every lane into a single scalar or finding the maximum of the lane). Reduction operations are the sanctioned way to collapse a vector down to one number.
ISPC exposes a couple of built-ins that describe the gang:
programCount - the number of instances in the gang (calculated from the vector width).programIndex - the index of the current instance within the gang, i.e. which lane "this" copy of the program is.Arrays behave exactly as they do in C and C++.
ISPC has pointers just like C and C++. Because both the pointer and the thing it points at can be uniform or varying, there are four combinations:
You can define structures just as in C/C++. They come in uniform and varying flavors, and ISPC can automatically convert between Array-of-Structures (AoS) and Structure-of-Arrays (SoA) layouts. SoA is generally friendlier to SIMD because it keeps each field contiguous, but the automatic conversion is not always the best-performing choice - layout still deserves the programmer's attention. Also sometimes the base choice of action is to do a Array of Structs of arrays with choice of parameter how many structs. The best way to find the optimum is usually to try many combinations with a valid benchmark.
ISPC supports the control-flow constructs you would expect, with a few SPMD-specific additions:
if, else, switch.for, while, do…while.foreach - the idiomatic parallel loop that distributes its iteration space across the gang.foreach_active - iterate over all active gang valuesforeach_unique - iterate over all unique gang valuesUnder the hood, divergent branches are handled by lane masking.
ISPC ships with a rich standard library covering logical operations, bitwise operations, math, random-number generation, reductions, masking/lane operations, and tasking (threads). It also provides several implementations of the math library so you can trade accuracy for speed:
The accompanying project implements the kernels a machine-learning model needs, in ISPC, and drives them from C++. The ISPC side is split into four source files - activations, backpropagation, linear algebra, and loss - and the C++ side wraps them into a tidy layer/model framework.
activation.ispc, linear.ispc, loss.ispc)The elementwise activations are textbook SPMD: a single foreach loop over the input, where every lane computes one element.
export void relu(uniform float input[], uniform float output[], uniform int count) {
foreach (i = 0 ... count) {
output[i] = max(input[i], 0.0f);
}
}
export void sigmoid(uniform float input[], uniform float output[], uniform int count) {
foreach (i = 0 ... count) {
output[i] = 1.0f / (1.0f + exp(-input[i]));
}
}
tanh_activation follows the same shape, computing the hyperbolic tangent from an exp-based formula. softmax is more interesting because it needs a normalization sum: it exponentiates each element into a varying accumulator, then uses reduce_add to collapse that vector into a single uniform total before dividing - a clean demonstration of the varying→uniform reduction rule:
export void softmax(uniform float input[], uniform float output[], uniform int n) {
varying float vsum = 0.0f;
foreach (i = 0 ... n) {
output[i] = exp(input[i]);
vsum += output[i];
}
uniform float total = reduce_add(vsum);
foreach (i = 0 ... n) {
output[i] /= total;
}
}
The linear-algebra kernels include a dot product (a multiply-accumulate into a varying partial, finished with reduce_add), elementwise helpers (add_vector, multiply_scalar), and two matrix-multiply implementations:
matmul - the naive triple loop, parallelized over the output's rows and columns with a two-dimensional foreach.matmul2 - a register-blocked, tiled version. It chooses a tile width at compile time from the target's vector width (TILE_WIDTH is derived from TARGET_WIDTH via the preprocessor) and a fixed TILE_HEIGHT of 8. It accumulates an entire tile of the output in registers while streaming over the shared dimension K, which dramatically improves data reuse and cache behavior. For inputs too small to tile, it falls back to the plain implementation.The mean-squared-error loss in loss.ispc accumulates squared differences into a partial and reduces:
export void mse_loss(uniform float predictions[], uniform float targets[],
uniform int n, uniform float result[]) {
float partial = 0.0f;
foreach (i = 0 ... n) {
float diff = predictions[i] - targets[i];
partial += diff * diff;
}
result[0] = reduce_add(partial) / (float)n;
}
backprop.ispc)Training requires gradients, so the project mirrors every forward kernel with a backward pass:
relu_backward - passes the upstream gradient through where the input was positive, zero elsewhere.sigmoid_backward and tanh_backward - use the cached forward output to avoid recomputing the expensive transcendental, applying σ(1−σ) and 1−tanh² respectively.softmax_backward - contracts the full Jacobian with the upstream gradient, again using a reduce_add to compute the shared dot product term.mse_backward - the simple 2·(pred − target)/n gradient.matmul_backward_A and matmul_backward_B - the two gradient matmuls (dA = dC·Bᵀ and dB = Aᵀ·dC) needed to backpropagate through a linear layer.Each backward kernel reuses the same SPMD patterns: elementwise foreach for the activations, nested foreach plus an inner uniform accumulation loop for the matrix gradients.
ISPC compiles each source file to an object file plus a generated header, which you then link against your C++ build. The workflow has three steps:
Compile the ISPC source to an object file and header, choosing a target ISA and optimization level:
ispc src/activation.ispc --target=avx -O3 -o target/activation.o -h target/activation.h
Compile your C++ against the generated headers:
g++ main.cpp -Itarget -O3 -Wall -c -o target/main.o
Link everything together:
g++ target/main.o target/activation.o -o target/main
On the C++ side, integration is minimal: include the generated header, and call the kernels. ISPC places the generated functions in the ispc namespace, so the C++ code calls ispc::relu(...), ispc::matmul(...), and so on.
The C++ layer turns these kernels into an ergonomic neural-network API:
ActivationFn is a small strategy object. Each activation bundles a forward kernel (which also caches whatever value its backward pass will need - the pre-activation for ReLU, the post-activation for sigmoid/tanh/softmax) and a backward kernel. Because the lambdas capture nothing, they decay to plain function pointers with no heap allocation.DenseLayer is a fully connected layer computing activation(W·x + b). Its forward calls ispc::matmul and ispc::add_vector; its backward calls the gradient matmuls and the activation's backward kernel; its update_weights applies SGD via ispc::multiply_scalar and ispc::add_vector.Model is a fluent builder that chains layers and validates that adjacent layer dimensions match.forward and train_step run inference and a full forward+backward+SGD step (returning the MSE loss) respectively.weight_io serializes and loads model weights in a small binary format with a magic header (MLWT), a version field, and per-layer shapes - with validation on load.run_model.cpp demonstrates the whole stack with a hand-tuned XOR network: a 2→8→1 architecture with sigmoid activations whose weights are baked in as constants, verifying that the ISPC kernels produce correct inference results across all four XOR inputs.Note on the numbers: the benchmark harness reports its "speedup" column as
scalar_time / ispc_time − 1rather than the conventionalscalar_time / ispc_time. That subtraction shifts every figure down by one and makes near-parity cases look like a fraction of a unit and slower cases appear negative. The discussion below describes the behavior the benchmarks reveal rather than reproducing the tabulated figures verbatim, because the printed column should be interpreted with that formula in mind.
The benchmark (benchmark.cpp) compares each ISPC kernel against an equivalent scalar C++ loop. The setup, as defined in the source:
[−1, 1].The qualitative picture that emerges is consistent and instructive:
exp) multiplies throughput across the lanes.tanh is the cautionary case: depending on the target and the math-library variant, the ISPC version can be slower than the scalar one. The scalar std::tanh is a heavily optimized library routine, while the project's ISPC tanh is built from an exp-based formula; on some targets that formulation costs more than it saves.matmul already gains substantially, and the tiled matmul2 - with its register blocking and target-aware tile width - pulls far ahead, with the advantage growing as the matrices get larger and data reuse matters more. The relative ranking of matmul2 over matmul is robust across targets, even though the absolute speedup figures vary widely by instruction set and vector width.The headline lesson is the one the programming-model section predicts: ISPC's benefit scales with the arithmetic intensity of the kernel. Workloads dominated by memory traffic see little benefit; workloads dominated by computation, especially those reusing data in cache like tiled GEMM, see the largest gains.
A second experiment sweeps the ReLU kernel across array sizes from a few kilobytes up to over a hundred megabytes, measuring effective memory throughput in GB/s. The result is a textbook memory hierarchy curve:
This reinforces why ReLU shows no SIMD speedup: once the array no longer fits in cache, the kernel is entirely bandwidth-bound, and the bottleneck is moving bytes, not computing on them.
ISPC is not a toy - it powers performance-critical code in well-known projects:
It is also a plausible drop-in for many other places where x86 efficiency matters.
The ecosystem continues to grow well beyond the basics in this article:
ISPC delivers SIMD performance through a programming model that stays close to plain C, sparing the programmer from intrinsics and assembly while giving far stronger guarantees than hope-it-vectorizes auto-vectorization. The neural-network case study makes the trade-offs concrete: the same foreach/uniform/varying/reduce_add patterns express activations, losses, and gradients cleanly, and the payoff is largest exactly where you would expect - compute-bound, cache-friendly kernels like tiled matrix multiplication - while memory-bound kernels remind us that no instruction set can outrun the memory bus. For data-parallel numerical code on the CPU, ISPC is a pragmatic and powerful tool.