What is Numeric Algorithms — in C++?

Complete reference for C++ numeric algorithms: accumulate, reduce, transform_reduce, prefix scans, gcd, lcm, midpoint, lerp, and parallel variants.

Which C++ standard introduced Numeric Algorithms — ?

Numeric Algorithms — was introduced in C++98.

What is the difficulty level of Numeric Algorithms — ?

Numeric Algorithms — is considered Basic-level C++ material.

Numeric Algorithms — <numeric>

<numeric>since C++98

The <numeric> header provides algorithms for numeric computation — sequential folds (accumulate, inner_product), prefix sums, and consecutive differences since C++98; parallel-capable counterparts and arithmetic utilities (gcd, lcm, midpoint, lerp) added through C++17 and C++20.

Overview

<numeric> complements <algorithm> for operations with numeric semantics. The C++98 set is entirely sequential and left-to-right ordered. C++17 introduced parallel-capable variants with deliberately relaxed ordering guarantees, and added arithmetic utilities. C++20 added midpoint and lerp. C++23 added std::ranges::fold_left / fold_right as the range-based fold story.

The central tension in this header: order matters. Every C++17 algorithm that can run in parallel sacrifices guaranteed left-to-right evaluation in exchange for that flexibility. That trade-off must be understood before substituting reduce for accumulate or inclusive_scan for partial_sum.

Examples

`accumulate` (C++98) vs `reduce` (C++17)

std::accumulate performs a strict left-to-right fold: result = init op v[0] op v[1] op .... This makes it correct for non-commutative operations, string concatenation, and anything with observable left-to-right order.

std::reduce makes no guarantee about evaluation order, enabling auto-vectorization and parallel execution. It is semantically equivalent to accumulate only when the binary operation is both associative and commutative.

cpp

Godbolt

#include <numeric>
#include <vector>
#include <functional>
#include <execution>  // C++17

std::vector<double> v = {1.0, 2.0, 3.0, 4.0, 5.0};

// accumulate — C++98, left-to-right, correct for any op
double sum  = std::accumulate(v.begin(), v.end(), 0.0);          // 15.0
double prod = std::accumulate(v.begin(), v.end(), 1.0,
                              std::multiplies<double>{});         // 120.0

// Wrong: accumulates as int until assignment
double bad  = std::accumulate(v.begin(), v.end(), 0);            // 15 (int arith)

// reduce — C++17, unspecified order; safe for associative+commutative ops
double rsum  = std::reduce(v.begin(), v.end());                  // init = 0.0
double rsum2 = std::reduce(v.begin(), v.end(), 100.0);           // 115.0
double rprod = std::reduce(v.begin(), v.end(), 1.0,
                           std::multiplies<>{});                 // 120.0

// Parallel reduce — requires -ltbb on GCC/Clang with libstdc++
double psum  = std::reduce(std::execution::par, v.begin(), v.end());  // C++17

The type of init determines the accumulation type. accumulate(v.begin(), v.end(), 0) on a vector<double> performs integer arithmetic for all intermediate steps.

`transform_reduce` (C++17)

Two overloads: a binary form (two ranges — transform pairs, then reduce) and a unary form (one range — unary transform, then reduce).

cpp

Godbolt

std::vector<double> a = {1.0, 2.0, 3.0};
std::vector<double> b = {4.0, 5.0, 6.0};

// Dot product (binary form): sum(a[i] * b[i])
double dot = std::transform_reduce(
    a.begin(), a.end(),
    b.begin(),
    0.0,                    // init
    std::plus<>{},          // reduce op
    std::multiplies<>{}     // element-wise transform
);  // 1*4 + 2*5 + 3*6 = 32.0

// L2 norm squared (unary form)
double norm_sq = std::transform_reduce(
    a.begin(), a.end(),
    0.0,
    std::plus<>{},
    [](double x) { return x * x; }
);  // 1 + 4 + 9 = 14.0

// Parallel, with a predicate — count elements meeting a condition
std::vector<int> data = {1, -2, 3, -4, 5};
int neg_count = std::transform_reduce(
    std::execution::par,
    data.begin(), data.end(),
    0,
    std::plus<>{},
    [](int x) { return x < 0 ? 1 : 0; }
);  // 2

inner_product (C++98) is the sequential predecessor to the binary form of transform_reduce. Prefer transform_reduce in new code; use inner_product only when the guaranteed left-to-right evaluation order matters.

Prefix Scans

partial_sum (C++98) computes out[i] = v[0] op v[1] op ... op v[i], left-to-right, no parallelism. The C++17 scan algorithms relax evaluation order for parallelism.

