NumPy versus a Python loop: a measured dot product¶

A single timing can be dominated by warm-up, scheduling, or library initialization. This notebook compares numpy.dot with an explicit Python loop across four input sizes and nine batched samples. Each batch runs for at least 50 ms where one call is shorter than that threshold. It checks that both implementations return the same value before reporting any timing.

This is a microbenchmark, not a universal claim about Python or NumPy. The result depends on the processor, BLAS implementation, thread configuration, array dtype, and the exact loop being compared.

In [1]:
import platform
import time

import numpy as np
import pandas as pd
from threadpoolctl import threadpool_info, threadpool_limits


def python_dot(left, right):
    return sum(float(x) * float(y) for x, y in zip(left, right))


def elapsed_batch(callable_, iterations):
    started = time.perf_counter()
    result = None
    for _ in range(iterations):
        result = callable_()
    return result, (time.perf_counter() - started) / iterations


def calibrate(callable_, target_seconds=0.05):
    iterations = 1
    while iterations < 1_048_576:
        started = time.perf_counter()
        for _ in range(iterations):
            callable_()
        elapsed = time.perf_counter() - started
        if elapsed >= target_seconds:
            break
        iterations *= max(2, int(target_seconds / max(elapsed, 1e-9)))
    return iterations


def benchmark(size, repeats=9, seed=42):
    rng = np.random.default_rng(seed)
    left = rng.random(size, dtype=np.float64)
    right = rng.random(size, dtype=np.float64)

    # Warm up both code paths before recording samples.
    expected = np.dot(left, right)
    observed = python_dot(left, right)
    np.testing.assert_allclose(observed, expected, rtol=1e-10)

    numpy_call = lambda: np.dot(left, right)
    loop_call = lambda: python_dot(left, right)
    numpy_iterations = calibrate(numpy_call)
    loop_iterations = calibrate(loop_call)

    numpy_samples = []
    loop_samples = []
    for sample_index in range(repeats):
        order = ('numpy', 'loop') if sample_index % 2 == 0 else ('loop', 'numpy')
        sample = {}
        for implementation in order:
            if implementation == 'numpy':
                sample['numpy'] = elapsed_batch(numpy_call, numpy_iterations)
            else:
                sample['loop'] = elapsed_batch(loop_call, loop_iterations)
        numpy_result, numpy_time = sample['numpy']
        loop_result, loop_time = sample['loop']
        np.testing.assert_allclose(loop_result, numpy_result, rtol=1e-10)
        numpy_samples.append(numpy_time)
        loop_samples.append(loop_time)

    numpy_median = float(np.median(numpy_samples))
    loop_median = float(np.median(loop_samples))
    return {
        'items': size,
        'samples': repeats,
        'numpy_calls_per_sample': numpy_iterations,
        'loop_calls_per_sample': loop_iterations,
        'numpy_median_ms': 1_000 * numpy_median,
        'numpy_iqr_ms': 1_000 * float(np.subtract(*np.percentile(numpy_samples, [75, 25]))),
        'loop_median_ms': 1_000 * loop_median,
        'loop_iqr_ms': 1_000 * float(np.subtract(*np.percentile(loop_samples, [75, 25]))),
        'median_speedup': loop_median / numpy_median,
    }


with threadpool_limits(limits=1):
    backend_info = threadpool_info()
    results = pd.DataFrame(
        [benchmark(size) for size in (1_000, 10_000, 100_000, 1_000_000)]
    )
results = results.round({
    'numpy_median_ms': 4, 'numpy_iqr_ms': 4,
    'loop_median_ms': 3, 'loop_iqr_ms': 3, 'median_speedup': 1,
})
results
Out[1]:
items samples numpy_calls_per_sample loop_calls_per_sample numpy_median_ms numpy_iqr_ms loop_median_ms loop_iqr_ms median_speedup
0 1000 9 69454 386 0.0007 0.0001 0.114 0.011 162.4
1 10000 9 22016 66 0.0033 0.0003 1.146 0.129 347.2
2 100000 9 3814 8 0.0228 0.0007 9.918 0.337 434.1
3 1000000 9 106 1 0.4519 0.0263 99.622 7.790 220.4
In [2]:
environment = {
    'python': platform.python_version(),
    'platform': platform.platform(),
    'numpy': np.__version__,
    'thread_limit_requested': 1,
    'thread_limit_observed': bool(backend_info) and all(
        pool.get('num_threads') == 1 for pool in backend_info
    ),
    'threadpools': backend_info or 'No compatible controller detected; backend thread count is uncontrolled.',
}
environment
Out[2]:
{'python': '3.12.13',
 'platform': 'macOS-26.5.2-arm64-arm-64bit',
 'numpy': '2.5.2',
 'thread_limit_requested': 1,
 'thread_limit_observed': False,
 'threadpools': 'No compatible controller detected; backend thread count is uncontrolled.'}

Reading the table¶

The median is the central measurement across nine batched samples; the interquartile range (IQR) shows the spread of the middle half. The two implementations alternate which one runs first. The speedup is calculated from their medians, and the displayed precision follows the scale of each measurement.

The benchmark requests one native-library thread through threadpoolctl and records whether a compatible controller was observed inside that context. If thread_limit_observed is false, the backend thread count was not verified and this run must be treated as uncontrolled on that dimension.

The explicit loop iterates over NumPy scalar values, so this benchmark measures Python interpreter overhead on top of NumPy-backed storage. numpy.dot dispatches the complete operation to compiled numerical code and may use an optimized, threaded BLAS library. Those are the mechanisms being compared.