UOJ Logo immortalCO的博客

博客

NOIP2018 D2T3 题解 + 关于动态 DP 等科技的一些总结

2018-11-12 11:29:08 By immortalCO

利益无关:猫锟不是本题的命题猫。

算法一

最小权覆盖集 = 全集 - 最大权独立集

强制取点、不取点可以使用把权值改成正无穷或负无穷实现。

接下来就是经典的“动态最大权独立集”了,这题的题解在我之前博客有,洛谷上也可以提交。

$O(n\log n)$。

算法二

考虑修改的两个点 $a,b$ 构成这条链。

可以把操作看作是:先在这条链伸出去的每棵子树上 DP,最后再在这条链上 DP。

伸出去的子树以及链的中间的点和修改无关,因而可以整合起来高效处理,因为它们的转移是不受影响的;但是 $a$ 和 $b$ 的转移是受影响的,需要单独处理。

因此可以先预处理出每个子树的 DP 值,然后用倍增或树剖维护“从一个点往上按照通用规则 DP 到另一个点的结果”,这样只需要特殊处理修改的点,其他地方只需要直接倍增。

$O(n\log n)$。

算法三

如果不想写倍增,也可以写这样一个算法:把所有的询问放在一起处理。

也就是说,由于每条链中间的 DP 转移全部相同,因而我们可以批量地对所有需要进行这一转移的 DP 进行转移,从而加快速度。

这个算法虽然实现较容易,但理解较为困难,因而不再叙述。

$O(n\log n)$。

总结

本题的三个做法对应着动态 DP、虚树计算和整体 DP 的思想,体现了三个算法(在一定条件下)的等价性。

动态 DP $\Rightarrow$ 虚树计算:动态维护被修改的点构成的虚树(可能会增加复杂度);

虚树计算 $\Rightarrow$ 动态 DP:一个一个加入虚树中的点;

动态 DP $\Rightarrow$ 整体 DP:每个修改看作时间的区间,用线段树批量应用修改(不一定可行);

整体 DP $\Rightarrow$ 动态 DP:一个一个询问进行处理,依次加入对称差(可能会增加复杂度);

虚树计算 $\Rightarrow$ 整体 DP:按询问分组批量进行;

整体 DP $\Rightarrow$ 虚树计算:通过 整体 DP $\Rightarrow$ 动态 DP $\Rightarrow$ 虚树计算 进行转化(不一定可行)。

可以看出,动态 DP 和整体 DP 具有很强的适用面,都有各自才能解决的问题;而虚树计算的能力则是最弱的,常常可以等复杂度地转化为另外两者(一些例外是需要在虚树上点分治的问题)。

基于变换合并的树上动态 DP 的链分治算法 & SDOI2017 切树游戏(cut)解题报告

2017-05-22 12:53:03 By immortalCO

最近听说 SDOI2017 R2 某题可以用我论文中的“基于链的分治优化树上动态 DP”的方法 A 掉……我很开心啊……

因此我把这道题出到了校内训练中……

在出题之前,肯定是要去把这题先做一遍的……然后,我发现我弄出了一些新的东西……!

题意简述

给出一棵无根树 $T$,节点从 $1$ 到 $n$ 编号,点 $i$ 的权值为 $v_i$。定义一棵树的权值为树中所有点的权值异或起来的结果。

你想要在树上玩一个游戏。你会做 $m$ 次以下两种操作:

  • C x y:将 $v_x$ 修改为 $y$。
  • Q k:询问 $T$ 有多少棵非空的连通子树,满足这棵子树的权值为 $k$。输出模 $10007$。

数据范围:$n,m\le 30000$,$0\le v_i,y,k < 128$。

时空限制:2s、128MB。

算法讨论

如果只有一次询问,非常容易想到暴力 DP。先转有根树。在全局记录答案数组 $ans(k)$ 表示权值为 $k$ 的子树个数。对每个点 $i$ 记录 $f(i, k)$ 表示子树中深度最小的点为 $i$ 且子树权值为 $k$ 的连通子树个数。记录 $g(i, j, k)$ 表示子树中深度最小的点为 $i$ 且所有其他的节点都在 $i$ 的前 $j$ 个子节点的子树中的连通子树个数。那么我们就有以下方程(设 $Ch(i)$ 为 $i$ 的子节点列表):

  • $g(i, 0, k) = [k=v_i]$
  • $g(i, j, k) = \sum_{t=0}^{127} g(i,j-1,t)\times (f(Ch(i)_j, k\oplus t) + [k\oplus t = 0])$
  • $f(i, k) = g(i, |Ch(i)|, k)$
  • $ans(k) = \sum_{i=1}^n f(i, k)$

总时间复杂度为 $O(nm\times 128^2)$。

接下来可以注意到第 2 个式子是一个“异或卷积”的形式,不难想到使用 FWT 可以优化到 $O(128\log 128)$。然后注意到 FWT 之后,加法和乘法都可以直接按位相加,因此可以在一开始就将所有数组 FWT,运算过程中全部使用 FWT 之后的数组,最后再将 $ans( * ) $ 数组 FWT 回去即可。这样就可以去掉一个 $\log 128$。时间复杂度为 $O((n + \log 128)\times 128)$。

再接下来就是优化修改复杂度了。看过我论文或做过 BZOJ 4712 的同学容易想到使用链分治维护树上动态 DP。首先将树进行轻重链剖分,然后按照重链构成的树的顺序进行 DP。如果这样以后每一条重链上的转移可以高效维护、支持修改,那么每次修改点 $p$ 之后,我们就可以高效地更新点 $p$ 到根的 $O(\log n)$ 条重链的信息即可。

首先 $ans(k)$ 是全局变量,不好维护。那么可以不记录 $ans( * )$,而是记录 $h(i, k)$ 表示 $i$ 子树中的 $f(i, k)$ 的和,那么这样整个 DP 就有子树的阶段性了。

可以发现 $f(i, k)$ 就是先将 $g(i, 0, * )$ 和所有子节点 $p$ 的 $f(p, k ) + [k = 0]$ 全部卷积起来的值。即如果设 $F_i(z)$ 表示 $f(i, * )$ 这一数组的生成函数,那么可以得出 $$F_i(z) = z^{v_i}\prod_{p\in Ch(i)} (F_p(z) + z^0) $$ 这里的卷积定义为异或卷积。那么对于一条重链上的每一个点 $i$,我们只需要将 $i$ 的所有轻儿子 $lp$ 的 $F_{lp}(z) + z^0$ 全部卷积起来,这样就考虑了所有轻儿子的子树中的贡献,设这个卷积的结果为 $LF_i(z)$。同样对于每个点我们记录 $LH_i(z)$ 表示这个点的每个轻儿子的 $H_{lp}(z)$ 之和(这里 $H_i(z)$ 的定义类似 $F_i(z)$,只不过是对 $h(i, * )$ 定义的)。每个点的轻边的信息和可以用线段树维护来支持高效修改。

Claris 大神说这里信息可减因此不用线段树,但我觉得这里的 $LF_i(z)$ 的信息相减需要做除法,如果出现 10007 的倍数则没有逆元,无法相除,因此我仍然采用线段树维护。

注意到上述算法只要求我们能求出 $F_{重链顶}(z)$ 和 $H_{重链顶}(z)$,就可以维护父亲重链的信息或答案了。因此现在只需要考虑所有过当前重链的子树。在这里我们有如下两种截然不同的思路。

基于维护序列信息的算法

论文中提到的方法是转化为序列上的问题,然后使用线段树维护。由于连通子树和重链的交一定也是一段连续的链,那么我们显然就可以像最大子段和问题那样,记录 $D_{a,b}(z)$ 表示 $a=[$左端点在连通子树中$]$、$b=[$右端点在连通子树中$]$ 的方案数。这个算法将修改复杂度优化为 $O(128\log^2 n + 128\log 128)$,已经可以通过本题了。

但是这个算法有一定的问题:首先它具有较大的常数因子,运行时间较慢。其次,这个算法仍然利用了具体题目的性质——连通子树和重链的交还是链。而并非所有的题都有这样的性质。最后,由于要不重不漏地计数,代码细节繁多,十分难写。

基于变换合并的算法

对于一条重链,设重链上的点按深度从小到大排序后为 $p_1,p_2,...,p_c$,那么我们可以得出以下方程:

  • $F_{p_c}(z) = H_{p_c}(z) = z^{v_{p_c}}$ (因为 $p_c$ 没有子节点)
  • $F_{p_i}(z) = LP_{p_i}(z) \times (F_{p_{i+1}}(z) + z^0) \times {z^{v_{p_i}}}$
  • $H_{p_i}(z) = H_{p_{i+1}}(z) + F_{p_{i}}(z)$

而我们所需要求的只有 $F_{p_1}(z)$ 和 $H_{p_1}(z)$。

可以观察到上面这个式子中,向量 $\left(F_{p_{i+1}}(z), H_{p_{i+1}}(z), z^0\right)$ 是通过一个线性变换得到向量 $\left(F_{p_i}(z), H_{p_i}(z), z^0\right)$,具体地来说是右乘上这样一个矩阵:

$$M_i=\begin{pmatrix} LF_{p_i}{z}\times {z^{v_{p_i}}} & LF_{p_i}{z}\times {z^{v_{p_i}}} & 0 \\ 0 & 1 & 0 \\ LF_{p_i}{z}\times {z^{v_{p_i}}} & LH_ {p_i} + LF_{p_i}{z}\times {z^{v_{p_i}}} & 1 \end{pmatrix}$$

而矩阵乘法是满足结合律的,也就是说,我们只需要用线段树支持单点修改某个矩阵 $M_i$、维护矩阵的积,我们就可以高效地求出我们所需要的向量 $(F_{p_1}(z), H_{p_1}(z), 1)$。而这是容易做到的,因此这个算法是完全可行的。这样,这个算法也将修改复杂度优化为了 $O(128\log^2 n + 128\log 128)$,可以通过本题。

简单优化这个算法的常数。注意到形如 $$\begin{pmatrix} \underline{a} & \underline{b} & 0 \\ 0 & 1 & 0 \\ \underline{c} & \underline{d} & 1 \end{pmatrix}$$ 的矩阵乘法对这个形式封闭,因为 $$\begin{pmatrix} \underline{a_1} & \underline{b_1} & 0 \\ 0 & 1 & 0 \\ \underline{c_1} & \underline{d_1} & 1 \end{pmatrix} \times \begin{pmatrix} \underline{a_2} & \underline{b_2} & 0 \\ 0 & 1 & 0 \\ \underline{c_2} & \underline{d_2} & 1 \end{pmatrix} = \begin{pmatrix} \underline{a_1 a_2} & \underline{b_1 + a_1 b_2} & 0 \\ 0 & 1 & 0 \\ \underline{a_2 c_1 + c_2} & \underline{b_2 c_1 + d_1 + d_2} & 1 \end{pmatrix}$$

因此我们只需要对每个矩阵维护 $a,b,c,d$ 四个变量即可。同时可以直接用等号右边的形式来计算矩阵乘法,这样就只需要做 $4$ 次而不是 $27$ 次生成函数乘法了,常数大大减小了。

比较与扩展

这两个算法的时间复杂度相同,并且都可以扩展到询问某一个有根树子树的信息——只需要对那一条重链询问一下信息和/变换和即可。

我们来审视一下后一个算法。首先,这个算法基于的是直接对这个 DP 本身进行优化这样一种思想,而不是通过分析具体题目的性质进行处理,因此这种算法具有更高的通用性。其次,由于这个算法是直接对这个 DP 本身进行优化,因此正确性显然,细节也要少于论文中介绍的在区间中维护 $D_{a,b}(z)$ 信息的方法(维护 $D_{a,b}(z)$ 这个方法必须严格分类,因此细节繁多,常数也较大)。因此这个算法比前一个的算法更加优秀。

然而,事实上这个算法同样利用了题目的一些性质——这题是计数类问题,而且转移是线性变换,因此可以用矩阵来维护,而矩阵恰恰是一种可以合并的变换。那么对于其他的题目,是否也能用这种基于变换合并的算法呢?

BZOJ 4712 加强版(不保证修改不降)

这是一道在论文中有提到的题目。论文中介绍的算法正是一种基于标记合并的算法,其标记为 $trans_{a,b}(x)=\min(a, x+b)$ 并且对加法封闭。这个算法的常数十分优秀。

树上动态最大权独立集问题

这同样是一道在论文中有提到的题目。论文中介绍的是使用线段树维护序列上信息的方法。

考虑新的思路。记录 $f(i)$ 表示所有点都在 $i$ 子树中的最大权独立集的权值(以下简称答案),$g(i)$ 表示所有点都在 $i$ 子树中且独立集不包括点 $i$ 的答案;同样记录 $Lf(i)$、$Lg(i)$ 表示所有轻儿子的 $f(lp)$ 或 $g(lp)$ 的和

