Back to Number Theory Dsa

Range Arithmetic Query Engine

September 23, 2026 Wasil Zafar 13 min read

Build a query engine around arithmetic monoids and precomputation.

Table of Contents

  1. Intuition, Vocabulary and Foundations
  2. Step-by-Step Derivations and Algorithms
  3. Implementation, Complexity and Engineering
  4. Applications, Case Studies and Exercises

Intuition, Vocabulary and Foundations

A range arithmetic query engine stores an array while answering questions about contiguous slices. The central trick is to retain only summaries that can be merged correctly. For example, in [12, 18, 30, 25], the gcd of indices 0 through 2 is $(6)$, while exactly three values in the full range are divisible by $(5)$. Gcd is naturally mergeable because $(gcd(gcd(A),gcd(B))=gcd(A\cup B))$; divisor-frequency queries instead need counts indexed by divisors.

Think of the engine as a warehouse. A segment tree combines compact labels from neighboring shelves, while sieve-derived metadata records which divisibility categories each item belongs to. Immutable workloads can be prepared aggressively with sparse tables or offline sorting. If updates occur, every changed value must update each summary that depends on it.

An lcm constraint requires particular care. Although lcm is associative, values can grow beyond the numeric type even when all inputs fit. Many applications need only a capped lcm, a divisibility predicate, or prime-exponent metadata rather than the exact integer. The API must state which interpretation it implements.

ConceptConcrete exampleNon-example or boundaryAssumption and question answered
Range gcd$(gcd(12,18,30)=6)$Adding child gcds does not produce a range gcdAssociative merge with identity $(0)$: what divides every value in the interval?
Range lcm$(operatorname{lcm}(4,6,10)=60)$Multiplying all values can overflow and is not the lcmUse gcd-based merging with a stated cap or wide integer: what is the least common multiple?
Divisor frequencyThree values in [10, 15, 20, 7] are divisible by $(5)$Counting only values equal to $(5)$For positive $(d)$, count indices satisfying $(d\mid a_i)$
Sieve metadataSPF data yields divisors of $(60)$ from $(2^2\cdot3\cdot5)$Trial-dividing every value for every queryValues lie within a declared preprocessing bound
Point updateReplacing 15 by 14 changes gcd and divisor countsUpdating the array but not its index structuresAssignment semantics: after set(i, x), all later queries observe $(x)$
Online workloadEach query may depend on earlier updatesReordering dependent operations for batchingPreserve operation order and support queries immediately
Offline workloadSort independent threshold queries by one endpointReordering updates without a time dimensionAll operations are known and may be reordered without changing meaning
Zero valueEvery positive divisor divides $(0)$, and $(gcd(0,x)=\lvert x\rvert)$Trying to enumerate all divisors of zeroTrack zeros separately in divisor-frequency summaries

Step-by-Step Derivations and Algorithms

Start with the baseline: scan each requested interval. A gcd or divisibility-count query then costs $(O(r-l+1))$, but this direct method is an excellent testing oracle. For static range gcd, idempotence permits a sparse table: precompute gcds of power-of-two intervals in $(O(n\log n))$ time and space, then answer a query from two overlapping blocks in $(O(1))$. Once point updates are allowed, a segment tree is the safer abstraction. Each node stores the gcd of its interval, so changing one leaf requires recomputing only its $(O(\log n))$ ancestors.

Divisor-frequency queries need a richer structure. Precompute smallest prime factors through a maximum value $(V)$, factor each nonzero array value, and enumerate its positive divisors. With square-root decomposition, each block stores a counter mapping divisor $(d)$ to the number of nonzero values divisible by $(d)$, plus a separate zero count. A query scans at most two partial blocks and reads counters for every complete block. If the block size is $(B)$, query time is $(O(B+n/B))$; choosing $(B)$ near $(sqrt n)$ balances the terms. A point update rebuilds one block, costing $(O(B))$ plus divisor-generation work.

