Algorithms

Fenwick tree

Prefix sums over values that change: where a row starts and which row is at an offset, in O(log n).

Behaviour

A prefix sum is the sum of the values before an index

fc.assert(
  fc.property(heights, (values) => {
    const tree = fenwickFrom(values);
    return values.every((_, index) => fenwickPrefix(tree, index) === sum(values, index));
  }),
);

After any changes, prefix sums are those of the changed values

fc.assert(
  fc.property(
    heights.filter((values) => values.length > 0),
    fc.array(fc.tuple(fc.nat(), fc.integer({ min: -50, max: 50 })), { maxLength: 40 }),
    (values, changes) => {
      const current = [...values];
      const tree = fenwickFrom(values);
      for (const [at, delta] of changes) {
        const index = at % current.length;
        current[index]! += delta;
        fenwickAdd(tree, index, delta);
      }
      return current.every((_, index) => fenwickPrefix(tree, index + 1) === sum(current, index + 1));
    },
  ),
);

An offset is found in the item whose span contains it

fc.assert(
  fc.property(
    heights.filter((values) => values.length > 0),
    fc.nat(),
    (values, at) => {
      const tree = fenwickFrom(values);
      const offset = at % sum(values, values.length);
      const index = fenwickFind(tree, offset);
      return sum(values, index) <= offset && offset < sum(values, index + 1);
    },
  ),
);

Items of no height are passed over: the item found is the last to start at or before the offset

expect(fenwickFind(fenwickFrom([10, 0, 10]), 10)).toBe(2);
expect(fenwickFind(fenwickFrom([0, 0, 10]), 0)).toBe(2);
expect(fenwickFind(fenwickFrom([10, 0, 0, 10]), 9)).toBe(0);
fc.assert(
  fc.property(
    fc.array(fc.integer({ min: 0, max: 3 }), { minLength: 1, maxLength: 300 }),
    fc.nat(),
    (values, at) => {
      const total = sum(values, values.length);
      fc.pre(total > 0);
      const offset = at % total;
      const index = fenwickFind(fenwickFrom(values), offset);
      return values[index]! > 0 && sum(values, index) <= offset && offset < sum(values, index + 1);
    },
  ),
);

An offset at or past the end finds the number of items

const tree = fenwickFrom([10, 20, 30]);
expect(fenwickFind(tree, 60)).toBe(3);
expect(fenwickFind(tree, 1000)).toBe(3);
expect(fenwickFind(tree, Infinity)).toBe(3);
expect(fenwickFind(fenwickFrom([]), 0)).toBe(0);
expect(fenwickFind(fenwickFrom([10, 0, 0]), 10)).toBe(3);
expect(fenwickFind(fenwickFrom([0, 0]), 0)).toBe(2);

An offset before the start, as a list scrolled past its top reports, finds the first item

expect(fenwickFind(fenwickFrom([10, 20, 30]), -40)).toBe(0);
expect(fenwickFind(fenwickFrom([0, 20, 30]), -1)).toBe(0);

An offset that is not a number is an error

expect(() => fenwickFind(fenwickFrom([10, 20, 30]), Number.NaN)).toThrow(RangeError);

Changing an item that is not there is an error

const tree = fenwickFrom([10, 20, 30]);
for (const index of [3, 4, 0.5, Number.NaN, Infinity, -1]) {
  expect(() => fenwickAdd(tree, index, 1)).toThrow(RangeError);
}
expect(() => fenwickAdd(fenwickFrom([]), 0, 1)).toThrow(RangeError);
expect([...tree]).toEqual([...fenwickFrom([10, 20, 30])]);

A prefix of fewer than no items, of more than there are or of part of one is an error

const tree = fenwickFrom([10, 20, 30]);
expect(fenwickPrefix(tree, 0)).toBe(0);
expect(fenwickPrefix(tree, 3)).toBe(60);
for (const count of [-1, 4, 1.5, Number.NaN, Infinity]) {
  expect(() => fenwickPrefix(tree, count)).toThrow(RangeError);
}

Heights in halves, quarters or any binary fraction of a pixel stay exact however often they change

const sixtyFourths = fc.integer({ min: 0, max: 64 * 200 }).map((n) => n / 64);
fc.assert(
  fc.property(
    fc.array(sixtyFourths, { minLength: 1, maxLength: 300 }),
    fc.array(fc.tuple(fc.nat(), sixtyFourths), { maxLength: 200 }),
    (values, changes) => {
      const current = [...values];
      const tree = fenwickFrom(values);
      for (const [at, height] of changes) {
        const index = at % current.length;
        fenwickAdd(tree, index, height - current[index]!);
        current[index] = height;
      }
      // The same cells a tree built afresh would have, so no error has built up.
      return fenwickFrom(current).every((cell, index) => cell === tree[index]);
    },
  ),
);

With heights of any fraction, the offset at which an item starts is found in that item, or past it when it has no height

// Thousandths of a pixel: of a like size, and almost none of them exact in binary.
const thousandths = fc.integer({ min: 0, max: 200_000 }).map((n) => n / 1000);
fc.assert(
  fc.property(fc.array(thousandths, { minLength: 1, maxLength: 300 }), (values) => {
    const tree = fenwickFrom(values);
    return values.every((_, index) => {
      // Both add the same cells in the same order, so rounding cannot part them.
      const start = fenwickPrefix(tree, index);
      const found = fenwickFind(tree, start);
      return found >= index && fenwickPrefix(tree, found) <= start;
    });
  }),
);

