标签：矩阵乘法   基本   坐标   usaco   org   his   move   解决   重载   

[ZJOI2007]捉迷藏

这是我最近写过最长的代码QAQ
码力太弱了QAQ
动态点分治模板题。
我们可以用三种堆来维护答案，这些堆要求支持删除非顶元素，以及查询次小值。我们把两个STL堆封装起来就可以实现。

三种堆：

d[x]表示以x为根的点分树中所有黑点到它分治爹的距离
c[x]表示以x为根的所有点分儿子d堆中的最大值
ans表示全局的最大值

我们从c中取出最大值和次大值就可以得到过这个点分根的最长链。我们不断用它来更新答案。

注意d的定义是到分治父亲的父亲的最大值
c数组要加入一个0

Code

#include<cstdio>
#include<queue>
#include<algorithm>
using namespace std;
namespace orz{
const int N=210000;
const int inf=2147483647;
struct heap{
    priority_queue<int>q,rm;
    heap(){
        q.push(-inf);
        rm.push(-inf+1);
    }
    inline void push(int x){
        q.push(x);
    }
    inline void remove(int x){
        if(q.top()==x)q.pop();
        else rm.push(x);
    }
    inline int top(){
        while(q.top()==rm.top())q.pop(),rm.pop();
        return q.top();
    }
    inline void pop(){
        while(q.top()==rm.top())q.pop(),rm.pop();
        q.pop();
    }
    inline int second(){
        int mx=top();
        if(mx==-inf)return -inf;
        q.pop();
        int ans=top();
        q.push(mx);
        return ans;
    }
    inline int size(){
        return q.size()-rm.size();
    }
    inline void show(){
        priority_queue<int>res;
        while(size()){
            res.push(top());
            printf("%d ",top());
            pop();
        }
        putchar('\n');
        while(res.size()){
            q.push(res.top());
            res.pop();
        }
    }
}ans,d[N],c[N];
//d表示所有的黑点到这个分治根的距离
//c是所有分治儿子的最长链
int head[N],ver[N],next[N],tot;
int size[N];
int father[N];
int dis[N][20];
int color[N],cnt;
int depth[N];
int n,m;
int root,all,rootmax;
bool rm[N];
inline void add(int x,int y){
    next[++tot]=head[x],head[x]=tot,ver[tot]=y;
    next[++tot]=head[y],head[y]=tot,ver[tot]=x;
}
inline int read(){
    int a=1,b=0;char t;
    do{t=getchar();if(t=='-')a=-1;}while(t>'9'||t<'0');
    do{b=b*10-'0'+t;t=getchar();}while(t>='0'&&t<='9');
    return a*b;
}
void write(int x){
    if(x)write(x/10),putchar('0'+x%10);
}
inline void print(int x){
    if(x<0){
        putchar('-');
        x=-x;
    }
    else if(!x){
        putchar('0');
    }
    write(x);
    putchar('\n');
}
inline bool getOrder(){
    char t;
    do{t=getchar();}while(t!='G'&&t!='C');
    return t=='G';
}
void getroot(int x,int fa){
    size[x]=1;
    int maxpart=-inf;
    for(int i=head[x];i;i=next[i]){
        if(rm[ver[i]]||ver[i]==fa)continue;
        getroot(ver[i],x);
        size[x]+=size[ver[i]];
        maxpart=max(maxpart,size[ver[i]]);
    }
    maxpart=max(maxpart,all-size[x]);
    if(maxpart<rootmax)rootmax=maxpart,root=x;
}
inline void addAns(heap &a){
    ans.push(a.top()+a.second());
}
inline void removeAns(heap &a){
    ans.remove(a.top()+a.second());
}
inline int ask(){
    if(cnt<=1)return cnt-1;
    return ans.top();
}
inline void BFS(int x,int dep){
    queue<int>q;
    q.push(x);
    int t;
    while(q.size()){
        t=q.front();
        q.pop();
        d[x].push(dis[t][depth[father[x]]]);
        for(int i=head[t];i;i=next[i]){
            if(rm[ver[i]])continue;
            if(dis[ver[i]][dep])continue;
            dis[ver[i]][dep]=dis[t][dep]+1;
            q.push(ver[i]);
        }
    }
}
int build(int x,int fa,int dep,int sum){
    //得到根
    rootmax=inf;
    all=sum;
    getroot(x,0);
    x=root;
    //保证size是对的
    getroot(x,0);
    //记录在点分树上的爸爸
    father[x]=fa;
    rm[x]=true;
    depth[x]=dep;
    //得到这一层所有点到分治根的距离
    BFS(x,dep);
    int son;
    //自己也算在所有链中
    c[x].push(0);
    for(int i=head[x];i;i=next[i]){
        if(rm[ver[i]])continue;
        son=build(ver[i],x,dep+1,size[ver[i]]);
        c[x].push(d[son].top());
    }
    addAns(c[x]);
    return x;
}
inline void turnOn(int x){
    //先更新对自己的影响
    removeAns(c[x]);
    c[x].remove(0);
    addAns(c[x]);
    //一层一层往上跳
    for(int i=x;father[i];i=father[i]){
        //先删去对答案为影响
        removeAns(c[father[i]]);
        c[father[i]].remove(d[i].top());
        d[i].remove(dis[x][depth[father[i]]]);
        c[father[i]].push(d[i].top());
        addAns(c[father[i]]);
    }
}
inline void turnOff(int x){
    removeAns(c[x]);
    c[x].push(0);
    addAns(c[x]);
    for(int i=x;father[i];i=father[i]){
        removeAns(c[father[i]]);
        c[father[i]].remove(d[i].top());
        d[i].push(dis[x][depth[father[i]]]);
        c[father[i]].push(d[i].top());
        addAns(c[father[i]]);
    }
}
inline void work(){
//  for(int i=1;i<=n;++i)
//      c[i].show();
//  for(int i=1;i<=3;++i){
//      for(int j=1;j<=n;++j)
//          printf("%d ",dis[j][i]);
//      putchar('\n');
//  }
    cnt=n;
    int x;
    m=read();
    for(int i=1;i<=m;++i){
        switch(getOrder()){
            case true:
                print(ask());
                break;
            case false:
                x=read();
                switch(color[x]){
                    case 0:
                        //开灯
                        turnOn(x);
                        color[x]=1;
                        --cnt;
                        break;
                    case 1:
                        //熄灯
                        turnOff(x);
                        color[x]=0;
                        ++cnt;
                        break;
                }
                break;
        }
    }
}
int QAQ(){
    n=read();
    for(int i=1;i<n;++i)
        add(read(),read());
    build(1,0,1,n);
    work();
    return false;
}
}
int main(){
    return orz::QAQ();
}

7.27爆零赛

这次考试题起名比较随意，题目难度不是按顺序排的，大概是倒序的。

随

给出\(n\)个正整数\(a_1,a_2…a_n\)和一个质数\(mod\).一个变量\(x\)初始为\(1\).进行\(m\)次操作.每次在\(n\)个数中随机选一个\(a_i\),然后\(x=xa_i\).问\(m\)次操作之后\(x\)的取值的期望。
答案输出a乘b的逆元的形式

NOIP模拟赛T1考O(n^2)矩阵乘法和原根，我佛了。
首先我们可以把问题转化为一个假期望，我们只需要求出每一种数值的方案数乘上数值大小除以\(n^m\)就可以了。
我们可以很容易的发现dp的转移是类似矩阵转移的形式。
所以就去打了一个\(O(mod^3\log m)\)的矩阵快速幂。然后因为常数太大被卡常卡死了。
题解中说这样可以得80分，然而一些人尝试去卡常只得了50分，而我卡了半天的才卡到了30分。
这个数据把mod设得太卡矩阵乘法导致了暴力反而能得更多的分，我觉得这很不合理QAQ。
我们很容易看出来这里需要一个\(O(n^2)\)的矩阵乘法，但我们现在只知道一个循环矩阵能做到。
所以正解肯定是一个循环矩阵。
但是我们发现因为原题中做的是乘法，所以转移是乱跳的，这不好。而出题人给了我们一些提示.

孙金宁教你学数学
质数P的原根g满足1<=rt<P,且rt的1次方,2次方…(P-1)次方在模P意义下可以取遍1到(P-1)的所有整数.
欧拉定理:对于质数P,1<=x<P的任意x的P-1次方在模P意义下都为1.
显然,原根的1次方,2次方…（P-2）次方在模P意义下都不为1,只有(P-1)次方在模P意义下为1.这也是一个数成为原根的充分必要条件.

因为原根的几次方可以取遍模p剩余系的所有数，反过来说就是每一个数都可以用原根的几次方的形式表示。那就很棒了，因为这样就可以将乘法变成加法，而加法是明显会形成一个循环矩阵的，所以这题就A了。

Code