cpp

Godbolt

std::vector<int> v = {1, 2, 3, 4, 5};
std::vector<int> out(v.size());

// partial_sum — C++98, sequential, inclusive (current element included)
std::partial_sum(v.begin(), v.end(), out.begin());
// out = {1, 3, 6, 10, 15}

// inclusive_scan — C++17, same semantics for associative ops, can parallelize
std::inclusive_scan(v.begin(), v.end(), out.begin());
// out = {1, 3, 6, 10, 15}

std::inclusive_scan(v.begin(), v.end(), out.begin(), std::multiplies<>{});
// out = {1, 2, 6, 24, 120}

// exclusive_scan — C++17, out[i] = reduce of v[0..i-1]; init is the first output
std::exclusive_scan(v.begin(), v.end(), out.begin(), 0);
// out = {0, 1, 3, 6, 10}

std::exclusive_scan(v.begin(), v.end(), out.begin(), 1, std::multiplies<>{});
// out = {1, 1, 2, 6, 24}

// transform_inclusive_scan — C++17, apply unary transform first, then scan
std::transform_inclusive_scan(
    v.begin(), v.end(), out.begin(),
    std::plus<>{},
    [](int n) { return n * n; }
);
// out = {1, 5, 14, 30, 55}  — prefix sums of squares

// transform_exclusive_scan — C++17, note: init comes before the binary op
std::transform_exclusive_scan(
    v.begin(), v.end(), out.begin(),
    0,               // init — before binary_op, unlike transform_inclusive_scan
    std::plus<>{},
    [](int n) { return n * n; }
);
// out = {0, 1, 5, 14, 30}

`adjacent_difference` (C++98)

Computes out[0] = v[0], out[i] = v[i] op v[i-1]. With default subtraction, it is the inverse of partial_sum with default addition.

cpp

Godbolt

// Financial: day-over-day price deltas
std::vector<double> prices = {100.0, 110.5, 98.0, 115.0, 112.0};
std::vector<double> deltas(prices.size());

std::adjacent_difference(prices.begin(), prices.end(), deltas.begin());
// deltas = {100.0, 10.5, -12.5, 17.0, -3.0}

// Ratios: using custom op
std::adjacent_difference(prices.begin(), prices.end(), deltas.begin(),
                         [](double cur, double prev) { return cur / prev; });
// deltas[0] = 100.0, deltas[1] = 110.5/100.0, ...

// In-place is safe: the algorithm reads each input before overwriting
std::adjacent_difference(prices.begin(), prices.end(), prices.begin());

`iota` (C++11) and `views::iota` (C++20)

cpp

Godbolt

std::vector<int> v(10);
std::iota(v.begin(), v.end(), 0);   // {0, 1, 2, ..., 9}
std::iota(v.begin(), v.end(), 1);   // {1, 2, 3, ..., 10}

// C++20: views::iota is lazy — no allocation, generates on demand
for (int i : std::views::iota(0, 10)) { /* ... */ }

// C++23: std::ranges::iota writes to a range (range-based iota algorithm)
std::ranges::iota(v, 0);  // C++23

Prefer std::views::iota (C++20) for iteration and range-based use; std::iota only when a materialized container is needed.

`gcd` and `lcm` (C++17)

cpp

Godbolt

std::gcd(12, 8);    // 4
std::gcd(0, 5);     // 5  — gcd(0, n) == |n|
std::gcd(0, 0);     // 0  — well-defined
std::lcm(4, 6);     // 12
std::lcm(0, 5);     // 0  — lcm(0, n) == 0

// Works with any integer types; result type is the common type
std::gcd(12LL, 8);  // long long

`midpoint` and `lerp` (C++20)

cpp

Godbolt

// midpoint — exact, no overflow; also works for pointers and iterators
std::midpoint(3, 7);                          // 5
std::midpoint(INT_MAX, INT_MAX - 2);          // INT_MAX - 1 — no overflow
std::midpoint(3.0, 7.0);                      // 5.0

// lerp — exact at t=0 and t=1; extrapolation allowed
std::lerp(0.0, 10.0, 0.0);   // 0.0  — exactly a
std::lerp(0.0, 10.0, 1.0);   // 10.0 — exactly b
std::lerp(0.0, 10.0, 0.3);   // 3.0
std::lerp(0.0, 10.0, 1.5);   // 15.0 — extrapolation

The naive (a + b) / 2 overflows for large integers. The naive a + t * (b - a) loses precision at the boundaries when a and b differ by many ULPs. std::midpoint and std::lerp avoid both issues.

Reference

