2D prefix sums and 2D Fenwick trees

Why this matters

Many problems store numbers in a grid. A grid is a table with rows and columns, like a spreadsheet. A game map stores a cost for every tile. A photo stores a brightness for every pixel. A seat chart stores which seats are taken. A common question is the total inside a rectangle of the grid, for example the total brightness of one part of a photo.

The direct way is to add up every cell in the rectangle. If the grid has 1000 rows and 1000 columns, a rectangle can hold a million cells, and a program that asks a thousand such questions does a billion additions. This lesson shows two faster methods. The first is a two-dimensional prefix table. It answers any rectangle in four lookups, but it cannot change after it is built. The second is a two-dimensional Fenwick tree. It answers a rectangle in at most about 324 steps and takes a change to one cell in at most 100 steps.

The one-dimensional versions of these ideas are the prefix sums and the Fenwick tree from the earlier lessons of this module. This lesson extends them to rows and columns at the same time.

Words you will need

• Grid: a table of values with a fixed number of rows and columns. In C++ it can be a vector of vectors of numbers.
• Cell: one value in the grid, found by its row and its column.
• Rectangle query: a question of the form "what is the total of the cells from row r1 to row r2 and from column c1 to column c2".
• Inclusive bounds: both ends of the range are part of the range. The rows 1 to 2 are row 1 and row 2. This lesson uses inclusive bounds everywhere.
• Index from 0: the first row is row 0 and the first column is column 0.
• Prefix table: a table where each cell holds the total of everything above and to the left of it.
• Inclusion and exclusion: a way to count a region by adding a big area, taking away the parts that should not be there, and adding back any part that was taken away twice.
• Fenwick tree: a table that stores partial totals. It can add a change and read a total in a few steps each.
• Lowest set bit: in the binary form of a number, the value of its rightmost 1 bit. For the number 6, which is 110 in binary, the lowest set bit is 2, which is 010. In C++ the expression i & -i gives it.
• Brute force: the slow, obvious method. It is used to check the fast method.

Rectangles and coordinates

Each cell is named by two numbers, its row and its column. The cell in row 1 and column 2 is written (1, 2). A rectangle is named by its top left corner and its bottom right corner. The rectangle (1, 0) to (2, 1) covers rows 1 and 2 and columns 0 and 1. That is four cells: row 1 column 0, row 1 column 1, row 2 column 0 and row 2 column 1.

Every example in this lesson uses the grid below, which has 3 rows and 4 columns. The numbers are chosen so that each sum is easy to check by hand.

     col 0  col 1  col 2  col 3
row 0   1      2      3      4
row 1   5      6      7      8
row 2   9     10     11     12

The rectangle from (1, 0) to (2, 1) holds 5, 6, 9 and 10. Their total is 30.

Part 1: the two-dimensional prefix table

The idea: a total from the corner

Start with the one-dimensional case. The prefix sum at position i is the total of all the values from position 0 to position i. Once the prefix sums are stored, the total from l to r is the prefix at r minus the prefix at l - 1.

In two dimensions the same idea uses a corner. Define P[i][j] as the total of all the cells with row less than i and column less than j. The table has one extra row and one extra column, so that P[0][j] and P[i][0] are both 0. These zero entries stand for "nothing above" and "nothing to the left". They let the formula work without special cases at the edges.

Here is the table for the grid above. Row 0 and column 0 of the table are all zeros, and the rest is filled in by the rule in the next part.

          0    1    2    3    4      <- column count, j
    0     0    0    0    0    0
    1     0    1    3    6   10
    2     0    6   14   24   36
    3     0   15   33   54   78

Read one cell to check it. The entry P[2][3] is 24. It is the total of rows 0 and 1 and columns 0, 1 and 2. That is 1 + 2 + 3 + 5 + 6 + 7, which is 24. The entry P[3][4] is 78, the total of the whole grid.

Building the table

Each cell of the table is built from three cells that are already known. Take the cell P[i + 1][j + 1]. Its region is every grid cell with row at most i and column at most j. That region is the region of P[i][j + 1] (rows above) plus the region of P[i + 1][j] (columns to the left). The two regions overlap in the region of P[i][j], so that overlap is counted twice and must be taken away once. The new grid value is then added in the corner.

P[i+1][j+1] = grid[i][j] + P[i][j+1] + P[i+1][j] - P[i][j]

Check it on the cell P[2][2]. The grid value at (1, 1) is 6. The table value P[1][2] is 3, and P[2][1] is 6, and P[1][1] is 1. So the entry is 6 + 3 + 6 - 1 = 14, which matches the table.

Building the table touches every cell once. For a grid with R rows and C columns, the cost is about R times C steps.

Answering a rectangle with four lookups