#include<cstdio>
#include<algorithm>
using namespace std;
namespace orz{
const int N=1100;
const int MOD=1000000007;
bool vis[N];
int b[N];
int c[N];
int n,m,mod;
int root;
int ans[N],base[N];
int res[N];
inline int read(){
    int a=1,b=0;char t;
    do{t=getchar();if(t=='-')a=-1;}while(t>'9'||t<'0');
    do{b=b*10-'0'+t;t=getchar();}while(t>='0'&&t<='9');
    return a*b;
}
inline bool check(long long x){
    long long res=x;
    for(int i=1;i<mod-1;++i){
        if(res==1)return false;
        res=res*x%mod;
    }   
    return true;
}
inline int pow(long long x,int y){
    long long ans=1;
    while(y){
        if(y&1)ans=ans*x%MOD;
        x=x*x%MOD;
        y>>=1;
    }
    return ans;
}
inline void mul(int *a,int *b){
    for(int i=0;i<mod-1;++i)
        res[i]=0;
    for(int i=0;i<mod-1;++i)
        for(int j=0;j<mod-1;++j)
            res[(i+j)%(mod-1)]=(res[(i+j)%(mod-1)]+1ll*a[i]*b[j]%MOD)%MOD;
    for(int i=0;i<mod-1;++i)
        a[i]=res[i];
}
int QAQ(){
    n=read(),m=read(),mod=read();
    for(int i=1;i<mod;++i)
        if(check(i)){
            root=i;
            break;
        }
    for(int i=1,x=root;i<mod-1;++i,x=x*root%mod)
        b[x]=i,c[i]=x;
    c[0]=1,b[1]=0;
    for(int i=1;i<=n;++i)
        ++base[b[read()]];
    ans[0]=1;
    int y=m;
    while(y){
        if(y&1)mul(ans,base);
        mul(base,base);
        y>>=1;  
    }
    long long fans=0;
    for(int i=0;i<mod;++i)
        fans=(fans+(long long)ans[i]*c[i]%MOD)%MOD;
    printf("%lld\n",fans*pow(pow(n,m),MOD-2)%MOD);
    return false;
}
}
int main(){
    return orz::QAQ();
}

后面两道题懒得转格式了

单

单车联通大街小巷.这就是出题人没有写题目背景的原因.
对于一棵树,认为每条边长度为1,每个点有一个权值a[i].dis(u,v)为点u到v的最短路径的边数.dis(u,u)=0.对每个点求出一个重要程度.点x的重要程度b[x]定义为其他点到这个点的距离乘上对应的点权再求和. 即:b[x]=a[1]dis(1,x)+a[2]dis(2,x)+....+a[n]*dis(n,x)
现在有很多树和对应的a数组，并求出了b数组.不幸的是，记录变得模糊不清了.幸运的是,树的形态完好地保存了下来,a数组和b数组至少有一个是完好无损的,但另一个数组完全看不清了.
希望你求出受损的数组.多组数据.

从a求b：
换根DP，秒了。
从b求a：
我们取1号节点为根。
sum表示所有节点a值之和，size表示子树a值之和。
我们考虑从a求b的过程，因为我们原来是换根求的b数组，所以相邻两个节点的b值相减后为\(sum-2size_y\)
如果我们能知道sum的话所有值就能求出来了，但是要怎么求sum呢，考场上没想出来QAQ。
总是想着把它们加起来可以搞出sum，结果不行。
我们发现一号节点没有这个式子，所以我们把它的式子暴力写出来
\(b_1=\sum\limits_{i=1}^{n}a_idep_i\)
我们惊奇的发现这个式子等于
\(\sum\limits_{i=2}^nsize_i\)
减一下就出来了。
考场上就差最后两行了，Orz。

Code

#include<cstdio>
#include<algorithm>
#define tkj 0
using namespace std;
namespace orz{
const int N=300000;
long long a[N],b[N],size[N];
int head[N],next[N],ver[N],tot;
long long qwq[N];
long long sum=0;
inline int read(){
    int a=1,b=0;char t;
    do{t=getchar();if(t=='-')a=-1;}while(t>'9'||t<'0');
    do{b=b*10-'0'+t;t=getchar();}while(t>='0'&&t<='9');
    return a*b;
}
inline void add(int x,int y){
    next[++tot]=head[x],head[x]=tot,ver[tot]=y;
    next[++tot]=head[y],head[y]=tot,ver[tot]=x;
}
inline void dky(int x,int father){
    size[x]=a[x];b[x]=0;
    for(int i=head[x];i;i=next[i]){
        if(ver[i]==father)continue;
        dky(ver[i],x);
        size[x]+=size[ver[i]];
        b[x]+=b[ver[i]]+size[ver[i]];
    }
}
inline void ykd(int x,int father){
    for(int i=head[x];i;i=next[i]){
        if(ver[i]==father)continue;
        qwq[ver[i]]=b[ver[i]]-b[x];
        ykd(ver[i],x);
    }
}
inline void sfd(int x,int father){
    a[x]=size[x]=(sum-qwq[x])/2;
    for(int i=head[x];i;i=next[i]){
        if(ver[i]==father)continue;
        sfd(ver[i],x);
        a[x]-=size[ver[i]];
    }
}
inline void dfs(int x,int father){
    for(int i=head[x];i;i=next[i]){
        if(ver[i]==father)continue;
        b[ver[i]]=b[x]+sum-2*size[ver[i]];
        dfs(ver[i],x);
    }
}
int QAQ(){
//  freopen("qaq.in","r",stdin);
    int t;
    t=read();
    while(t--){
        int n=read();
        for(int i=1;i<=n;++i)
            head[i]=0;
        tot=0;
        for(int i=1;i<n;++i)
            add(read(),read());
        switch(read()){
            case 0:
                sum=0;
                for(int i=1;i<=n;++i)
                    sum+=(a[i]=read());
                dky(1,tkj);
                dfs(1,tkj);
                for(int i=1;i<=n;++i)
                    printf("%lld ",b[i]);
                break;
            case 1:
                sum=0;
                for(int i=1;i<=n;++i)
                    b[i]=read();
                ykd(1,tkj);
                for(int i=2;i<=n;++i)
                    sum+=qwq[i];
                sum=(sum+2*b[1])/(long long)(n-1);
                sfd(1,tkj);
                a[1]=sum;
                for(int i=2;i<=n;++i)
                    a[1]-=a[i];
                for(int i=1;i<=n;++i)
                    printf("%lld ",a[i]);
                break;
        }
        putchar('\n');
    }
    return false;
}
}
int main(){
    return orz::QAQ();
}

题

出个题就好了.这就是出题人没有写题目背景的原因.
你在平面直角坐标系上.
你一开始位于(0,0).
每次可以在上/下/左/右四个方向中选一个走一步.
即:从(x,y)走到(x,y+1),(x,y-1),(x-1,y),(x+1,y)四个位置中的其中一个.
允许你走的步数已经确定为n.现在你想走n步之后回到(0,0).但这太简单了.你希望知道有多少种不同的方案能够使你在n步之后回到(0,0).当且仅当两种方案至少有一步走的方向不同,这两种方案被认为是不同的.
答案可能很大所以只需要输出答案对10^9+7取模后的结果.(10^9+7=1000000007,1和7之间有8个0)
这还是太简单了,所以你给能够到达的格点加上了一些限制.一共有三种限制,加上没有限制的情况,一共有四种情况,用0,1,2,3标号:
0.没有任何限制,可以到达坐标系上所有的点,即能到达的点集为{(x,y)|x,y为整数}
1.只允许到达x轴非负半轴上的点.即能到达的点集为{(x,y)|x为非负数,y=0}
2.只允许到达坐标轴上的点.即能到达的点集为{(x,y)|x=0或y=0}
3.只允许到达x轴非负半轴上的点,y轴非负半轴上的点以及第1象限的点.即能到达的点集为{(x,y)|x>=0,y>=0}

type0,1,3 组合数，秒了
type2的数据较小，明显是一个DP，可是我不会。
有10分数据是<=100的，复杂度大一些也能过，暴力建图后跑矩阵快速幂。复杂度\(O((2n)^3log n)\),成功水到分。
正解是DP，我们把问题转化，变成了从（0，0）出发沿一个方向瞎走，回到（0，0），换（或不换）个方向再次瞎走。于是就变成了一个类似背包的东西，原行走序列被分成了一段一段的，考场上想到这里就没往下想，因为如果在半路走回零会导致统计重复，所以我们的一个单位应该是从零出去再回来中间不经过（0，0），这是（i-2）/2的卡特兰数，DP就完了。

Code