Function	Description	Since
`accumulate(b, e, init[, op])`	Left-to-right fold	C++98
`inner_product(b, e, b2, init[, op1, op2])`	Sequential generalized dot product	C++98
`partial_sum(b, e, out[, op])`	Sequential inclusive prefix fold	C++98
`adjacent_difference(b, e, out[, op])`	Consecutive differences	C++98
`iota(b, e, val)`	Fill with incrementing values	C++11
`reduce(b, e[, init[, op]])`	Unordered fold, parallelizable	C++17
`transform_reduce(b, e, b2, init, rop, top)`	Pairwise transform then fold	C++17
`transform_reduce(b, e, init, rop, uop)`	Unary transform then fold	C++17
`inclusive_scan(b, e, out[, op[, init]])`	Parallel-capable inclusive scan	C++17
`exclusive_scan(b, e, out, init[, op])`	Parallel-capable exclusive scan	C++17
`transform_inclusive_scan(b, e, out, op, uop)`	Transform then inclusive scan	C++17
`transform_exclusive_scan(b, e, out, init, op, uop)`	Transform then exclusive scan	C++17
`gcd(m, n)`	Greatest common divisor	C++17
`lcm(m, n)`	Least common multiple	C++17
`midpoint(a, b)`	Overflow-safe exact midpoint	C++20
`lerp(a, b, t)`	Linear interpolation with exact endpoints	C++20
`std::ranges::fold_left(r, init, op)`	Range-based left fold	C++23

Best Practices

Use reduce over accumulate for large numeric ranges. For integer sums and floating-point reductions where reordering is acceptable, reduce enables auto-vectorization and execution::par parallelism. Reserve accumulate when left-to-right order is semantically load-bearing — custom state machines, string building, observable side effects.

Prefer transform_reduce over accumulate with a transform. accumulate(b, e, init, [](auto acc, auto x){ return acc + f(x); }) forces sequential evaluation of both transform and reduction. transform_reduce separates the two operations, opening the door to parallelism and making intent explicit.

Enable parallel policies deliberately. On GCC and Clang with libstdc++, std::execution::par requires linking -ltbb. The code compiles without it, but the parallel policy may silently serialize or abort at runtime. Benchmark before assuming parallelism wins — for small ranges, the scheduling overhead dominates.

Use views::iota (C++20) for iteration. std::iota fills an already-allocated container. std::views::iota(0, n) generates values lazily at zero cost and composes with other range adaptors.

Common Pitfalls

reduce with floating-point gives non-reproducible results. Floating-point addition is not strictly associative — (a + b) + c ≠ a + (b + c) in general due to rounding. reduce may reorder operations even without a parallel policy. Use accumulate when reproducibility matters.

inclusive_scan and exclusive_scan have the same ordering constraints as reduce. The C++17 parallel scan algorithms require an associative (and for parallel policies, commutative) binary op. Using subtraction, division, or any non-associative operation produces unspecified results — not UB, but non-deterministic output.

transform_inclusive_scan and transform_exclusive_scan have different init positions. In transform_exclusive_scan, init is the third output argument (before binary_op). In transform_inclusive_scan, init is the last and optional argument. Swapping them silently compiles, producing wrong results.

lcm overflows for large inputs. If lcm(m, n) is not representable in the common type, the behavior is undefined — even if neither m nor n is large. lcm(46341, 46341) overflows int. Use long long or validate inputs.

gcd with INT_MIN is undefined behavior. The standard requires that |m| be representable as type M. Since |INT_MIN| exceeds INT_MAX, std::gcd(INT_MIN, x) is UB. Cast to unsigned or long long first.

std::reduce vs std::accumulate

std::accumulate is strictly sequential and left-to-right; std::reduce may reorder operations. For large numeric ranges where the operation is associative and commutative, reduce enables auto-vectorisation and explicit parallelism. When evaluation order is observable — string building, custom state, subtraction — only accumulate is correct.

cpp

Godbolt

#include <numeric>
#include <execution>
#include <vector>

std::vector<double> v = {1.0, 2.0, 3.0, 4.0, 5.0};

// accumulate — C++11, guaranteed left-to-right
double seq_sum = std::accumulate(v.begin(), v.end(), 0.0);   // 15.0

// reduce — C++17, unspecified evaluation order; correct only for assoc+comm ops
double red_sum = std::reduce(v.begin(), v.end(), 0.0);       // C++17

// Parallel sum — requires linking -ltbb on GCC/Clang with libstdc++
double par_sum = std::reduce(                                 // C++17
    std::execution::par_unseq, v.begin(), v.end(), 0.0);