首先写出 $p_i$ 关于 $p_{i+1}$ 的转移:

  • $g(p_i) = Lf(p_i) + f(p_{i+1})$
  • $f(p_i) = \max(g(p_i), v_i + Lg(p_i) + g_(p_{i+1})$

考虑对这些 $g_{p_i}$、$f_{p_i}$ 定义新的加法和乘法运算:其中“加法” $a\oplus b$ 定义为 $\max(a, b)$、“乘法” $a\otimes b$ 定义为 $a+b$——根据加法分类、乘法分步的原理定义。那么转移方程可以写成这样:

  • $g(p_i) = Lf(p_i) \otimes f(p_{i+1})$
  • $f(p_i) = Lf(p_i) \otimes f(p_{i+1}) \oplus v_i \otimes Lg(p_i) \otimes g_(p_{i+1})$ // 这里将 $g(p_i)$ 代入了

这样,我们的转移就同样是“线性变换”,可以用矩阵来维护矩阵乘法了。算法的正确性可以用将这种新定义运算的矩阵乘法展开后证明,也可以将这种算法类比为 Floyd 求最短路来理解——它和 Floyd 的实现其实是一样的。

神奇的子图

这就是论文题……

首先考虑树上的算法。记录 $f(i, k)$ 表示所有点都在 $i$ 中且 $i$ 的点度为 $k$ 时的连通子图最大权值和,记录 $g(i)$ 表示 $i$ 子树中所有点 $j$ 的 $\max f(j, * )$。转移的时候,$f(i, * )$ 是对子节点跑重量为 $1$、价值为 $\max f(p, k < K )$ 的背包,同时更新 $g(i)$。

可以注意到,背包可以写成卷积的形式,而卷积可以写成矩阵乘法的形式,因此使用矩阵来表示变换是完全可行的。维护向量 $(f(p_i, 0), f(p_i,1), ..., f(p_i, k), g(p_i), 1)$,这样不难写出转移矩阵,就可以进行转移了。使用线段树维护即可做到和论文中相同的复杂度。同样,结合圆方树可以得到本题的另一个满分算法,复杂度相同。

这个算法的细节比原来的算法要少非常多,常数也小一些。

UOJ 268

终于不是论文中提到的题目了……

题意经过简化后变成这样:给出一棵树,每个点有 $a_i,b_i$ 两个权值,定义一条链 $(x, y)$ 的权值为链上所有点的 $a_i$ 加上 $b_{\mathrm{LCA}(x, y)}$。要求单点修改 $a_i$、链上修改 $b_i$(修改都是加一个数或减一个数),动态维护权值最大的链。

可以类比《神奇的子图》记录状态 $f(i, k\in\{0,1,2\})$、$g(i)$,意义和转移类似《神奇的子图》,只是 $f(i, * )$ 更新 $g(i)$ 时要加上 $b_i$。单点修改 $a_i$ 容易支持,考虑如何支持链上修改 $b_i$。不妨将 $b_i$ 差分,即设 $b'_i = b_i - b_{fa(i)}$,那么我们同样可以用 $b'_i$ 写出转移方程,并且链上修改变成了单点修改,易于维护。

转化为单点修改后,虽然要进行 $4$ 个节点的单点修改,但单点修改的常数小于链上修改,因此这个算法常数优于之前的算法。

神奇的简单路径

又是论文中的题了……

题意简述:给出一张图,保证每个点双 $\le 8$,每次修改点权,或询问两点之间最长简单路径,或询问全局最长简单路径。

如果只有全局最长简单路径,那么就是《神奇的子图》的弱化版。现在考虑第一个询问。这种询问相当于在圆方树上做链上询问。如果是暴力的话,那么就是沿着链一个一个点 DP 上来。这种链上的 DP 显然也可以写成变换的形式,那么我们也只需要维护这条链上的变换和即可,具体实现同样用矩阵,复杂度和原算法相同。

要注意的是一个一个 DP 上来时有向上的一段和向下的一段,因此线段树中维护矩阵的积时也要维护一正一反两种信息。

总结

基于变换合并的树上动态 DP 算法,完全地采用了对 DP 本身进行优化的思想,思路较于一般的树上动态 DP 算法思路更加直观、正确性更加显然。它不仅可以通用地解决之前所有的问题,还可以获得更小的思维难度、更小的代码难度和更快的运行速度,是论文中介绍的原算法的一个很大的改进。

当然,这个算法还有一些局限性。希望大家能和我分享关于这方面的研究,希望在我们的共同努力下能有更多更好的算法出现。

用无旋转 Treap 来写 LCT

2017-02-17 11:18:51 By immortalCO

之前本猫做 http://uoj.ac/problem/173 的时候,发现如果将 Treap 的 size 记为整棵子树(而非当前层)的大小,复杂度好像是对的。之前过不了好像是建树时的 push_up 次数常数过大,我就把初始化改成先不 push_up,建好之后再 push_up,直接用优化成线性;后面改成每次修改完之后才 push_up 后,操作次数就只剩三百多万了。

然后就想用 Treap 来写 LCT,看一看能不能达到同样是 $O(n\log n)$ 的复杂度。本猫用的是非旋转、无权值的 Treap,利用 splitmerge 来维护序列,在 merge 时以 size 作为概率决定哪一个点在上面。

以下以 http://uoj.ac/problem/3 题为例。

结构体

struct Node
{
    Node *p, *l, *r;    // p:为了方便,记录节点父亲;l、r:左右儿子
    int val, max;    // val:自己的权值;max:当前重链上的最大权值
    uint cnt, sum;    // cnt:这个点的所有轻边的子树大小 + 1;sum:当前子树的 cnt 之和
    bool rev;  // 翻转标记
}

为了保证复杂度,本猫采用了“LCT 维护子树大小”的做法,将 cntsum 维护成子树大小,并据此作为 Treap 的随机种子。

access

这是 LCT 的最基本的操作。伪代码如下:

access(node):
    p = node 所在的重链
    q = 空链
    while p != null:
        # 设 b 为 node 的轻儿子
        将 p 所在的重链从 node 处拆分为 a、b 两段(node 是 a 段的末尾)

        # 维护轻边子树大小之和
        node.cnt += b.sum - q.sum
        更新 p 使得其正确维护信息

        # 设 q 为 node 的重儿子
        p = 连接 a、q 
        更新 p 使得其正确维护信息

        # 下一步循环
        q = p
        node = q 的重链顶端的父亲
        p = node 所在的重链

那么我们就要支持对 Treap 进行从某一个点开始拆分以及合并两个 Treap。显然合并就是 merge 操作。我们只需要实现一个自下而上的 split 操作,循环即可实现。注意这些操作中都要维护好每个节点的父亲。

其他 make_rootlinkcut 操作都能利用上述操作实现(link 只需要设置重链父亲并维护 cntcut 需要一次 split)。

效果

提交为 http://uoj.ac/submission/129275 。感觉上速度挺正常的,不比 Splay 版的慢多少(吐槽:同样的代码连续提交两次使用时间差了 0.7s!)。

为了验证我们维护子树大小确实能欧化复杂度,修改 cnt 的维护方式(修改后为 $O(n\log^2n)$)再次提交为 http://uoj.ac/submission/129302 。慢了一倍以上。

但仍然无法查明其复杂度是否低于 $O(n\log^2n)$。

有啥用

可能可以持久化?(然而还要解决 access 拼接次数是均摊的问题……)

可能可以用其他重量平衡树(如旋转版的 Treap)代替这里的非旋转 Treap,能解决更多 SAM + LCT 的题,或更好的解决 BZOJ 4234,甚至出到子树询问和仙人掌上?(之前有过块链 LCT 做法)

一种高效处理无修改区间或树上询问的数据结构(附代码)

2016-11-02 19:45:47 By immortalCO

问题描述

给出一个某种元素的序列 $a_1,a_2,...,a_n$,要求进行 $m$ 次询问,每一次是询问一段区间 $[l,r]$ 的某种支持结合律和快速合并的信息,要求在线。

这类问题比较通用,比如在 DP 的优化中就常常见到。

要求复杂度尽量优秀,假设合并复杂度为 $O(1)$。

以下的时间复杂度 $O(A) + O(B) + O(C)$ 表示预处理复杂度、单次询问的复杂度和空间复杂度。

各种性质的具体问题

可减信息,如区间和、区间异或和

直接用前缀和实现,复杂度 $O(n)+O(1)+O(n)$。

可重复贡献信息,如区间最值

如果序列很长,使用线段树即可,复杂度 $O(n)+O(\log n)+O(n)$。

如果询问很多,但序列不是特别长,可以用倍增的 RMQ,复杂度 $O(n\log n)+O(1) + O(n\log n)$。

如果序列很长,询问也很多,可以对序列线性建立笛卡尔树转为 LCA 问题,然后转为正负 1 RMQ,每 $\log n$ 分一段打表预处理,复杂度 $O(n)+O(1)+O(n)$。

仅满足支持结合律和快速合并,如最大子段和、区间的串哈希值

一般使用线段树实现,复杂度 $O(n) + O(\log n) + O(n)$。

一种新的做法

上面一般做法的问题是:询问复杂度不够优秀。

如果问题没有任何性质,一般只能使用线段树。然而在实际的问题(或所优化的问题)中,由于询问不受空间大小的限制,一般询问的次数是大大多于序列的长度的,比如有的莫队需要查询 $O(n\sqrt n)$ 次的区间最值,这样倍增 RMQ 就有了用武之地。因此,我们考虑能否像 RMQ 那样,在一般的问题上,以预处理的时间和空间,换取快速的询问

考虑一个询问 $[l,r]$。如果 $l=r$,则我们可以直接得出答案。否则,考虑 $[l,r]$ 在线段树上定位时,是在哪个线段树节点 $p$ 中同时递归到左右两个儿子(如果 $[l,r]$ 只定位为线段树上的一个节点,则视为同时递归进入两个儿子)。设该节点的区间 $[s, t)$,其中点为 $m$(两个一半的区间是 $[s,m)$ 和 $[m,t)$),那么区间一定是满足 $s\le l < m \le r < t$。如果我们对于每个 $l\le i < m$ 都预处理了 $[i,m)$ 的信息和 $prel(p,i)$,对于每个 $m \le i < r$ 都预处理了 $[m,i]$ 的信息和 $prer(p,i)$,那么所求的信息和就是 $prel(p,l)+prer(p,r)$。这一步是 $O(1)$ 的。现在唯一的问题是我们能否快速求出 $p$。

一般的线段树没有什么好的性质,但如果一棵线段树的值域是 $[0, 2^k-1)$,并且编号遵从堆式存储($i$ 的儿子为 $2i$ 和 $2i+1$),那么 $[l,r]$ 所定位的区间,就是 $[l,r]$ 所对应的两个节点的 LCA,也就是 $l,r$ 两个位置对应节点的编号的二进制数的 LCP!获得 $x,y$ 的二进制数的 LCP 只需要使用 x >> Log2[x ^ y](找出并丢弃第一个不同的二进制位及后面的位),其中 Log2 数组可以线性预处理,因此求出 $p$ 也只需要 $O(1)$ 的时间。这样,我们就得到了一种 $O(n\log n)+O(1)+O(n\log n)$ 的优秀算法。

我们称这一数据结构为猫树

应用

经典问题:区间最大子段和

给出一个序列,多次询问区间的最大子段和。

这是一个经典的模型。不同于经典做法,我们只需要记录 $prel$、$prer$ 为对应前(后)缀的最大子段和、最大前(后)缀和即可。

复杂度:$O(n\log n)+O(1)+O(n\log n)$。实测在 $n=m=200000,a_i\le 10^9$ 的情况下,此做法的运行时间接近经典做法(非递归线段树实现)的 $2\over 3$。

经典问题:NAND

给出一个序列,多次询问一个 $x$ 对一个区间中所有数按顺序依次 NAND 的结果。

NAND 没有结合律,因此我们要用一些经典的处理方式。NAND 是按位独立的,因此我们可以对每一位维护信息:如果这一位刚开始是 0,那么按顺序 NAND 了这个区间中的所有数后会变成什么;如果这一位是 1 那么会变成多少。用位运算可以优化为 $O(1)$ 的信息合并。

使用猫树,即可直接支持 $O(n\log n)+O(1)+O(n\log n)$。

经典问题:区间凸包

现在来看看不那么傻逼的题。

给出一个平面上的点的序列,每次询问一个区间的凸包的面积,保证横坐标单调递增,强制在线。

这题的凸包不是贡献独立的,因此我们必须求出这个凸包。经典做法是用线段树套可持久化平衡树维护凸包,在线段树的节点直接保存凸包,用可持久化线段树的合并保证复杂度,代码非常繁琐,复杂度为 $O(n \log^2n)+O(\log^2n)+O(n\log^2n)$。

结合猫树,只要我们可持久化地保存 $prel,prer$ 为区间的凸包,这样我们就只需要资磁对左右两边预处理的凸包二分出一个公切线,然后进行计算即可,并不需要直接合并两个平衡树维护的凸包。而这里的凸包由于只有 $x$ 单调的插入,用简单的可持久化数组(线段树或平衡树)即可直接维护,代码量很小。复杂度为 $O(n\log^2n)+O(\log n)+O(n\log^2n)$。

BZOJ4540

给出一个区间,每次询问一个区间 $[l,r]$ 的所有子区间的区间最小值之和。

根据猫树的思想,我们考虑分开来计算代价。其中,$[l,m)$ 和 $[m,r]$ 的答案可以用单调栈线性预处理,就只需要考虑跨越部分了。对 $[l,m)$ 中记录每个数到 $m$ 的区间最小值,考虑这个最小值贡献答案的区间的右端点 $r$,这个最右的右端点可以通过左右的归并得出。这样预处理之后,就能做到 $O(1)$ 的询问复杂度。这里的细节比较繁琐,不过多说明。

这个做法实际运行时间十分优秀,在 $n=100000,m=2\times 10^7$ 的数据可以 $0.3s$ 出解。

扩展?

优化 DP

如果某种 DP 需要按顺序决定,但有很多次的区间查询,我们可以用猫树结合二进制分组做到相同的复杂度。

在查询次数很多的情况下,此做法比直接线段树维护具有优势。

树上无修改询问

考虑把猫树扩展到树上,以资磁快速的链上信息询问。我们对树进行点分治,预处理所有点到重心的信息(需要保存“包括中心”和“不包括重心”两种信息,如果有方向要求要专门考虑),这样同样能拆成两个已经预处理的信息的和。

询问的时候,我们需要找出这条链上最上面的重心。对点分治建立点分树(把重心连成树),然后我们要找到的重心就是询问的链端点的 LCA,这可以用倍增 RMQ 来预处理。

复杂度和区间问题相同,$O(n\log n)+O(1)+O(n\log n)$。

其实 BZOJ4568 是这种做法的离线版本。

UPD:代码参考

区间最大子段和
    const int MaxN = 200010;
    const ll inf = 1000000000000000000ll;
    int val[MaxN], pos[MaxN];
    struct Node
    {
        struct Data {ll pre, sum;} *dl, *dr;
        int mov_l, mov_r;
        IL ll query(RG int l, RG int r)
        {
            return dl[l -= mov_l].sum > dr[r -= mov_r].sum 
                ? dmax(dl[l].sum, dl[l].pre + dr[r].pre)
                : dmax(dr[r].sum, dl[l].pre + dr[r].pre);
        }
    } T[524288];
    int N;
    void build(RG int i, RG int l, RG int r)
    {
        if(l + 1 == r || N <= l) return;
        RG int m = (l + r) >> 1;
        build(i << 1, l, m);
        build(i << 1 | 1, m, r);

        RG Node *o = T + i;
        o->dl = new Node::Data[m - l], o->mov_l = l;
        o->dr = new Node::Data[r - m], o->mov_r = m;
        RG ll pre, max_pre, sum, max_sum;
        pre = sum = 0, max_pre = max_sum = -inf;
        for(RG int k = m - 1; k >= l; --k)
        {
            pre += val[k],                     cmax(max_pre, pre);
            sum = dmax(sum, 0) + val[k],    cmax(max_sum, sum);
            o->dl[k - l] = (Node::Data) {max_pre, max_sum};
        }
        pre = sum = 0, max_pre = max_sum = -inf;
        for(RG int k = m; k <= r - 1; ++k)
        {
            pre += val[k],                     cmax(max_pre, pre);
            sum = dmax(sum, 0) + val[k],    cmax(max_sum, sum);
            o->dr[k - m] = (Node::Data) {max_pre, max_sum};
        }
    }
    int pre[1024];

    IL void main()
    {
        RG int (*F)() = read<int>;
        RG int n = N = F();
        Rep(i, 0, n) val[i] = F();
        RG int D = 1; while(D < n) D <<= 1;
        build(1, 0, D);
        RG int l, r, d, v;
        pre[0] = 0;
        Rep(i, 1, 1024) pre[i] = pre[i >> 1] + 1;
        Rep(_, 0, F())
        {
            l = F() - 1, r = F() - 1;
            if(l == r) Print(val[l]);
            else
            {
                l += D, r += D;
                v = l ^ r;
                d = (v < 1024 ? pre[v] : 10 + pre[v >> 10]);
                Print(T[l >> d].query(l - D, r - D));
            }
        }
        io::flush();
    }
BZOJ4540
    const int MaxN = 131072;
    int a[MaxN], N;
    // 对单调栈中的数进行归并,并确定每一位对应的单调栈中元素
    struct Node
    {
        int m;
        struct Data
        {
            int pos;            // 在左右单调栈中的排名
            ll sum, ans;        // 答案的和
            ll pre;                // 到中点这一段的和
        } *pre;
        IL ll query(RG int l, RG int r)
        {
            RG int pl = pre[l].pos, pr = pre[r].pos;
            return pre[l].pre + pre[r].pre + ((pl < pr)
                ? pre[l].ans + (m - l)       * (pre[r].sum - (pl == m - l        ? 0 : pre[l + pl - 1].sum))
                : pre[r].ans + (r - m + 1) * (pre[l].sum - (pr == r - m + 1 ? 0 : pre[r - pr + 1].sum)));
        }
    } T[262144];
    void build(RG const int i, RG const int l, RG const int r)
    {
        if(l + 1 == r || N <= l) return;
        RG const int m = T[i].m = (l + r) >> 1;
        build(i << 1, l, m);
        build(i << 1 | 1, m, r);

        RG Node::Data *f = T[i].pre = (new Node::Data[r - l]) - l;

        static int stack[MaxN], b[MaxN];
        RG int top, cur, V; RG ll pre, ran, sum;
        pre = sum = top = ran = 0; cur = 2000000000; *stack = m;
        for(RG int i = m - 1; i >= l; --i)
        {
            V = a[i];
            for(; top && a[stack[top]] >= V; --top)
                ran -= (ll) a[stack[top]] * (stack[top - 1] - stack[top]);
            ran += (ll) V * (stack[top] - i);
            stack[++top] = i;
            cmin(cur, V); 
            f[i].sum = (sum += (b[i] = cur));
            f[i].pre = (pre += ran);
        }
        pre = sum = top = ran = 0; cur = 2000000000; *stack = m - 1;
        for(RG int i = m; i <= r - 1; ++i)
        {
            V = a[i];
            for(; top && a[stack[top]] >= V; --top)
                ran -= (ll) a[stack[top]] * (stack[top] - stack[top - 1]);
            ran += (ll) V * (i - stack[top]);
            stack[++top] = i;
            cmin(cur, V);
            f[i].sum = (sum += b[i] = cur);
            f[i].pre = (pre += ran);
        }
        sum = 0;
        RG int pl = m - 1, pr = m, N = 0, tmp;
        while(l <= pl || pr < r)
            (pr >= r || (l <= pl && b[pl] > b[pr]))
                ? (f[pl].pos = ++N, (tmp = (N - (m - pl)))         ? sum += (ll) b[pl] * tmp : 0, f[pl--].ans = sum)
                : (f[pr].pos = ++N, (tmp = (N - (pr - m + 1)))     ? sum += (ll) b[pr] * tmp : 0, f[pr++].ans = sum);
    }

    IL void main()
    {
        RG int (*F)() = read<int>;
        RG int n = N = F(), m = F(), D = 1;
        while(D < n) D <<= 1;
        Rep(i, 0, n) a[i] = F();
        build(1, 0, D);
        RG int l, r, v, d;
        static int pre[1024];
        Rep(i, 1, 1024) pre[i] = pre[i >> 1] + 1;
        while(m--)
        {
            l = F() - 1, r = F() - 1;
            if(l == r) Print(a[l]);
            else
            {
                l += D, r += D;
                v = l ^ r;
                d = (v < 1024 ? pre[v] : 10 + pre[v >> 10]);
                Print(T[l >> d].query(l - D, r - D));
            }
        }
        io::flush();
    }
代码的模板
//pb_ds 20160711
#include <cstdio>
#include <algorithm>
#include <cstring>
#include <cassert>
#include <cmath>
#include <iostream>
using namespace std;
#define RG register
#define IL __inline__ __attribute__((always_inline))
#define For(i, a, b) for(RG int i = a, ___u = b; i <= ___u; ++i)
#define Dor(i, a, b) for(RG int i = b, ___d = a; i >= ___d; --i)
#define Rep(i, a, b) for(RG int i = a, ___u = b; i != ___u; ++i)
#define dmin(a, b) ((a) < (b) ? (a) : (b))
#define dmax(a, b) ((a) > (b) ? (a) : (b))
#define cmin(a, b) ((a) > (b) ? (a) = (b) : 1)
#define cmax(a, b) ((a) < (b) ? (a) = (b) : 1)
#define ddel(a, b) ((a) > (b) ? (a) - (b) : (b) - (a))
#define dabs(i) ((i) > 0 ? (i) : -(i))
typedef long long ll;
typedef unsigned uint;
typedef unsigned long long ull;
typedef long double ld;

#include <queue>
#include <vector>

namespace pb_ds
{   
    // 输入输出优化模板从此开始
    namespace io
    {
        const int MaxBuff = 1 << 15;
        const int Output = 1 << 23;
        char B[MaxBuff], *S = B, *T = B;
        //#define getc() getchar()
        #define getc() ((S == T) && (T = (S = B) + fread(B, 1, MaxBuff, stdin), S == T) ? 0 : *S++)
        char Out[Output], *iter = Out;
        IL void flush()
        {
            fwrite(Out, 1, iter - Out, stdout);
            iter = Out;
        }
    }

    template<class Type> IL Type read()
    {
        using namespace io;
        RG char ch; RG Type ans = 0; RG bool neg = 0;
        while(ch = getc(), (ch < '0' || ch > '9') && ch != '-')     ;
        ch == '-' ? neg = 1 : ans = ch - '0';
        while(ch = getc(), '0' <= ch && ch <= '9') ans = ans * 10 + ch - '0';
        return neg ? -ans : ans;
    }
    template<class Type> IL void Print(RG Type x, RG char ch = '\n')
    {
        using namespace io;
        if(!x) *iter++ = '0';
        else
        {
            if(x < 0) *iter++ = '-', x = -x;
            static int s[100]; RG int t = 0;
            while(x) s[++t] = x % 10, x /= 10;
            while(t) *iter++ = '0' + s[t--];
        }
        *iter++ = ch;
    }
    // 输入输出优化模板到此结束

    /** 把上面的代码拷贝到这里即可运行 */
    /** 把上面的代码拷贝到这里即可运行 */
    /** 把上面的代码拷贝到这里即可运行 */
    /** 把上面的代码拷贝到这里即可运行 */
    /** 把上面的代码拷贝到这里即可运行 */
    /** 把上面的代码拷贝到这里即可运行 */
}

int main()
{
    #define FILE "dev"
    //freopen(FILE ".in", "r", stdin), freopen(FILE ".out", "w", stdout);
    pb_ds::main();
    fclose(stdin), fclose(stdout);  
}

完结撒花

指针和数组的代码哲学和一些数组的小技巧(雾)

2016-08-07 14:18:20 By immortalCO

一直很喜欢在 struct 时用指针。除了在速度上有优势(数组要做额外的加法),还有一点就是它的代码是符合中英文的语序的。

举个例子:

k->fa->val->max = 1;

如果要解读这段代码,我们可以这样:

设置 
    k 
        的 fa 
            的 val 
                的 max 
为
    1

因为指针的代码是符合我们常用语言的语序的。如果用数组,我们就必须这样写:

max[val[fa[k]]] = 1

强行解读的话,就变成:

设置
    在所有 max 中,是
        在所有 val 中,是
            在所有 fa 中,是
                k
            的那一个 fa
        的那一个 val
    的那一个 max
为
    1

简直佶屈聱牙,难以成颂!然而,一种特殊的 C++ 数组写法可以解决这个问题。在 C++ 中,a[b] 等价于 *(a + b),因此

a[b] = *(a + b) = *(b + a) = b[a]

也就是说,我们可以倒转数组名和下标的位置。比如我们可以

a[i] = i[a]
b[1] = 1[b]
c[2][2] = 2[c[2]] = 2[2[c]] // 每一个 [] 都有两种写法
d[1][2][3][4][5] = 5[4[3[2[1[d]]]]]
A[B[C[D[E[i]]]]] = i[E][D][C][B][A] // [] 是从左往右计算的

这个的第一个用途是混乱代码,防止在 CF 上被 Hack。试想如果你有个 5 维 DP,那么每一维你都有 2 种写法,总共就有 32 种写法,然后你在 DP 中混用这几种写法……(嘿嘿嘿!!!)

观察最后一行左边的那个形式,不正是我们 max[val[fa[k]]] 的形式吗?因此我们可以写成

max[val[fa[k]]] = k[fa][val][max]

嘿嘿嘿!熟悉的语序回来了!而且没有方括号的嵌套,减小了漏打、错打括号导致 CE 的可能性。那么在不方便使用指针(比如 64 位机卡内存)时,我们同样可以用数组做到优美的代码风格。

以下是数组版《树上两点距离》,用树链剖分求 LCA。在 query 时,语序十分自然。

#include <cstdio>
const int MaxN = 1000010;
int to[MaxN << 1], len[MaxN << 1], next[MaxN << 1], fir[MaxN];
void link(int x, int y, int v)
{
    static int tot = 1;
    ++tot, tot[next] = x[fir], tot[to] = y, tot[len] = v, x[fir] = tot;
    ++tot, tot[next] = y[fir], tot[to] = x, tot[len] = v, y[fir] = tot; 
}
int pre[MaxN], dep[MaxN], cnt[MaxN], cho[MaxN], top[MaxN];
void dfs_init(int i)
{
    i[cnt] = 1;
    for(int k = i[fir]; k; k = k[next])
        if(!k[to][cnt])
        {
            k[to][pre] = i;
            k[to][dep] = i[dep] + k[len];
            dfs_init(k[to]);
            i[cnt] += k[to][cnt];
            if(k[to][cnt] > i[cho][cnt])
                i[cho] = k[to];
        }
}
void dfs_make(int i)
{
    i[top] = (i == i[pre][cho]) ? i[pre][top] : i;
    for(int k = i[fir]; k; k = k[next])
        if(!k[to][top]) dfs_make(k[to]);
}
int query(int x, int y)
{
    int ans = x[dep] + y[dep];
    while(x[top] != y[top])
        x[top][dep] > y[top][dep]
            ? x = x[top][pre]
            : y = y[top][pre];
    ans -= (x[dep] < y[dep] ? x[dep] : y[dep]) << 1;
    return ans;
}

int main()
{
    int n, m, x, y, v; scanf("%d%d", &n, &m);
    for(int i = n; --i; ) scanf("%d%d%d", &x, &y, &v), link(x, y, v);
    dfs_init(1), dfs_make(1);
    while(m--) scanf("%d%d", &x, &y), printf("%d\n", query(x, y));
}

圆方树——处理仙人掌的利器

2016-08-01 09:11:36 By immortalCO

营员交流课件

链接: http://pan.baidu.com/s/1hrM1Ic8 密码: dc4j

概述

随着业界毒瘤的推动,我们开始见到越来越多的仙人掌题。对于这种题,我们总是想办法把仙人掌转化为树来做,从而将熟悉的树的算法修改后用于仙人掌上。

圆方树是一种由我们提出的树的结构,它能轻松的将仙人掌转化为树,并能资磁仙人掌上点分治、虚仙人掌、仙人掌剖分等经典的树上问题,在 UOJ 上若干道仙人掌题目中,均有极高的效率。

构造

首先我们要构造原仙人掌 $G$ 的圆方树 $T_G$。我们定义原仙人掌上的点为圆点,现在我们考虑转化。

我们对仙人掌进行点双的 Tarjan 算法,因为仙人掌上每个环都是一个点双,而且在栈中的顺序就是环上点的顺序。那么考虑一个点 $i$ 的一条出边 $(i,j)$ 满足 $dfn(i) < low(j)$,那么说明 $(i, j)$ 是一条树边,我们在 $T_G$ 中直接连上这两个点即可;如果 $dfn(i) = low(j)$,那么我们找到了一个环(可能是重边造成的二元环),则从栈中取出点直到取出 $j$ 为止,设这样从栈中取出的点集为 $R$,则 $R \cup \{i\}$ 是一个环;对于其他情况,我们按照 Tarjan 的步骤执行即可,无需特殊处理。

对于一个环 $R \cup \{i\}$,设 $i$ 为环的 。我们新建一个方点 $P$,在 $T_G$ 中连接所有的 $(P, q), q\in R \cup \{i\}$。这样连的边数等于环上点数,而点又多了一个,因此最后连出的 $T_G$ 是树。

性质

圆方树有优美的性质。

  1. 无论如何换根,圆方树形态不变 //因为你是把环连成菊花而不是别的什么

  2. 圆方树的子树 = 原仙人掌的子仙人掌 //分类讨论各种情况可以证明

  3. 方点不会和方点相连

还有其他神奇的性质,要具体题目具体分析。

应用

最短路

类比树上最短路,我们在 LCA 处考虑。然而圆方树的树上距离并不一定就等于仙人掌上的最短路长度。

注意到圆方树中的边权都还没有确定,那么我们就可以现在确定了。如果一条边是圆圆边(圆点指向圆点的边,这里已经任意选择一圆点为根),那么它是一个树边,边权等于原仙人掌上边权;如果是方圆边,那么边权等于环的根到那个圆点的最短路径长度;如果是圆方边,那么边权为 0。

这样,如果一对询问点在圆方树上的 LCA 是圆点,那么其圆方树上长度就是原仙人掌上长度,因为路径上唯一的不同就在于对于每个环是走哪一侧(这个决策对每个环是独立的),而前面边权的确定已经解决了这一问题;如果 LCA 是方点,则我们只需要考虑在 LCA 这个环中是否要走经过根的那一侧。我们需要找出这两个询问点 $(x, y)$ 分别是在这个方点的哪两个儿子 $(p_x, p_y)$ 的子树中(这就是个 jump 操作,可以用倍增或树链剖分轻松实现)。则 $x \rightarrow p_x$ 和 $p_y \rightarrow y$ 的树上距离都是最短路径,现在就只需要考虑同一个环上的点 $(p_x, p_y)$ 的最短路径(其实也就是走哪一侧的问题)。 前面先记录环上边权的前缀和,这样就能 $O(1)$ 计算出不经过根那条路径的长度,再通过记录整个环的边权和比较出走哪一侧更优。

例题:http://www.lydsy.com/JudgeOnline/problem.php?id=2125

DP

仙人掌有若干种 DP,比如最大独立集和直径。

最大独立集比较简单,只需要 DFS 时,树边按照树上独立集的方式转移,环上用一个 $f(i, 0..1, 0..1)$ 的简单 DP 去转移即可。

例题:http://www.lydsy.com/JudgeOnline/problem.php?id=4316

求直径的话,我们类比求树直径的方式,在 LCA 处考虑。设 $f(i, j \in {1,2})$ 为 $i$ 子树中的第 $j$ 深的点,那么圆点可以轻松转移、更新答案。在方点转移也容易,但更新答案不能直接更新(因为有哪一侧的问题,直接取不一定是最短路,会导致答案偏大)。我们可以把环复制一份用单调队列来解决,也可以正反各做一次前缀和。

例题:http://www.lydsy.com/JudgeOnline/problem.php?id=1023

虚仙人掌

对一个点集的子集进行询问。分例题进行讲解。

UOJ 87

http://uoj.ac/problem/87

题意:给出一个点集,求点集中两两最短路径的最大值。

首先对原仙人掌建立 BZOJ 2125 那种带边权的圆方树。然后每次对圆方树建立虚树,在虚树上按照 BZOJ 1023 的方法进行简单 DP 即可。

UOJ 189

http://uoj.ac/problem/189

题意:给出仙人掌,要求资磁:(1)询问若干条最短、最长路径的并的边权和 (2)修改边权

先对所有端点的点建立圆方树的虚树。现在我们要考虑两部分的贡献:(1)对于虚树上一条边,对应圆方树上一条链,我们要一起考虑这条链的贡献;(2)对于虚树上的一个方点,对应一个环,我们要考虑它有哪几段计入答案。

我们对于每条边这么操作:首先找到 LCA,判断出是哪一段符合要求(最长或最短),然后像求圆的周长(某些计算几何题)那样记录两个括号用来扫描线;然后找出不在环上的一段(如果 LCA 是方点,要利用 jump 操作求出 $p_x, p_y$),利用子树前缀和的思想记录一下是否存在最长、最短路的询问(也就是在 $x$ 处 $+1$,在 $p_x$ 处 $-1$,然后一个点的权值为其子树中的权值之和)。

那么对于一条虚树上的边,我们就有 4 种情况,表示最长路和最短路分别是否有。我们可以通过一种特殊的“深度”来计算这些值。设 $f(i)$ 为从圆方树根到 $i$ 的树边长度之和,$g(i)$ 为从根到 $i$ 的最短路径,$h(i)$ 为从根到 $i$ 的最长路径。那么如果经过边$(i, j)$的又有最长路径又有最短路径,则贡献为 $g(j) + h(j) - f(j) - (g(i) + h(i) - f(i))$;如果只有最短路径则贡献为 $g(i) - g(j)$,只有最长路径贡献为 $h(i) - h(j)$。

环上的情况,我们像求圆上周长那样,对之前添加的括号排个序即可。

如果要修改,那么既有可能影响环上距离,又有可能影响那 3 种深度。我们分别用树状数组维护均可。注意如果修改的是环上的边,那么影响的深度数组到底是哪几段需要分类讨论。

点分治

UOJ 23

http://uoj.ac/problem/23

题意:求出仙人掌上从 1 出发长度为 $[1..n]$ 的简单路径的条数。

暴力是用背包的思想,设 $f(i, j)$ 表示从 $i$ 出发往下走 $j$ 步的路径条数。注意到转移是卷积的形式,那么我们可以用 FFT 来解决这一问题。

进行点分治。设当前子树中重心为 $G$,考虑求出 $f(G)$ 为这棵子树的背包数组。首先 $G$ 的包含当前连通块最高点(记为 $top(G)$)这一分支可以递归去做。现在考虑 $G$ 的儿子们,如果 $G$ 是圆点,那么直接把出边的 $f$ 加起来求和即可;否则如果是方点,那么每个儿子都有两种选择:向左或向右,两种情况都要加上去。现在还需要考虑 $top(G)$ 到 $G$ 这一段对 $G$ 出边贡献的 $f$ 的影响。显然我们可以再套一个分治 FFT 把 $top(G)$ 到 $G$ 这一条链上的数组算出来,但这样太笨了。我们维护 $g(G)$ 表示 $top(G)$ 到 $G$ 这一段的背包数组,那么我们发现从 $fa(G)$ 开始一直跳 $top$,一段一段跳上去,不会超过 $\log(top(G)\ 到\ G\ 的链长)$ 段,而且后一段的 $g$ 数组长度是前一段的两倍。这样暴力 FFT 卷积起来,复杂度只有 $top(G)\ 到\ G\ 的链长 \times \log(top(G)\ 到\ G\ 的链长)$,那么这个问题就在 $O(n \log^2 n)$ 的时间内完美解决。

一个细节是设方点点权为环的大小,点分治要找带权重心。

这种分治就是我们熟悉的树分治,代码难度减小了不少。

链剖分

UOJ 158

http://uoj.ac/problem/158

题意:给出仙人掌要求资磁 3 种操作:(1)一个点到根的最短路径上取反(2)一个点到根的最长路径上取反(3)子树内黑点个数询问。

树上肯定是树剖。仙人掌的话,我们可以对圆方树进行树剖。

考虑怎么支持我们的操作。首先由于要支持跳重链,我们必须能高效修改一整条重链上最长或最短路上的点。首先必须经过的点(即割点)肯定是要取反的,然而环上的点还得分类考虑,有时只要修改一部分点。

那么我们就可以把点进行分类。考虑每一条重链,将点分成 3 类:(1)割点(2)在重链上作为最短路径出现的点(3)在重链上作为最长路径出现的点。那么显然对于第 1 类点,直接在树上修改即可。然而对于后两种点,树上路径不一定会扫过需要更改的点。

我们考虑这样一个事情:在 DFS 序中,如果一个点是方点,那么我们先访问它的所有儿子,然后再按照重边先行的顺序访问儿子的子树(不包括那些儿子了),这样的话,我们修改一段区间的时候,就会扫过这个环上所有点;而我们已经把点分类了,因此也可以选择只修改最短路径或最长路径上的点。

那么这个问题就完美解决了。写起来比较码农,在进入重链时需要特殊处理。

【UNR #1】果冻运输 改 checker 记

2016-07-17 21:29:56 By immortalCO

题目地址

http://uoj.ac/problem/215

考场上的做法

数据这么小,显然要么手玩要么爆搜咯。

爆搜模拟好难写啊!!!

然后就就纯手玩了 3.5h,玩的欲仙欲死才 47 分。

再次乱搞

试试爆搜。

这次出题人良心福利,checker 居然是发源代码,而且代码挺可读的。那就研究一下。

move(x, y, z) 是执行移动命令,print() 是输出地图,check() 是判断是否完成。

改改改! move 改成返回值为 bool 表示是否移动成功,check() 改成返回里面的变量 cnt 表示连通块个数。只有 a[N][N] 是地图内容,那就弄个 struct 包起来。

check() 的返回值不错,可以当估价,那就写个 A* 吧。实现 hashcode() 表示哈希函数,以及 operator ==operator < 作为哈希和堆的操作符。哈希表不知道开多少,那就分块分配,每次 new 1024 个。

开跑!绝大多数点在 2min 内都跑了出来。16 一直 RE,把估价函数改成返回 0/1 就出来了。只有 17 跑了巨久,跑了 40min 才出解,用掉了我 5G 内存。

交一发,发现只有 99 分,第 19 个点少了 1 分。然而我的步数是 19,而 AC 提交的步数是 29,不科学啊。 把估价改成 0,就跑出了 18,于是就 AC 啦!

提交记录

http://uoj.ac/submission/84887

修改版 checker

用 clang-format 给代码加了空格,不喜勿喷

FILE 里面的数字改掉就可以跑不同的测试点了。

#include <bits/stdc++.h>
using namespace std;
#define X first
#define Y second
#define mp make_pair
#define pb push_back
#define rep(i, a, n) for (int i = (a); i <= (n); ++i)
#define SZ(x) ((int)(x).size())
typedef pair<int, int> pii;
const int DX[] = {-1, 1, 0, 0}, DY[] = {0, 0, 1, -1};
const int N = 32;
const int Base = 233;
const int mod = 3333331;
int n, m;
struct State;
struct Record{State *last; int k0, k1, k2;};

int tmp[N][N], sv[N][N];
set<int> S[10];
queue<pii> Q;
int bo;
void get(int &p) 
{
    int x;
    for (x = getchar();
             (x < '0' || x > '9') && x != '|' && x != '-' && x != '.' && x != ' ';
             x = getchar())
        ;
    if (x == '|' || x == '-')
        p = 1;
    else if (x == '.' || x == ' ')
        p = 0;
    else if (x == '0')
        p = -1;
    else
        p = x - '0';
}
struct State
{
    int a[N][N], p, e;
    Record rec;
    void dfs1(int x, int y, int z) 
    {
        if (x < 1 || x > 2 * n - 1 || y < 1 || y > 2 * m - 1) return;
        if (!a[x][y]) return;
        if (a[x][y] == -1) bo = 0;
        tmp[x + DX[z] * 2][y] = a[x][y];
        a[x][y] = 0;
        rep(i, 0, 3) dfs1(x + DX[i], y + DY[i], z);
        if (a[x + DX[z] * 2][y]) Q.push(mp(x + DX[z] * 2, y));
    }
    void dfs2(int x, int y) 
    {
        if (x < 1 || x > 2 * n - 1 || y < 1 || y > 2 * m - 1) return;
        if (!a[x][y]) return;
        if (tmp[x][y]) return;
        tmp[x][y] = 1;
        rep(i, 0, 3) dfs2(x + DX[i], y + DY[i]);
        dfs2(x, y + 2);
    }
    bool move(State *last, int x, int y, int z) 
    {
        rec = (Record) {last, x, y, z}; ++p;
        if (!a[2 * x - 1][2 * y - 1]) return 0;
        bo = 1;
        rep(i, 1, 2 * n - 1) rep(j, 1, 2 * m - 1) tmp[i][j] = 0, sv[i][j] = a[i][j];
        while (Q.size()) Q.pop();
        dfs1(x * 2 - 1, y * 2 - 1, z);
        while (Q.size()) 
        {
            pii k = Q.front();
            Q.pop();
            dfs1(k.X, k.Y, z);
        }
        if (!bo) 
        {
            rep(i, 1, 2 * n - 1) rep(j, 1, 2 * m - 1) a[i][j] = sv[i][j];
            return 0;
        }
        rep(i, 1, 2 * n - 1) rep(j, 1, 2 * m - 1) if (tmp[i][j]) a[i][j] = tmp[i][j];

        while (1) 
        {
            rep(i, 1, 2 * n - 1) rep(j, 1, 2 * m - 1) tmp[i][j] = 0;
            rep(i, 1, 2 * n - 1) rep(j, 1, 2 * m - 1) if (a[i][j] == -1) dfs2(i, j);
            int bo = 0;
            rep(j, 1, 2 * m - 1) rep(i, 1, 2 * n - 1) if (a[i][j] && !tmp[i][j])
                    bo = 1,
                    a[i][j - 1] = a[i][j], a[i][j] = 0;
            if (!bo) break;
        }

        rep(i, 1, 2 * n - 1) rep(j, 1, 2 * m - 1) if ((i & 1) && (j & 1))
                rep(k, 0, 3) if (a[i][j] >= 1 && a[i][j] <= 5 && a[i][j] == a[i + DX[k] * 2][j + DY[k] * 2]) 
                    a[i + DX[k]][j + DY[k]] = 1;
        return 1;
    }
    void print(FILE *x = stdout) 
    {
        fprintf(x, "p = %d e = %d\n", p, estimate());
        rep(i, 1, 2 * m - 1) 
        {
            rep(j, 1, 2 * n - 1) if ((i & 1) && (j & 1)) 
            {
                if (a[j][2 * m - 1 - i + 1] == 0)
                    fprintf(x, ".");
                else if (a[j][2 * m - 1 - i + 1] == -1)
                    fprintf(x, "0");
                else
                    fprintf(x, "%c", a[j][2 * m - 1 - i + 1] + '0');
            }
            else 
            {
                fprintf(x, "%c", a[j][2 * m - 1 - i + 1]
                    ? i & 1 ? '-' : '|'
                    : (i & 1) && (j & 1) ? '.' : ' ');
            }
            fprintf(x, "\n");
        }
    }
    void dfs3(int x, int y, int z) const
    {
        if (x < 1 || x > 2 * n - 1 || y < 1 || y > 2 * m - 1) return;
        if (a[x][y] < 1 || a[x][y] > 5) return;
        if (tmp[x][y]) return;
        tmp[x][y] = 1;
        S[a[x][y]].insert(z);
        rep(i, 0, 3) if (a[x + DX[i] * 2][y + DY[i] * 2] == a[x][y])
                dfs3(x + DX[i] * 2, y + DY[i] * 2, z);
    }
    int check()
    {
        rep(i, 1, 5) S[i].clear();
        int cnt = 0;
        rep(i, 1, 2 * n - 1) rep(j, 1, 2 * m - 1) tmp[i][j] = 0;
        rep(i, 1, 2 * n - 1) rep(
                j, 1, 2 * m - 1) if ((i & 1) && (j & 1) && a[i][j] >= 1 && a[i][j] <= 5)
                dfs3(i, j, ++cnt);
        cnt = 0;
        rep(i, 1, 5) if (SZ(S[i])) cnt += SZ(S[i]) - 1;
        return e = cnt;
    }
    void read_map()
    {
        p = 0;
        scanf("%d%d", &m, &n);
        rep(i, 1, 2 * m - 1) rep(j, 1, 2 * n - 1) get(a[j][2 * m - 1 - i + 1]);
        rep(i, 1, 2 * n - 1) rep(j, 1, 2 * m - 1) if ((i & 1) && (j & 1))
            rep(k, 0, 3) if (a[i][j] >= 1 && a[i][j] <= 5 && a[i][j] == a[i + DX[k] * 2][j + DY[k] * 2]) 
                a[i + DX[k]][j + DY[k]] = 1;
    }
    int estimate() const {return p + e;} 
    //第 16 和第 19 个点会出问题,要把 + e 删掉。
    bool operator == (const State& S) const
    {
        rep(i, 1, 2 * n - 1) rep(j, 1, 2 * m - 1)
            if(a[i][j] != S.a[i][j]) return 0;
        return 1;
    }
    int hashcode() const
    {
        int ans = 0;
        rep(i, 1, 2 * n - 1) rep(j, 1, 2 * m - 1)
            ans = (ans * Base + a[i][j] + 23) % mod;
        return ans;
    }
    void output()
    {
        if(rec.last) 
        {
            rec.last->output();
            printf("%d %d %d\n", rec.k0, rec.k1, rec.k2);
        }
    }
} Init, Tmp;

struct Hash{State S; Hash *next;} *fir[mod];
State* find(const State& S)
{
    int h = S.hashcode();
    for(Hash *k = fir[h]; k; k = k->next)
        if(k->S == S) return 0;
    static Hash *mem, *top;
    if(mem == top) top = (mem = new Hash[1024]) + 1024;
    *--top = (Hash) {S, fir[h]}, fir[h] = top;
    return &top->S;
}

struct Comp
{
    bool operator() (State *a, State *b) const
    {
        return a->estimate() > b->estimate();
    }
};
priority_queue<State*, vector<State*>, Comp> heap;

int main() 
{
    #define FILE "C20"
    freopen(FILE ".in", "r", stdin);
    Init.read_map();
    heap.push(find(Init));
    int counter = 0;
    while(heap.size())
    {
        State *S = heap.top(); heap.pop();
        Tmp = *S;
        printf("new State %d S = %d p = %d e = %d\n", ++counter, S->hashcode(), S->p, S->estimate());
        for(int i = 1; i <= n; ++i)
            for(int j = 1; j <= m; ++j)
                for(int k = 0; k != 2; ++k)
                {
                    Tmp = *S;
                    if(!Tmp.move(S, i, j, k)) continue;
                    State *T = find(Tmp);
                    if(!T) continue;
                    if(!T->check())
                    {
                        puts("Found!");
                        T->print();
                        freopen(FILE ".out", "w", stdout);
                        T->output();
                        return 0;
                    }
                    heap.push(T);
                }
    }
    return 0;
}

UOIP十合一乱搞记+详细题解附代码

2016-07-06 22:33:03 By immortalCO

题目地址

http://uoj.ac/problem/208

题意简述

给出一张有向图,求有多少边的子集是DAG。

开头

为啥会做这一题呢?

几周前,机房高一神犇 WrongAnswer 在这题刚出来不久就秒掉了它。我也看了看,觉得比较interesting,而且我的数数题一向很弱,就决定开始搞了。从此入坑,一发不可收拾。

真是好题啊!搞了我整整一周。不过,其中的一些做法还是非常有趣的。

分点题解+日志

Case 1

打开1。卧槽,n和m都这么大,暴力肯定跑不出来。观察了一下,发现似乎前面几条边的左端点都不大于右端点?写个程序assert一发,确实如此。

那就很容易了,把自环扔掉,剩下k条边,答案就是$2^k$。跑了一下,过check了,就扔一边了。

#include <bits/stdc++.h>
using namespace std;
const int mod = 1000000007;

int main()
{
    assert(freopen("uoip1.in", "r", stdin));
    assert(freopen("uoip1.out", "w", stdout));
    int n, m;
    cin >> n >> m;
    int ans = 1;
    while(m--)
    {
        int x, y;
        cin >> x >> y;
        if(x != y)
        {
            assert(x < y);
            ans = 2ll * ans % mod;
        }
    }
    cout << ans;
    cerr << "check = " << (ans ^ 597855228) % 181 << endl;
}

Case 3

看了一下2,感觉没啥性质,就扔一边先去看3了。

一眼n=m。环套树即视感。仔细看一下,每条边左端点都是i。果然是环套树。写个程序assert一下,发现果然是环套树森林。

那么显然找环。不在环上的边,删不删没啥关系。环上的边就至少得删一条。于是就求出所有环大小,答案就是$2^{不在环上的边}\prod(2^{环的大小}-1) $咯。刚开始没想清楚WA了好几次,改清楚就过了。

#include <bits/stdc++.h>
using namespace std;
#define IL __inline__ __attribute__((always_inline))
#define RG register
const int mod = 1000000007;
const int MaxN = 100010;

namespace pb_ds
{
    int v[MaxN], timer;
    int to[MaxN], dep[MaxN];
    int ring;
    void dfs(int i)
    {
        v[i] = timer; if(!v[to[i]]) return (dep[to[i]] = dep[i] + 1, dfs(to[i]));
        if(v[to[i]] == timer) ring = (dep[i] - dep[to[i]] + 1); 
        else ring = -1;
    }
}

int main()
{
    using namespace pb_ds;
    assert(freopen("uoip3.in", "r", stdin));
    assert(freopen("uoip3.out", "w", stdout));

    int n, m;
    cin >> n >> m; assert(m == n);
    for(int i = 1; i <= n; ++i)
    {
        int x; cin >> x >> to[i];
        assert(x == i);
        assert(to[i] != i);
    }
    int ans = 1, other = n;
    for(int i = 1; i <= n; ++i) if(!v[i]) 
    {
        ++timer;
        dfs(i);
        if(ring == -1) continue;
        //printf("i = %d ring = %d\n", i, ring);
        other -= ring;
        int tmp = 1; while(ring--) tmp = (tmp << 1) % mod;
        ans = (long long) ans * (tmp - 1) % mod;
    }
    //printf("other = %d\n", other);
    while(other--) ans = (ans << 1) % mod;

    cout << ans;
    cerr << "check = " << (ans ^ 597855228) % 181 << endl;
}

Case 2 & 4

看4。一眼还是没啥性质。拉到下面去,似乎左右端点的宽度都差不多(我承认我是这么看出性质的)?写个程序验证一下,发现左右端点差不超过7。 想到2,发现左右端点同样不超过7。原来两个测试点能用一样的方法做。

这性质怎么用呢?貌似可以按顺序DP。考虑一条边(x,y),它能被加入当且仅当原图中不存在一条从y到x的路径;而由于端点距离不超过7,我们只需要考虑这个点往前7个点(包括本身,共8个点)的连通性。那我们就可以以最后8个点的传递闭包为状态进行DP。

我们可以用unsigned long long恰好64位整数来状压。然而这样状态很多,而我们知道合法的状态很有限(比如如果(1,2)和(2,3)为1,则(1,3)就必须为1)。我懒得预处理所有状态和转移,就直接用map<ull, int>来记录状态,开两个做滚动数组。每个位置还要编号前移,不过这也只是唯一的细节。

刚开始一直WA,后来发现1 << x的时候没用1ll。我很弱。

然后就把2和4都跑过了。2闪出,4开了O2只跑了8秒左右。

#include <bits/stdc++.h>
using namespace std;
#define IL __inline__ __attribute__((always_inline))
#define RG register
const int mod = 1000000007;
const int MaxN = 4010;
vector< pair<int, int> > e[MaxN];
#define pb push_back
#define mkp make_pair
#define fir first
#define sec second
typedef unsigned long long ull;

namespace pb_ds
{
    IL int inc(RG int x, RG int v) {x += v; return x >= mod ? x - mod : x;}
    IL void add(RG int& x, RG int v) {x = inc(x, v);}
    //number from 0 to 7
    //chuan di bi bao
    const ull empty = 9241421688590303745ull;
    IL bool test(RG const ull& s, RG int x, RG int y) {return s >> (x << 3 | y) & 1;}
    IL void set(RG ull& s, RG int x, RG int y) {s |= 1ull << (x << 3 | y);}
    IL int out_set(RG const ull& s, RG int i) {return (s >> (i << 3)) & 255;}
    IL void print(RG const ull& s, RG int n = 8)
    {
        for(RG int i = 0; i != n; ++i)
            for(RG int j = 0; j != n; ++j)
                printf("%d%c", test(s, i, j), " \n"[j + 1 == n]);
    }
    IL ull add_edge(RG ull s, RG int x, RG int y)
    {
        if(test(s, y, x)) return 0;
        for(RG int i = 0; i != 8; ++i) if(test(s, i, x))
            s |= (ull) out_set(s, y) << (i << 3);
        return s;
    }
    IL ull move_next(RG const ull& s)
    {
        //(x, y) -> (x + 1, y + 1)
        RG ull ret = empty;
        for(RG int i = 0; i + 1 != 8; ++i)    
            for(RG int j = 0; j + 1 != 8; ++j)
                if(test(s, i, j)) set(ret, i + 1, j + 1);
        return ret;
    }

    map<ull, int> Q[2];
    template<class T> void reinit(RG T& t) {t.~T(); new (&t)T; }
}

int main()
{
    using namespace pb_ds;
    assert(freopen("uoip4.in", "r", stdin));
    assert(freopen("uoip4.out", "w", stdout));
    int n, m;
    cin >> n >> m;
    while(m--)
    {
        int x, y;
        cin >> x >> y;
        assert(fabs(x - y) <= 7);
        if(x == y) continue;
        e[max(x, y)].pb(mkp(x, y));
    }
    RG bool d = 0;
    Q[d][empty] = 1;
    RG map<ull, int>::iterator iter, end;
    #define formap(M) for(iter = M.begin(), end = M.end(); iter != end; ++iter)
    RG ull next;
    for(RG int p = 1; p <= n; ++p)
    {
        reinit(Q[d ^= 1]);
        formap(Q[!d]) add(Q[d][move_next(iter->fir)], iter->sec);
        for(RG unsigned i = 0; i != e[p].size(); ++i)
        {
            RG int x = p - e[p][i].fir, y = p - e[p][i].sec;
            reinit(Q[d ^= 1]);
            formap(Q[!d]) 
            {
                add(Q[d][iter->fir], iter->sec);
                next = add_edge(iter->fir, x, y);
                if(next) add(Q[d][next], iter->sec);
            }
        }
    }
    RG int ans = 0;
    formap(Q[d]) add(ans, iter->sec);
    cout << ans;// << endl;
    cerr << "check = " << (ans ^ 597855228) % 181 << endl;
}

Case 5

看起来并没有啥规律……

我相信出题人一定隐藏了什么奇怪的性质,毕竟这才第5个点,不可能用啥高端算法吧。用Mathematica画了个图,没看出啥性质。写了个找连通块的程序,每个连通块也不小。

最后我写了一个找强连通的程序,发现大多数强连通分量都是1,最大一个只有4。于是就可以每个强连通分量暴力枚举边,每个强连通分量方案相乘,强连通分量之间的边直接用$2^{个数}$乘在答案上即可。

这个点的性质真神奇…………………

#include <bits/stdc++.h>
using namespace std;
#define IL __inline__ __attribute__((always_inline))
#define RG register
const int mod = 1000000007;
const int MaxN = 1010;
#define pb push_back
#define mkp make_pair
#define fir first
#define sec second

namespace pb_ds
{
    vector<int> e[MaxN];
    int dfn[MaxN], low[MaxN], timer;
    int stack[MaxN], top;
    int bel[MaxN], bid[MaxN], total;
    vector<int> blo[MaxN];

    void dfs(RG int i)
    {
        dfn[i] = low[i] = ++timer;
        stack[++top] = i;
        for(RG unsigned j = 0; j != e[i].size(); ++j)
            if(!dfn[e[i][j]])
            {
                dfs(e[i][j]);
                low[i] = min(low[i], low[e[i][j]]);
            }
            else if(!bel[e[i][j]])
                low[i] = min(low[i], dfn[e[i][j]]);
        if(low[i] == dfn[i])
        {
            RG int tmp;
            ++total;
            do 
            {
                tmp = stack[top--];
                bel[tmp] = total;
                bid[tmp] = blo[total].size();
                blo[total].pb(tmp);
            }
            while(tmp != i);
        }
    }

    struct Edge{int x, y;};
}

int main()
{
    using namespace pb_ds;
    assert(freopen("uoip5.in", "r", stdin));
    assert(freopen("uoip5.out", "w", stdout));

    int n, m;
    cin >> n >> m;
    for(int i = 1; i <= m; ++i)
    {
        int x, y; cin >> x >> y;
        if(x == y) {--i, --m; continue;}
        e[x].pb(y);
    }
    for(int i = 1; i <= n; ++i) if(!dfn[i]) dfs(i);
    int ans = 1;
    for(int B = 1; B <= total; ++B)
    {
        //cout << blo[B].size() << endl;
        static Edge edge[MaxN];
        int M = 0, size = blo[B].size();
        for(int i = 0; i != size; ++i)
        {
            int p = blo[B][i];
            for(int j = 0; j != (int) e[p].size(); ++j)
                if(bel[e[p][j]] == B)
                {
                    edge[M++] = (Edge) {i, bid[e[p][j]]};
                }
        }
        m -= M;
        int cur = 0;
        for(int s = 0; s != (1 << M); ++s)
        {
            static bool v[20][20];
            for(int i = 0; i != size; ++i)
                memset(v[i], 0, size), v[i][i] = 1;
            for(int i = 0; i != M; ++i)
                if(s >> i & 1) v[edge[i].x][edge[i].y] = 1;
            for(int k = 0; k != size; ++k)
                for(int i = 0; i != size; ++i)
                    for(int j = 0; j != size; ++j)
                        v[i][j] |= v[i][k] && v[k][j];
            for(int i = 0; i != size; ++i)
                for(int j = 0; j != i; ++j)
                    if(v[i][j] && v[j][i]) goto failed;
            ++cur;
            failed: ;
        }
        ans = (long long) ans * cur % mod;
    }
    while(m--) ans = (ans << 1) % mod;
    cout << ans;
    cerr << "check = " << (ans ^ 597855228) % 181 << endl;
}

Case 6 & 7 & 8

开6,发现是100个点的完全图。据说是个挺难推的式子。先放着去看7、8。

7和8中,点数很少。想起陈老师的《主旋律》那题。似乎可以容斥?找了些资料,发现可以对于每个子集S枚举源点集合s,则我们就有这些规则:

  • 任何以s里面的点为右端点的边都不能选
  • s指向S-s中的点可以任意选择连不连
  • S-s必须也是DAG。

于是就可以状压DP。记f(S)表示S集合为DAG的方案数,枚举源点集合s,则贡献为$2^{s连向S-s的边数}f(S-s)$。注意到这里会算重,但是符合容斥的性质,因此乘上$(-1)^{|S|+1}$即可。但是直接这么做是$n3^n$的(求s向S-s连边数),8可能要跑很久。我们可以利用无重边的性质,利用位运算加速,优化到$3^n$,这样8就只要跑80s左右了。

再回去做6,现在6就很简单了。由于结构相同,相同大小的子集的f都相同。于是就转移时枚举|s|,乘组合数即可。似乎许多本来要暴力的容斥,有了特殊的性质后,就可以高效解决了。

Case 7 & 8 程序

#include <bits/stdc++.h>
using namespace std;
#define IL __inline__ __attribute__((always_inline))
#define RG register
typedef long long ll;
const int mod = 1000000007;
IL int F() 
{
    RG char ch;
    while(ch = getchar(), ch < '0' || ch > '9')    ;
    RG int ans = ch - '0';
    while(ch = getchar(), '0' <= ch && ch <= '9') ans = ans * 10 + ch - '0';
    return ans;
}
IL int inc(RG int x, RG int v) {x += v; return x >= mod ? x - mod : x;}
IL int dec(RG int x, RG int v) {x -= v; return x < 0 ? x + mod : x;}

const int MaxN = 22;
int f[1 << MaxN], g[1 << MaxN];
int e_i[MaxN], e_o[MaxN];
int Pow[MaxN * MaxN], Pop[1 << MaxN], Fir[1 << MaxN];
#define Bit(i) (1 << (i))

int main()
{
    assert(freopen("uoip8.in", "r", stdin));
    assert(freopen("uoip8.out", "w", stdout));

    RG int n = F(), m = F();
    Pow[0] = 1;
    for(RG int i = 1; i <= m; ++i)
    {
        RG int x = F() - 1, y = F() - 1;
        e_o[x] |= 1 << y;
        e_i[y] |= 1 << x;
        Pow[i] = inc(Pow[i - 1], Pow[i - 1]);
    }
    f[0] = 1;
    for(RG int set = 1, sub, tmp, rev, pos, val; set != Bit(n); ++set)
    {
        cerr << set << '\n';
        Pop[set] = __builtin_popcount(set);
        Fir[set] = __builtin_ctz(set);
        tmp = (Pop[set] & 1) ? 1 : mod - 1;
        for(sub = set & (set - 1); sub; sub = set & (sub - 1))
        {
            pos = Fir[rev = set - sub];
            val = g[rev] = g[rev - Bit(pos)] + Pop[e_o[pos] & sub] - Pop[e_i[pos] & rev];
            tmp = (Pop[sub] & 1)
                ? inc(tmp, (ll) Pow[val] * f[rev] % mod)
                : dec(tmp, (ll) Pow[val] * f[rev] % mod);
        }
        f[set] = tmp;
    }
    RG int ans = f[Bit(n) - 1];
    cout << ans;
    cerr << "check = " << (ans ^ 597855228) % 181 << endl;    
}

Case 6 程序

#include <bits/stdc++.h>
using namespace std;
#define IL __inline__ __attribute__((always_inline))
#define RG register
typedef long long ll;
const int mod = 1000000007;
IL int F() 
{
    RG char ch;
    while(ch = getchar(), ch < '0' || ch > '9')    ;
    RG int ans = ch - '0';
    while(ch = getchar(), '0' <= ch && ch <= '9') ans = ans * 10 + ch - '0';
    return ans;
}
const int MaxN = 3010;
int comb[MaxN][MaxN], p[MaxN * MaxN], f[MaxN];

int main()
{
    assert(freopen("uoip6.in", "r", stdin));
    assert(freopen("uoip6.out", "w", stdout));

    RG int n = F();
    comb[0][0] = p[0] = 1;
    for(RG int i = 1, u = n * n; i <= u; ++i)
        p[i] = (p[i - 1] << 1) % mod;
    f[0] = 1;
    for(RG int i = 1, j, tmp; i <= n; ++i)
    {
        tmp = 0;
        for(j = comb[i][0] = 1; j <= i; ++j)
        {
            comb[i][j] = (comb[i - 1][j] + comb[i - 1][j - 1]) % mod;
            tmp = (j & 1)
                ? int((tmp + (ll) comb[i][j] * p[j * (i - j)] % mod * f[i - j]) % mod)
                : int((tmp - (ll) comb[i][j] * p[j * (i - j)] % mod * f[i - j]) % mod);
        }
        f[i] = tmp;
    }
    RG int ans = (f[n] + mod) % mod;    
    cout << ans;
    cerr << "check = " << (ans ^ 597855228) % 181 << endl;    
}

Case 9

一眼还是没啥性质…

发现除了前面两三个点,接下来这么多行中,可以两行两行分成一组,每一组左端点或右端点中都有一个相同的。也就是相当于一个一个点加入图中,每个点朝已经加入的任意两个点连任意方向的边(似乎是一种广义的树?)。

然而发现了这个似乎并没有什么帮助。先拿5的程序来跑,最大的双连通分量非常大,这个方法不可行。然后想起某模拟赛的一题“图染色计数”,发现外向和内向DAG都可以直接删掉,答案是一个乘在外面的常数。试着写了一些代码,删掉了四五十个点,然而剩下的点还是很多。接下来就是那些度为1的点了,它们构成了一些链。注意到第8个点中那个DP之所以可行,是因为能计算出一条边连或不连的方案数。把链缩成边,一条由$k$条边组成的新边,连的方案数是1,不连的方案数是$2^k-1$。这样还有些重边,那么也可以用类似的方法求出重边合并后连和不连的方案数。

写代码中出现了各种逗比错误,我真是太弱。最后缩完链后,点只有23个了。然而我不太会继续往下缩,也不是特别会处理有具体方案数的边,就直接写了一个$m3^{23}$的算法……半夜去跑……跑了14h终于跑了出来……

要是在考场上肯定跑不出来了……据 WrongAnswer 说他手动缩掉了两个点,然后又强行优化成$3^{21}$,真是令人膜拜不已。

(千万不要尝试用这个程序去跑……)

#include <bits/stdc++.h>
using namespace std;
#define IL __inline__ __attribute__((always_inline))
#define RG register
typedef long long ll;
const int mod = 1000000007;
IL int F() 
{
    RG char ch;
    while(ch = getchar(), ch < '0' || ch > '9')    ;
    RG int ans = ch - '0';
    while(ch = getchar(), '0' <= ch && ch <= '9') ans = ans * 10 + ch - '0';
    return ans;
}

#define pb push_back
#define mkp make_pair
#define fir first
#define sec second

namespace pb_ds
{
    IL int inc(RG int x, RG int v) {x += v; return x >= mod ? x - mod : x;}
    IL int dec(RG int x, RG int v) {x -= v; return x < 0 ? x + mod : x;}

    const int MaxN = 10010;
    vector<int> e_i[MaxN], e_o[MaxN];
    int Pow[MaxN], d_i[MaxN], d_o[MaxN];
    int pid[MaxN];
    bool del[MaxN];

    vector< pair<int, pair<int, int> > > e[MaxN];
    #define ForSTL(it, mp) for(RG typeof(mp.begin()) it = mp.begin(), ___e = mp.end(); it != ___e; ++it)

    const int MaxT = 25;
    int f[1 << MaxT], Pop[1 << MaxT], Fir[1 << MaxT];
}

int main()
{
    using namespace pb_ds;
    assert(freopen("uoip9.in", "r", stdin));
    assert(freopen("uoip9.out", "w", stdout));

    RG int n = F(), m = F();
    Pow[0] = 1;
    for(RG int i = 1; i <= m; ++i)
    {
        RG int x = F(), y = F();
        if(x == y) {--i, --m; continue;}
        assert(!e[x][y]);
        e_o[x].pb(y); ++d_o[x];
        e_i[y].pb(x); ++d_i[y];
        Pow[i] = inc(Pow[i - 1], Pow[i - 1]);
    }
    RG int ans = 1;

    static int q[MaxN];
    RG int l = 1, r = 0;
    for(RG int i = 1; i <= n; ++i) if(!d_i[i]) q[++r] = i;
    while(l <= r) 
    {
        RG int i = q[l++]; del[i] = 1;
        ans = (ll) ans * Pow[d_o[i]] % mod;
        for(RG unsigned j = 0; j != e_o[i].size(); ++j)
            if(!--d_i[e_o[i][j]]) q[++r] = e_o[i][j];
    }
    for(RG int i = 1; i <= n; ++i) if(d_i[i] && !d_o[i]) q[++r] = i;
    while(l <= r) 
    {
        RG int i = q[l++]; del[i] = 1;
        ans = (ll) ans * Pow[d_i[i]] % mod;
        for(RG unsigned j = 0; j != e_i[i].size(); ++j)
            if(!--d_o[e_i[i][j]]) q[++r] = e_i[i][j];
    }
//for(RG int i = 1; i <= n; ++i) if(d_i[i] * d_o[i]) printf("%d\n", i);
    //step 1
    //if(0)
    {
        RG int N = n; n = 0;
        memset(pid + 1, -1, N << 2);
        for(RG int i = 1; i <= N; ++i) if(d_i[i] * d_o[i] > 1) pid[i] = n++;
        //for(RG int i = 1; i <= N; ++i) if(d_i[i] * d_o[i]) pid[i] = n++;
        map< pair<int, int>, pair<int, int> > edge;
        for(RG int i = 1; i <= N; ++i) if(pid[i] >= 0)
        {
            for(RG unsigned j = 0; j != e_o[i].size(); ++j)
            {
                if(del[e_o[i][j]]) continue;
                RG int p = e_o[i][j], v = 1;
                while(pid[p] == -1) 
                {
                    assert(!del[p]);
                    del[p] = 1;
                    ++v;
                    for(RG unsigned k = 0; ; ++k)
                        if(!del[e_o[p][k]])
                        {
                            p = e_o[p][k];
                            break;
                        }
                }
                edge[mkp(pid[i], pid[p])].fir += (v);
                if(!edge[mkp(pid[i], pid[p])].sec) edge[mkp(pid[i], pid[p])].sec = 1;
                edge[mkp(pid[i], pid[p])].sec = (ll) edge[mkp(pid[i], pid[p])].sec * dec(Pow[v], 1) % mod;
            }
            //del[i] = 1;
        }
        for(RG int i = 1; i <= N; ++i) if(pid[i] == -1 && !del[i]) 
        {
            //there is a ring
            RG int p = i, v = 1;
            for( ; ; )
            {
                del[p] = 1;
                ForSTL(k, e_o[p]) if(pid[*k] == -1 && !del[*k])
                    {p = *k, ++v; goto next_step;}
                break;
                next_step: ;
            }
            ans = (ll) ans * dec(Pow[v], 1) % mod;
        }
        ForSTL(k, edge)
        {
            k->sec.fir = Pow[k->sec.fir];
            e[k->fir.fir].pb(mkp(k->fir.sec, k->sec));
        }
    }
    f[0] = 1;
    #define Bit(i) (1 << (i))
    for(RG int set = 1, tmp, sub, rev, pos, sum, gay, j, s; set != Bit(n); ++set)
    {
        cerr << set << '\n';
        Pop[set] = __builtin_popcount(set);
        Fir[set] = __builtin_ctz(set);
        tmp = 0;
        for(sub = set; sub; sub = set & (sub - 1))
        {
            rev = set - sub;
            sum = 1;
            for(gay = rev; gay; gay -= 1 << pos)
            {
                pos = Fir[gay];
                for(j = 0, s = e[pos].size(); j != s; ++j)
                    if(Bit(e[pos][j].fir) & sub) sum = (ll) sum * e[pos][j].sec.sec % mod;
            }
            for(gay = sub; gay; gay -= 1 << pos)
            {
                pos = Fir[gay];
                for(j = 0, s = e[pos].size(); j != s; ++j)
                    if(Bit(e[pos][j].fir) & rev) sum = (ll) sum * e[pos][j].sec.fir % mod;
                    else if(Bit(e[pos][j].fir) & sub) sum = (ll) sum * e[pos][j].sec.sec % mod;
            }
            tmp = (Pop[sub] & 1)
                ? inc(tmp, (ll) sum * f[rev] % mod)
                : dec(tmp, (ll) sum * f[rev] % mod);
        }
        f[set] = tmp;
    }
    ans = (ll) ans * f[Bit(n) - 1] % mod;
    cout << ans;
    cerr << "check = " << (ans ^ 597855228) % 181 << endl;    
}

Case 10

这幅图的性质和9类似,但是是一个点向之前选3个点连边,而且点数多了10倍。这下就不容易缩边了……

听说貌似树分解能做这个测试点?然而我并没有去WC2015,并不是很熟悉这一套理论。于是b*i*u一发,然而只找到一篇收费的文章,而且讲的内容寥寥无几……G*o*l*了一发也没啥收获。最后我灵机一动,在b*i*u搜Tree Decomposition,就搜到了一篇不错的论文:

http://dept-info.labri.fr/~dicky/PUF/M1Projects/TreeDecomposition-2011-2012.pdf

然后就研究了一段时间。

一个无向图$G = (V, E)$的树分解$(T, X)$表示一棵树$T$,其每个节点都有一个点集$X_i$满足以下条件

  1. $\bigcup_{i} X_i = V$,即任意一个点至少属于一个集合
  2. $\forall (x, y) \in E, \exists i \in T, s.t. x \in X_i, y \in Y_i$,即任意一条边的两个端点,都存在一个树上点集同时包含它们
  3. 对于任意$x, y\in T$,若存在$p \in X_x, p\in X_y$,则对于$x$到$y$唯一路径上的所有点$z$,都有$p \in X_z$。即包含一个特定图上点的树上点构成一个连通块。

对有向图的基图执行树分解。如果每个$X_i$都不大,那么我们就可以用4的方法进行传递闭包状压DP。由于性质3成立,我们也可以像测试点4那样,只需要记录需要的点,因为如果一个点不再需要被记录,它在后面也不会影响转移。似乎可以设计一种树形DP来求方案数。然而我们要先构造一个树分解。

论文中提供了一种方法。每次选择一个点度最小的点p,将p相邻的点连成团,然后将二元组(p, 与p相连的点集合S)压入栈中,最后从原图上删除p及其所有边,一直重复到只剩下两个点。建树时,令剩下的两个点为根,每次取栈顶,寻找编号最小的且包含S的点作为其父亲,并设当前点的$X$集合为{p}$\cup$与p相连的点集合S。这样构造出来的树分解,每个儿子的$X$集合只和其父亲差1,给转移带来了莫大的方便。

试着实现了这个算法,发现所有点的$X$集合大小均不超过4。真棒!那就可以DP了。然而我们该如何DP呢?刚开始我想用map<ull, int>作为DP数组进行记忆化搜索DP,然而后来却发现这样的转移具有后效性——子树中加的点可能会影响连通性,因此这里的转移不能直接用乘法原理!

感觉怎么做复杂度都会炸的样子。我突发奇想,干脆不记方案数,而是用map<ull, map<ull, int> >记录一个状态队列,表示进入子树时状态是$s$,离开子树时会有多少种状态、它们的方案数分别是多少。这样就可以解决这个问题,在儿子与儿子之间轻松转移了。

代码要注意从父亲切换到儿子和从儿子回到父亲时,状压中点编号的变动,十分繁琐……拍了好久终于拍过了……虽然看起来复杂度很糟糕,然而也只跑了2s就出解了。

提示:下面这个代码可以直接用来跑9。

#include <bits/stdc++.h>
using namespace std;
#define IL __inline__ __attribute__((always_inline))
#define RG register
typedef long long ll;
typedef unsigned long long ull;
const int mod = 1000000007;
IL int F() 
{
    RG char ch;
    while(ch = getchar(), ch < '0' || ch > '9')    ;
    RG int ans = ch - '0';
    while(ch = getchar(), '0' <= ch && ch <= '9') ans = ans * 10 + ch - '0';
    return ans;
}

#define pb push_back
#define mkp make_pair
#define fir first
#define sec second

namespace pb_ds
{
    IL int inc(RG int x, RG int v) {x += v; return x >= mod ? x - mod : x;}
    IL int dec(RG int x, RG int v) {x -= v; return x < 0 ? x + mod : x;}
    IL void add(RG int& x, RG int v) {x = inc(x, v);}

    const int MaxN = 10010;
    bool e[MaxN][MaxN];
    bool E[MaxN][MaxN];
    int d[MaxN];
    bool v[MaxN];

    struct Piece {int p; vector<int> s;} S[MaxN];

    const ull empty = 9241421688590303745ull;
    IL bool test(RG const ull& s, RG int x, RG int y) {return s >> (x << 3 | y) & 1;}
    IL void set(RG ull& s, RG int x, RG int y) {s |= 1ull << (x << 3 | y);}
    IL void flip(RG ull& s, RG int x, RG int y) {s ^= 1ull << (x << 3 | y);}
    IL int out_set(RG const ull& s, RG int i) {return (s >> (i << 3)) & 255;}
    IL void print(RG const ull& s, RG int n = 8)
    {
        for(RG int i = 0; i != n; ++i)
            for(RG int j = 0; j != n; ++j)
                printf("%d%c", test(s, i, j), " \n"[j + 1 == n]);
    }
    IL ull add_edge(RG ull s, RG int x, RG int y)
    {
        if(test(s, y, x)) return 0;
        for(RG int i = 0; i != 8; ++i) if(test(s, i, x))
            s |= (ull) out_set(s, y) << (i << 3);
        return s;
    }

    struct Edge{int x, y;};
    template<class T> void reinit(RG T& t) {t.~T(); new (&t)T; }
    #define formap(III, EEE, M) for(III = M.begin(), EEE = M.end(); III != EEE; ++III)

    IL int value(RG map<ull, int>& M)
    {
        RG int ans = 0;
        RG map<ull, int>::iterator iter, end;
        formap(iter, end, M) add(ans, iter->sec);
        return ans;
    }

    struct Node *_T;
    struct Node
    {
        int p;
        vector<int> vec, son;
        map<ull, map<ull, int> > mem;
        IL bool include(RG const vector<int> &_V) const
        {
            RG int s = vec.size(), j = 0, S = _V.size();
            RG const int *v = &*vec.begin(), *V = &*_V.begin();
            for(RG int i = 0; i != S; ++i)
            {
                while(j != s && v[j] < V[i]) ++j;
                if(j == s || v[j] != V[i]) return 0;
            }
            return 1;
        }
        map<ull, int>& dp(RG const ull& s)
        {
            if(mem.find(s) != mem.end())
                return mem.find(s)->sec;

            RG map<ull, int> Q[2]; RG bool d = 0;
            RG map<ull, int>::iterator iter, end;
            Q[d][s] = 1;
            for(RG unsigned k = 0; k != son.size(); ++k)
            {
                RG Node &T = _T[son[k]];
                vec.pb(T.p);
                RG int n = vec.size(), np = n - 1, m = 0;

                RG Edge e[10];
                for(RG int i = 0; i != np; ++i)
                {
                    if(pb_ds::e[vec[i]][T.p]) e[m++] = (Edge) {i, np};
                    if(pb_ds::e[T.p][vec[i]]) e[m++] = (Edge) {np, i};
                }
                for(RG int i = 0; i != m; ++i)
                {
                    reinit(Q[d ^= 1]);
                    formap(iter, end, Q[!d])
                    {
                        add(Q[d][iter->fir], iter->sec);
                        RG ull next = add_edge(iter->fir, e[i].x, e[i].y);
                        if(next) add(Q[d][next], iter->sec);
                    }
                }

                RG int id[10];
                for(RG int i = 0; i != n; ++i)
                {
                    RG const int *k = find(&*T.vec.begin(), &*T.vec.end(), vec[i]);
                    if(k == &*T.vec.end()) id[i] = -1;
                    else id[i] = k - &*T.vec.begin();
                }

                reinit(Q[d ^= 1]);
                formap(iter, end, Q[!d])
                {
                    RG ull set = iter->fir, conv = empty;
                    for(RG int i = 0; i != n; ++i) if(~id[i])
                        for(RG int j = 0; j != n; ++j) if(~id[j])
                            if(test(set, i, j))
                                pb_ds::set(conv, id[i], id[j]);
                    RG map<ull, int>& ret = T.dp(conv);
                    RG map<ull, int>::iterator Iter, End;
                    formap(Iter, End, ret)
                    {
                        RG ull set = Iter->fir, conv = iter->fir;
                        for(RG int i = 0; i != n; ++i) if(~id[i])
                            for(RG int j = 0; j != n; ++j) if(i != j && ~id[j])
                                if(test(set, id[i], id[j]))
                                    conv = add_edge(conv, i, j);
                        for(RG int i = 0; i != np; ++i)
                        {
                            if(test(conv, i, np)) flip(conv, i, np);
                            if(test(conv, np, i)) flip(conv, np, i);
                        }
                        add(Q[d][conv], (ll) iter->sec * Iter->sec % mod);
                    }
                }
                vec.pop_back();
            }
            mem.insert(mkp(s, Q[d]));
            return mem.find(s)->sec;
        }
    } T[MaxN];
}

int main()
{
    using namespace pb_ds;
    assert(freopen("uoip10.in", "r", stdin));
    assert(freopen("uoip10.out", "w", stdout));
    RG int n = F(), m = F();
    for(RG int i = 1; i <= m; ++i)
    {
        RG int x = F(), y = F();
        if(x == y) {--i, --m; continue;}
        ++d[x]; 
        ++d[y];
        E[x][y] = E[y][x] = 1;
        e[x][y] = 1;
        //e[y][x] = -1;
    }
    RG int top = 0;
    for(RG int lim = 1, k = n; k > 1; ++lim)
    {
        for(RG int p; k > 1; )
        {
            p = 0;
            for(RG int i = 1; i <= n; ++i)
                if(!v[i] && d[i] < lim) {p = i; break;}
            if(!p) break;
            --k;
            v[p] = 1;
            S[++top].p = p;
            for(RG int i = 1; i <= n; ++i)
                if(!v[i] && E[p][i]) 
                {
                    S[top].s.pb(i);
                    --d[i];
                }
            RG int D = d[p], *s = &*S[top].s.begin();
            for(RG int i = 0; i != D; ++i)
                for(RG int j = 0; j != i; ++j)
                    if(!E[s[i]][s[j]])
                    {
                        ++d[s[i]], ++d[s[j]];
                        E[s[i]][s[j]] = E[s[j]][s[i]] = 1;
                    }
        }
    }
    RG int N = 0;
    while(top)
    {
        RG Piece& P = S[top--];
        sort(&*P.s.begin(), &*P.s.end());
        RG Node cur = (Node) {P.p, P.s};
        cur.vec.pb(P.p);
        inplace_merge(&*cur.vec.begin(), &*cur.vec.end() - 1, &*cur.vec.end());
        if(!N) T[++N] = cur;
        else for(RG int i = 1; i <= N; ++i)
            if(T[i].include(P.s)) 
            {
                T[++N] = cur, T[i].son.pb(N);
                break;
            }
    }

    static Edge edge[10]; RG int tot = 0;
    for(RG unsigned i = 0; i != T[1].vec.size(); ++i)
        for(RG unsigned j = 0; j != T[1].vec.size(); ++j)
            if(e[T[1].vec[i]][T[1].vec[j]]) 
                edge[tot++] = (Edge) {(int) i, (int) j};
    _T = T;
    RG map<ull, int> Q[2];
    RG map<ull, int>::iterator iter, end;
    RG bool d = 0;
    Q[d][empty] = 1;
    for(RG int i = 0; i != tot; ++i)
    {
        reinit(Q[d ^= 1]);
        formap(iter, end, Q[!d])
        {
            add(Q[d][iter->fir], iter->sec);
            RG ull next = add_edge(iter->fir, edge[i].x, edge[i].y);
            if(next) add(Q[d][next], iter->sec);
        }
    }
    RG int ans = 0;
    formap(iter, end, Q[d]) 
        ans = (ans + (ll) iter->sec * value(T[1].dp(iter->fir))) % mod;

    cout << ans;// << endl;
    cerr << "check = " << (ans ^ 597855228) % 181 << endl;
}

/*
5 7
2 1
1 3
2 3
4 1
3 4
5 2
4 5

^/\*lgg\*/

总结

这真是一道神题啊!

orz matthew99

orz WrongAnswer

我真是太弱了!

完结撒花!

【广告向】求更优算法:带Link、Cut树上路径k小值

2015-08-04 19:13:08 By immortalCO

本萌萌哒弱菜出了一题(其实思路是集训队论文给的),但~算法不是很优美,求各位神犇去虐一虐~~~

提交地址

http://www.lydsy.com/JudgeOnline/problem.php?id=4234

测试数据

测试数据: http://pan.baidu.com/s/1qWDV9la 密码: c7nj

命题报告:http://pan.baidu.com/s/1o6qtgJC (写的很傻逼,神犇轻D)

如果有了更优的算法,一定记得告诉本弱菜一声哦~~~

代码什么的

写的很长,因为前面有2k的输入输出优化……

lgg语句是debug语句

具体实现和命题报告里的有一些差别

//pb_ds NOI template 20150710
#include <cstdio>
#include <algorithm>
#include <iostream>
#include <cstring>
#include <cassert>
#include <cmath>
#include <vector>
#include <set>
#include <utility>
#include <queue>
#define For(i, a, b) for(register int i = a, ___u = b; i <= ___u; ++i)
#define ForDown(i, a, b) for(register int i = b, ___d = a; i >= ___d; --i)
#define cmax(i, j) ((i) < (j) ? (i) = (j) : (i))
#define cmin(i, j) ((i) > (j) ? (i) = (j) : (i))
#define dmax(i, j) ((i) > (j) ? (i) : (j))
#define dmin(i, j) ((i) < (j) ? (i) : (j))
#define ddel(i, j) ((i) > (j) ? (i) - (j) : (j) - (i))

#include <queue>
#include <ext/pb_ds/assoc_container.hpp>
#include <ext/pb_ds/priority_queue.hpp>
#include <ext/pb_ds/tree_policy.hpp>

namespace io
{
    const int MAXBUF = 1 << 18;
    char B[MAXBUF], *S = B, *T = B;
    #define getc() (S == T && (T = (S = B) + fread(B, 1, MAXBUF, stdin), S == T) ? 0 : *S++)
    inline int F()
    {
        register int ch, aa = 0, bb = 0;
        register char *S = io::S, *T = io::T;
        for(ch = getc(); (ch < '0' || ch > '9') && ch != '-'; ch = getc())
            ;
        for(ch == '-' ? bb = 1 : aa = ch - '0', ch = getc(); '0' <= ch && ch <= '9'; ch = getc())
            aa += (aa << 3) + aa + ch - '0';
        io::S = S, io::T = T;
        return bb ? -aa : aa; 
    }

    inline double Ff()
    {
        register double aa, bb;
        register char ch;
        register char *S = io::S, *T = io::T;
        while(ch=getc(),(ch<'0'||ch>'9'))
            ;aa=ch-'0';
        while(ch=getc(),(ch>='0'&&ch<='9'))aa=aa*10+ch-'0';
        if(ch=='.'){bb=1;while(ch=getc(),ch>='0'&&ch<='9')bb*=0.1,aa+=bb*(ch-'0');}
        io::S = S, io::T = T;
        return aa;
    }

    char buff[MAXBUF], *iter = buff;
    template<class T>inline void P(register T x, register char ch = '\n')
    {
        static int stack[110];
        register int O = 0;
        register char *iter = io::iter;
        if(!x)*iter++ = '0';
        else
        {
            (x < 0) ? x = -x, *iter++ = '-' : 1;
            for(; x; x /= 10)
                stack[++O] = x % 10;
            for(; O; *iter++ = '0' + stack[O--])
                ;
        }
        *iter++ = ch, io::iter = iter;
    }

    inline int gets(register char *o)
    {
        register char *s = o, ch = getc();
        for(; ch && (ch == ' ' || ch == '\n'); ch = getc())
            ;
        for(; ch && (ch != ' ' && ch != '\n'); ch = getc())
            *s++ = ch;
        *s = 0;
        return s - o;
    }

    inline void output() {fwrite(buff, 1, iter - buff, stdout);}
}

#include <ext/pb_ds/priority_queue.hpp>

namespace pb_ds
{
    using io::F;
    using io::P;
    #define RG register
    int n;
    const int MAXN = 200010, MAXM = 400010;
    namespace lct
    {
        const int MAXN = pb_ds::MAXN << 1;
        int tot;
        int fa[MAXN], son[MAXN][2], time[MAXN], cho[MAXN];
        bool rev[MAXN], mark[MAXN], have[MAXN];
        #define type(i) (son[fa[i]][1] == (i))
        #define not_root(i) (son[fa[i]][type(i)] == (i))
        inline void push_up(register int i)
        {
            static int A, B;
            have[i] = have[son[i][0]] || mark[i] || have[son[i][1]];
            time[i] < time[cho[i] = ((time[A = cho[son[i][0]]] < time[B = cho[son[i][1]]]) ? A : B)]
                ? cho[i] = i : 1;
///*lgg*/assert(time[cho[i]] <= time[i] && time[cho[i]] <= time[cho[son[i][0]]] && time[cho[i]] <= time[cho[son[i][1]]]);
        }
        inline void push_down(register int i)
        {
            if(rev[i])
            {
                rev[i] = 0, rev[son[i][0]] ^= 1, rev[son[i][1]] ^= 1;
                register int t = son[i][0]; son[i][0] = son[i][1], son[i][1] = t;
            }
        }
        inline void rotate(register int i)
        {
            register int f = fa[i], t = type(i);
            fa[i] = fa[f], not_root(f) ? son[fa[f]][type(f)] = i : 1;
            (son[f][t] = son[i][!t]) ? fa[son[f][t]] = f : 1;
            push_up(son[fa[f] = i][!t] = f);
        }
        inline void preview(register int i)
        {
            if(not_root(i)) preview(fa[i]);
            push_down(i);
        }
        inline void splay(register int i, register bool need = 1)
        {
            if(need) preview(i);
            for( ; not_root(i); rotate(i))
                if(not_root(fa[i])) rotate(type(i) ^ type(fa[i]) ? i : fa[i]);
            push_up(i);
        }
        inline void dfs_print(register int i)
        {
            if(!i)return;
            dfs_print(son[i][0]);
            printf("%d ", i);
            dfs_print(son[i][1]);
        }
        inline void debug(register int root)
        {
            splay(root);
            printf("root = %d chain = ", root);
            dfs_print(root);
            puts("");
        }
        inline int access(register int i)
        {
            register int j = 0;
            for( ; i; splay(i), son[i][1] = j, j = i, i = fa[i])
                ;
            return j;
        }
        inline int query(RG int x, RG int y) {return access(x), access(y);}
        inline void make_root(register int i) {access(i), splay(i, 0), rev[i] ^= 1;}
        inline void link(register int i, register int j) 
        {
///*lgg*/printf("link %d ~ %d\n", i, j);
            make_root(i), fa[i] = j;
///*lgg*///printf("linked "); debug(j);
        }
        inline void cut(register int i, register int j) 
        {
            make_root(i), access(j), splay(i, 0), son[i][1] = fa[j] = 0;
///*lgg*/ assert(son[i][1] == j);
///*lgg*/// puts("cut"); printf("i: "); debug(i); printf("j: "); debug(j);
        }
        inline int find_root(register int i)
        {
            for(access(i), splay(i, 0); son[i][0]; i = son[i][0])
                ;
            return i;
        }
        inline int expose(register int i, register int j) 
        {
            if(i == j) return j;
            make_root(i), access(j), splay(j, 0);
            return j;
        }
        inline void set_mark(RG int i, RG bool b)
        {
            splay(i), mark[i] = b, push_up(i);
        }    
    }
    namespace graph
    {
        int fir[MAXN], tot = 1, next[MAXM], to[MAXM];
        inline void reset() {memset(fir, 0, sizeof fir); tot = 1;}
        inline void link(RG int x, RG int y)
        {
            next[++tot] = fir[x], to[fir[x] = tot] = y;
            next[++tot] = fir[y], to[fir[y] = tot] = x;
        }    
    }
    namespace segment
    {
        const int MAXNODE = 10000000;
        const int LL = 1, RR = 1000000;
        int tot, lef[MAXNODE], rig[MAXNODE], size[MAXNODE];
        inline int insert(RG int val, RG int old, RG int l = LL, RG int r = RR)
        {
            RG int at = ++tot;
            size[at] = size[old] + 1;
///*lgg*/if(l==LL&&r==RR)printf("insert val = %d at = %d l = %d r = %d size = %d\n", val, at, l, r, size[at]);
            if(l < r)
            {
                RG int m = (l + r) >> 1;
                (val <= m)
                    ? (lef[at] = insert(val, lef[old], l, m), rig[at] = rig[old])
                    : (lef[at] = lef[old], rig[at] = insert(val, rig[old], m + 1, r));
///*lgg*/assert(size[at] == size[lef[at]] + size[rig[at]]);
            }
            return at;
        }
    }
    int val[MAXN];
    struct E {int x, y;} e[MAXN];
    int root[MAXN], top[MAXN], fa[MAXN], dep[MAXN] = {0, 1};
    bool chosen[MAXN], vis[MAXN];
    inline int dfs_init(RG int at)
    {
///*lgg*///printf("at = %d fa = %d\n", at, fa[at]);
        using namespace graph;
        root[at] = segment::insert(val[at], root[fa[at]]);
        vis[at] = 1;
        RG int size = 1, cho = 0, max_cho = 0, d = dep[at] + 1, temp;
        for(RG int i = fir[at]; i; i = next[i])
            !vis[to[i]]
                ?    fa[to[i]] = at, dep[to[i]] = d,
                     size += (temp = dfs_init(to[i])), 
                    temp > max_cho ? cho = to[i], max_cho = temp : 1 : 1;
        chosen[cho] = 1;
///*lgg*///printf("dfs_init at = %d val = %d dep = %d size = %d cho = %d\n", at, val[at], dep[at], size, cho);
        return size;
    }
    inline void dfs_make(RG int at)
    {
///*lgg*/printf("at = %d fa = %d\n", at, fa[at]);
        using namespace graph;
        vis[at] = 0,
        top[at] = (chosen[at] ? top[fa[at]] : at);
///*lgg*///printf("at = %d top = %d\n", at, top[at]);
        for(RG int i = fir[at]; i; i = next[i])
            if(vis[to[i]]) dfs_make(to[i]);
    }
    inline int lca(RG int x, RG int y)
    {
        while(top[x] != top[y])
            dep[top[x]] > dep[top[y]]
                ? x = fa[top[x]] : y = fa[top[y]];
        return dep[x] < dep[y] ? x : y;
    }
    struct Change {int at, val, cnt; } change[MAXN];
    struct Segment {int root, cnt; } seg[MAXN];
    int change_cnt, seg_cnt;
    int marked[MAXN], marked_cnt;
    namespace lct
    {
        inline void dfs(RG int at)
        {
            if(!have[at]) return;
            push_down(at);
            dfs(son[at][0]);
///*lgg*///printf("at = %d\n", at);
            if(mark[at])
            {
///*lgg*/printf("found mark %d\n", at);
                marked[++marked_cnt] = at;
            }
            dfs(son[at][1]);
        }
    }
    inline void build()
    {
        RG int x, y;
        graph::reset();
        segment::tot = change_cnt = 0;
        memset(chosen, 0, sizeof chosen);
        For(i, 1, n) 
            if(lct::mark[i]) lct::set_mark(i, 0);
        For(i, n + 1, n + n - 1)
        {
            graph::link(x = e[i].x, y = e[i].y);
            if(lct::mark[i])
                lct::set_mark(i, 0);
///*lgg*/assert(!lct::mark[i]);
///*lgg*/printf("link %d ~ %d time = %d\n", x, y, lct::time[i]+(1<<30));
        }
///*lgg*/For(i,0,n*2)assert(!lct::mark[i]);
        fa[1] = 0;
        dfs_init(1);
        dfs_make(1);
    }


    inline void main()
    {
        n = F();
        RG int q = F(), H = 2 * sqrt(n), timer = -1 << 30;
        For(i, 1, n) val[i] = F();
        RG int x, y;
        For(i, n + 1, n + n - 1)
        {
            e[i] = (E) {x = F(), y = F()};
            lct::time[i] = ++timer;
            lct::link(x, i);
            lct::link(i, y);
        }
        build();
        RG int k, modify_tot = 0, cho;
///*lgg*/printf("H = %d\n", H);
        For(aaaa, 1, q)
        {
            k = F(), x = F(), y = F();
///*lgg*/if(k==-2)build();else
            if(k == -1)
            {
                cho = lct::cho[lct::expose(x, y)];
                lct::cut(e[cho].x, cho);
                lct::cut(cho, e[cho].y);
///*lgg*/printf("replace %d (%d ~ %d) time = %d with (%d ~ %d) time = %d\n", cho, e[cho].x, e[cho].y,lct::time[cho]+(1<<30),x,y,timer+1+(1<<30));
                e[cho] = (E) {x, y};
                lct::time[cho] = ++timer;
                lct::mark[cho] = lct::mark[x] = lct::mark[y] = 1;
                lct::link(x, cho);
                lct::link(cho, y);
                if(!(++modify_tot % H))
                    build();
                lct::expose(x, y);
///*lgg*///lct::debug(lct::expose(x,y));
            }
            else if(!k)
            {
                change[++change_cnt] = (Change) {x, val[x], -1};
///*lgg*/printf("change at = %d cnt[%d] -1 cnt[%d] +1\n", x,val[x],y);
                change[++change_cnt] = (Change) {x, val[x] = y, 1};
                if(!(++modify_tot % H))
                    build();
            }
            else
            {
///*lgg*/printf("query %d %d\n", x, y);
                if(x == y)
                {
                    P(val[x]);
                    continue;
                }
                seg_cnt = marked_cnt = 0;
                if(!lct::mark[x]) marked[++marked_cnt] = x;
                lct::dfs(lct::expose(x, y));
                if(!lct::mark[y]) marked[++marked_cnt] = y;

                //现在是处理mark的点的时候
                marked[++marked_cnt] = MAXN;
                RG int a, b, c;
                For(i, 1, marked_cnt)
                {
                    a = marked[i];
                    while(marked[i] <= n) ++i;
                    b = marked[i - 1];
                    if(a == b)
                    {
                        seg[++seg_cnt] = (Segment) {root[a], 1};
                        if(fa[a])seg[++seg_cnt] = (Segment) {root[fa[a]], -1};    
///*lgg*/printf("add_point %d fa = %d\n", a, fa[a]);                    
                    }
                    else
                    {
                        c = lca(a, b);
                        seg[++seg_cnt] = (Segment) {root[a], 1};
                        seg[++seg_cnt] = (Segment) {root[b], 1};
                        seg[++seg_cnt] = (Segment) {root[c], -1};
                        if(fa[c])seg[++seg_cnt] = (Segment) {root[fa[c]], -1};
///*lgg*/printf("add_chain %d ~ %d lca = %d flca = %d\n", a, b, c, fa[c]);
                    }
                }

                static int temp[2][MAXN];
                RG int *f = temp[0], *g = temp[1], *t, l = segment::LL, r = segment::RR, left, m;
                f[0] = 0, lct::make_root(x);
                For(i, 1, change_cnt)
                    if(lct::query(y, change[i].at) == change[i].at)
                        f[++f[0]] = i;
///*lgg*/For(i, 1,seg_cnt)printf("root = %d cnt = %d size = %d\n", seg[i].root, seg[i].cnt, segment::size[seg[i].root]); 
                while(l < r)
                {
///*lgg*/printf("l = %d r = %d k = %d m = %d\n", l, r, k, (l+r)/2);
///*lgg*/For(i,1,f[0])printf("assoc change at = %d val = %d cnt = %d\n", change[f[i]].at, change[f[i]].val, change[f[i]].cnt);
                    left = 0;
                    m = (l + r) >> 1;
                    For(i, 1, seg_cnt)
                        left += segment::size[segment::lef[seg[i].root]] * seg[i].cnt;
///*lgg*/printf("segment left = %d\n", left);
                    For(i, 1, f[0])
                        if(change[f[i]].val <= m)
                            left += change[f[i]].cnt;
///*lgg*/printf("segment + change left = %d\n", left);    
                    if(k <= left)
                    {
                        r = m;
                        For(i, 1, seg_cnt)
                            seg[i].root = segment::lef[seg[i].root];
                        g[0] = 0;
                        For(i, 1, f[0])
                            if(change[f[i]].val <= m)
                                g[++g[0]] = f[i];
                    }                
                    else
                    {
                        l = m + 1;
                        k -= left;
                        For(i, 1, seg_cnt)
                            seg[i].root = segment::rig[seg[i].root];
                        g[0] = 0;
                        For(i, 1, f[0])
                            if(m < change[f[i]].val)
                                g[++g[0]] = f[i];
                    }
                    t = g, g = f, f = t;
///*lgg*/getchar();    
                }
///*lgg*/printf("ans = %d\n", l);
                P(l);

            }
//if(!(modify_tot % H)) fprintf(stderr, "rebuild i = %d\n", aaaa);
        }
        io::output();
    }
}

int main()
{
    pb_ds::main();
}
共 9 篇博客