#include<cstdio>
#include<algorithm>
#include<cstring>
using namespace std;
namespace orz{
const int N=210000;
const int MOD=1000000007;
long long fac[N];
long long inv[N];
long long f[N];
int tot;
inline int read(){
    int a=1,b=0;char t;
    do{t=getchar();if(t=='-')a=-1;}while(t>'9'||t<'0');
    do{b=b*10-'0'+t;t=getchar();}while(t>='0'&&t<='9');
    return a*b;
}
inline int pow(long long x,int y){
    long long ans=1;
    while(y){
        if(y&1)ans=ans*x%MOD;
        x=x*x%MOD;
        y>>=1;
    }
    return ans;
}
inline long long C(int n,int m){
    return fac[n]*inv[m]%MOD*inv[n-m]%MOD;
}
inline long long Catalan(int n){
    return (C(2*n,n)-C(2*n,n-1)+MOD)%MOD;
}
int QAQ(){
    int n=read();
    fac[0]=1;
    for(int i=1;i<=n;++i)
        fac[i]=fac[i-1]*i%MOD;
    inv[n]=pow(fac[n],MOD-2);
    for(int i=n-1;i>=0;--i)
        inv[i]=inv[i+1]*(i+1)%MOD;
    long long ans;
    switch(read()){
        case 0:
            ans=0;
            for(int i=0;i<=n;i+=2){
                ans=(ans+C(n,i)*C(i,i/2)%MOD*C(n-i,(n-i)/2)%MOD)%MOD;
            }
            printf("%lld\n",ans);
            break;
        case 1:
            printf("%lld\n",Catalan(n/2));
            break;
        case 2:
            f[0]=1;
            for(int i=0;i<=n;i+=2)
                for(int j=0;j<i;j+=2)
                    f[i]=(f[i]+f[j]*4*Catalan((i-j)/2-1))%MOD;
            printf("%lld\n",f[n]);
            break;
        case 3:
            ans=0;
            for(int i=0;i<=n;i+=2)
                ans=(ans+C(n,i)*Catalan(i/2)%MOD*Catalan((n-i)/2)%MOD)%MOD;
            printf("%lld\n",ans);
            break;
    }
    return false;
}
}
int main(){
    return orz::QAQ();
}

循环矩阵学习笔记

循环矩阵是一个方阵，大概形式是这样的：
12345
51234
45123
34512
23451
它有一些很好的性质，比如循环矩阵加循环矩阵还是循环矩阵，循环矩阵乘循环矩阵还是一个循环矩阵。
所以我们的矩阵乘法就可以由\(O(n^3)\)优化到\(O(n^2)\)我们只需要把矩阵\(n^2\)的乘出来一行，之后就都是循环的了，最暴力的方法就是\(O(n^2)\)赋值，这样的复杂度是对的，可是它不够优美。我们可以用一个\(n^2\)的循环来代替这个过程，我们只需要计算好每一个位置对哪里有贡献就可以了。
如果矩阵是从0开始的，贡献的位置是\((i+j)\mod n\)
如果矩阵是从1开始的，贡献的位置是\((i+j-2)\mod n+1\)

[P5056][模板]插头DP

插头DP模板题，大力分类讨论就完了。

Code

#include<cstdio>
#include<algorithm>
#include<cstring>
using namespace std;
namespace orz{
const int MOD=233333,N=1000000;
int map[14][14];
int n,m;
int ex,ey;
long long ans;
struct HashMap{
    int head[MOD],val[N],next[N],tot;
    long long cnt[N];
    inline void clear(){
        memset(head,0,sizeof(head));
        tot=0;
    }
    inline long long& operator[](const int &x){
        int pos=x%MOD;
        for(int i=head[pos];i;i=next[i])
            if(val[i]==x)
                return cnt[i];
        next[++tot]=head[pos],head[pos]=tot,val[tot]=x,cnt[tot]=0;
        return cnt[tot];
    }
}f[2];
inline int find(int &s,int id){
    return s>>((id-1)<<1)&3;
}
inline void set(int &s,int id,int val){
    id=(id-1)<<1;
    s&=~(3<<id);
    s|=val<<id;
}
inline int link(int &s,int id){
    int delta=(find(s,id)==1?1:-1);
    int t,cnt=0;
    for(int pos=id;pos<=m+1;pos+=delta){
        t=find(s,pos);
        if(t==1)++cnt;
        else if(t==2)--cnt;
        if(!cnt)return pos;
    }
    return -1;
}
inline void solve(int x,int y){
    HashMap &now=f[((x-1)*m+y)&1];
    HashMap &last=f[(((x-1)*m+y)&1)^1];
    now.clear();
    int tot=last.tot;
    int t1,t2;
    long long val;
    int state;
    for(int i=1;i<=tot;++i){
        state=last.val[i];
        val=last.cnt[i];
        t1=find(state,y);
        t2=find(state,y+1);
        if(!t1&&!t2){
            if(map[x][y]&&map[x+1][y]&&map[x][y+1]){
                set(state,y,1);
                set(state,y+1,2);
                now[state]+=val;
            }
            else if(!map[x][y]){
                now[state]+=val;
            }
        }
        else if(!t1&&t2){
            if(map[x][y+1])now[state]+=val;
            if(map[x+1][y]){
                set(state,y,t2);
                set(state,y+1,0);
                now[state]+=val;
            }   
        }
        else if(t1&&!t2){
            if(map[x+1][y])now[state]+=val;
            if(map[x][y+1]){
                set(state,y+1,t1);
                set(state,y,0);
                now[state]+=val;    
            }
        }
        else if(t1==1&&t2==1){
            set(state,link(state,y+1),1);
            set(state,y,0);
            set(state,y+1,0);
            now[state]+=val;
        }
        else if(t1==1&&t2==2){
            if(x==ex&&y==ey)
                ans+=val;
        }
        else if(t1==2&&t2==1){
            set(state,y,0);
            set(state,y+1,0);
            now[state]+=val;
        }
        else if(t1==2&&t2==2){
            set(state,link(state,y),2);
            set(state,y,0);
            set(state,y+1,0);
            now[state]+=val;
        }
    }
}
inline int read(){
    int a=1,b=0;char t;
    do{t=getchar();if(t=='-')a=-1;}while(t>'9'||t<'0');
    do{b=b*10-'0'+t;t=getchar();}while(t>='0'&&t<='9');
    return a*b;
}
inline int getMap(){
    char t;
    do{t=getchar();}while(t!='*'&&t!='.');
    return t=='.';
}
int QAQ(){
//  freopen("qaq.in","r",stdin);
    n=read(),m=read();
    for(int i=1;i<=n;++i)
        for(int j=1;j<=m;++j)
            map[i][j]=getMap();
    for(int i=1;i<=n;++i)
        for(int j=1;j<=m;++j)
            if(map[i][j])
                ex=i,ey=j;
    f[0][0]=1;
    for(int i=1;i<=n;++i){
        for(int j=1;j<=m;++j)
            solve(i,j);
        if(i!=n){
            HashMap &now=f[(i*m)&1];
            for(int i=1;i<=now.tot;++i)
                now.val[i]<<=2;
        }
    }
    printf("%lld\n",ans);
    return false;
}
}
int main(){
    return orz::QAQ();
}

BigInt类

Code

class BigInt{
#define LEN 6
#define base 100000000
    private:
        int len,s[LEN];
    public:
        BigInt(){
            len=0;
            memset(s,0,sizeof(s));
        }
        BigInt(int x){
            len=0;
            memset(s,0,sizeof(s));
            while(x){
                s[++len]=x%base;
                x/=base;
            }
        }
        BigInt operator+(const BigInt &b)const{
            BigInt res;
            res.len=max(this->len,b.len);
            for(int i=1;i<=res.len;++i)
                res.s[i]+=this->s[i]+b.s[i],res.s[i+1]+=res.s[i]/base,res.s[i]%=base;
            if(res.s[res.len+1])++res.len;
            return res;
        }
        void operator+=(const BigInt &b){
            *this=*this+b;
        }
        inline void print(){
            printf("%d",s[len]);
            for(int i=len-1;i>=0;--i)
                printf("%08d",s[i]);
        }
};

[SDOI2011]地板

这一天，OIer们终于想起了，被插头DP支配的恐惧

毒瘤插头DP题又写了170行。
不过好像是因为我没有去压行。
如果把转移封装起来的话应该能短不少。

在这一题中我们设2中插头，转过向的和没转过向的，然后就可以大力分类讨论了。
最后输出答案可以输出最后一次的状态0。

Code