// PITFALL: subtraction is neither commutative nor associative
// accumulate gives 1-(2-(3-(4-5))) evaluated strictly left-to-right: -13
double safe = std::accumulate(v.begin(), v.end(), 0.0,
                              std::minus<>{});   // C++11, deterministic
// reduce with minus: implementation-defined result — do NOT use
// double bad = std::reduce(v.begin(), v.end(), 0.0, std::minus<>{});  // C++17 UB-adjacent

std::reduce without an explicit execution policy may still reorder operations compared to accumulate. For floating-point sums where reproducibility matters, always use accumulate.

std::transform_reduce for Dot Products

std::transform_reduce fuses a per-element transform with an unordered reduction in a single pass, enabling both vectorisation and parallelism. The binary overload operates on two ranges (classic dot product); the unary overload applies a single transform before reducing.

cpp

Godbolt

#include <numeric>
#include <execution>
#include <vector>
#include <string>

std::vector<double> a = {1.0, 2.0, 3.0};
std::vector<double> b = {4.0, 5.0, 6.0};

// Dot product — binary form                                              // C++17
double dot = std::transform_reduce(
    a.begin(), a.end(),
    b.begin(),
    0.0,
    std::plus<>{},
    std::multiplies<>{}
);  // 1*4 + 2*5 + 3*6 = 32.0

// Parallel dot product                                                   // C++17
double par_dot = std::transform_reduce(
    std::execution::par_unseq,
    a.begin(), a.end(),
    b.begin(),
    0.0,
    std::plus<>{},
    std::multiplies<>{}
);

// Word-count example — unary form summing string lengths               // C++17
std::vector<std::string> words = {"hello", "world", "cpp"};
std::size_t total_chars = std::transform_reduce(
    words.begin(), words.end(),
    std::size_t{0},
    std::plus<>{},
    [](const std::string& s) { return s.size(); }
);  // 5 + 5 + 3 = 13

The manual equivalent — std::transform into a temporary buffer then std::accumulate — allocates memory and forces two passes. transform_reduce does both in a single pass with no intermediate allocation.

std::prefix_sum: inclusive_scan and exclusive_scan

std::inclusive_scan produces a prefix sum where output element i includes input element i; std::exclusive_scan excludes it, using an explicit init value as the first output. Both are parallelisable when the binary operation is associative.

cpp

Godbolt

#include <numeric>
#include <execution>
#include <vector>

std::vector<int> v = {1, 2, 3, 4};
std::vector<int> out(v.size());

// inclusive_scan — C++17: out[i] = v[0] + v[1] + ... + v[i]
std::inclusive_scan(v.begin(), v.end(), out.begin());
// out = {1, 3, 6, 10}

// exclusive_scan — C++17: out[i] = init + v[0] + ... + v[i-1]
std::exclusive_scan(v.begin(), v.end(), out.begin(), 0);
// out = {0, 1, 3, 6}

// Use case: offset table for variable-length records
// exclusive_scan gives the start offset of each record; last element + size gives end
std::vector<int> lengths = {3, 5, 2, 8};
std::vector<int> offsets(lengths.size());
std::exclusive_scan(lengths.begin(), lengths.end(), offsets.begin(), 0);
// offsets = {0, 3, 8, 10}  — byte offsets into a packed buffer

// Parallel inclusive scan                                               // C++17
std::inclusive_scan(std::execution::par,
                    v.begin(), v.end(), out.begin());

The key distinction: for input {1, 2, 3, 4}, inclusive scan yields {1, 3, 6, 10} while exclusive scan with init 0 yields {0, 1, 3, 6}. Exclusive scan is the natural choice when the output value at index i represents "everything before i" — offset tables, histogram cumulative counts, etc.

std::iota and Arithmetic Sequences

std::iota fills a range with consecutively incrementing values, using operator++ on the start value. It is the standard way to initialise index arrays or generate monotone test data without a hand-written loop.

cpp

Godbolt

#include <numeric>
#include <ranges>    // C++20
#include <vector>
#include <algorithm>

std::vector<int> idx(8);
std::iota(idx.begin(), idx.end(), 0);   // C++11: {0, 1, 2, 3, 4, 5, 6, 7}
std::iota(idx.begin(), idx.end(), 1);   // C++11: {1, 2, 3, 4, 5, 6, 7, 8}

// Lazy alternative — no allocation, composes with range adaptors          // C++20
for (int i : std::views::iota(0, 8)) { /* generates 0..7 on demand */ }

// Generate even numbers lazily                                             // C++20
auto evens = std::views::iota(0, 8)
           | std::views::transform([](int i) { return i * 2; });
