[Data Structure]Binary Indexed Tree/Fenwick Tree

本文主要参考topcoder讲解BIT(Binary Index Tree)的文章，图片均来自链接的文章。

Binary Index Tree/Fenwick Tree是一种应用在区间查询的数据结构。其基本的思想是任何整数都可以分解成2的指数为基的线性组合；以此类推，任何一个区间也可以分解成不同的长度为2的指数的区间的并集。那么BIT是如何运作的呢？我们以数组arr(长度为N)为例，所求的为区间和，对应的BIT也是一个长度为N的数组tree，对于index = i，tree[i]存的是：

• sum(arr[i - 2^j + 1]: arr[i])，其中j为i的二进制表示从右到左第一个1的位置(从0开始)
• 换句话说，如果i的二进制表示为a1b，其中1为从右到左第一个1，b全为0。那么tree[i]存的是从arr[a0b + 1]到arr[a1b]的和
• 以12为例，其二进制为1100；去掉最低位的1为1000，十进制为8。那么tree[12] = arr[12] + arr[11] + arr[10] + arr[9]

如图所示，一个长度为16的数组对应的BIT数组每个索引所存的区间和如下图所示。注意BIT的数组下标从1开始，所以长度比原数组多1：

如果我们用多叉树的形式表示，如下图所示：

这个树具有以下特点：

• 每个节点所表示的区间依然如上述，即节点i表示[i - 2^j + 1, i]的区间, j为LSB(lowest significant bit)所在的位置
• 每个节点n的所有子节点表示的区间的是连续相邻的区间，并且这些区间的并集为和n表示的区间相等
• 每个节点n的子节点为n - 2^(j - 1), n - 2^(j - 2)..., n - 2^0, j为LSB所在的位置
• 每个节点n的父节点为n + 2^j, j为LSB所在的位置

对应更新和查询的操作为：

• 查询区间[0, i]的和: 每次查询tree[i]并且抹掉i的最后一个1直到1等于0。累加所有对应的tree[i]的值。
• 拿11(1011)举例，抹掉最后一个1之后为10(1010)。第一次我们查询tree[11] = arr[11] 
• 抹掉10(1010)最后一个1之后为8(1000)。第二次查询tree[10] = arr[10] + arr[9]
• 抹掉8(1000)最后一个1之后为0。最后一次查询tree[8] = tree[8] + tree[7] + ... + tree[1]
• 所以tree[11] + tree[10] + tree[8] = arr[11] + arr[10] + ... + arr[1]
• 查询区间[i, j]的和 = 区间和[0, i] - 区间和[0, j - 1]
• 单点更新arr[i]的值：我们要更新所有包含i的区间。所以我们每一次更新 i += 1 << j，j为i的最低位的1。
• 假设输入数组长度为16
• 拿i = 10(1010)举例，第一次我们更新tree[10] = arr[10] + arr[9]
• 更新i之后i = 10 + 2 = 12(1100)。第二次更新tree[12] = arr[12] + ... + arr[9]
• 更新i之后i = 12 + 4 = 16(10000)。第三次更新tree[16] = arr[16] + ... + arr[1]
• i之后会超过16，停止

实现的话有一个小技巧，我们可以通过i & -i来取i的最低为1的位。具体的原理可以参考这篇文章。实现起来代码非常的简洁。假设区间长度为N，单点更新和区间查询的复杂度为O(log N)。

假设我们要实现如下功能：

• 更新处于[0, N]区间的数, [0, N]区间初始为0
• 统计[i, j]的区间的区间和，其中i <= j <= N

那么实现代码如下：

	class BIT
	{
	public:
	BIT(int n)
	{
	N = n + 1;
	tree = vector<int>(N, 0);
	}
	void update(int key, int val)
	{
	int i = key + 1;
	while(i < N)
	{
	tree[i] += val;
	i += i & -i;
	}
	}
	int query(int s, int e)
	{
	return queryUntil(e) - queryUntil(s - 1);
	}
	private:
	int N;
	vector<int> tree;
	int queryUntil(int key)
	{
	int i = key + 1, res = 0;
	while(i)
	{
	res += tree[i];
	i -= i & -i;
	}
	return res;
	}
	};

view raw BIT.cpp hosted with ❤ by GitHub