#include<cstdio>
#include<algorithm>
#include<cstring>
using namespace std;
namespace orz{
const int N=2000000;
const int MOD=499211;
const int P=20110520;
struct HashMap{
    int head[MOD],next[N],val[N],cnt[N],tot;
    HashMap(){
        memset(head,0,sizeof(head));
        tot=0;
    }
    inline void clear(){
        memset(head,0,sizeof(head));
        tot=0;
    }
    inline int& operator[](const int &x){
        int pos=x%MOD;
        for(int i=head[pos];i;i=next[i])
            if(val[i]==x)
                return cnt[i];
        next[++tot]=head[pos],head[pos]=tot,val[tot]=x,cnt[tot]=0;
        return cnt[tot];
    }
}f[2];
inline int find(int &s,int id){
    return (s>>((id-1)<<1))&3;
}
inline void set(int &s,int id,int val){
    id=(id-1)<<1;
    s&=~(3<<id);
    s|=(val<<id);
}
bool qwq[200][200];
bool map[200][200];
int n,m;
inline int read(){
    int a=1,b=0;char t;
    do{t=getchar();if(t=='-')a=-1;}while(t>'9'||t<'0');
    do{b=b*10-'0'+t;t=getchar();}while(t>='0'&&t<='9');
    return a*b;
}
inline bool getMap(){
    char t;
    do{t=getchar();}while(t!='*'&&t!='_');
    return t=='_';
}
inline void show(int s){
    for(int i=1;i<=m+1;++i)
        printf("%d",find(s,i));
}
inline void solve(int x,int y){
    HashMap &now=f[((x-1)*m+y)&1];
    HashMap &last=f[(((x-1)*m+y)&1)^1];
    now.clear();
    int tot=last.tot;
    int t1,t2,s,cnt;
    for(int i=1;i<=tot;++i){
        s=last.val[i];
        cnt=last.cnt[i];
        t1=find(s,y);
        t2=find(s,y+1);
        if(!t1&&!t2){
            if(!map[x][y]){
                set(s,y,0);
                set(s,y+1,0);
                now[s]=(now[s]+cnt)%P;
                continue;
            }
            if(map[x][y]&&map[x+1][y]&&map[x][y+1]){
                set(s,y,2);
                set(s,y+1,2);
                now[s]=(now[s]+cnt)%P;
            }
            if(map[x][y]&&map[x+1][y]){
                set(s,y,1);
                set(s,y+1,0);
                now[s]=(now[s]+cnt)%P;
            }
            if(map[x][y]&&map[x][y+1]){
                set(s,y,0);
                set(s,y+1,1);
                now[s]=(now[s]+cnt)%P;
            }
        }
        else if(!t1&&t2==1){
            if(map[x][y+1]){
                set(s,y,0);
                set(s,y+1,2);
                now[s]=(now[s]+cnt)%P;
            }
            if(map[x+1][y]){
                set(s,y,1);
                set(s,y+1,0);
                now[s]=(now[s]+cnt)%P;
            }
        }
        else if(!t1&&t2==2){
            if(map[x+1][y]){
                set(s,y,2);
                set(s,y+1,0);
                now[s]=(now[s]+cnt)%P;
            }
            set(s,y,0);
            set(s,y+1,0);
            now[s]=(now[s]+cnt)%P;
        }
        else if(t1==1&&!t2){
            if(map[x][y+1]){
                set(s,y,0);
                set(s,y+1,1);
                now[s]=(now[s]+cnt)%P;
            }
            if(map[x+1][y]){
                set(s,y,2);
                set(s,y+1,0);
                now[s]=(now[s]+cnt)%P;
            }
        }
        else if(t1==1&&t2==1){
            set(s,y,0);
            set(s,y+1,0);
            now[s]=(now[s]+cnt)%P;
        }
        else if(t1==2&&!t2){
            if(map[x][y+1]){
                set(s,y,0);
                set(s,y+1,2);
                now[s]=(now[s]+cnt)%P;
            }
            set(s,y,0);
            set(s,y+1,0);
            now[s]=(now[s]+cnt)%P;
        }
    }
}
int QAQ(){
    n=read(),m=read();
    for(int i=1;i<=n;++i)
        for(int j=1;j<=m;++j)
            qwq[i][j]=getMap();
    if(n<m){
        for(int i=1;i<=n;++i)
            for(int j=1;j<=m;++j)
                map[j][i]=qwq[i][j];
        swap(n,m);
    }
    else {
        for(int i=1;i<=n;++i)
            for(int j=1;j<=m;++j)
                map[i][j]=qwq[i][j];
    }
    f[0][0]=1;
    for(int i=1;i<=n;++i){
        for(int j=1;j<=m;++j)
            solve(i,j);
        if(i!=n){
            HashMap &last=f[(i*m)&1];
            for(int i=1;i<=last.tot;++i)
                last.val[i]<<=2;
        }
    }
    HashMap &ans=f[(n*m)&1];
    printf("%d\n",ans[0]);
    return false;
}
}
int main(){
    return orz::QAQ();
}

[COGS775]山海经

听说是线段树恶心题，结果一遍过了。
就是区间最长子段和，但是要求询问具体方案。而且不设SpecialJudge，要求输出字典序最小的解。
区间最长子段和很容易维护，细节全在答案的维护。
因为要求字典序最小，所以我们在答案相等是时候也尽量要取靠左的，所以很容易就A了。

Code

#include<cstdio>
#include<algorithm>
using namespace std;
namespace orz{
#define lc x<<1
#define rc x<<1|1
const int N=1100000;
int a[N];
int sum[N];
struct node{
    int l,r;
    int ans,lans,rans;
    int lp,rp,llp,rrp;
}t[(N<<2)+233];
inline int calc(const int &l,const int &r){
    return sum[r]-sum[l-1]; 
}
inline node merge(node a,node b){
    node res;
    //先维护无关紧要的信息
    res.l=a.l;
    res.r=b.r;
    //=================
    //从左往右更新值
    //更新左答案，左答案位置
    //左答案的右端点要尽量小
    res.lans=a.lans;
    res.llp=a.llp;
    if(calc(a.l,a.r)+b.lans>res.lans){
        res.llp=b.llp;
        res.lans=calc(a.l,a.r)+b.lans;
    }
    //更新右答案，右答案位置
    res.rans=b.rans;
    res.rrp=b.rrp;
    if(calc(b.l,b.r)+a.rans>=res.rans){
        res.rrp=a.rrp;
        res.rans=calc(b.l,b.r)+a.rans;
    }
    //更新答案，更新答案位置
    if(a.ans>=b.ans){
        res.ans=a.ans;
        res.lp=a.lp;
        res.rp=a.rp;
    }
    else{
        res.ans=b.ans;
        res.lp=b.lp;
        res.rp=b.rp;
    }
    if(a.rans+b.lans>res.ans){
        res.ans=a.rans+b.lans;
        res.lp=a.rrp;
        res.rp=b.llp;
    }
    else if(a.rans+b.lans==res.ans){
        if(a.rrp<res.lp){
            res.ans=a.rans+b.lans;
            res.lp=a.rrp;
            res.rp=b.llp;
        }
    }
    return res;
}
void build(int x,int l,int r){
    if(l==r){
        t[x].ans=t[x].lans=t[x].rans=a[l];
        t[x].l=t[x].r=t[x].lp=t[x].rp=t[x].llp=t[x].rrp=l;
        return ;
    }
    int mid=(l+r)>>1;
    build(lc,l,mid);
    build(rc,mid+1,r);
    t[x]=merge(t[lc],t[rc]);
}
inline node query(int x,int l,int r){
    if(t[x].l==l&&t[x].r==r)return t[x];
    int mid=(t[x].l+t[x].r)>>1;
    if(r<=mid)return query(lc,l,r);
    else if(l>mid)return query(rc,l,r);
    else return merge(query(lc,l,mid),query(rc,mid+1,r));
}
inline int read(){
    int a=1,b=0;char t;
    do{t=getchar();if(t=='-')a=-1;}while(t>'9'||t<'0');
    do{b=b*10-'0'+t;t=getchar();}while(t>='0'&&t<='9');
    return a*b;
}
int QAQ(){
    int n,m;
    int l,r;
    node res;
    n=read(),m=read();
    for(int i=1;i<=n;++i)
        a[i]=read();
    for(int i=1;i<=n;++i)
        sum[i]=sum[i-1]+a[i];
    build(1,1,n);
    while(m--){
        l=read(),r=read();
        res=query(1,l,r);
        printf("%d %d %d\n",res.lp,res.rp,res.ans);
    }
    return false;
}
}
int main(){
    return orz::QAQ();
}

[LOJ10222]佳佳的Fibonacci

矩阵乘法模板题
我们把nfn的式子拆开就可以求出矩阵
乘就完事了

Code

#include<cstdio>
#include<algorithm>
#include<cstring>
using namespace std;
namespace orz{
const int N=6;
int MOD;
struct Matrix{
    int a[N][N];
    Matrix(){
        memset(a,0,sizeof(a));
    }
    Matrix operator*(const Matrix &b)const{
        Matrix res;
        for(int i=1;i<=5;++i)
            for(int j=1;j<=5;++j)
                for(int k=1;k<=5;++k)
                    res.a[i][j]=(res.a[i][j]+1ll*a[i][k]*b.a[k][j]%MOD)%MOD;
        return res;
    }
}ans,base;
inline int read(){
    int a=1,b=0;char t;
    do{t=getchar();if(t=='-')a=-1;}while(t>'9'||t<'0');
    do{b=b*10-'0'+t;t=getchar();}while(t>='0'&&t<='9');
    return a*b;
}
int QAQ(){
    int n=read();
    MOD=read();
    Matrix ans,base;
    base.a[1][1]=1;
    base.a[2][1]=1;
    base.a[2][2]=1;
    base.a[3][2]=1;
    base.a[4][2]=1;
    base.a[5][2]=2;
    base.a[2][3]=1;
    base.a[4][4]=1;
    base.a[5][4]=1;
    base.a[4][5]=1;
    ans.a[1][1]=1;
    ans.a[1][2]=2;
    ans.a[1][3]=1;
    ans.a[1][4]=1;
    ans.a[1][5]=1;
    int y=n-1;
    while(y){
        if(y&1)ans=ans*base;
        y>>=1;
        base=base*base;
    }
    printf("%d\n",ans.a[1][1]);
    return false;
}
}
int main(){
    return orz::QAQ();
}