// evens produces {0, 2, 4, 6, 8, 10, 12, 14} without allocation

// C++23: std::ranges::iota — range-based overload writing into a range
std::ranges::iota(idx, 0);              // C++23

// Sort-by-index idiom — sort indices by value in another array
std::vector<double> data = {3.1, 1.4, 2.7, 0.5};
std::iota(idx.begin(), idx.begin() + 4, 0);   // C++11
std::sort(idx.begin(), idx.begin() + 4,
          [&](int a, int b) { return data[a] < data[b]; });
// idx now holds the permutation that sorts data

Prefer std::views::iota (C++20) when only iteration is needed. Use std::iota when a materialised container of sequential values is explicitly required.

Parallel Numeric Algorithms (C++17)

All C++17 numeric algorithms that accept an ExecutionPolicy can run in parallel, but the binary operation must be both associative and commutative for par_unseq — violating this constraint gives implementation-defined (not undefined) results. Link against Intel TBB (-ltbb) on GCC/Clang to activate parallelism; without it the policy may silently degrade to sequential execution.

cpp

Godbolt

#include <numeric>
#include <execution>
#include <vector>

std::vector<double> v(1'000'000, 1.0);
std::vector<double> w(1'000'000, 2.0);

// Parallel sum                                                           // C++17
double sum = std::reduce(
    std::execution::par_unseq, v.begin(), v.end(), 0.0);

// Parallel dot product                                                   // C++17
double dot = std::transform_reduce(
    std::execution::par_unseq,
    v.begin(), v.end(),
    w.begin(),
    0.0,
    std::plus<>{},
    std::multiplies<>{}
);

// Parallel prefix sum                                                    // C++17
std::vector<double> prefix(v.size());
std::inclusive_scan(std::execution::par, v.begin(), v.end(), prefix.begin());

// PITFALL: non-commutative op with par_unseq — implementation-defined
// std::reduce(std::execution::par_unseq, v.begin(), v.end(), 0.0, std::minus<>{});

Parallel policies pay off at roughly N > 10,000 for CPU-bound operations on modern hardware. For smaller ranges the thread-scheduling and SIMD-setup overhead exceeds the computation cost. Always measure with realistic data before enabling par or par_unseq.

std::midpoint and std::lerp (C++20)

std::midpoint(a, b) computes the exact arithmetic midpoint without integer overflow; std::lerp(a, b, t) returns the exact linear interpolation with guaranteed exact results at t == 0 and t == 1. Both were introduced to replace common idioms that contain subtle numerical errors.

cpp

Godbolt

#include <numeric>  // midpoint, lerp

// midpoint — integer, no overflow                                        // C++20
int m = std::midpoint(INT_MAX - 1, INT_MAX);  // INT_MAX - 1 + 1 rounds down = INT_MAX - 1... wait
// Correct: (INT_MAX-1 + INT_MAX) / 2 via the safe formula, no intermediate overflow
int safe_mid = std::midpoint(2'000'000'000, 2'000'000'001);  // 2'000'000'000

// PITFALL: naive formula overflows
// int bad = (2'000'000'000 + 2'000'000'001) / 2;  // signed overflow — UB

// midpoint on floating-point                                             // C++20
double fmid = std::midpoint(0.1, 0.3);   // correctly rounded midpoint

// midpoint on pointers (iterators into the same array)                  // C++20
std::vector<int> v = {1, 2, 3, 4, 5};
auto* mid_ptr = std::midpoint(v.data(), v.data() + v.size());  // points to v[2]

// lerp — exact at endpoints, monotone                                   // C++20
double a = 0.0, b = 10.0;
std::lerp(a, b, 0.0);   // exactly 0.0
std::lerp(a, b, 1.0);   // exactly 10.0
std::lerp(a, b, 0.3);   // 3.0
std::lerp(a, b, 1.5);   // 15.0 — extrapolation is supported

// PITFALL: naive lerp has cancellation error near t=0 and t=1
// a + t * (b - a)  — when a and b differ by many ULPs, t*(b-a) loses low bits
// std::lerp avoids this via a branch: returns a when t==0, b when t==1

// Animation use case: smooth transition over 60 frames
for (int frame = 0; frame <= 60; ++frame) {
    double pos = std::lerp(0.0, 100.0, frame / 60.0);  // C++20, monotone, exact endpoints
}

std::lerp also handles extrapolation (t outside [0,1]) and NaN propagation correctly. The naive formula a + t * (b - a) is not monotone in general and does not guarantee lerp(a, b, 1.0) == b due to floating-point rounding, making std::lerp essential for animation, easing functions, and binary search midpoints.