An offset is found among any power of two of items, and one more or fewer

for (let power = 0; power <= 20; power += 1) {
  for (const count of [2 ** power - 1, 2 ** power, 2 ** power + 1]) {
    const tree = fenwickFrom(new Float64Array(count).fill(24));
    expect(fenwickFind(tree, count * 24 - 1)).toBe(Math.max(count - 1, 0));
    expect(fenwickFind(tree, count * 24)).toBe(count);
  }
}

A tree holds up to 2^31 − 1 items; building a larger one is an error

// Only the length is read before the refusal, so nothing this large is made.
expect(() => fenwickFrom({ length: 2 ** 31 })).toThrow("at most 2^31 − 1 values");
expect(() => fenwickFrom({ length: 2 ** 32 + 2 })).toThrow("at most 2^31 − 1 values");

Among a million items, each operation touches at most twenty-one cells

const tree = counted(fenwickFrom(new Float64Array(1_000_000).fill(24)));
for (const index of [0, 1, 524_287, 524_288, 999_999]) {
  tree.reset();
  fenwickAdd(tree.array, index, 8);
  // Each cell on the way up is read and written.
  expect(tree.accesses() / 2).toBeLessThanOrEqual(21);
  tree.reset();
  fenwickPrefix(tree.array, index + 1);
  expect(tree.accesses()).toBeLessThanOrEqual(21);
  tree.reset();
  fenwickFind(tree.array, index * 24);
  expect(tree.accesses()).toBeLessThanOrEqual(21);
}

React

Install

npx shadcn@latest add https://ui.publictheta.com/r/fenwick.json

Source

src/algorithms/fenwick.ts

/**
 * A Fenwick tree: prefix sums over a sequence of numbers that change.
 *
 * It answers the two questions a list of rows with differing heights asks on
 * every scroll — where does row i start, and which row is at offset y — and
 * absorbs a row changing height, each in O(log n), with nothing allocated.
 *
 * The tree is a plain `Float64Array` one longer than the sequence; these
 * functions are all there is to it. It holds up to 2^31 − 1 values, so that
 * every position is a 32-bit integer; sums are 64-bit floats and have no such
 * limit.
 *
 * Sums are exact while the values are integers, or multiples of one binary
 * fraction — halves, quarters, sixty-fourths of a pixel — and the total stays
 * under 2^53. Other fractions are rounded at each addition, so after many
 * changes a sum may differ in its last bits from one made afresh: round sizes
 * to such a fraction before they come here, or rebuild the tree now and then.
 * Either way `fenwickFind` and `fenwickPrefix` add the same cells in the same
 * order, and so agree with each other to the bit.
 *
 * A position that is not in the tree is a mistake in the caller, and throws a
 * `RangeError`.
 */

const MAX_COUNT = 2 ** 31 - 1;

/**
 * The largest power of two not above `count`, by its leading zeros: exact,
 * where a logarithm is rounded. For 0 the shift wraps to a negative number,
 * which ends the loops below before they begin.
 */
const highestBit = (count: number) => 1 << (31 - Math.clz32(count));

/** The tree for `values`, built in O(n). */
export function fenwickFrom(values: ArrayLike<number>): Float64Array {
  const count = values.length;
  if (count > MAX_COUNT) throw new RangeError(`A Fenwick tree holds at most 2^31 − 1 values, not ${count}.`);
  const tree = new Float64Array(count + 1);
  for (let i = 1; i <= count; i += 1) {
    tree[i]! += values[i - 1]!;
    const parent = i + (i & -i);
    if (parent <= count) tree[parent]! += tree[i]!;
  }
  return tree;
}

/** Add `delta` to the value at `index`. */
export function fenwickAdd(tree: Float64Array, index: number, delta: number): void {
  // One comparison refuses negatives, fractions and what is not a number.
  if (index >>> 0 !== index || index >= tree.length - 1) {
    throw new RangeError(`No value at index ${index} among ${tree.length - 1}.`);
  }
  for (let i = index + 1; i < tree.length; i += i & -i) tree[i]! += delta;
}

/** The sum of the first `count` values: where the item at index `count` starts. */
export function fenwickPrefix(tree: Float64Array, count: number): number {
  if (count >>> 0 !== count || count >= tree.length) {
    throw new RangeError(`No prefix of ${count} values among ${tree.length - 1}.`);
  }
  let sum = 0;
  let position = 0;
  // From the largest block down, as `fenwickFind` goes, so the two round alike.
  for (let step = highestBit(count); step > 0; step >>= 1) {
    if ((count & step) !== 0) {
      position += step;
      sum += tree[position]!;
    }
  }
  return sum;
}

/**
 * The index of the item that `offset` falls in, when no value is negative: the
 * largest count whose prefix sum does not exceed `offset`. So an item of no
 * size is passed over for the one after it, an offset at or past the total
 * gives the number of values, and an offset before the start gives 0.
 */
export function fenwickFind(tree: Float64Array, offset: number): number {
  if (offset !== offset) throw new RangeError("An offset must be a number.");
  const count = tree.length - 1;
  let position = 0;
  let sum = 0;
  // Descend from the largest power of two, taking each block that still fits.
  for (let step = highestBit(count); step > 0; step >>= 1) {
    const next = position + step;
    if (next > count) continue;
    const taken = sum + tree[next]!;
    if (taken <= offset) {
      position = next;
      sum = taken;
    }
  }
  return position;
}