7.29爆零赛

爆零了QAQ
上来一看，T1扫描线，T2启发式合并平衡树，T3期望不会。
自闭了，打了暴力分。

这三道题目名字也挺不正经的说。

辣鸡

这题其实你如果要扫描线拧干的话估计也能干出来，可是我并不会扫描线。
这题\(n^2\)就完事了，循环的时候加一个小优化，先把它排了序，如果以后再也不能对答案贡献了就跳出。
然后就可以见证\(n^2\)过\(100000\)了。

Code

#include<cstdio>
#include<algorithm>
using namespace std;
namespace orz{
const int N=200000;
struct point{
    long long x1,y1,x2,y2;
    bool operator<(const point &b)const{
        return this->x1^b.x1?this->x1<b.x1:this->y1<b.y1;
    }
}a[N];
inline int read(){
    int a=1,b=0;char t;
    do{t=getchar();if(t=='-')a=-1;}while(t>'9'||t<'0');
    do{b=b*10-'0'+t;t=getchar();}while(t>='0'&&t<='9');
    return a*b;
}
inline long long calc(long long l1,long long r1,long long l2,long long r2){
    long long res=0,tl=max(l1-1,l2),tr=min(r1-1,r2);
    res+=max(0ll,tr-tl+1);
    tl=max(l1+1,l2),tr=min(r1+1,r2);
    res+=max(0ll,tr-tl+1);
    return res;
}
int QAQ(){
    int n=read();
    long long ans=0;
    for(register int i=1;i<=n;++i)
        a[i].x1=read()+1,a[i].y1=read()+1,a[i].x2=read()+1,a[i].y2=read()+1;
    for(register int i=1;i<=n;++i)
        ans+=(a[i].x2-a[i].x1)*(a[i].y2-a[i].y1)*2ll;
    sort(a+1,a+n+1);
    for(register int i=1;i<=n;++i){
        for(register int j=i+1;j<=n;++j){
            if(a[j].x1>a[i].x2+1)break;
            if(a[j].x1==a[i].x2+1){ans+=calc(a[i].y1,a[i].y2,a[j].y1,a[j].y2);}
            else if(a[j].y1==a[i].y2+1){ans+=calc(a[i].x1,a[i].x2,a[j].x1,a[j].x2);}
            else if(a[j].x2==a[i].x1-1){ans+=calc(a[i].y1,a[i].y2,a[j].y1,a[j].y2);}
            else if(a[j].y2==a[i].y1-1){ans+=calc(a[i].x1,a[i].x2,a[j].x1,a[j].x2);}
        }
    }
    printf("%lld\n",ans);
    return false;
}
}
int main(){
    return orz::QAQ();
}

模板

题解是启发式合并动态开点线段树，因为我比较懒，所以还是去写启发式合并Splay了。
平衡树启发式合并跑的还挺快的。
过了样例之后出了两个问题，一个是没判断k=0的情况T了，第二个是边数组没开2倍T了。
Orz。

Code

#include<cstdio>
#include<algorithm>
#include<set>
using namespace std;
namespace orz{
#define IT set<WSL>::iterator
const int N=110000;
struct SplayNode{
    int son[2],father,size,cnt,val,time,color;
}t[N*18];
int tot=0;
int res;
int head[N],ver[N*2],next[N*2],cwy;
inline void add(int x,int y){
    next[++cwy]=head[x],head[x]=cwy,ver[cwy]=y;
    next[++cwy]=head[y],head[y]=cwy,ver[cwy]=x;
}
struct WSL{
    int color;
    mutable int time;
    WSL(int c,int t){
        color=c;
        time=t;
    }
    bool operator<(const WSL&b)const{
        return this->color<b.color;
    }
};
struct Splay{
    int root;
    inline void update(int x){
        t[x].size=t[t[x].son[0]].size+t[t[x].son[1]].size+1;
        t[x].cnt=t[t[x].son[0]].cnt+t[t[x].son[1]].cnt+t[x].val;
    }
    inline bool get(int x){
        return t[t[x].father].son[1]==x;
    }
    inline void rotate(int x){
        int y=t[x].father,z=t[y].father,k=get(x);
        t[z].son[get(y)]=x;
        t[x].father=z;
        t[y].son[k]=t[x].son[k^1];
        t[t[y].son[k]].father=y;
        t[x].son[k^1]=y;
        t[y].father=x;
        update(y);
    }
    inline void splay(int x){
        int y,z;
        while(t[x].father){
            y=t[x].father;z=t[y].father;
            if(z)
                rotate(get(x)^get(y)?x:y);
            rotate(x);
        }
        update(x);
        root=x;
    }
    inline int findAns(int k){
        int x=root;
        int ans=0;
        if(k>=t[x].size){
            return t[x].cnt;
        }
        if(k==0){
            return 0;
        }
        while(true){
            if(k<=t[t[x].son[0]].size)x=t[x].son[0];
            else if(k<=t[t[x].son[0]].size+1)return ans+t[t[x].son[0]].cnt+t[x].val;
            else ans+=t[t[x].son[0]].cnt+t[x].val,k-=t[t[x].son[0]].size+1,x=t[x].son[1];
        }
    }
    inline void insert(int time,int color,int val){
        int father=0,x=root;
        while(x&&t[x].time!=time){
            father=x;
            x=t[x].son[time>t[x].time];
        }
        x=++tot;
        if(father)t[father].son[time>t[father].time]=x;
        t[x].time=time;t[x].size=1;t[x].cnt=t[x].val=val;
        t[x].father=father;t[x].color=color;
        splay(x);
    }
    inline void change(int time){
        int x=root;
        while(x&&time!=t[x].time){
            x=t[x].son[time>t[x].time];
        }
        t[x].val=0;
        splay(x);
    }
    inline int size(){
        return t[root].size;
    }
};
struct QWQ{
    set<WSL>dky;
    Splay s;
    inline int size(){
        return s.size();
    }
    inline void insert(int time,int color){
        IT it=dky.find(WSL(color,0));
        if(it==dky.end()){
            dky.insert(WSL(color,time));
            s.insert(time,color,1);
        }
        else if(time>it->time){
            s.insert(time,color,0);
        }
        else{
            s.change(it->time);
            it->time=time;
            s.insert(time,color,1);
        }
    }
}a[N];
int p[N],tong[N],ans[N];
inline int read(){
    int a=1,b=0;char t;
    do{t=getchar();if(t=='-')a=-1;}while(t>'9'||t<'0');
    do{b=b*10-'0'+t;t=getchar();}while(t>='0'&&t<='9');
    return a*b;
}
inline void merge(int x){
    if(!x)return;
    merge(t[x].son[0]);
    a[res].insert(t[x].time,t[x].color);
    merge(t[x].son[1]);
}
inline void dfs(int x,int father){
    for(int i=head[x];i;i=next[i]){
        if(ver[i]==father)continue;
        dfs(ver[i],x);
        if(a[p[x]].size()<a[p[ver[i]]].size()){
            res=p[ver[i]];
            merge(a[p[x]].s.root);
            p[x]=p[ver[i]];
        }
        else{
            res=p[x];
            merge(a[p[ver[i]]].s.root);
        }
    }
    ans[x]=a[p[x]].s.findAns(tong[x]);
}
int QAQ(){
    int n,m,q,x,c;
    n=read();
    for(int i=1;i<n;++i)
        add(read(),read());
    for(int i=1;i<=n;++i)
        tong[i]=read(),p[i]=i;
    m=read();
    for(int i=1;i<=m;++i)
        x=read(),c=read(),a[p[x]].insert(i,c);  
    dfs(1,0);
    q=read();
    for(int i=1;i<=q;++i)
        printf("%d\n",ans[read()]);
    return false;
}
}
int main(){
    return orz::QAQ();
}

大佬

题解有看不懂的DP做法，实际上有很优秀的快速幂做法，时间复杂度是\(m\log k\)的，比题解的\(n^2m\)不知道高到哪里去了。
开始的时候我设的dp状态一维是天数，因为前面对后面有影响导致不可转移。实际上我们的所有情况包括了m取值的所有情况，而所有长度为k的区间的所有情况也都存在，所以我们只需要考虑长度为k的区间，乘上\(n-k+1\)就可以了。
而对于长度为k的区间，最大值为i的方案数为\(i^k-(i-1)^k\)

Code