Now the rectangle sum. Suppose the question is the rectangle with rows r1 to r2 and columns c1 to c2. The table gives the total of everything from row 0 down to row r2, and from column 0 across to column c2. That region is the big block P[r2 + 1][c2 + 1].

The rectangle is not the whole block. It is missing three parts. The top strip, rows 0 to r1 - 1, is above it. The left strip, columns 0 to c1 - 1, is to its left. The top-left corner, rows 0 to r1 - 1 and columns 0 to c1 - 1, lies inside both strips, so it has been counted twice.

So the answer is the big block, minus the top strip, minus the left strip, plus the corner that was taken away twice. In table entries the rule is:

sum = P[r2 + 1][c2 + 1] - P[r1][c2 + 1] - P[r2 + 1][c1] + P[r1][c1]

Use it on the rectangle from (1, 0) to (2, 1). Here r1 = 1, c1 = 0, r2 = 2, c2 = 1. The four lookups are P[3][2] = 33, P[1][2] = 3, P[3][0] = 0, and P[1][0] = 0. The sum is 33 - 3 - 0 + 0 = 30, which matches the hand count.

The formula always uses four lookups, whatever the size of the rectangle. The rectangle that covers the whole grid uses P[R][C] and the three zero entries, so it returns the grand total.

Worked example 1: a prefix table with three queries

This program reads a grid, builds the table, prints it, and then answers a list of rectangle queries. The class stores the table, and the method rectSum applies the four-lookup formula.

Input:

3 4
1 2 3 4
5 6 7 8
9 10 11 12
3
0 0 2 3
1 1 2 2
2 0 2 0

Output:

table:
0 0 0 0 0
0 1 3 6 10
0 6 14 24 36
0 15 33 54 78
sum(0,0)-(2,3) = 78
sum(1,1)-(2,2) = 34
sum(2,0)-(2,0) = 9

The table printed by the program matches the hand-built table above. The first query covers the whole grid and returns 78. The second covers the four cells 6, 7, 10 and 11, and returns 34. The third covers one cell, the value 9 in row 2 and column 0, and returns 9. A single cell is a rectangle with the same top left and bottom right corner, so the formula handles it without a special case.

The limit of the prefix table

The table answers queries fast, but it is fixed. If the cell at (1, 2) changes from 7 to 100, every table entry at or below and to the right of that cell is wrong. In the grid of this example, that is 4 of the 20 table entries. A bigger grid is worse. A change near the top left corner of a 1000 by 1000 grid changes about a million table entries. So the prefix table suits a grid that is built once and queried many times. A grid with frequent changes needs the Fenwick tree in Part 2.

Part 2: the two-dimensional Fenwick tree

The one-dimensional idea, briefly

A one-dimensional Fenwick tree stores a partial total at each index. Index i stores the total of the values from i - lowbit(i) + 1 up to i, where lowbit(i) is the lowest set bit of i. The lowest set bit of 6 is 2, so index 6 stores the total of positions 5 and 6. The lowest set bit of 8 is 8, so index 8 stores the total of positions 1 to 8.

Two loops over these indexes do all the work. To read a prefix total up to index i, start at i, add the stored value, then remove the lowest set bit with i -= i & -i. To add a change at position i, start at i, add the change, then move up with i += i & -i. Each loop takes about log2 of the size steps, because each step clears or carries one bit.

Two loops, one inside the other

A two-dimensional Fenwick tree does the same thing for rows and for columns. The row loop runs over one set of indexes and the column loop runs over another. For each row index it runs the full column loop. The stored value at (i, j) is the total of a block of cells whose size comes from the lowest set bits of i and j. The program does not need to know the block sizes, because the two loops visit exactly the blocks it needs.

Indexes start at 1 inside the tree, because the rule i & -i gives 0 for index 0 and the loop would never end. So a cell at row r and column c is stored at tree index (r + 1, c + 1). The class below hides this shift. Outside the class, all indexes are 0-based.

The update method is:

void add(long long r, long long c, long long delta) {
    for (long long i = r + 1; i <= rows_; i += i & -i) {
        for (long long j = c + 1; j <= cols_; j += j & -j) {
            tree_[i][j] += delta;
        }
    }
}

The prefix method, which gives the total of all cells with row at most r and column at most c, is:

long long prefix(long long r, long long c) const {
    long long total = 0;
    for (long long i = r + 1; i > 0; i -= i & -i) {
        for (long long j = c + 1; j > 0; j -= j & -j) {
            total += tree_[i][j];
        }
    }
    return total;
}

If r or c is -1, the loop for that index never starts, and prefix returns 0. This is what makes the rectangle formula safe at the top and left edges.

The rectangle formula stays the same

The rectangle sum uses the same inclusion and exclusion as the prefix table. The four prefix calls replace the four table lookups:

long long rectSum(long long r1, long long c1, long long r2, long long c2) const {
    return prefix(r2, c2) - prefix(r1 - 1, c2) - prefix(r2, c1 - 1) + prefix(r1 - 1, c1 - 1);
}

Each prefix call reads at most 81 cells for a grid of 1000 by 1000. So a rectangle sum costs about four times that, and an update costs about the same. Neither one copies the grid.

Worked example 2: updates and queries

This program keeps a 2D Fenwick tree for the grid from before. It reads commands. A line a r c d adds d to the cell at row r and column c. A line q r1 c1 r2 c2 prints the rectangle sum. The grid starts with all zeros, and the first twelve commands fill it in.

Input. The first twelve lines set the cells to the grid above. The rest are queries and changes:

3 4
20
a 0 0 1
a 0 1 2
a 0 2 3
a 0 3 4
a 1 0 5
a 1 1 6
a 1 2 7
a 1 3 8
a 2 0 9
a 2 1 10
a 2 2 11
a 2 3 12
q 0 0 2 3
q 1 1 2 2
a 1 2 -7
q 1 1 2 2
q 0 2 2 2
a 0 0 10
q 0 0 0 0
q 0 0 2 3

Output:

sum(0,0)-(2,3) = 78
sum(1,1)-(2,2) = 34
sum(1,1)-(2,2) = 27
sum(0,2)-(2,2) = 14
sum(0,0)-(0,0) = 11
sum(0,0)-(2,3) = 81

Check the lines by hand. The first two match the prefix table example. Then the command a 1 2 -7 takes 7 away from the cell at row 1 and column 2, which leaves 0 there. The rectangle from (1, 1) to (2, 2) now holds 6, 0, 10 and 11, a total of 27. The column 2 query covers 3, 0 and 11, which is 14. The command a 0 0 10 raises the top left cell from 1 to 11, and the whole grid total goes up by 10 and down by 7, so 78 + 10 - 7 = 81.

Notice the input has 20 commands, and the second line of input says so. The program reads the count first, as the earlier examples do.

Worked example 3: a check against brute force, and the cost

A new method needs a test. This example uses the same class, with main replaced. It keeps a plain copy of the grid, runs random changes and random queries against both, and counts how often the two answers differ. The random numbers come from a formula with a fixed seed, so every run is the same.

int main() {
    const long long R = 40, C = 50;
    Fenwick2D fen(R, C);
    vector<vector<long long>> grid(R, vector<long long>(C, 0));
    unsigned long long seed = 7;
    auto next = [&seed](long long mod) {
        seed = seed * 6364136223846793005ULL + 1442695040888963407ULL;
        return static_cast<long long>((seed >> 33) % static_cast<unsigned long long>(mod));
    };
    int checks = 0, mismatches = 0;
    for (int step = 0; step < 3000; step++) {
        if (next(2) == 0) {
            long long r = next(R), c = next(C), d = next(201) - 100;
            fen.add(r, c, d);
            grid[r][c] += d;
        } else {
            long long r1 = next(R), r2 = next(R), c1 = next(C), c2 = next(C);
            if (r1 > r2) swap(r1, r2);
            if (c1 > c2) swap(c1, c2);
            long long brute = 0;
            for (long long r = r1; r <= r2; r++)
                for (long long c = c1; c <= c2; c++) brute += grid[r][c];
            checks++;
            if (brute != fen.rectSum(r1, c1, r2, c2)) mismatches++;
        }
    }
    cout << "checks: " << checks << ", mismatches: " << mismatches << "\n";
    return 0;
}

Output:

checks: 1454, mismatches: 0

Zero mismatches in 1454 checks is good evidence that the four-prefix formula and the two loops are right. It is not a proof. A test like this should use grids of several shapes, including one row and one column, before the code is trusted.

The cost can be measured too. A variant of the same loops that counts how many cells one call visits prints these numbers for a 1000 by 1000 grid:

add at (0,0) on 1000x1000 touches 100 cells
add at (500,500) on 1000x1000 touches 16 cells
prefix at (999,999) reads 36 cells
table of 1001x1001 longs: 8016008 bytes

An update at the top left corner touches 100 cells, which is 10 row indexes times 10 column indexes. An update at the bottom right corner touches only one cell, because the first index it reaches is already past 1000. The prefix read at the bottom right corner reads 36 cells. The table of the same grid takes about 8 MB, and the Fenwick tree stores the same number of entries.

Common mistakes

Mistake 1: the wrong corner in the formula

The most common error is to subtract the wrong prefix. The formula needs r1 - 1 and c1 - 1, not r1 and c1. If the program uses r1 where it should use r1 - 1, it includes or removes an extra row or column.