2D Binary Index Tree

2D BIT和2D Segment Tree的思想是类似的。对于输入为m x n的数组，对应的2D BIT就是一个(m + 1) x (n + 1)的数组。每一行就是对应输入数组那一行的1D BIT。更新和查询依然是很类似的，我们要读/写对应行列的cell的数据。以下图为例，输入数组是一个16 * 8的2D矩阵，如果我们要update input矩阵的(4, 2)，也就是BIT中的(5, 3)：

• 要update的行为(5, 6, 8 16)
• 要update的列为(3, 4, 8)
• 要update的cell如下图蓝色所示

时间复杂度的话单点更新和查询就变为O(log m) * O(log n)。具体的实现可以参考这道求2D区间和的题目。

BIT与RMQ

我们可以看出，区间和的运算是一个可逆的过程。比如sum(0, 4) = 5, sum(5, 7) = 6, 那么sum(0, 7) = sum(0, 4) + sum(5, 7) = 11。同理我们可以反过来推出sum(5,7) = sum(0, 7) - sum(0, 4) = 11 - 5 = 6。但是并不是所有的区间查询都拥有这样的性质。一个很明显的例子，RMQ(Range Minimum Query)问题，我们可以推导出min(0, 7) = min(min(0, 4), min(5, 7))，但是这个过程不是可逆的，也就是说已知min(0, 7)和min(0, 4)的情况下，我么不一定能推导出min(5, 7)的值。我们之前实现区间和BIT的时候利用了区间和的这种性质，比如查询的时候我们查询两个从0开始的区间和然后相减。更新的时候我们直接在每个父节点的区间和上累加。但是对于RMQ问题，我们都无法进行这种操作，查询不必说，更新的话如果我们把arr中位于i位置的值oldVal更新为val，BIT的节点对应的值如何改变有以下几种情况：

1. 如果oldVal > val，那么我们从i对应的节点向上一直更新，即currNode.val = min(currNode.val, val)
2. 如果oldVal == val，我们不需要做任何事
3. 如果oldVal < val，情况会比较复杂：
1. 首先如果currNode.val != oldVal，代表oldVal不是当前区间的最小值，所以对于当前节点和父节点我们不需要做任何事
2. 如果currNode.val != oldVal，代表oldVal是当前区间的最小值，我们把oldVal更新为更大的val，那么当前区间的最小值是可能变的。那么我们需要去当前节点的所有子节点查询新的最小值，再更新当前区间的最小值。同理，对父节点我们也要做同样的事情，知道我们到达根节点或者发现不需要继续更新了。

更新的思路我们清楚了，但是如果按照Region Sum的方式存，我们还是没有办法解决查询的问题。为了解决这个问题，我们需要改变存储的方式，用两个BIT来存数据，具体的详情可以参照这篇paper。我们在这里简单介绍：

• 第一个BIT为我们之前介绍的BIT，节点i存的是[i - 2^j + 1, i]区间的最小值，j为LSB所在的位置
• 第二个BIT为反向的BIT，节点i存的是[i, i + 2^j - 1]区间的最小值，j为LSB所在的位置。节点n的父亲为i - (i & i), 儿子为n + 2^(j - 1), n + 2^(j - 2)..., n + 2^0, j为LSB所在的位置
• 此外我们还需要一个数组val来存原始的输入数组，val[i] = arr[i]。我们总共需要这三个数组存储。

两个BIT如下图所示(图片来自上面的论文，忽略每个节点存的值，我们只看结构)，我们实现的时候BIT2还是当做从1开始，所以第二个图会有些许不同:

对于区间[s, e]的查询，输入查询区间为[s - 1, e - 1]，因为我们存BIT从0开始:

1. 在BIT1中我们从s对应的节点开始向上爬，保证s不大于e。
2. 在BIT2中我们从e对应的节点开始向上爬，保证e不小于s。
3. 在两次攀爬的过程中我们最终停在的节点的index是一样的。证明如下: 我们可以把s和e写成二进制表达式，并且假设两个数第一个不一样的位从p + 1(从左到右)开始:
• s = a1a2a3...ap0bp+2...bn
• e = a1a2a3...ap1cp+2...cn
• 我们每次更新s += s & -s，我们最终会到达k = a1a2a3....ap10...0
• 我们每次更新e -= s & -s，我们最终也会到达k
4. 在1中我们经过的所有节点的index，我们去BIT2中取index对应节点对应区间的最小值
5. 在2中我们经过的所有节点的index，我们去BIT1中取index对应节点对应区间的最小值
6. 4和5对应的所有区间在加上val[k]的最小值，就是查询区间的最小值