#include<cstdio>
#include<algorithm>
using namespace std;
namespace orz{
const long long MOD=1000000007;
inline long long read(){
    long long a=1,b=0;char t;
    do{t=getchar();if(t=='-')a=-1;}while(t>'9'||t<'0');
    do{b=b*10-'0'+t;t=getchar();}while(t>='0'&&t<='9');
    return a*b;
}
inline long long pow(long long x,long long y){
    long long ans=1;
    while(y){
        if(y&1)ans=ans*x%MOD;
        y>>=1;
        x=x*x%MOD;
    }
    return ans%MOD;
}
int QAQ(){
    long long n=read(),m=read(),k=read();
    long long ans=0;
    if(n<k){printf("0\n");return false;}
    for(int i=1;i<=m;++i)
        ans=(ans+read()*((pow(i,k)-pow(i-1,k))%MOD+MOD)%MOD*(n-k+1)%MOD)%MOD;
    printf("%lld\n",ans*pow(pow(m,k),MOD-2)%MOD);
    return false;
}
}
int main(){
    return orz::QAQ();
}

[SCOI2014]方伯伯的玉米田

这是个DP
看起来好像是一个LIS，我们考虑如何处理区间同时加一。
因为我们求的是最长不下降子序列，所以一次区间加一定是从一个点开始加到最后，否则会导致后面的变小。
所以我们设\(f_{i,j}\)表示考虑了前i个玉米，进行了j次拔高的最大长度，之后就是一个LIS了。
注意第二层循环要反着来，否则会自己加自己。

Code

#include<cstdio>
#include<algorithm>
using namespace std;
namespace orz{
const int N=20000;
const int mx=5550;
const int mxk=550;
int a[N];
int c[5600][600];
int n,k;
inline int ask(int x,int y){
    int ans=0;
    for(int i=x;i;i-=i&-i)
        for(int j=y;j;j-=j&-j)
            ans=max(ans,c[i][j]);
    return ans;
}
inline void add(int x,int y,int val){
    for(int i=x;i<=mx;i+=i&-i)
        for(int j=y;j<=mxk;j+=j&-j)
            c[i][j]=max(c[i][j],val);
}
inline int read(){
    int a=1,b=0;char t;
    do{t=getchar();if(t=='-')a=-1;}while(t>'9'||t<'0');
    do{b=b*10-'0'+t;t=getchar();}while(t>='0'&&t<='9');
    return a*b;
}
int QAQ(){
    n=read(),k=read();
    for(int i=1;i<=n;++i)
        a[i]=read();
    int ans=0;
    int fans=0;
    for(int i=1;i<=n;++i)
        for(int j=k+1;j;--j){
            ans=ask(a[i]+j,j)+1;
            add(a[i]+j,j,ans);
            fans=max(fans,ans);
        }
    printf("%d\n",fans);
    return false;
}
}
int main(){
    return orz::QAQ();
}

[SDOI2011]拦截导弹

第一问很简单，裸的CDQ优化DP。
第二问就很难计算了。
它大概是一个假的概率，等于经过某个导弹的方案数除以总方案数。
我们考虑如何判断一个点是否在x到y的最短路上，我们可以从x跑一遍最短路，从y跑一遍最短路，如果一个点两个dis的值加起来等于x到y的最短路长度加一，那么这个点就是好的。这个显然。
我们用相同的思路处理这一道题。
我们先正着跑一遍最长不上升子序列,dp数组为\(f_i\),方案数为\(fc_i\)，再倒着跑一遍最长不下降子序列,dp数组为\(g_i\)，方案数为\(gc_i\)如果\(f_i+g_i=ans\)那么经过它的方案数为\(fc_i*gc_i\)。

之后我们发现我们还是不会求，CDQ怎么统计方案数？
我们发现我们的树状数组只能求出最优方案，而不能求出方案数。
这时候就需要我们去魔改一下我们的树状数组了。
我们再记一个cnt数组，我们在更新最大值的时候把它更新，如果最大值和插入值相等就累加。
然后我们发现跑出来的全是0，因为我们的计数没有初值,把所有数初值赋为0吗？那你必WA，不用想也知道那样一定会出错。
我们发现只有一个数作为第一个数出现时才会从0转移，所以当右区间的一个数从0转移时我们将它跳过，最后跑到它本身如果还是0那说明它是从0转移的，我们为它赋上初值。
这题方案数还挺多的，需要开double，比如10 10 9 9 8 8 7 7 6 6 5 5 4 4 3 3 2 2 1 1就有\(2^{10}\)种总方案。

Code

#include<cstdio>
#include<algorithm>
using namespace std;
namespace orz{
const int N=200000;
struct node{
    int t,h,v;
    bool operator<(const node&b)const{
        return this->h<b.h;
    }
}a[N],t[N];
int hc[N],vc[N],htot,vtot,c[N];
double cnt[N];
int n,f[N],g[N];
double fc[N],gc[N];
inline void add(int x,int ans,double tot){
    for(;x<=n+2;x+=x&-x){
        if(c[x]<ans)c[x]=ans,cnt[x]=tot;
        else if(ans==c[x])cnt[x]+=tot;
    }
}
inline int ask(int x,double &res){
    res=0;
    int ans=0;
    for(;x;x-=x&-x){
        if(ans<c[x])ans=c[x],res=cnt[x];
        else if(ans==c[x])res+=cnt[x];
    }
    return ans;
}
inline void remove(int x){
    for(;x<=n;x+=x&-x)
        c[x]=cnt[x]=0;
}
inline int read(){
    int a=1,b=0;char t;
    do{t=getchar();if(t=='-')a=-1;}while(t>'9'||t<'0');
    do{b=b*10-'0'+t;t=getchar();}while(t>='0'&&t<='9');
    return a*b;
}
void cdq(int l,int r){
    if(l==r){
        if(!f[a[l].t]){
            f[a[l].t]=1;
            fc[a[l].t]=1;
        }
        return;
    }
    int mid=(l+r)>>1;
    cdq(l,mid);
    int tl=l,tr=mid+1,cwy;
    double res;
    for(int i=l;i<=r;++i)t[i]=a[i];
    sort(t+l,t+mid+1);
    sort(t+mid+1,t+r+1);
    for(int i=l;i<=r;++i){
        if((tl<=mid&&t[tl].h<=t[tr].h)||tr>r){
            add(t[tl].v,f[t[tl].t],fc[t[tl].t]);
            ++tl;
        }
        else{
            cwy=ask(t[tr].v,res)+1;
            if(cwy==1){++tr;continue;}
            if(f[t[tr].t]<cwy)f[t[tr].t]=cwy,fc[t[tr].t]=res;
            else if(f[t[tr].t]==cwy)fc[t[tr].t]+=res;
            ++tr;
        }
    }
    for(int i=l;i<=mid;++i)
        remove(t[i].v);
    cdq(mid+1,r);
}
int QAQ(){
    n=read();
    for(int i=1;i<=n;++i)
        a[i].h=hc[i]=read(),a[i].v=vc[i]=read(),a[i].t=i;
    sort(hc+1,hc+n+1);
    sort(vc+1,vc+n+1);
    htot=unique(hc+1,hc+n+1)-hc-1;
    vtot=unique(vc+1,vc+n+1)-vc-1;
    for(int i=1;i<=n;++i){
        a[i].h=lower_bound(hc+1,hc+htot+1,a[i].h)-hc;
        a[i].v=lower_bound(vc+1,vc+vtot+1,a[i].v)-vc;
    }
    for(int i=1;i<=n;++i)a[i].h=n+1-a[i].h;
    for(int i=1;i<=n;++i)a[i].v=n+1-a[i].v;
    cdq(1,n);
    for(int i=1;i<=n;++i)a[i].h=n+1-a[i].h;
    for(int i=1;i<=n;++i)a[i].v=n+1-a[i].v;
    for(int i=1;i<=n;++i)g[i]=f[i],gc[i]=fc[i];
    for(int i=1;i<=n;++i)f[i]=fc[i]=0;
    for(int i=1,j=n;i<j;++i,--j)swap(a[i],a[j]);
    cdq(1,n);
    int ans=0;
    double tot=0;
    for(int i=1;i<=n;++i)
        ans=max(ans,g[i]);
    for(int i=1;i<=n;++i)
        if(g[i]==ans)
            tot+=gc[i];
    printf("%d\n",ans);
    for(int i=1;i<=n;++i){
        if(f[i]+g[i]==ans+1)printf("%lf ",(fc[i]*gc[i])/tot);
        else printf("0.0000000 ");
    }
    return false;
}
}
int main(){
    return orz::QAQ();
}

后缀数组(SA)学习笔记

由于HZ学长讲的过于清晰，所以蒟蒻并没有听懂。
以下内容是我从网上各位大神的博客学习后总结的，感谢他们的无私付出。
后半部分基本来自2009年罗穗骞的论文，感谢他写出了如此清晰的国家集训队论文，这是少数几个蒟蒻也能看懂的国集论文。