This is the bad version of rectSum:

return prefix(r2, c2) - prefix(r1, c2) - prefix(r2, c1 - 1) + prefix(r1, c1 - 1);

Run on the input of worked example 2, the bad version prints this.

sum(0,0)-(2,3) = 68
sum(1,1)-(2,2) = 21
sum(1,1)-(2,2) = 21
sum(0,2)-(2,2) = 11
sum(0,0)-(0,0) = 0
sum(0,0)-(2,3) = 61

The first answer is 68, not 78. The bad formula subtracts row 0, whose total is 10, from the whole grid, so the answer is too small by 10. No error message is printed. A test with a small grid and a hand-counted answer catches this at once.

Mistake 2: an index outside the grid is ignored

The tree has no check on its indexes. If a change is given a row past the last row, the first loop never starts, so the change is dropped without a message. Take a grid with 3 rows and add a value to row 3.

3 4
2
a 3 0 5
q 0 0 2 3

Output:

sum(0,0)-(2,3) = 0

The grid has no row 3, so the value 5 went nowhere. The program gives no error, and the answer is 0. A safe program checks 0 <= r < rows and 0 <= c < cols before calling add.

Mistake 3: the table is not updated

The prefix table from Part 1 is built once. If the program changes one grid cell and does not rebuild the table, every later query reads the old totals. The answer is wrong and the program gives no warning. Either rebuild the table after each change, which costs R times C steps each time, or use the Fenwick tree, which costs about a few dozen steps per change.

Mistake 4: passing the wrong number of arguments

A query needs four numbers, the two corners. Passing only three gives a compile error. This is the error g++ prints for a class with one rectSum method that takes four arguments.

l13_err.cpp:10:22: error: no matching function for call to 'Fenwick2D::rectSum(int, int, int)'
   10 |     cout << f.rectSum(0, 0, 2) << "\n";
      |             ~~~~~~~~~^~~~~~~~~
l13_err.cpp:6:15: note: candidate 1: 'long long int Fenwick2D::rectSum(long long int, long long int, long long int, long long int) const'

The note tells you the method takes four long long values. The fix is to pass the column of the bottom right corner as well, as in f.rectSum(0, 0, 2, 3).

Practice

1. On the grid of this lesson, the rectangle from (1, 0) to (2, 1) sums to 30. Compute the same sum from the prefix table using the four-lookup formula. Write down the four table values you use.

Show the answer

The four values are P[3][2] = 33, P[1][2] = 3, P[3][0] = 0 and P[1][0] = 0. The sum is 33 - 3 - 0 + 0 = 30, which matches the hand count of 5 + 6 + 9 + 10.

2. A query asks for prefix(-1, 5) on a Fenwick tree. What does it return, and which part of the code makes that happen?

Show the answer

It returns 0. The outer loop starts with i = r + 1, which is 0 when r is -1. The condition i > 0 is false from the start, so the loop body never runs and total stays 0.

3. A grid has 1000 rows. For an update at row 0, how many row indexes does the outer loop visit? Use the loop for (i = r + 1; i <= rows; i += i & -i).

Show the answer

The loop visits the indexes 1, 2, 4, 8, 16, 32, 64, 128, 256 and 512. That is 10 indexes. The next index would be 1024, which is above 1000, so the loop stops.

4. Why does the prefix table have one extra row and one extra column? Explain in two sentences what would go wrong without them.

Show the answer

The extra row and column hold zeros for "nothing above" and "nothing to the left". Without them, a rectangle that starts at row 0 or column 0 would need the table at index -1, which does not exist, so the formula would need a special case at every edge.

5. A 1000 by 1000 grid receives one change between queries, and the program runs 100000 changes and queries in turn. Compare the cost of rebuilding a prefix table after each change with the cost of a Fenwick update. Use the cell counts from worked example 3 and give a rough total for each approach.

Show the answer

Rebuilding the table after each change writes about 1000000 cells, so 100000 changes write about 100000000000 cells. A Fenwick update touches at most 100 cells on a 1000 by 1000 grid, so 100000 updates touch at most about 10000000 cells. The Fenwick approach is much cheaper for a grid that changes often.

Summary

• A prefix table stores the total of everything above and to the left of each cell, with an extra row and column of zeros.
• A rectangle sum is the big block minus the strip above, minus the strip to the left, plus the corner that was taken away twice.
• Building a prefix table costs one step per cell, and it must be rebuilt after any change.
• A 2D Fenwick tree uses two nested loops over the lowest-bit steps, so a change takes at most 100 steps and a query at most about 324 for a 1000 by 1000 grid.
• Use inclusive bounds and check that the rectangle formula handles row 0 and column 0.
• Test a fast method against a plain brute force on random inputs before you trust it.