树上差分
树上差分,顾名思义,意思就是在树上做差分。
至于什么是差分呢?好吧,那就先讲讲差分是什么东西吧,作为树上差分的前置技能。。。。
差分:
先从差分的定义讲起。
定义:差分(difference)又名差分函数或差分运算,差分的结果反映了离散量之间的一种变化,是研究离散数学的一种工具。它将原函数f(x) 映射到f(x+a)-f(x+b) 。差分运算,相应于微分运算,是微积分中重要的一个概念。总而言之,差分对应离散,微分对应连续。差分又分为前向差分、向后差分及中心差分三种。(一般套路,没什么实际意义,只是用来看看的)。
相信大家都知道前缀和是个什么玩意儿(后面引出差分)。
前缀和:定义一个数组f[i], 使得sum[i] = sum[i - 1] + sum[i]; 即sum[i] 是所有 小于等于i 的 k 所表示的 a[k] 的累加和。(注:这里讨论的所有字母均为大于0的整型,不存在浮点型)
即 sum[i] = ∑a[i] (i 属于[1,i]),sum[i] 即为 a[i] 的 前缀和数组.
用前缀和可以有效应对多次查询,我们可以得到,属于区间[i, j ]的k的所有a[k]的累加和 为 sum[j] - sum[i - 1].
应用这个公式,先用一遍求出所有前缀和sum[i], 然后面对所有k次的查询区间[a, b],我们均可以用O(1)的算法:sum[b] - sum[a - 1]得出区间累加和。
总时间复杂度为(O(N + K));
关于前缀和还有许多应用,这里不一一赘述,只是顺带提及一下这个概念。
好了下面进入差分环节:
前面普及了一下前缀和,由前缀和数组,我们知道: sum[i] = sum[i - 1] + a[i];
移项得:a[i] = sum[i] - sum[i - 1];...........................(一)
此时,我们变换一下数组,定义一个数组f[i] , 它表示a[i] 的差分数组,它使得f[i] = a[i] - a[i - 1];
所以:f[i] = a[i] - a[i - 1];...........................(二)
对比(一)(二)两式子,我们发现,原数组a[i] 也可表示为前缀和数组sum[i] 的差分数组.同时原数组a[i]也可以表示为差分数组f[i]的前缀和数组。所以,我们发现,差分和前缀和其实是一对互逆运算。我们同样也可以用前缀和的思路逆推来求出原数组的差分数组。
即 a[i] = ∑f[i] (i 属于[1,i]).............................(三)
仔细思考一下,若我们对属于区间[a,b]的 i 的每一个a[i] 都加1, 该肿么办?每一个数都模拟,一个一个加吗?
显然一个一个加效率太低。其实,我们想一想,是不是也可以用类似前缀和的方法,用O(1)的算法解决?
好,以下进入讲解时刻:
由公式(三),因为属于区间[a, b]中的 i 的a[i]都同时加了1,所以相对于区间[a + 1,b]中,差分数组f[i]不变。(i 属于[a + 1,b])
但是对于f[a] ,它应该 +1, 使得后面所有的a[i](i >= a) 都 +1. 由公式a[i] = ∑f[i] (i 属于[1,i])得出。
而a[b + 1] 及以后的数组不需要增加,而前面的类似链式反应的操作使它 +1了,所以我们要把它减回去,即f[b + 1] - 1;
所以,我们通过对f[i]数组操作,令f[a] + 1, f[b + 1] - 1,就能表示出这一操作了,至于具体的a[i] 的值,则可以用前缀和a[i] = ∑f[i](i 属于[1,i])求出.下面给出一个例子:
【差分】 糖块
题目描述
现在有n(1 <= N <= 1,000,000, N 是奇数)个盒子,编号是1..n。
数学老师为了惩罚他,决定让他做一个难题,他让小x会对这些盒子做k(1 <=k <= 25,000)次放糖块的操作(这得多少糖块呀).
数学老师每次会给小x一个区间[a,b],这时小x就会从编号是a的盒子到编号是b的盒子每个盒子都放一个糖块。
做完k次操作后,数学老师会问小x,在所有盒子的糖块个数中,这些糖块个数的中位数是多少。(最中间的值)。
因为n是奇数,所以这个答案肯定是唯一的。
输入格式
第一行:两个整数 n k。 接下来k行,每行一个区间 ai bi ,表示要放糖块的区间。
输出格式
一个整数,表示中位数的值。
样例数据
input
7 4
5 5
2 4
4 6
3 5
output
1
【样例解释】一共有7个盒子,4个操作,第一次在5号盒子里放1个,第二次在2 3 4号盒子里放。
放完糖块之后,盒子里的糖块依次是0,1,2,3,3,1,0.排过序后,数字1是中位数。
数据规模与约定:
时间限制:1s
空间限制:256MB
---------------------------------------------------我是完美的分割线------------------------------------------------------------------
这道题就是一道典型的差分模板题,对于每一个操作,按照之前描述的f[i]差分,用前缀和求出原数组,再排一下序,用求出中位数就可以了。
代码如下:
1 #include<bits/stdc++.h>
2 using namespace std;
3 int n, m;
4 int f[1000010], a[1000100];
5 int main()
6 {
7 //freopen("candy.in","r",stdin);
8 //freopen("candy.out","w",stdout);
9 memset(f, 0,sizeof(f));
10 scanf("%d%d",&n, &m);
11 for(int i = 1;i <= m;i++)
12 {
13 int ai, bi;
14 scanf("%d%d",&ai, &bi);
15 f[ai]++, f[bi+1]--;
16 }
17 for(int i = 1;i <= n;i++) a[i] = a[i - 1] + f[i];
18 sort(a + 1, a + n + 1);
19 printf("%d
", a[n/2+1]);
20 return 0;
21 }
差分差不多讲完了。下面来说说树的lca和倍增优化吧,也作为树上差分的前置技能。
LCA(最近公共祖先)
先说一下lca 的定义吧:LCA(Lowest Common Ancestors),即最近公共祖先,是指在有根树中,找出某两个结点u和v最近的公共祖先。
LCA的实现有多种算法,有朴素(暴力)算法,倍增,线段树,tarjan等等。
这里只讨论朴素(暴力)算法和倍增优化。
暴力枚举(朴素算法)
对于有根树T的两个结点u、v,首先将u,v中深度较深的那一个点向上蹦到和深度较浅的点,然后两个点一起向上蹦,直到蹦到同一个点,这个点就是u,v的最近公共祖先,记作LCA(u,v)。但是这种方法的时间复杂度在极端情况下会达到O(N)。特别是有多组数据求解时,时间复杂度将会达到O(N * M)。
考虑一下倍增优化。
先说说倍增是个什么东西吧。
倍增:
算法时间复杂度为log 级的,相比一步一步走到终点的复杂度(O(N))确实是快了不少。
倍增的主要应用就是LCA和RMQ(区间最值查询)。
顺带讲一下RMQ吧:
RMQ定义:对于长度为n的数列A,回答若干询问RMQ(A,i,j)(i,j<=n),返回数列A中下标在i,j里的最小(大)值,也就是说,RMQ问题是指求区间最值的问题。
面对这个问题,显然,暴力求区间最小(大)值是行不通的。
那么,我们能否像差分那样,用倍增的方法先预处理,用O(1)的时间复杂度回答每一个查询呢?
显然可以做到。
这就要用到动态规划(DP)和倍增的思想了(ST算法)。
例如:
a数列为:3 2 4 5 6 8 1 2 9 7
首先:设a[i]是要求区间最值的数列,f[i, j]表示从第i个数起连续2 ^ j个数中的最大值。(DP的状态)
f[1,0]表示第1个数起,长度为2 ^ 0=1的最大值,其实就是3这个数。同理 f[1,1] = max(3,2) = 3,f[1,2]=max(3,2,4,5) = 5,f[1,3] = max(3,2,4,5,6,8,1,2) = 8;
并且我们可以容易的看出f[i,0]就等于A[i]。(DP的初始值)
我们把f[i,j]平均分成两段(因为f[i,j]一定是偶数个数字),从 i 到i + 2 ^ (j - 1) - 1为一段,i + 2 ^ (j - 1)到i + 2 ^ j - 1为一段(长度都为2 ^ (j - 1))。
于是我们得到了状态转移方程f[i, j]=max(f[i,j - 1], f[i + 2 ^ (j - 1),j - 1])。
预处理完成!
接下来是查询:
假如我们需要查询的区间为(i,j),那么我们需要找到覆盖这个闭区间(左边界取i,右边界取j)的最小幂(可以重复,比如查询1,2,3,4,5,我们可以查询1234和2345)。
因为这个区间的长度为j - i + 1,所以我们可以取k=log2( j - i + 1),则有:RMQ(i, j)=max{f[i , k],,f[ j - 2 ^ k + 1, k]}。
举例说明,要求区间[1,5]的最大值,k = log2(5 - 1 + 1)= 2,即求max(f[1,2],f[5 - 2 ^ 2 + 1, 2])=max(f[1, 2],f[2, 2]);
为了方便理解,以下给出(ST)预处理代码:
1 void ST_prework()
2 {
3 memset(f, 10, sizeof(f));
4 for (int i = 1; i <= n;i++) f[i][0] = a[i];
5 int t = log(n) / log(2) + 1;
6 for (int j = 1; j < t;j++)
7 for (int i = 1; i <= n - (1 << j) + 1;i++)
8 f[i][j] = min(f[i][j-1], f[i + (1 << (j-1))][j-1]);//区间的最小值,最大值同理
9 }
查询代码:
1 int ST_query(int l, int r)
2 {
3 int k = log(r - l + 1) / log(2);
4 return min(f[l][k], f[r - (1 << k) + 1][k]);
5 }
以上就是RMQ大概了,好像讲的有点偏。
给一道例题加深理解:
[倍增]忠诚
题目描述
老管家是一个聪明能干的人。他为财主工作了整整10年,财主为了让自已账目更加清楚。
要求管家每天记k次账,由于管家聪明能干,因而管家总是让财主十分满意。但是由于一些人的挑拨,财主还是对管家产生了怀疑。
于是他决定用一种特别的方法来判断管家的忠诚,他把每次的账目按1,2,3…编号,然后不定时的问管家问题,
问题是这样的:在a到b号账中最少的一笔是多少?为了让管家没时间作假他总是一次问多个问题。
输入格式
输入中第一行有两个数m,n表示有m(m<=100000)笔账,n表示有n个问题,n<=100000。
第二行为m个数,分别是账目的钱数
后面n行分别是n个问题,每行有2个数字说明开始结束的账目编号。
输出格式
输出文件中为每个问题的答案。具体查看样例。
样例数据
input
10 3
1 2 3 4 5 6 7 8 9 10
2 7
3 9
1 10
output
2 3 1
数据规模与约定
时间限制:1s
空间限制:256MB
--------------------------------我是完美的分割线--------------------------------------------------------
这道题是RMQ的模板题,对每次进行ST算法就行了,做法及细节不多说,比较基础。
代码如下:
1 #include <bits/stdc++.h>
2 using namespace std;
3 int f[100000][50];
4 int n, m;
5 int a[1000001];
6 void ST_prework()
7 {
8 memset(f, 10, sizeof(f));
9 for (int i = 1; i <= n;i++) f[i][0] = a[i];
10 int t = log(n) / log(2) + 1;
11 for (int j = 1; j < t;j++)
12 for (int i = 1; i <= n - (1 << j) + 1;i++)
13 f[i][j] = min(f[i][j-1], f[i + (1 << (j-1))][j-1]);
14 }
15 int ST_query(int l, int r)
16 {
17 int k = log(r - l + 1) / log(2);
18 return min(f[l][k], f[r - (1 << k) + 1][k]);
19 }
20 int main()
21 {
22 //freopen("input.in", "r", stdin);
23 //freopen("output.out" , "w", stdout);
24 scanf("%d%d",&n ,&m);
25 for (int i = 1; i <= n;i++) scanf("%d", &a[i]);
26 ST_prework();
27 for (int i = 1; i <= m;i++)
28 {
29 int x, y;
30 scanf("%d%d", &x, &y);
31 printf("%d", ST_query(x, y));
32 if (i != m) printf(" ");
33 }
34 return 0;
35 }
以上就是这样了。
通过RMQ这个例子,我们也可以体会到倍增的强大力量!!!
下面进入LCA的环节:
同样运用倍增的思想优化,
注意到u,v走到最近公共祖先w之前,u,v所在结点不相同。而到达最近公共祖先w后,再往上走仍是u,v的公共祖先,即u,v走到同一个结点,这具有二分性质。
同样,类似于解决RMQ的方法,我们预处理fa[ i ][ j ]数组表示u往上走2 ^ k步走到的结点.
不妨假设depth[u] < depth[v]
①将v往上走d = depth[v] - depth[u]步,此时u,v所在结点深度相同,该过程可用二进制(倍增)优化。由于d是确定值,将d看成2的次方的和值,d = 2k ^ 1 + 2k ^ 2 + ... + 2k ^ m,利用fa数组,如v = fa[k1][v],v = fa[k2][v]加速。
②若此时u = v,说明lca(u,v)已找到。
③利用fa数组加速u,v一起往上走到最近公共祖先w的过程。令d = depth[u] - depth[w],虽然d是个未知值,但依然可以看成2的次方的和。从高位到低位枚举d的二进制位,设最低位为第0位,若枚举到第k位,有fa[k][u] != fa[k][v],则令u = fa[k][u],v = fa[k][v]。最后最近公共祖先w = fa[0][u] = fa[0][v],即u和v的父亲。
如何预处理?
k=0时,fa[k][u]为u在有根树中的父亲,令根结点fa[k][root]=-1。
k>0时,fa[k][u]=fa[k-1][fa[k-1][u]]。树的高度最多为n,k是logn级别。
时间复杂度:
预处理O(n log n)
单次查询O(log n)
这里给出广搜(bfs)的预处理代码:(用的是邻接表存储和STL的队列,也有深搜的写法)
1 void bfs()
2 {
3 q.push(1);
4 d[1] = 1;
5 while(!q.empty())
6 {
7 int x = q.front();
8 q.pop();
9 for(int i = head[x];i;i = e[i].next)
10 {
11 int y = e[i].to;
12 if(d[y]) continue;
13 d[y] = d[x] + 1;
14 dist[y] = dist[x] + e[i].to;
15 f[y][0] = x;
16 for(int j = 1;j <= t;j++)
17 f[y][j] = f[f[y][j - 1]][j - 1];
18 q.push(y);
19 }
20 }
21 }
接下来是LCA的主代码:
1 int lca(int x, int y)
2 {
3 if(d[x] > d[y]) swap(x, y);
4 for(int i = t;i >= 0;i--)
5 {
6 if(d[f[y][i]] >= d[x])
7 y = f[y][i];
8 }
9 if(x == y) return x;
10 for(int i = t;i >= 0;i--)
11 if(f[x][i] != f[y][i])
12 x = f[x][i], y = f[y][i];
13 return f[x][0];
14 }
以上就是LCA的倍增优化算法,至于其他的算法,以后再说。
好了,终于到正题了!!前面的都是铺垫(前置技能),接下来的才是好戏!
树上差分:
再来一波普及前置知识:
需要知道的树的性质:
(一).树上任意两个点的路径唯一.
(二).任何子节点的父亲节点唯一.(可以认为根节点是没有父亲的)
树上差分的形式有两种:
1. 点差分 2.边差分
先讲点差分:
令cnt_i为节点i被经过的次数.
例如,我们从 s-->t ,求这条路径上的点被经过的次数.
很明显的,我们需要找到他们的LCA(这个点时中转点)
我们需要让cnt[s]++,让cnt[t]++,而让他们的cnt[lca]--,cnt[faher(lca)]--;(差分)
考虑:我们搜索到s,向上回溯.
下面以u表示当前节点,son_i代表i的儿子节点.(如果一些son不给出下标,即代表当前节点u的儿子
每个u统计它的子树大小,顺着路径标起来.(即cnt[u]+=cnt[son])
我们会发现第一次从s回溯到它们的LCA时候,cnt[LCA]+=cnt[son[LCA]]cnt[LCA]=0!
"不是LCA会被经过一次嘛,为什么是0!"
别急,我们继续搜另一边.
继续:我们搜索到t,向上回溯.依旧统计每个u的子树大小cnt[u]+=cnt[son]再度回到LCA 依旧 是cnt[LCA]+=cnt[son[LCA]]
这个时候 cnt[LCA]=1 这就达到了我们要的效果。
担忧: 万一我们再从LCA向上回溯的时候使得其父亲节点的子树和为1怎么办?
这样我们不就使得其父亲节点被经过了一次?
因此我们需要在cnt[faher(lca)]--
这样就达到了标记我们路径上的点的要求!
边差分:
原理和点差分差不多,只要把边塞给深度较深的点即可。
其实本质就是差分的方法,应用到树上,找到lca ,进行差分。最后统计的时候用一下前缀和就可以了。
如果lca ,树 , 差分 ,前缀和 等知识都掌握得很牢的话,树上差分就是小ks。
以下给一道例题体会:
【bzoj 4390】 [Usaco2015 dec]Max Flow
题目描述
Farmer John has installed a new system of N-1 pipes to transport milk between the N stalls in his barn (2≤N≤50,0002 leq N leq 50,0002≤N≤50,000), conveniently numbered 1…N1…N. Each pipe connects a pair of stalls, and all stalls are connected to each-other via paths of pipes.
FJ is pumping milk between KKK pairs of stalls (1≤K≤100,0001 leq K leq 100,0001≤K≤100,000). For the iiith such pair, you are told two stalls sis_isi and tit_iti, endpoints of a path along which milk is being pumped at a unit rate. FJ is concerned that some stalls might end up overwhelmed with all the milk being pumped through them, since a stall can serve as a waypoint along many of the KKK paths along which milk is being pumped. Please help him determine the maximum amount of milk being pumped through any stall. If milk is being pumped along a path from sis_isi to tit_iti, then it counts as being pumped through the endpoint stalls sis_isi and tit_iti, as well as through every stall along the path between them.
FJ给他的牛棚的N(2≤N≤50,000)个隔间之间安装了N-1根管道,隔间编号从1到N。所有隔间都被管道连通了。
FJ有K(1≤K≤100,000)条运输牛奶的路线,第i条路线从隔间si运输到隔间ti。
一条运输路线会给它的两个端点处的隔间以及中间途径的所有隔间带来一个单位的运输压力,你需要计算压力最大的隔间的压力是多少。
输入格式
The first line of the input contains N and K.
The next N?1 lines each contain two integers x and y (x≠y) describing a pipe between stalls x and y.
The next K lines each contain two integers s and t describing the endpoint stalls of a path through which milk is being pumped.
输出格式
An integer specifying the maximum amount of milk pumped through any stall in the barn.
样例数据
input
5 10
3 4
1 5
4 2
5 4
5 4
5 4
3 5
4 3
4 3
1 3
3 5
5 4
1 5
3 4
output
9
数据规模与约定
时间限制:1s
空间限制:256MB
------------------------------我是完美的分割线------------------------------------
这道题留给读者思考。
以下为参考程序(若非穷途末路,请勿查看):
1 #include<bits/stdc++.h>
2 using namespace std;
3 int d[100010], f[100010][50], dist[100010];
4 int head[1000100], cnt = 0, pows[100010];
5 int n, m, ans = 0;
6 int t;
7 queue <int > q;
8 struct node
9 {
10 int to, next;
11 }e[1000100];
12 void add(int x, int y)
13 {
14 cnt++;
15 e[cnt].to = y;
16 e[cnt].next = head[x];
17 head[x] = cnt;
18 }
19 void bfs()
20 {
21 q.push(1);
22 d[1] = 1;
23 while(!q.empty())
24 {
25 int x = q.front();
26 q.pop();
27 for(int i = head[x];i;i = e[i].next)
28 {
29 int y = e[i].to;
30 if(d[y]) continue;
31 d[y] = d[x] + 1;
32 dist[y] = dist[x] + e[i].to;
33 f[y][0] = x;
34 for(int j = 1;j <= t;j++)
35 f[y][j] = f[f[y][j - 1]][j - 1];
36 q.push(y);
37 }
38 }
39 }
40 int lca(int x, int y)
41 {
42 if(d[x] > d[y]) swap(x, y);
43 for(int i = t;i >= 0;i--)
44 {
45 if(d[f[y][i]] >= d[x])
46 y = f[y][i];
47 }
48 if(x == y) return x;
49 for(int i = t;i >= 0;i--)
50 if(f[x][i] != f[y][i])
51 x = f[x][i], y = f[y][i];
52 return f[x][0];
53 }
54 void dfs(int x, int fa)
55 {
56 for(int i = head[x];i;i = e[i].next)
57 {
58 int y = e[i].to;
59 if(y == fa) continue;
60 dfs(y, x);
61 pows[x] += pows[y];
62 }
63 ans = max(ans, pows[x]);
64 }
65 int main()
66 {
67 //freopen("maxflow.in","r",stdin);
68 //freopen("maxflow.out","w",stdout);
69 scanf("%d%d", &n, &m);
70 t = log2(n);
71 for(int i = 1;i < n;i++)
72 {
73 int xx, yy;
74 scanf("%d%d",&xx, &yy);
75 add(xx, yy);
76 add(yy, xx);
77 }
78 bfs();
79 for(int i = 1;i <= m;i++)
80 {
81 int xx, yy;
82 scanf("%d%d",&xx, &yy);
83 int l = lca(xx, yy);
84 pows[xx]++, pows[yy]++;
85 pows[l]--, pows[f[l][0]]--;
86 }
87 dfs(1, 0);
88 printf("%d
", ans);
89 return 0;
90 }
树上差分的内容差不多就到这儿,本蒟蒻的手也打累了,休息一会儿qwq。