没听懂QAQ
后缀数组是处理字符串的一个好东西。
它是一个数组。
我们先把所有的后缀排个序，然后我们就可以得到后缀数组。
以下所有第i个指的是原字符串i到len的子串。
排名第i的所有后缀排过序后的第i个后缀。
这里的排序指的是按照字典序排序。
a<ab
一些数组定义(以下都是数组)：
rank表示第i个后缀的排名
sa表示排名第i的后缀是第几个后缀
还有一个常用的height数组之后再定义吧。

我们可以发现rank数组和sa数组是互逆的，求一个就可以了。
但是怎样优秀的求就是一个问题。
重载运算符?
复杂度为O(n^2logn)

有一些优秀的算法可以解决这个问题。
其中倍增法最为常用。
时间复杂度为\(O(n\log n)\)

倍增算法是基于后缀的一些性质。
第i+1的后缀的第一个字符是第i个后缀的第二个字符。
所以我们在算出第一个字符顺序的同时也得到了第二个字符的排序，我们算出1,2个字符的排序后也得到了第3,4个字符的排序。
所以我们就可以倍增的用双关键字的排序来做了。
直接用是sort是\(O(n\log^2n)\)的。
用基数排序会更优秀。

写了很多注释的代码

#include<cstdio>
#include<algorithm>
#include<cstring>
using namespace std;
namespace orz{
const int N=1100000;
char s[N];
int n,m;
int x[N],y[N],c[N],sa[N];
inline void SA(){
    //c是一个桶，现在统计了对应元素的个数
    for(int i=1;i<=n;++i)++c[x[i]=s[i]];
    //求了一个前缀和
    for(int i=1;i<=m;++i)c[i]+=c[i-1];
    //求出只用第一个字母的SA值
    //如果所有字符都不相同的话，对应的排名就是前缀和
    //为了保证sa互不相同要每次减一
    //这里正着倒着都对
    //其实好像乱着也对
    for(int i=n;i>=1;--i)sa[c[x[i]]--]=i;
    //开始倍增
    //k的意义为目前已经处理出了关键字为前的SA，扩展到2×k的SA
    for(int k=1;k<=n;k<<=1){
        //x表示的是第i个后缀当前的rank，相同rank会重复。
        //p只是一个计数器
        int p=0;
        //y表示以第二关键字排序的SA值
        //先把长度不够的给排到前面
        //这里正着倒着都对，因为之前已经排好了
        for(int i=n-k+1;i<=n;++i)y[++p]=i;
        //倍增法求SA的核心操作
        //第二关键字可以从第一关键字得到
        for(int i=1;i<=n;++i)if(sa[i]>k)y[++p]=sa[i]-k;
        //清空桶
        for(int i=1;i<=m;++i)c[i]=0;
        //用第一关键字信息更新桶
        for(int i=1;i<=n;++i)++c[x[i]];
        //做前缀和
        for(int i=1;i<=m;++i)c[i]+=c[i-1];
        //与上面的排序差不多，只是按y的顺序来
        for(int i=n;i>=1;--i)sa[c[x[y[i]]]--]=y[i];
        //更新出下一次的x
        //现在的SA数组已经是按2k排序的SA数组了
        //先把原来的y废了
        //现在的y是上一次的x
        swap(x,y);
        //计数器清空
        p=1;
        //设置第一位
        x[sa[1]]=1;
        //如果相同的话拥有相同rank
        for(int i=2;i<=n;++i)
            x[sa[i]]=y[sa[i]]==y[sa[i-1]]&&y[sa[i]+k]==y[sa[i-1]+k]?p:++p;
        //如果p=n说明所有后缀都排好了
        if(p==n)break;
        m=p;
    }
}
inline int read(){
    int a=1,b=0;char t;
    do{t=getchar();if(t=='-')a=-1;}while(t>'9'||t<'0');
    do{b=b*10-'0'+t;t=getchar();}while(t>='0'&&t<='9');
    return a*b;
}
int QAQ(){
    scanf("%s",s+1);
    n=strlen(s+1),m='z';
    SA();
    for(int i=1;i<=n;++i)
        printf("%d ",sa[i]);
    return false;
}
}
int main(){
    return orz::QAQ();
}

没有注释的板子

#include<cstdio>
#include<algorithm>
#include<cstring>
using namespace std;
namespace orz{
const int N=1100000;
char s[N];
int x[N],y[N],c[N],sa[N];
int n,m;
inline int read(){
    int a=1,b=0;char t;
    do{t=getchar();if(t=='-')a=-1;}while(t>'9'||t<'0');
    do{b=b*10-'0'+t;t=getchar();}while(t>='0'&&t<='9');
    return a*b;
}
inline void SA(){
    for(int i=1;i<=n;++i)++c[x[i]=s[i]];
    for(int i=1;i<=m;++i)c[i]+=c[i-1];
    for(int i=n;i;--i)sa[c[x[i]]--]=i;
    for(int k=1;k<=n;k<<=1){
        int p=0;
        for(int i=n-k+1;i<=n;++i)y[++p]=i;
        for(int i=1;i<=n;++i)if(sa[i]>k)y[++p]=sa[i]-k;
        for(int i=1;i<=m;++i)c[i]=0;
        for(int i=1;i<=n;++i)++c[x[i]];
        for(int i=1;i<=m;++i)c[i]+=c[i-1];
        for(int i=n;i;--i)sa[c[x[y[i]]]--]=y[i];
        swap(x,y);
        p=1;
        x[sa[1]]=1;
        for(int i=2;i<=n;++i)
            x[sa[i]]=y[sa[i]]==y[sa[i-1]]&&y[sa[i]+k]==y[sa[i-1]+k]?p:++p;
        if(p==n)break;
        m=p;
    }
}
int QAQ(){
    scanf(" %s",s+1);
    n=strlen(s+1);
    m='z';
    SA();
    for(int i=1;i<=n;++i)
        printf("%d ",sa[i]);
    return false;
}
}
int main(){
    return orz::QAQ();
}

我们终于出了SA数组，可是仍然不知道该怎么用。
我们一般还需要height数组。
\(height_i\)表示排名为i的后缀与排名为i-1的后缀的LCP
如果直接求就\(O(n^2)\)了。
那就很不好
设\(h_i=height_{rank_i}\)
即h数组表示第i个后缀的height值
我们有一个性质:\(h_i>=h_{i-1}-1\)

证明

我们设k为第i-1个后缀排名前一个后缀
则\(lcp(k,i-1)=h_{i-1}\)
把这两个后缀都去掉首字符
变成了k+1和i
而它们的lcp少了1
这个只是方便理解的证明，实际上好像不是佷严谨。

带求height的板子

inline void SA(){
    for(int i=1;i<=n;++i)++c[x[i]=a[i]];
    for(int i=1;i<=m;++i)c[i]+=c[i-1];
    for(int i=1;i<=n;++i)sa[c[x[i]]--]=i;
    for(int k=1;k<=n;k<<=1){
        int p=0;
        for(int i=n-k+1;i<=n;++i)y[++p]=i;
        for(int i=1;i<=n;++i)if(sa[i]>k)y[++p]=sa[i]-k;
        for(int i=1;i<=m;++i)c[i]=0;
        for(int i=1;i<=n;++i)++c[x[i]];
        for(int i=1;i<=m;++i)c[i]+=c[i-1];
        for(int i=n;i;--i)sa[c[x[y[i]]]--]=y[i];
        swap(x,y);
        p=1;
        x[sa[1]]=1;
        for(int i=2;i<=n;++i)
            x[sa[i]]=y[sa[i]]==y[sa[i-1]]&&y[sa[i]+k]==y[sa[i-1]+k]?p:++p;
        if(p==n)break;
        m=p;
    }
    for(int i=1;i<=n;++i)
        rank[sa[i]]=i;
    for(int i=1,p=0;i<=n;++i){
        if(rank[i]==1)continue;
        if(p)--p;
        int j=sa[rank[i]-1];
        while(a[i+p]==a[j+p])++p;
        height[rank[i]]=p;
    }
}

SA一些常用套路的总结

最长可重叠重复子串

重复子串的定义：在字符串中至少出现两次的子串
输出height的最大值就可以了

最长不可重叠重复子串

最长可重叠k重复子串

我们依旧是把问题转化为判定，使用二分来解决.
二分长度。如果存在一段连续的height都大于等于这个长度且这些数的个数大于等于k-1，那么这个长度就是可行的。

USACO06DEC牛奶模式Milk Patterns

Code