For a static batch asking how many values are divisible by each $(d)$, build a global value-frequency array and compute $(F(d)=\sum_{k\ge1}\operatorname{freq}(kd))$. This divisor transform costs $(O(V\log V))$ and answers each whole-array query in $(O(1))$, but it cannot represent arbitrary subranges or interleaved updates without additional structure.

Range Arithmetic Query Engine algorithm map
flowchart TD
    A[Classify workload and query semantics] --> B{Are updates present?}
    B -->|No| C{Query family}
    C -->|Range gcd| D[Sparse table]
    C -->|Whole-array divisor counts| E[Offline divisor transform]
    C -->|Mixed subranges| F[Block metadata]
    B -->|Yes| G[Segment tree for gcd]
    G --> H[SPF-derived divisors]
    H --> I[Per-block divisor counters]
    I --> J[Update leaf and rebuild one block]

The segment-tree invariant is that every node equals the gcd of exactly its covered interval; the block invariant is that each counter equals the number of nonzero block values divisible by its key, with zeros counted separately. Construction uses $(O(n))$ segment-tree space plus stored divisor incidences, gcd operations take $(O(\log n))$, and block divisibility queries and updates take $(O(B+n/B))$ and $(O(B))$, respectively.

Implementation, Complexity and Engineering

The reference engine supports assignment updates, half-open ranges [left, right), range gcd, and counts of values divisible by a positive integer. It fixes a maximum absolute value at construction: updates outside that bound are rejected rather than silently using incomplete sieve metadata. Negative values are interpreted through their absolute values, and zero is tracked separately because it is divisible by every positive query divisor.

The segment tree and block index deliberately use independent update paths. This makes the contract visible and permits differential tests against brute force. In a performance-oriented implementation, divisor lists should be cached, counters should use compact arrays when the value bound is modest, and block size should be benchmarked rather than assumed to equal an exact square root.

from collections import Counter
from math import gcd, isqrt
from random import Random


def smallest_prime_factors(limit: int) -> list[int]:
    if limit < 1:
        raise ValueError("maximum absolute value must be positive")
    spf = list(range(limit + 1))
    spf[0] = 0
    spf[1] = 1
    for prime in range(2, isqrt(limit) + 1):
        if spf[prime] == prime:
            for multiple in range(prime * prime, limit + 1, prime):
                if spf[multiple] == multiple:
                    spf[multiple] = prime
    return spf