以查询区间[4, 12]距离，对应BIT的index的查询区间为[5, 13]:

1. 我们首先在BIT1中从5开始爬，会经过5, 6最后停在8，在BIT2中对应的区间为[5, 5], [6, 7]
2. 同理在BIT2中从13开始爬，会经过13, 12最后停在8，在BIT1中对应的区间为[13, 13], [9, 12]
3. 那么val[8]加上1中[5, 5], [6, 7]和2中[13, 13], [9, 12]的区间正好对应了[5, 13]的区间，那么查询的结果就是[5, 13]的最小值

查询的时间复杂度为O(log n)，更新的话因为对于每个节点我们可能需要查一下所有的子节点，所以时间复杂度为O(log^2 n)，代码如下：

	//Binay Indexed Tree for RMQ
	class BIT
	{
	public:
	BIT(vector<int>& input)
	{
	int n = input.size();
	N = n + 1;
	tree1 = vector<int>(N, INT_MAX);
	tree2 = vector<int>(N, INT_MAX);
	vals = vector<int>(N, INT_MAX);
	for (int i = 0; i < n; ++i)add(i, input[i]);
	}
	void update(int key, int val)
	{
	int i = key + 1, original = vals[i], cache = val;
	vals[i] = val;
	//update BIT1
	while (i < N)
	{
	//if new value is greater than original value
	//there is chance that tree1[i] will increase
	//since orignal can be minval for current node
	//we need to find out what should be the new
	//minimum val for curr node. Since all children
	//plus original[i] covers the range that current
	//node represents
	if (val > original)
	{
	if(tree1[i] == original)
	{
	val = min(val, vals[i]);
	//check all children
	int r = 1;
	while ((i & -i) >> r)
	{
	int j = (i & -i) >> r;
	val = min(val, tree1[i - j]);//update val based on infromation in children nodes
	++r;
	}
	}
	else
	break;
	}
	if (val == tree1[i])break;
	tree1[i] = val;//update BIT1 node with new value
	i += (i & -i);//go to parent
	}
	//update BIT2, similar to code above
	i = key + 1;
	val = cache;
	while (i)
	{
	if (val > original)
	{
	if(tree2[i] == original)
	{
	val = min(val, vals[i]);
	int r = 1;
	while ((i & -i) >> r)
	{
	int j = (i & -i) >> r;
	val = min(val, tree2[i + j]);
	++r;
	}
	}
	}
	else
	break;
	if (val == tree2[i])break;
	tree2[i] = val;
	i -= (i & -i);
	}
	}
	int query(int start, int end)
	{
	++start;
	++end;
	int res = INT_MAX;
	int prev = start, curr = prev + (prev & -prev);
	while (curr <= end)
	{
	res = min(res, tree2[prev]);
	prev = curr;
	curr += curr & -curr;
	}
	prev = end;
	curr = prev - (prev & -prev);
	while (curr >= start)
	{
	res = min(res, tree1[prev]);
	prev = curr;
	curr -= curr & -curr;
	}
	//prev stays at the same place after each of
	//these two iteration
	res = min(vals[prev], res);
	return res;
	}
	private:
	int N;
	//tree1: normal BIT, bigger value as parent, node i covers [i - 2^r + 1, i],
	//where r is LSB(least significant bit) of 1 in i.
	//tree2: reversed BIT, smaller value as parent, node i covers [i, i + 2^r - 1]
	//where r is LSB of 1 in i.
	//vals: vals store original value at i
	vector<int> tree1, tree2, vals;
	void add(int key, int val)
	{
	int i = key + 1;
	vals[i] = val;
	while (i < N)
	{
	tree1[i] = min(tree1[i], val);
	i += i & -i;
	}
	i = key + 1;
	while (i)
	{
	tree2[i] = min(tree1[i], val);
	i -= i & -i;
	}
	}
	};

view raw BIT4RMQ.cpp hosted with ❤ by GitHub