#include<cstdio>
#include<algorithm>
using namespace std;
namespace orz{
const int N=40000;
int x[N],y[N],c[1100000],sa[N],rank[N],height[N];
int a[N];
int n,k,m;
inline int read(){
    int a=1,b=0;char t;
    do{t=getchar();if(t=='-')a=-1;}while(t>'9'||t<'0');
    do{b=b*10-'0'+t;t=getchar();}while(t>='0'&&t<='9');
    return a*b;
}
inline void SA(){
    for(int i=1;i<=n;++i)++c[x[i]=a[i]];
    for(int i=1;i<=m;++i)c[i]+=c[i-1];
    for(int i=1;i<=n;++i)sa[c[x[i]]--]=i;
    for(int k=1;k<=n;k<<=1){
        int p=0;
        for(int i=n-k+1;i<=n;++i)y[++p]=i;
        for(int i=1;i<=n;++i)if(sa[i]>k)y[++p]=sa[i]-k;
        for(int i=1;i<=m;++i)c[i]=0;
        for(int i=1;i<=n;++i)++c[x[i]];
        for(int i=1;i<=m;++i)c[i]+=c[i-1];
        for(int i=n;i;--i)sa[c[x[y[i]]]--]=y[i];
        swap(x,y);
        p=1;
        x[sa[1]]=1;
        for(int i=2;i<=n;++i)
            x[sa[i]]=y[sa[i]]==y[sa[i-1]]&&y[sa[i]+k]==y[sa[i-1]+k]?p:++p;
        if(p==n)break;
        m=p;
    }
    for(int i=1;i<=n;++i)
        rank[sa[i]]=i;
    for(int i=1,p=0;i<=n;++i){
        if(rank[i]==1)continue;
        if(p)--p;
        int j=sa[rank[i]-1];
        while(a[i+p]==a[j+p])++p;
        height[rank[i]]=p;
    }
}
inline bool check(int x){
    int res=0;
    for(int i=2;i<=n;++i){
        if(height[i]>=x)++res;
        else res=0;
        if(res>=k-1)return true;
    }   
    return false;
}
int QAQ(){
    n=read(),k=read();
    for(int i=1;i<=n;++i)
        m=max(m,a[i]=read());
    SA();
    int l=0,r=n,mid;
    while(l<r){
        mid=(l+r+1)>>1;
        if(check(mid))l=mid;
        else r=mid-1;
    }
    printf("%d\n",l);
    return false;
}
}
int main(){
    return orz::QAQ();
}

不同子串个数

一个很显然的转化是将子串变成后缀的前缀。
所以问题也就变成了求不同的后缀的前缀的个数
我们按照排完的顺序依次插入后缀，每插入一个后缀，都会多\(n-sa_i-height_i+1\)个本质不同的子串。把它们加起来就可以得到个数，注意第一个后缀是\(n-sa_1+1\)

LuoguP2408不同子串个数

Code

#include<cstdio>
#include<algorithm>
#include<cstring>
using namespace std;
namespace orz{
const int N=200000;
int x[N],y[N],c[N],sa[N],height[N],rank[N];
char s[N];
int n,m;
inline int read(){
    int a=1,b=0;char t;
    do{t=getchar();if(t=='-')a=-1;}while(t>'9'||t<'0');
    do{b=b*10-'0'+t;t=getchar();}while(t>='0'&&t<='9');
    return a*b;
}
inline void SA(){
    for(int i=1;i<=n;++i)++c[x[i]=s[i]];
    for(int i=1;i<=m;++i)c[i]+=c[i-1];
    for(int i=n;i;--i)sa[c[x[i]]--]=i;
    for(int k=1;k<=n;k<<=1){
        int p=0;
        for(int i=n-k+1;i<=n;++i)y[++p]=i;
        for(int i=1;i<=n;++i)if(sa[i]>k)y[++p]=sa[i]-k;
        for(int i=1;i<=m;++i)c[i]=0;
        for(int i=1;i<=n;++i)++c[x[i]];
        for(int i=1;i<=m;++i)c[i]+=c[i-1];
        for(int i=n;i;--i)sa[c[x[y[i]]]--]=y[i];
        swap(x,y);
        p=1;
        x[sa[1]]=1;
        for(int i=2;i<=n;++i)
            x[sa[i]]=y[sa[i]]==y[sa[i-1]]&&y[sa[i]+k]==y[sa[i-1]+k]?p:++p;
        if(p==n)break;
        m=p;
    }
    for(int i=1;i<=n;++i)
        rank[sa[i]]=i;
    for(int i=1,p=0;i<=n;++i){
        if(rank[i]==1)continue;
        if(p)--p;
        int j=sa[rank[i]-1];
        while(s[i+p]==s[j+p])++p;
        height[rank[i]]=p;
    }
}
int QAQ(){
    n=read();
    scanf("%s",s+1);
    m='z';
    SA();
//  for(int i=1;i<=n;++i)
//      printf("%d %d\n",sa[i],height[i]);
//  putchar('\n');
    long long ans=0;
    ans+=n-sa[1]+1;
    for(int i=2;i<=n;++i)
        ans+=(long long)n-sa[i]+1-height[i];
    printf("%lld\n",ans);
    return false;
}
}
int main(){
    return orz::QAQ();
}

对于这个思路，我们还可以发现这样插入的子串是按照字典序从小到大的，所以我们还可以二分取出第k大的本质不同的字串。

Code

#include<cstdio>
#include<algorithm>
#include<cstring>
using namespace std;
namespace orz{
const int N=100;
int x[N],y[N],c[N],sa[N],rank[N],height[N];
long long sum[N];
char s[N];
int f[N][20];
int Log[N];
int n,m;
inline long long read(){
    long long a=1,b=0;char t;
    do{t=getchar();if(t=='-')a=-1;}while(t>'9'||t<'0');
    do{b=b*10-'0'+t;t=getchar();}while(t>='0'&&t<='9');
    return a*b;
}
inline void SA(){
    for(int i=1;i<=n;++i)++c[x[i]=s[i]];
    for(int i=1;i<=m;++i)c[i]+=c[i-1];
    for(int i=n;i;--i)sa[c[x[i]]--]=i;
    for(int k=1;k<=n;k<<=1){
        int p=0;
        for(int i=n-k+1;i<=n;++i)y[++p]=i;
        for(int i=1;i<=n;++i)if(sa[i]>k)y[++p]=sa[i]-k;
        for(int i=1;i<=m;++i)c[i]=0;
        for(int i=1;i<=n;++i)++c[x[i]];
        for(int i=1;i<=m;++i)c[i]+=c[i-1];
        for(int i=n;i;--i)sa[c[x[y[i]]]--]=y[i];
        swap(x,y);
        p=1;
        x[sa[1]]=1;
        for(int i=2;i<=n;++i)
            x[sa[i]]=y[sa[i]]==y[sa[i-1]]&&y[sa[i]+k]==y[sa[i-1]+k]?p:++p;
        if(p==n)break;
        m=p;
    }
    for(int i=1;i<=n;++i)
        sa[rank[i]]=i;
    for(int i=1,p=0;i<=n;++i){
        if(rank[i]==1)continue;
        if(p)--p;
        int j=sa[rank[i]-1];
        while(s[i+p]==s[j+p])++p;
        height[rank[i]]=i;
    }
}
inline void preLog(){
    for(int i=0;(1<<i)<=n;++i)
        Log[1<<i]=i;
    for(int i=0,res=0;i<=n;++i){
        if(Log[i])res=Log[i];
        else Log[i]=res;
    }
}
inline void findString(long long x,int &tl,int &tr){
    int l=1,r=n,mid;
    while(l<r){
        mid=(l+r+1)>>1;
        if(sum[mid]>=x)r=mid-1;
        else l=mid;
    }
    tl=sa[l];
    tr=sa[l]+x-sum[l]+height[l];
}
int QAQ(){
    long long tot=0;
    int q;
    long long x,y;
    int xl,xr,yl,yr;
    n=read(),q=read();
    m='z';
    sum[1]=n-sa[1]+1;
    for(int i=2;i<=n;++i)
        sum[i]=n-sa[i]-height[i]+1;
    for(int i=1;i<=n;++i)
        sum[i]+=sum[i-1];
    while(m--){
        x=read(),y=read();
        findString(x,xl,xr);
        findString(y,yl,yr);
        printf("%d %d %d %d\n",xl,xr,yl,yr);
    }
    return false;
}
}
int main(){
    return orz::QAQ();
}

没写完2

#include<cstdio>
#include<algorithm>
using namespace std;
namespace orz{
const int N=100;
int n,m;
int c[N],x[N],y[N];
int Log[N];
inline void log2(int x){

}
inline void logPre(){
    
}
struct String{
    int f[N][20];
    inline void getMin(){
    
    }
    inline void SA(){
    
    }
}
inline int read(){
    int a=1,b=0;char t;
    do{t=getchar();if(t=='-')a=-1;}while(t>'9'||t<'0');
    do{b=b*10-'0'+t;t=getchar();}while(t>='0'&&t<='9');
    return a*b;
}
int QAQ(){
    return false;
}
}
int main(){
    return orz::QAQ();
}

2019年7月博客汇总下

标签：矩阵乘法   基本   坐标   usaco   org   his   move   解决   重载   

原文地址：https://www.cnblogs.com/oiertkj/p/12203542.html