class RangeArithmeticEngine:
    def __init__(self, values: list[int], max_abs_value: int | None = None):
        if not values:
            raise ValueError("values must be nonempty")
        observed = max(abs(value) for value in values)
        self.max_value = observed if max_abs_value is None else max_abs_value
        if self.max_value < max(1, observed):
            raise ValueError("maximum bound does not cover the input")
        self.values = list(values)
        self.n = len(values)
        self.spf = smallest_prime_factors(self.max_value)
        self.divisor_cache: dict[int, tuple[int, ...]] = {1: (1,)}
        self.size = 1
        while self.size < self.n:
            self.size *= 2
        self.tree = [0] * (2 * self.size)
        for index, value in enumerate(self.values):
            self.tree[self.size + index] = abs(value)
        for node in range(self.size - 1, 0, -1):
            self.tree[node] = gcd(self.tree[2 * node], self.tree[2 * node + 1])
        self.block_size = max(1, isqrt(self.n))
        self.block_count = (self.n + self.block_size - 1) // self.block_size
        self.block_divisors = [Counter() for _ in range(self.block_count)]
        self.block_zeros = [0] * self.block_count
        for block in range(self.block_count):
            self._rebuild_block(block)

    def _check_range(self, left: int, right: int) -> None:
        if not (0 <= left < right <= self.n):
            raise ValueError("range must satisfy 0 <= left < right <= n")

    def _divisors(self, value: int) -> tuple[int, ...]:
        value = abs(value)
        if value == 0:
            raise ValueError("zero has no finite divisor list")
        cached = self.divisor_cache.get(value)
        if cached is not None:
            return cached
        remaining = value
        factors: list[tuple[int, int]] = []
        while remaining > 1:
            prime = self.spf[remaining]
            exponent = 0
            while remaining % prime == 0:
                remaining //= prime
                exponent += 1
            factors.append((prime, exponent))
        result = [1]
        for prime, exponent in factors:
            existing = tuple(result)
            power = 1
            for _ in range(exponent):
                power *= prime
                result.extend(divisor * power for divisor in existing)
        answer = tuple(result)
        self.divisor_cache[value] = answer
        return answer

    def _rebuild_block(self, block: int) -> None:
        start = block * self.block_size
        stop = min(self.n, start + self.block_size)
        counts = Counter()
        zeros = 0
        for value in self.values[start:stop]:
            if value == 0:
                zeros += 1
            else:
                counts.update(self._divisors(value))
        self.block_divisors[block] = counts
        self.block_zeros[block] = zeros

    def update(self, index: int, value: int) -> None:
        if not 0 <= index < self.n:
            raise IndexError("update index out of range")
        if abs(value) > self.max_value:
            raise ValueError("updated value exceeds the preprocessing bound")
        self.values[index] = value
        node = self.size + index
        self.tree[node] = abs(value)
        node //= 2
        while node:
            self.tree[node] = gcd(self.tree[2 * node], self.tree[2 * node + 1])
            node //= 2
        self._rebuild_block(index // self.block_size)

    def range_gcd(self, left: int, right: int) -> int:
        self._check_range(left, right)
        left += self.size
        right += self.size
        result_left = result_right = 0
        while left < right:
            if left & 1:
                result_left = gcd(result_left, self.tree[left])
                left += 1
            if right & 1:
                right -= 1
                result_right = gcd(self.tree[right], result_right)
            left //= 2
            right //= 2
        return gcd(result_left, result_right)

    def count_divisible(self, left: int, right: int, divisor: int) -> int:
        self._check_range(left, right)
        if divisor <= 0:
            raise ValueError("divisor must be positive")
        if divisor > self.max_value:
            return sum(value == 0 for value in self.values[left:right])
        total = 0
        while left < right and left % self.block_size:
            total += self.values[left] % divisor == 0
            left += 1
        while left + self.block_size <= right:
            block = left // self.block_size
            total += self.block_divisors[block][divisor]
            total += self.block_zeros[block]
            left += self.block_size
        while left < right:
            total += self.values[left] % divisor == 0
            left += 1
        return total


engine = RangeArithmeticEngine([12, 18, 30, 25, 0], max_abs_value=100)
assert engine.range_gcd(0, 3) == 6
assert engine.count_divisible(0, 5, 5) == 3
engine.update(1, 20)
assert engine.range_gcd(0, 3) == 2
assert engine.count_divisible(0, 5, 5) == 4
assert engine.count_divisible(4, 5, 97) == 1

rng = Random(7)
values = [rng.randint(-50, 50) for _ in range(40)]
tested = RangeArithmeticEngine(values, max_abs_value=100)
for _ in range(500):
    if rng.randrange(3) == 0:
        index = rng.randrange(len(values))
        value = rng.randint(-100, 100)
        values[index] = value
        tested.update(index, value)
    else:
        left = rng.randrange(len(values))
        right = rng.randrange(left + 1, len(values) + 1)
        divisor = rng.randint(1, 120)
        expected_gcd = 0
        for value in values[left:right]:
            expected_gcd = gcd(expected_gcd, abs(value))
        expected_count = sum(value % divisor == 0 for value in values[left:right])
        assert tested.range_gcd(left, right) == expected_gcd
        assert tested.count_divisible(left, right, divisor) == expected_count

for operation in (
    lambda: tested.range_gcd(3, 3),
    lambda: tested.count_divisible(0, 1, 0),
    lambda: tested.update(0, 101),
):
    try:
        operation()
        raise AssertionError("invalid input was accepted")
    except ValueError:
        pass

print("trace:", engine.range_gcd(0, 5), engine.count_divisible(0, 5, 5))
print("all tests passed")
MethodBuild, query, and update costsBest workloadFailure mode and mitigation
Brute-force scanBuild $(O(1))$; query $(O(k))$; update $(O(1))$Small inputs and testing oracleToo slow on long ranges; retain it for differential tests
Gcd sparse tableBuild $(O(n\log n))$; query $(O(1))$; no efficient updatesImmutable arrays with many gcd queriesOverlapping-block logic works because gcd is idempotent; do not generalize it blindly
Gcd segment treeBuild $(O(n))$; query and update $(O(\log n))$Online mixed gcd workloadsWrong identity or stale ancestors; use identity $(0)$ and recompute to the root
Block divisor countersBuild proportional to stored divisor incidences; query $(O(B+n/B))$; update $(O(B))$Moderate online divisibility workloadsZero has infinitely many divisors; track zero counts separately
Offline divisor transform$(O(V\log V))$ build; $(O(1))$ whole-array queryStatic global frequency questionsDoes not answer arbitrary subranges; select only when the query shape matches
Exact lcm treeQuery and update $(O(\log n))$, but integers may grow rapidlySmall values or arbitrary-precision environmentsOverflow; use checked multiplication or cap at a documented sentinel
Mo-style orderingRoughly $(O((n+q)\sqrt n))$ add/remove workStatic offline queries with reversible frequency stateInterleaved updates need a time dimension; otherwise preserve online order

Applications, Case Studies and Exercises

In an approachable analytics task, an array represents batch sizes. Users ask for the gcd of a range and how many entries are divisible by a requested machine capacity. A segment tree plus block divisor counters answers both while supporting corrections to individual entries. The brute-force version remains small enough to serve as a trusted oracle.

A professionally realistic workload might contain millions of ledger-derived quantities, a mostly static daily snapshot, occasional corrections, and several query families. The engine should classify queries before choosing structures: sparse tables for immutable gcd snapshots, offline divisor transforms for whole-dataset reports, and an update-aware tree or block index only where freshness requires it. Memory, preprocessing latency, update latency, and numeric limits are part of the public service contract.

Exercises can progress from recalling the gcd identity and half-open range convention, through proving segment-tree and block-counter invariants, to implementing capped lcm queries. Debugging tasks should target forgotten zero counts, inclusive-versus-exclusive endpoint errors, stale metadata after updates, and sieve-bound violations. An extension task can compare square-root blocks, a divisor-indexed Fenwick design for restricted queries, and Mo's algorithm on the same recorded workload.

Deliverable

Apply Range Arithmetic Query Engine

Case study 1 — trace a mixed workload. Start with [12, 18, 30, 25, 0] and preprocessing bound 100. Query gcd on [0, 3), count multiples of $(5)$ on [0, 5), assign index 1 the value 20, and repeat both queries. Checkpoints: the gcd changes from $(6)$ to $(2)$, while the divisibility count changes from $(3)$ to $(4)$. Confirm that zero contributes to every positive-divisor query.

Case study 2 — select a production design. Generate a reproducible array of 200,000 values bounded by $(10^6)$ and a trace containing 70 percent gcd queries, 25 percent divisibility queries, and 5 percent updates. Benchmark brute force on a sampled prefix, then compare block sizes $(256)$, $(448)$, and $(700)$. Record build time, median and p95 latency by operation, peak memory, and oracle mismatches. Success requires zero mismatches and a documented structure choice justified by measurements.

Hands-on deliverable. Expose gcd left right, count-divisible left right d, and set index value commands using half-open ranges and assignment updates. Add tests for single-element ranges, all-zero ranges, negative values, divisors above the sieve bound, repeated updates, and invalid indices. The engine passes when 10,000 randomized operations match brute force, its memory fits a stated budget, and its README documents semantics and asymptotic costs.

Choose by workload: scan for tiny instances, use a sparse table for immutable gcd queries, a segment tree for updates, divisor transforms for static global counts, and block or offline techniques for richer frequency state. Common mistakes are treating every arithmetic statistic as a monoid, ignoring zero and overflow semantics, or updating only one index; the prerequisite material on arithmetic functions and Möbius-based counting provides the next tools for more advanced divisor summaries.