「ROI 2019 Day1」运输 20/19

Posted on 2023-08-12 Edited on 2024-04-19 In CP 题解

「ROI 2019 Day1」运输 20/19

题目大意：

给定一个带权 \(DAG\) ， \(1\) 为起始点，给定小常数 \(p\)

每次查询一个点 \(u\) ，一个权值 \(r\) ，问是否存在一条路径 \(1\ldots u\) ，其长度 \(x\) 满足 \(r\leq x\leq \frac{p}{p-1}\cdot r\)

转换一下， \(dp\) 每个点是否存在 \(r\) ，那么对于路径的权值 \(x\) ，合法的 \(r\) 即为 \([\frac{p-1}{p}\cdot x,x]\)

对于任意两个区间，如果其相交，则可以合并，并且用两个区间中最小和最大的 \(x\) 来表示这个区间

而不相交的区间最多只有 \(\log_{\frac{p}{p-1}} w\) 段，大概 \(700\) 段

任意时刻，每个点的 \(dp\) 情况可以用 \(700\) 段不交的区间表示，转移可以归并数组进行

因此维护 \(dp\) 复杂度为 \(O(m\log_{\frac{p}{p-1}} w)\) 常数极小，单次查询复杂度为 \(O(\log \log_{\frac{p}{p-1}} w)\)

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
#define pb push_back
#define rep(i,a,b) for(int i=a,i##end=b;i<=i##end;++i)
char IO;
template <class T=int> T rd(){
	T s=0; int f=0;
	while(!isdigit(IO=getchar())) if(IO=='-') f=1;
	do s=(s<<1)+(s<<3)+(IO^'0');
	while(isdigit(IO=getchar()));
	return f?-s:s;
}

const int N=5e5+10;

int n,m,q,p;
struct Interval{
	ll l,r,x,y;
	Interval(ll a,ll b){ 
		x=a,y=b;
		l=(x*(p-1)+p-1)/p,r=y;
	}
	void Add(ll c){
		x+=c,y+=c;
		l=(x*(p-1)+p-1)/p,r=y;
	}
	bool operator & (const Interval &__) const { return min(r,__.r)>=max(l,__.l)-1; }
	Interval operator + (const Interval &__) const { return Interval(min(x,__.x),max(y,__.y)); }
	bool operator < (const ll &x) const {
		return r<x;
	}
};
vector <Interval> dp[N];

struct Edge{
	int to,nxt; ll w;
} e[N];
int head[N],ecnt;
void AddEdge(int u,int v,ll w){
	e[++ecnt]=(Edge){v,head[u],w};
	head[u]=ecnt;
}
void Trans(vector <Interval> &x,vector <Interval> y,ll c){
	for(auto &i:y) i.Add(c);
	vector <Interval> res;
	int p1=0,p2=0,s1=x.size(),s2=y.size();
	while(p1<s1 || p2<s2) {
		if(p1<s1 && (p2==s2 || x[p1].l<=y[p2].l)) {
			if(res.size() && *res.rbegin()&x[p1]) res[res.size()-1]=*res.rbegin()+x[p1++];
			else res.pb(x[p1++]);
		} else {
			if(res.size() && *res.rbegin()&y[p2]) res[res.size()-1]=*res.rbegin()+y[p2++];
			else res.pb(y[p2++]);
		}
	}
	x=res;
}

int main(){
	rep(kase,1,rd()) {
		n=rd(),m=rd(),q=rd(),p=rd();
		rep(i,1,n) head[i]=ecnt=0;
		rep(i,1,m){
			int u=rd(),v=rd(); ll w=rd<ll>();
			AddEdge(u,v,w);
		}
		rep(i,1,n) dp[i].clear();
		dp[1].pb(Interval(0,0));
		rep(u,1,n) {
			for(int i=head[u];i;i=e[i].nxt) {
				Trans(dp[e[i].to],dp[u],e[i].w);
			}
		}
		while(q--){
			int x=rd(); ll y=rd<ll>();
			auto p=lower_bound(dp[x].begin(),dp[x].end(),y);
			if(p!=dp[x].end() && p->l<=y) putchar('1');
			else putchar('0');
		}
		puts("");
	}
}

「ROI 2018 Day 2」无进位加法

Posted on 2023-08-12 Edited on 2024-04-19 In CP 题解

「ROI 2018 Day 2」无进位加法

题目大意：

给出二进制数 \(a_1,\ldots a_n\) ，对于 \(b_1\ldots b_n\)

满足 \(a_i\leq b_i\) ， \(\bigoplus b_i=\sum b_i\) ，其中 \(\bigoplus\) 为异或和

求 \(\sum b_i\) 最小值

设长度量级为 \(N=\sum len(a_i)\)

\(O(N^2-N^3)\) ，从高到低确定答案的每一个位

枚举当前位为0，下面的位为1，贪心确定是否存在方案

检查一个答案是否合法：

动态维护一个倒序的 \(a_i\) 集合，从高到低考虑每一个位置

1.如果当前位为0：

如果 \(a_i\) 中存在大于等于这一位的数，非法

2.如果当前位为1：

2-1.如果 \(a_i\) 中存在2个当前位为1的数，非法

2-2.如果 \(a_i\) 中存在恰好一个，则将这个1用于这个 \(a_i\) ，并将 \(a_i\) 去掉最高位后放回集合

2-3.不存在，用这个 \(1\) 删除最大的一个 \(a_i\)

实际看来，这个贪心本身效率并不高

\[\ \]

优化1：快速确定答案最高位的可能范围

令 \(B=\max\{ len(a_i)+i-1\}\)

则 \(len(Ans)\in[B,B+1]\)

上下界均可以由上面的贪心模拟得到

\[ \ \]

优化2：快速维护 \(a_i\) 倒序

显然在不断更改的过程中，当前的 \(a_i\) 一定是原先的某一个 \(a_i\) 的一段后缀

考虑将所有这样的后缀排序，为了方便，用每一个最高的1来表示一个合法的后缀

显然可以先按照后缀长度分类，同长度的后缀，按照后缀中下一个1出现的位置排序

也就是一个类似基数排序个过程，额外维护每一个后缀中下一个出现的 \(1\) 所对应的后缀即可

预处理复杂度为 \(O(N\log N)\)

同时，也可以用线段树快速维护插入/删除的排名，得到 \(B\) 的值，单次操作复杂度 \(O(\log N)\)

\[ \ \]

优化3

称满足 \(len(a_i)+i-1=B\) 的 \(i\) 为 \(\text{critical number}\)

令 \(p\) 为最小的 \(\text{critical number}\) ，也就是在贪心过程中第一个出现情况2-1./2-2.的位置

决策答案为 \(B\) 还是为 \(B+1\) ，也就是决策

是用 \(len(a_p)\) 这个位置删除 \(a_p\) 的最高位，还是用 \(len(a_p)+1\) 的位置删除 \(a_p\)

( \([1,p-1]\) 的部分一定会被删掉)

\(\text{intended solution}\) 采用暴力递归来完成确定每一位的这个操作

Function Solve(Limit) Limit为当前可以使用的最高位的1
	求得 B,p
	删除 a[1,p-1]
	删除 a[p]最高位
	if B<=Limit and Solve(p-1) then
		ans[len(a[p]),B]=1
		return True
	删除a[p]
	if B+1<=Limit and Solve(p) then 
		ans[len(a[p])+1,B+1]=1
		return True
	else return False
end

至于复杂度，官方题解给出为 \(O(N)\) 次递归和删除/加入操作，最终复杂度为 \(O(N\log N)\)

#include<bits/stdc++.h>
using namespace std;
#define pb push_back
#define rep(i,a,b) for(int i=a,i##end=b;i<=i##end;++i)
#define drep(i,a,b) for(int i=a,i##end=b;i>=i##end;--i)
template <class T> inline void cmin(T &a,const T &b){ ((a>b)&&(a=b)); }
template <class T> inline void cmax(T &a,const T &b){ ((a<b)&&(a=b)); }
char IO;
int rd(){
	int s=0;
	while(!isdigit(IO=getchar()));
	do s=(s<<1)+(s<<3)+(IO^'0');
	while(isdigit(IO=getchar()));
	return s;
}

typedef vector <int> V;
const int N=3e5+10,INF=1e9+10;

int n,m,I[N],L;
char s[N];
int fir[N],nxt[N],rk[N],len[N],id[N];
V A[N];

struct Affirmation_Of_My_Existence{
	int s[N<<2],t[N<<2];
	void Down(int p){
		rep(v,p<<1,p<<1|1) t[v]+=t[p],s[v]+=t[p];
		t[p]=0;
	}
	void Upd(int p,int l,int r,int ql,int qr,int x) {
		if(ql>qr) return;
		if(ql<=l && r<=qr) {
			s[p]+=x,t[p]+=x;
			return;
		}
		Down(p);
		int mid=(l+r)>>1;
		if(ql<=mid) Upd(p<<1,l,mid,ql,qr,x);
		if(qr>mid) Upd(p<<1|1,mid+1,r,ql,qr,x);
		s[p]=max(s[p<<1],s[p<<1|1]);
	}
	void Build(int p,int l,int r){ 
		s[p]=len[id[l]]-INF;
		if(l==r) return;
		int mid=(l+r)>>1;
		Build(p<<1,l,mid),Build(p<<1|1,mid+1,r);
	}
	void Add(int x,int k){
		x=rk[x];
		Upd(1,1,m,x,x,INF*k),Upd(1,1,m,x+1,m,k);
	}
	// Find the first critical position "p", and return all the bits in [1,p]
	void Get(int p,int l,int r,int x,V &R){
		if(s[p]<0) return;
		if(l==r) return R.pb(id[l]);
		Down(p);
		int mid=(l+r)>>1;
		Get(p<<1,l,mid,x,R);
		if(s[p<<1]!=x) Get(p<<1|1,mid+1,r,x,R);
	}
} T;

int Solve(int L){
	// L denotes the maxmium bit we can use
	int B=T.s[1];
	if(B<0) return 1;
	if(B>L) return 0;
	V R; T.Get(1,1,m,B,R);
	int p=*R.rbegin(),l=len[p];
	for(int i:R) T.Add(i,-1);
	
	// Try ans B , so we use bit [nxt,B] to delete the number [1,p-1] 
    // and the number a[p] will be set to a[p]-2^l
	if(nxt[p]) T.Add(nxt[p],1);
	if(Solve(l-1)) {
		rep(i,l,B) s[i]=1;
		return B+1;
	}

	// Try ans B+1 , so we use bit [nxt+1,B+1] to delete the [1,p]
	if(nxt[p]) T.Add(nxt[p],-1);
	if(B<L && Solve(l)) {
		rep(i,l+1,B+1) s[i]=1;
		return B+2;
	}
	for(int i:R) T.Add(i,1);
	return 0;
}

int main(){
	rep(i,1,n=rd()) {
		scanf("%s",s); int l=strlen(s);
		cmax(L,l);
		drep(j,l-1,0) if(s[j]=='1') {
			nxt[++m]=fir[i];
			A[len[m]=l-j-1].pb(m);
			fir[i]=m;
		}
	}
	rk[0]=1e9;
	int k=m;
	rep(i,0,L-1) {
		k-=A[i].size();
		sort(A[i].begin(),A[i].end(),[&](int x,int y){ return rk[nxt[x]]<rk[nxt[y]]; });
		for(int j:A[i]) id[rk[j]=++k]=j;
		k-=A[i].size();
	}
	T.Build(1,1,m);
	rep(i,1,n) T.Add(fir[i],1);
	memset(s,0,sizeof s);
	drep(i,Solve(INF)-1,0) putchar(s[i]^48);
}

「USACO 2020 US Open Platinum」Exercise

Posted on 2023-08-12 Edited on 2024-04-19 In CP 题解

「USACO 2020 US Open Platinum」Exercise

做法与模数是否是质数无关

问题可能比较复杂，需要多步分析

1.对于一个已知的排列

显然这样的置换会构成若干个环，设每个环长度为 \(a_i,i\in [1,m]\) ，显然答案就是 \(lcm(a_i)\)

2.对于已知的 \(a_i\) 序列(注意这里说的是有序的)，计算其方案数

考虑已经排列的个数为 \(i\) ，加入一个环大小为 \(j\)

为了避免重复，应当固定这个环的初始位置为1号点，其余位置按照原先顺序插入

则方案数可以分为两部分考虑:

2-1.环内排列，固定的环首不可排列，即 \((j-1)!\)

2-2.剩下的 \(j-1\) 一个点位置未知，从未固定的 \(i+j-1\) 个点中选择

即 \(C(i+j-1,j-1)\)

所以就是 \(C(i+j-1,j-1)(j-1)!=\frac{(i+j-1)!}{i!}\)

归纳一下，发现更形象的描述就是 \(\begin{aligned}\frac{n!}{\prod_{i=1}^{m} (\sum_{j=1}^i a_j)}\end{aligned}\)

也就是每次除掉转移时的大小，将 \(n!\) 分成若干段，这似乎有利于理解下面的dp优化

3.计算 \(lcm\) 之积

考虑对于每个质因数计算其出现的幂次，注意这个幂次是对于 \(\varphi\) 取模的

原先是求恰好包含 \(x^k\) 的方案数，得到的 \(dp\) 不好优化，考虑转换为:

求质因数 \(x\) 出现在答案里的幂次 \(\ge k\) 的方案数 \(F_k\) ，答案就是 \(x^{\sum F_k}\)

Solution 1

那么反向求解，令 \(dp_i\) 表示当前已经确定了 \(i\) 个点，没有出现 \(x^k\) 倍数大小的联通块

暴力转移，枚举 \(i\) 从所有 \(x^k\not|j\) 转移过来即可，单次求解复杂度为 \(O(n^2)\) ，不可行

优化1:

考虑分解系数 \(\frac{(i+j-1)!}{i!}\) ，累前缀和，对于 \(j\) 为 \(x^k\) 倍数的情况枚举减掉

这样单次求解复杂度为 \(O(\frac{n^2}{x^k})\) ，总复杂度为 \(O(n^2\ln n)\) ，且不好处理阶乘逆元

优化2:

不枚举 \(x^i\) 的倍数，直接再用一个前缀和数组 \(s_i\) 记录下来，让 \(s_i\) 从 \(s_{i-x^k}\) 转移过来即可

如何将系数 \(\frac{(i+j-1)!}{i!}\) 分解？

每次 \(i+j\) 增大1，就多乘上一个 \(i+j\) 即可

当 \(s_{i}\) 从 \(s_{i-x^k}\) 转移过来时，需要补上 \(\begin{aligned}\prod_{j=i-x^k+1}^{i-1}j\end{aligned}\)

也就是模拟了上面提到的把 \(n!\) 分段的过程

这样就去掉了阶乘逆元的求解

Tips:发现需要预处理 \(T_{i,j}=\prod_{k=i}^j k\) ，可以滚动一下会快一点，内存为 \(O(n)\)

\[ \ \]

不同的 \(x^i\) 上限为 \(O(n)\) 种，实际大概可能是 \(O(\pi(n)\log \log n=\frac{n\log \log n}{\log n})\) ?

因此复杂度为 \(O(n^2)\)

可以看到代码还是很简单的

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
#define Mod2(x) ((x<0)&&(x+=P2))
#define rep(i,a,b) for(int i=a,i##end=b;i<=i##end;++i)
#define drep(i,a,b) for(int i=a,i##end=b;i>=i##end;--i)
 
const int N=7510;
 
int n,P,P2;
int mk[N];
int s[N],T[N],dp[N];
 
ll qpow(ll x,ll k=P-2) {
    ll res=1;
    for(;k;k>>=1,x=x*x%P) if(k&1) res=res*x%P;
    return res;
}
 
 
int main(){
    scanf("%d%d",&n,&P),P2=P-1;
    int ans=1;
    rep(i,2,n) if(!mk[i]) {
        for(int j=i+i;j<=n;j+=i) mk[j]=1;
        for(int x=i;x<=n;x*=i) mk[x]=i;
    }
    rep(i,1,n) T[i]=1;
    int r=1;
    rep(i,1,n) r=1ll*r*i%P2;
    rep(x,2,n) {
        rep(j,1,n-x+1) T[j]=1ll*T[j]*(j+x-2)%P2;
        // 滚动求解区间乘积
        if(mk[x]<=1) continue;
        dp[0]=1,s[0]=1;
        int sum=1;
 
        rep(i,1,n) {
            s[i]=0;
            if(i>=x) s[i]=1ll*s[i-x]*T[i-x+1]%P2;
            dp[i]=sum-s[i]; Mod2(dp[i]);
            s[i]=(1ll*s[i]*i+dp[i])%P2;
            sum=(1ll*sum*i+dp[i])%P2;
        }
        ans=1ll*ans*qpow(mk[x],P2+r-dp[n])%P;
    }
    printf("%d\n",ans);
}

Solution 2

为了便于表达，设满足条件为至少出现一个 \(x\) 的倍数

实际用min-max容斥确实比较好理解，设对于集合 \(S\) ，求其最大值

\(\begin{aligned} \max \lbrace S\rbrace =\sum_{T\sube S} (-1)^{|T|+1}\min\lbrace T \rbrace\end{aligned}\)

简要证明的话:

把 \(S\) 中的元素倒序排成一排分别为 \(S_i\)

对于 \(S_1\) 即最大值，显然被计算一次

对于剩下的值 \(S_i(i>1)\) ，则它作为最小值产生贡献意味着选的数都在 \(1-i\) 内，显然有 \(2^{i-1}\) 次为奇数集合大小， \(2^{i-1}\) 为偶数集合大小，两部分抵消

\[ \ \]

要计算最大值为1的方案数，那么就要计算最小值为1的子集方案数

考虑强制一个子集中每一个环大小均为 \(x\) 的倍数，设选出了 \(i\) 个这样的环，总大小为 \(j\) 的方案数为 \(dp_{i,j}\)

则实际对答案的贡献还要考虑这样的子集出现的次数

考虑选择子集的位置，以及剩下的 \(n-j\) 个点任意排布，方案数应该为 \(C(n,j)\cdot (n-j)!=\frac{n!}{j!}\)

如果真的用 \(dp_{i,j}\) ，复杂度显然太高，考虑 \(i\) 这一维的影响只在于系数 \((-1)^{|T|+1}\) ，可以直接在转移过程中解决

因此可以直接记录大小 \(i\) ，从前面转移过来

(可以看到依然需要访问上面提到的 \(T_{i,j}\) ，要滚动的话还会更难处理)

这样的 \(dp\) 状态有 \(\frac{n}{x}\) 种，转移为 \((\frac{n}{x})^2\) ，最后统计复杂度为 \(O(\frac{n}{x})\)

实际上这样的复杂度已经足够了，是 \(O(\sum \frac{n^2}{i^2}\approx n^2)\)

以下是 \(O(n^2)\) 内存的代码

#include<bits/stdc++.h>
using namespace std;
#define rep(i,a,b) for(int i=a;i<=b;++i)
#define drep(i,a,b) for(int i=a;i>=b;--i)
typedef long long ll;

const int N=7510;

int n,P,P2,T[N][N],mk[N],dp[N];
ll qpow(ll x,ll k) {
	ll res=1;
	for(;k;k>>=1,x=x*x%P) if(k&1) res=res*x%P;
	return res;
}
int Calc(int x){
	int m=n/x,s=P2-1;
	rep(i,1,m) {
		s=1ll*s*T[(i-1)*x+1][i*x-1]%P2;
		dp[i]=P2-s;
		s=(1ll*s*i*x+dp[i])%P2;
	}
	int ans=0;
	rep(i,1,m) ans=(ans+1ll*dp[i]*T[i*x+1][n])%P2;
	return ans;
}

int main(){
	scanf("%d%d",&n,&P),P2=P-1;
	rep(i,2,n) if(!mk[i]) {
		for(int j=i;j<=n;j+=i) mk[j]=1;
		for(int j=i;j<=n;j*=i) mk[j]=i;
	}
	rep(i,1,n+1){
		T[i][i-1]=1;
		rep(j,i,n) T[i][j]=1ll*T[i][j-1]*j%P2;
	}
	int ans=1;
	rep(i,2,n) if(mk[i]>1) ans=ans*qpow(mk[i],Calc(i))%P;
	printf("%d\n",ans);
}

\[ \ \]

最终优化

可以看到，这个做法和Sol1的转移有十分的相同之处，因此考虑用同样的方法优化掉

但是预处理系数的部分依然是 \(O(n^2)\) 的，如何解决呢?

1.线段树大法

预处理复杂度为 \(O(n)\) ，查询复杂度为 \(O(\log n)\) ，总复杂度 \(O(n\log^2 n)\)

空间复杂度为 \(O(n)\)

Oh这个代码是ZKW线段树

#include<bits/stdc++.h>
using namespace std;
#define rep(i,a,b) for(int i=a;i<=b;++i)
typedef long long ll;

const int N=7510;

int n,P,P2;
struct Tree {
	int s[N<<2],bit;
	void Build() {
		for(bit=1;bit<=n;bit<<=1);
		rep(i,1,n) s[i+bit]=i;
		for(int i=bit;i>=1;--i) s[i]=1ll*s[i<<1]*s[i<<1|1]%P2;
	}
	int Que(int l,int r){
		if(l>r) return 1;
		if(l==r) return l;
		int res=1;
		for(l+=bit-1,r+=bit+1;l^r^1;l>>=1,r>>=1){
			if(~l&1) res=1ll*res*s[l^1]%P2;
			if(r&1) res=1ll*res*s[r^1]%P2;
		}
		return res;
	}
} T;

int mk[N];
ll qpow(ll x,ll k) {
	ll res=1;
	for(;k;k>>=1,x=x*x%P) if(k&1) res=res*x%P;
	return res;
}

int dp[N];
int Calc(int x){
	int m=n/x,s=P2-1;
	rep(i,1,m) {
		s=1ll*s*T.Que((i-1)*x+1,i*x-1)%P2;
		dp[i]=P2-s;
		s=(1ll*s*i*x+dp[i])%P2;
	}
	int ans=0;
	rep(i,1,m) ans=(ans+1ll*dp[i]*T.Que(i*x+1,n))%P2;
	return ans;
}

int main(){
	scanf("%d%d",&n,&P),P2=P-1;
	T.Build();
	rep(i,2,n) if(!mk[i]) {
		for(int j=i;j<=n;j+=i) mk[j]=1;
		for(int j=i;j<=n;j*=i) mk[j]=i;
	}
	int ans=1;
	rep(i,2,n) if(mk[i]>1) ans=ans*qpow(mk[i],Calc(i))%P;
	printf("%d\n",ans);
}

2.模数分解+CRT(Chinese Remainder Theory)

分解后，可以用模逆元处理，然后就直接做，最后CRT合并一下，其实我并不会实现。。。

3.猫树~~(嘿嘿)~~

这是一个 \(O(n\log n)\) 预处理， \(O(1)\) 查询的数据结构，空间复杂度为 \(O(n\log n)\)

因此复杂度为 \(O(n\log n)\)

#include<bits/stdc++.h>
using namespace std;
#define rep(i,a,b) for(int i=a;i<=b;++i)
#define drep(i,a,b) for(int i=a;i>=b;--i)
typedef long long ll;

const int N=7510;

int n,P,P2;
struct SuckCat{
	int s[14][8200],Log[8200],bit;
	void Build() {
		for(bit=1;bit<=n;bit<<=1);
		rep(i,2,bit) Log[i]=Log[i>>1]+1;
		for(int l=1,d=0;l*2<=bit;l<<=1,d++){
			for(int i=l;i<=bit;i+=l*2){
				s[d][i-1]=i-1,s[d][i]=i;
				drep(j,i-2,i-l) s[d][j]=1ll*s[d][j+1]*j%P2;
				rep(j,i+1,i+l-1) s[d][j]=1ll*s[d][j-1]*j%P2;
			}
		}
	}
	int Que(int l,int r){
		if(l>r) return 1;
		if(l==r) return l;
		int d=Log[l^r];
		return 1ll*s[d][l]*s[d][r]%P2;
	}
} T;

int mk[N];
ll qpow(ll x,ll k) {
	ll res=1;
	for(;k;k>>=1,x=x*x%P) if(k&1) res=res*x%P;
	return res;
}

int dp[N];
int Calc(int x){
	int m=n/x,s=P2-1;
	rep(i,1,m) {
		s=1ll*s*T.Que((i-1)*x+1,i*x-1)%P2;
		dp[i]=P2-s;
		s=(1ll*s*i*x+dp[i])%P2;
	}
	int ans=0;
	rep(i,1,m) ans=(ans+1ll*dp[i]*T.Que(i*x+1,n))%P2;
	return ans;
}

int main(){
	scanf("%d%d",&n,&P),P2=P-1;
	T.Build();
	rep(i,2,n) if(!mk[i]) {
		for(int j=i;j<=n;j+=i) mk[j]=1;
		for(int j=i;j<=n;j*=i) mk[j]=i;
	}
	int ans=1;
	rep(i,2,n) if(mk[i]>1) ans=ans*qpow(mk[i],Calc(i))%P;
	printf("%d\n",ans);
}

「ROI 2019 Day2」模式串查找 (口胡)

Posted on 2023-08-12 Edited on 2024-04-19 In CP 题解

「ROI 2019 Day2」模式串查找 (口胡)

设 \(S=\sum |w_i|\)

显然我们需要一个树形数据结构来维护题目中添加字符的操作

归纳一下，需要实现的操作就是：

1.添加一个新串

2.在当前串中分裂一段区间 \([L,R]\)

3.将一个串复制 \(k\) 次

将每一个单字符视为一个节点，考虑用一个可持久化的平衡树来维护上述操作，比如可持久化非旋 \(\text{Treap}\)

题解中给出的2-3-Tree不知道效果怎么样

2,3操作增加的数量为 \(\log\) 级别，最终的节点总数为 \(O(S+n\log S)\)

接下来需要实现的操作是合并两个子串的信息，显然在合并时计算跨过两个节点的串匹配模板串的次数

容易想到记录后缀匹配的 \(kmp\) 指针，但是这样的指针难以完成合并操作

为了计算这个，我们需要维护这个串的前缀/后缀在原串上最长的匹配子段

此时一个节点的信息可以用 \((l1,r1),(l2,r2),k\) 来表示，其中 \(k\) 标记了这个串是否整个出现在模板串中

合并子串 \(A,B\) 信息：

1.求出前缀/后缀匹配，以前缀为例：

1-1.如果 \(A\) 无法完全匹配，则匹配前缀为 \(A\) 的匹配前缀

1-2. \(A\) 能够完全匹配，设 \(A\) 在原串对应 \(L,R\) ， \(B\) 的最长匹配前缀在模板串匹配起始位置为 \(P\)

我们需要在 \([L,R]\) 之后借上 \(P\) 开始的一段后缀，而 \([L,R]\) 出现在模板串中的位置对应着后缀数组上一段 \(rank\) 区间

取模板串反向后缀数组，求出与 \(R\) 匹配长度超过 \(R-L+1\) 的 \(rank\) 区间 \([l,r]\)

则我们需要找到 \([l,r]\) 中 \(sa[i]+1\) 与 \(P\) 最长的 \(LCP\) ，显然是一个求临近 \(rank\) 的问题，可以在线主席树二分解决

由此得到最长前缀为 \(A\) 串再加上额外匹配得到的部分

\(O(\log m)\) 完成

\[ \ \]

2.求出跨过两个串的完美匹配个数

容易得到前串后缀对应的 \(\text{kmp}\) 指针，后串前缀对应反向的 \(\text{kmp}\) 指针

实际上就是求出这两个指针不断失配时相加恰好完全匹配的个数

注意到，实际上这个询问是完全可以离线的

可以通过在线得到的匹配情况，离线得到询问得到询问答案，再从叶子开始重新计算每个节点的答案

比较暴力的做法是：

建立两棵 \(\text{kmp}\) 树，在第一棵树上 \(\text{dfs}\) ，加入祖先对应的位置，然后再第二棵树上查询祖先匹配的个数

可以用 \(\text{dfs}\) 序+差分树状数组维护，复杂度为 \(O(\log m)\)

\[ \ \]

因此总体复杂度为 \(O((S+n\log S)\log m)\) ~~实际上常数非常大？~~

「ROI 2018 Day 1」量子隐形传态

Posted on 2023-08-12 Edited on 2024-04-19 In CP 题解

「ROI 2018 Day 1」量子隐形传态

题目大意：

在 \(N\times M\) 的网格上给定 \(K\) 个点 \(1\ldots K\) ，定义两点间的距离为 \(\displaystyle 2^{\max\{|x_i-x_j|,|y_i-y_j|\}}\)

\(N,M,K\leq 10^4\) ，求 \(1\) 到 \(k\) 的最短路，下文认为 \(N,M\) 同阶

如何存储距离

显然距离是一个不超过 \(10^4\) 位的二进制数，用 \(\text{bitset}\) 存下来

每一次转移需要维护一个位+1操作，比较大小操作，都可以 \(O(\frac{N}{w})\) 实现，其中 \(\text{w}\) 为压位数

\[ \ \]

Lemma:

对于点 \(A(x,y)\) ，将平面分为 \(8\) 个部分

注意对于 \(x'=x\) 或者 \(y'=y\) 的区域一定要分离

则有：在任意一个平面区域中，有效的转移点一定是距离 \((x,y)\) 最近的点

简要证明：

对于 \(A\) 来说，切比雪夫距离相同的的点构成一条带

设最近的点为 \(B\) ，那么对于任意一个其它点 \(C\) ，显然有 \(dis(A,C)>dis(B,C),dis(A,C)>dis(A,B)\)

故走 \(A\rightarrow B\rightarrow C\) 不劣

快速完成转移

这样的 \(B\) 显然不唯一存在，每次转移需要的是 \(\text{L}\) 形的段

故可以对于每行每列用线段树优化区间连边

故得到一个 \(O(K)\) 点数， \(O(K\log K)\) 边数的图

用 \(\text{Dijkstra}\) 完成最短路，复杂度为 \(O(K\log ^2 K\frac{N}{W})\)

ps: 这里没有考虑找到最近点的过程，下面的代码是直接暴力找的。。。

#include<bits/stdc++.h>
using namespace std;
typedef unsigned long long ull;
typedef pair <int,int> Pii;
#define reg register
#define pb push_back
#define mp make_pair
#define rep(i,a,b) for(int i=a,i##end=b;i<=i##end;++i)
#define drep(i,a,b) for(int i=a,i##end=b;i>=i##end;--i)
template <class T> inline void cmin(T &a,const T &b){ ((a>b)&&(a=b)); }
template <class T> inline void cmax(T &a,const T &b){ ((a<b)&&(a=b)); }

char IO;
int rd(){
	int s=0;
	while(!isdigit(IO=getchar()));
	do s=(s<<1)+(s<<3)+(IO^'0');
	while(isdigit(IO=getchar()));
	return s;
}

const int N=30011,INF=1e9+10,U=10000,D=6;
bool Mbe;

int n,m,k;
int X[N],Y[N],pre[N],vis[N];
vector <Pii> E[N];
typedef unsigned long long ull;
struct Bitset{
	ull a[N/64+10];
	int l;
	void Add(int x){
		while(1) {
			ull t=a[x>>D];
			a[x>>D]+=1ull<<(x&63);
			if(t>a[x>>D]) x=((x>>D)+1)<<D;
			else break;
		}
		cmax(l,((x>>D)<<D)+63-__builtin_clzll(a[x>>D]));
	}
	bool operator < (const Bitset &__) const {
		if(__.l>10009) return 1;
		if(l!=__.l) return l<__.l;
		drep(i,(l>>D)+1,0) if(a[i]!=__.a[i]) return a[i]<__.a[i];
		return 0;
	}
} dis[N];
struct Queue{
	int s[N<<2],bit;
	void Up(int p) { s[p]=dis[s[p<<1]]<dis[s[p<<1|1]]?s[p<<1]:s[p<<1|1]; }
	void Build(){
		for(bit=1;bit<=n+1;bit<<=1);
		s[bit+1]=1;
		for(int p=bit+1;p;p>>=1) s[p]=1;
	}
	void push(int p) { for(s[p+bit]=p,p+=bit;p>>=1;) Up(p); }
	int top(){
		int res=s[1],p=s[1];
		for(s[p+=bit]=0;p>>=1;) Up(p);
		return res;
	}
} que;

int Dis(int x,int y) { return max(abs(X[x]-X[y]),abs(Y[x]-Y[y])); }
int Ans[N],Ac;
int Min[N][8];
int Dir(int u,int v){
	int x=X[v]-X[u],y=Y[v]-Y[u];
	if(y==0) return x>0?0:4;
	if(y>0) return x==0?2:(x>0?1:3);
	return x==0?6:(x>0?7:5);
}

typedef vector <int> V;
V A[N],B[N];
int rtx[N],rty[N],ls[N],rs[N];
int Build(const V &vec,int l,int r){
	if(l==r) return vec[l];
	int mid=(l+r)>>1,u=++n;
	ls[u]=Build(vec,l,mid),rs[u]=Build(vec,mid+1,r);
	E[u].pb(mp(ls[u],-1)),E[u].pb(mp(rs[u],-1));
	return u;
}
V Res;
void Que(const V &vec,int p,int l,int r,int ql,int qr){
	if(!p) return;
	if(ql<=vec[l] && vec[r]<=qr) return Res.pb(p);
	int mid=(l+r)>>1;
	if(ql<=vec[mid]) Que(vec,ls[p],l,mid,ql,qr);
	if(qr>=vec[mid+1]) Que(vec,rs[p],mid+1,r,ql,qr);
}
void AddX(int u,int x,int l,int r){
	int d=abs(x-X[u]);
	if(rtx[x]) Que(A[x],rtx[x],0,A[x].size()-1,l,r);
	for(int v:Res) {
		E[u].pb(mp(v,d));
	}
	Res.clear();
}
void AddY(int u,int y,int l,int r){
	int d=abs(y-Y[u]);
	if(rty[y]) Que(B[y],rty[y],0,B[y].size()-1,l,r);
	for(int v:Res) {
		E[u].pb(mp(v,d));
	}
	Res.clear();
}

void Init(){
	rep(i,1,k) A[X[i]].pb(i),B[Y[i]].pb(i);
	rep(i,1,U) {
		if(A[i].size()) {
			sort(A[i].begin(),A[i].end(),[&](int x,int y){ return Y[x]<Y[y]; });
			rtx[i]=Build(A[i],0,A[i].size()-1);
			for(int &j:A[i]) j=Y[j];
		}
		if(B[i].size()) {
			sort(B[i].begin(),B[i].end(),[&](int x,int y){ return X[x]<X[y]; });
			rty[i]=Build(B[i],0,B[i].size()-1);
			for(int &j:B[i]) j=X[j];
		}
	}
	rep(i,1,k) rep(j,0,7) Min[i][j]=INF;
	rep(i,1,k) rep(j,i+1,k){
		int d=Dir(i,j),dis=Dis(i,j);
		cmin(Min[i][d],dis);
		cmin(Min[j][(d+4)&7],dis);
	}
	rep(i,1,k) {
		if(Min[i][0]!=INF) AddX(i,X[i]+Min[i][0],Y[i],Y[i]);
		if(Min[i][1]!=INF) {
			AddX(i,X[i]+Min[i][1],Y[i]+1,Y[i]+Min[i][1]);
			AddY(i,Y[i]+Min[i][1],X[i]+1,X[i]+Min[i][1]);
		}
		if(Min[i][2]!=INF) AddY(i,Y[i]+Min[i][2],X[i],X[i]);
		if(Min[i][3]!=INF) {
			AddX(i,X[i]-Min[i][3],Y[i]+1,Y[i]+Min[i][3]);
			AddY(i,Y[i]+Min[i][3],X[i]-Min[i][3],X[i]-1);
		}
		if(Min[i][4]!=INF) AddX(i,X[i]-Min[i][4],Y[i],Y[i]);
		if(Min[i][5]!=INF) {
			AddX(i,X[i]-Min[i][5],Y[i]-Min[i][5],Y[i]-1);
			AddY(i,Y[i]-Min[i][5],X[i]-Min[i][5],X[i]-1);
		}
		if(Min[i][6]!=INF) AddY(i,Y[i]-Min[i][6],X[i],X[i]);
		if(Min[i][7]!=INF) {
			AddX(i,X[i]+Min[i][7],Y[i]-Min[i][7],Y[i]-1);
			AddY(i,Y[i]-Min[i][7],X[i]+1,X[i]+Min[i][7]);
		}
	}
}

bool Med;
int main(){
	fprintf(stderr,"%.2lf\n",(&Med-&Mbe)/1024.0/1024.0);
	n=rd(),m=rd(),k=rd(),n=k;
	rep(i,1,k) X[i]=rd(),Y[i]=rd();
	Init();
	que.Build();
	dis[0].Add(10111);
	rep(i,2,n) dis[i].Add(10110);
	while(que.s[1]) {
		int u=que.top();
		vis[u]=1;
		for(auto t:E[u]) {
			int v=t.first;
			Bitset w=dis[u]; if(~t.second) w.Add(t.second);
			if(w<dis[v]) dis[v]=w,que.push(v),pre[v]=u;
		}
	}

	for(int u=k;u;u=pre[u]) if(u<=k) Ans[++Ac]=u;
	printf("%d\n",Ac);
	drep(i,Ac,1) printf("%d ",Ans[i]);
}

「USACO 2020.12 Platinum」Spaceship

Posted on 2023-08-12 Edited on 2024-04-19 In CP 题解

「USACO 2020.12 Platinum」Spaceship

看到题目第一个想到的是：按照路径长度可以确定按钮次数和路径次数

然而路径长度是 \(2^k\) 级别的。。

下文认为 \(n,q,k\) 同阶

既然无法考虑长度，那么就直接在 \(dp\) 时将路径作为状态压入

令 \(dp_{i,s,t}\) 表示前面 \(i\) 个按钮未按过，从 \(s\) 走到 \(t\) 的路径数

显然 \(dp_{i}\) 可以看做一个类似矩阵的转移，设邻接矩阵为 \(E\)

那么得到 \(dp_{i}=dp_{i-1}+dp_{i-1}\cdot E\cdot dp_{i-1}\)

( \(dp_i\) 不一定按了 \(i\) 这个按钮，所以考虑按或者不按)

那么对于 \(b_s,b_t\) 的情况，考虑把操作序列分成两段

显然操作序列有唯一一个按过恰好一次的最大的按钮 \(max\) ，可以在这里将操作序列分成两段

处理出 \(b_s\) 到 \(max\) 和 \(max\) 到 \(b_t\) 的两部分方案数合并即可

以处理 \(b_s\rightarrow max\) 为例

类似上面的操作，令 \(F_{i,j}\) 表示最大按钮为 \(i\) ，且按过 \(i\) 之后没有按过其他按钮，停留在 \(j\) 的方案数 ( \(F_i\) 是一个向量)

初始显然有 \(F_{b_s,s}=1\)

\(\displaystyle i>b_s,F_{i}=\sum_{j<i} F_{j}\cdot dp_{j-1}\cdot E\) (按过 \(j\) 之后可以继续按 \([1,j-1]\) )

前缀和即可，复杂度为 \(O(n^3)\)

同理处理出 \(max\rightarrow b_t\) 的权值合并即可

#include<bits/stdc++.h>
using namespace std;
#define Mod1(x) ((x>=P)&&(x-=P))
#define rep(i,a,b) for(int i=a,i##end=b;i<=i##end;++i)
#define drep(i,a,b) for(int i=a,i##end=b;i>=i##end;--i)

const int N=62,P=1e9+7,INF=1e9+10;

int n,m,q;
struct Mat{
	int a[N][N];
	void clear(){ memset(a,0,sizeof a); }
	Mat operator * (const Mat &x) const {
		Mat res; res.clear();
		rep(i,1,n) rep(j,1,n) if(a[i][j]) rep(k,1,n) res.a[i][k]=(res.a[i][k]+1ll*a[i][j]*x.a[j][k])%P;
		return res;
	}
	Mat operator + (const Mat &x) const {
		Mat res;
		rep(i,1,n) rep(j,1,n) res.a[i][j]=a[i][j]+x.a[i][j],Mod1(res.a[i][j]);
		return res;
	}
} I,E,dp[N];
int F[N][N],G[N][N];
int S[N][N];

int main(){
	freopen("spaceship.in","r",stdin),freopen("spaceship.out","w",stdout);
	scanf("%d%d%d",&n,&m,&q);
	rep(i,1,n) I.a[i][i]=1;
	rep(i,1,n) rep(j,1,n) scanf("%1d",&E.a[i][j]);
	dp[0]=I;
	rep(i,1,m) dp[i]=dp[i-1]+dp[i-1]*E*dp[i-1];
	while(q--) {
		int bs,s,bt,t; scanf("%d%d%d%d",&bs,&s,&bt,&t);
		memset(F,0,sizeof F),memset(G,0,sizeof G);
		memset(S,0,sizeof S);
		F[bs][s]=1;
		rep(i,1,n) S[bs][i]=dp[bs-1].a[s][i];
		rep(i,bs+1,m) {
			rep(a,1,n) if(S[i-1][a]) rep(b,1,n) F[i][b]=(F[i][b]+1ll*S[i-1][a]*E.a[a][b])%P;
			rep(a,1,n) S[i][a]=S[i-1][a];
			rep(a,1,n) if(F[i][a]) rep(b,1,n) S[i][b]=(S[i][b]+1ll*F[i][a]*dp[i-1].a[a][b])%P;
		}

		memset(S,0,sizeof S);
		G[bt][t]=1;
		rep(i,1,n) S[bt][i]=dp[bt-1].a[i][t];
		rep(i,bt+1,m) {
			rep(a,1,n) rep(b,1,n) G[i][a]=(G[i][a]+1ll*S[i-1][b]*E.a[a][b])%P;
			rep(a,1,n) S[i][a]=S[i-1][a];
			rep(a,1,n) rep(b,1,n) S[i][a]=(S[i][a]+1ll*G[i][b]*dp[i-1].a[a][b])%P;
		}

		int ans=0;
		rep(i,1,m) rep(j,1,n) ans=(ans+1ll*F[i][j]*G[i][j])%P;
		printf("%d\n",ans);
	}
}

「USACO 2020.12 Platinum」Cowmistry

Posted on 2023-08-12 Edited on 2024-04-19 In CP 题解

「USACO 2020.12 Platinum」Cowmistry

看到样例解释突然就有了思路.jpg

令 \(m=\min\{2^t|2^t>k\}\) ，也就是 \(k\) 最高位+1

对于数轴，按照 \([im,(i+1)m)\) 分组，显然跨过分组的数之间异或 \(\ge m>k\) ，不合法，扔掉

对于每组，直接看做是 \([0,m-1]\) ，此时令 \(d=\frac{m}{2}\)

将 \([0,m-1]\) 分为 \([0,d-1],[d,m-1]\) ，显然两段之内的数相互异或的结果 \(<d\) ，一定合法

也就是说，长度为 \(d\) 的段内随意选3个都合法

下面考虑跨过 \(d\) 的贡献，显然是一边选一个，一边选两个，此时这些数之间的异或和 \(\ge d\)

以左边选一个为例，令 \(k'=k-d\)

\(O(k\log k)\)

考虑异或和中一定有 \(d\) 这一位，下面只需要对于 \(i\in[0,d-1]\) 暴力统计出 \([d,m]\) 中有多少个数 \(j\) 满足 \(i\oplus j\leq k'\)

可以前缀和之后 \(\log k\) 做到，从高到低依次考虑 \(k\) 每一个为1的二进制位 \(i\) ，如果查询的数这一位和 \(x\) 相同，那么下面可以任意选择

否则将 \(x\rightarrow x\oplus 2^i\)

实现如下

int Que(int x,int *A,int k){
	int ans=0;
	drep(i,19,0) if(k&(1<<i)) {
		int l=x&(~((1<<i)-1)),r=l+(1<<i)-1;
		ans+=A[r]-A[l-1];
		x^=1<<i;
	}
	return ans+A[x]-A[x-1];
}

得到个数 \(c\) 之后，接下来就是为 \(x\) 在 \(c\) 里选择两个匹配，就是 \(\sum \frac{c(c-1)}{2}\)

由此得到一个 \(O(k\log k)\) 做法，可以通过前三档数据

namespace part_klogk{
	const int N=1<<20|10;
	int m,A[N],B[N],ans;
	int Que(int x,int *A,int k){
		int ans=0;
		drep(i,19,0) if(k&(1<<i)) {
			int l=x&(~((1<<i)-1)),r=l+(1<<i)-1;
			ans+=A[r]-A[l-1];
			x^=1<<i;
		}
		return ans+A[x]-A[x-1];
	}
	void Solve(int l,int r) {
		int mid=(l+r)>>1;
		int c1=0,c2=0;
		rep(i,l,mid) c1+=A[i];
		rep(i,mid+1,r) c2+=A[i];
		ans=(ans+D3(c1))%P,ans=(ans+D3(c2))%P;
		rep(i,l,mid) if(A[i]) {
			int c=Que(i-l,B+mid+2,k-(m/2));
			ans=(ans+1ll*c*(c-1)/2)%P;
		}
		rep(i,mid+1,r) if(A[i]) {
			int c=Que(i-mid-1,B+l+1,k-(m/2));
			ans=(ans+1ll*c*(c-1)/2)%P;
		}
	}
	void Solve() {
		for(m=1;m<=k;) m<<=1;
		rep(i,1,n) rep(j,L[i],R[i]) A[j]++;
		rep(i,1,N-1) B[i]=B[i-1]+A[i-1];
		rep(i,0,R[n]/m) Solve(i*m,(i+1)*m-1);
		printf("%d\n",ans);
	}
}

\[\ \]

\(O(n\log k-n\log ^2 k)\)

考虑对于 \(k'\) ，问题降阶为求区间 \(a\) 中的每一个数，在区间 \(b\) 中查询合法的个数 \(cnt_a\)

其中 \(a,b\) 区间可以认为对应相同的 \([0,L-1]\) ，但是出现元素不同

此时，继续采用上面的方法进行分组，分组对象变为两组区间

令 \(m'=\min\{2^t|2^t>k'\},d'=\frac{m}{2}\) ，分组之后

对于 \([0,d'-1],[d',m-1]\) ，显然两段之内异或 \(\leq d\) ，一定合法，加入答案 \(cnt_a\) 中

对于跨过区间的贡献，用下面的方法处理

左边区间对于右边区间继续递归进行查询，将得到的结果加上左边区间中数的个数

右边区间继续对于左边区间递归进行查询，将得到的结果加上右边区间中数的个数

问题不断降阶为子问题，会有 \(\log k\) 层子问题

每层子问题所有的区间可以分为 \(O(n)\) 段特殊的段

其他的部分要么完全覆盖要么为空，这些部分可以快速求出

发现答案如果用 \((num,cnt)\) 表示当前查询结果为 \(num\) 的个数为 \(cnt\)

每层可以用 \(O(n)\) 个不同的二元对表示结果

如此求得后，再自底向上合并得到答案

\[ \ \]

关于实现

如果真的按照上面的方法，一层层向下分裂区间，会面临着常数大，难以实现的问题

考虑将每个区间 \([L_i,R_i]\) 插入线段树 \([0,2^{30}-1]\)

显然得到 \(O(n\log k)\) 个节点，在底层完全覆盖的节点打上标记

先递归进行第一层的分裂区间操作，对于打上标记的节点可以分成 \(\frac{r-l+1}{m}\) 个完全覆盖区间

这些完全覆盖区间的答案可以 \(O(1)\) 求出

\[ \ \]

分裂完成后，每次查询的区间对用 \((a,b,L,k,f1,f2)\) 表示

其中 \(a\) 表示查询区间对应节点， \(b\) 表示被查询区间对应节点

\(L\) 表示区间长度， \(f1,f2\) 表示 \(a,b\) 是否继承到上层的完全覆盖标记

如此每次合并 \((ls_a,rs_b),(rs_a,ls_b)\) 的查询结果，加上另一边的贡献即可

当 \(b\) 为空或者为完全覆盖时，答案可以 \(O(1)\) 得到

当 \(a\) 为空时，可以得到空答案

从这样的底层向上合并，每个元素被合并次数 \(O(\log k)\)

\[ \ \]

没有仔细分析复杂度，应该就是 \(O(n\log k-n\log ^2 k)\)

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
typedef pair <int,int> Pii;
#define mp make_pair
#define pb push_back
#define rep(i,a,b) for(int i=a,i##end=b;i<=i##end;++i)
#define drep(i,a,b) for(int i=a,i##end=b;i>=i##end;--i)
char IO;
int rd(){
    int s=0;
    while(!isdigit(IO=getchar()));
    do s=(s<<1)+(s<<3)+(IO^'0');
    while(isdigit(IO=getchar()));
    return s;
}

const int N=2e4+10,P=1e9+7,I2=(P+1)/2,I3=(P+1)/3,U=1<<30;

int n,m,k;
int L[N],R[N];
int Check(){
    for(int i=k;i;i>>=1) if(~i&1) return 0;
    return 1;
}
int D3(int n){ return 1ll*n*(n-1)/2%P*(n-2)%P*I3%P; }
int D2(int n){ return 1ll*n*(n-1)/2%P; }

const int M=N*60;
int s[M],ls[M],rs[M],t[M],cnt,rt;
//在线段树上插入节点，打上标记，同时处理出size便于下面计算
void Ins(int &p,int l,int r,int ql,int qr){
    if(!p) p=++cnt;
    s[p]+=min(qr,r)-max(ql,l)+1;
    if(ql<=l && r<=qr) { t[p]=1; return; }
    int mid=(l+r)>>1;
    if(ql<=mid) Ins(ls[p],l,mid,ql,qr);
    if(qr>mid) Ins(rs[p],mid+1,r,ql,qr);
}

int ans;
// O(1)求出m
int Up(int k){ return 1<<(32-__builtin_clz(k)); }
typedef vector <Pii> V;
V Union(const V &A,const V &B){
    int p1=0,p2=0,s1=A.size(),s2=B.size();
    V R;
    // 这里是否归并并不影响复杂度
    while(p1<s1 || p2<s2) {
        if(p1<s1 && (p2==s2||A[p1]<B[p2])) {
            if(R.size() && R.back().first==A[p1].first) R.back().second+=A[p1].second;
            else R.pb(A[p1]);
            p1++;
        } else {
            if(R.size() && R.back().first==B[p2].first) R.back().second+=B[p2].second;
            else R.pb(B[p2]);
            p2++;
        }
    }
    return R;
}
V Calc(int a,int b,int L,int k,int f1,int f2){
    f1|=t[a],f2|=t[b];
    V Res;
    // 底层情况O(1)解决
    if(!f1 && !a) return Res;
    if(f2) return Res.pb(mp(k+1,f1?L:s[a])),Res;
    if(!b) return Res.pb(mp(0,f1?L:s[a])),Res;
    int m=Up(k);
    // L>m说明还要继续进行分裂
    if(L>m) return Union(Calc(ls[a],ls[b],L/2,k,f1,f2),Calc(rs[a],rs[b],L/2,k,f1,f2));
    // 进入降阶子问题，左查右，右查左
    k-=m/2;
    V A=Calc(ls[a],rs[b],L/2,k,f1,f2);
    for(auto &i:A) i.first+=f2?L/2:s[ls[b]];
    V B=Calc(rs[a],ls[b],L/2,k,f1,f2);
    for(auto &i:B) i.first+=f2?L/2:s[rs[b]];
    return Union(A,B);
}

int cm;
// 第一次分裂
void Solve(int p,int L){
    if(!p) return;
    // 完全覆盖的部分答案可以快速求出
    if(t[p]) { cm+=L/m; return; }
    // 已经完成分裂，此时进入第一层子问题
    if(L==m) {
        // 贡献分为 
        // 左选3 ， 右选 3
        ans=(ans+D3(s[ls[p]]))%P,ans=(ans+D3(s[rs[p]]))%P;
        // 左1，右2
        V A=Calc(ls[p],rs[p],m/2,k-m/2,0,0);
        // 左2，右1
        V B=Calc(rs[p],ls[p],m/2,k-m/2,0,0);
        for(auto i:A) ans=(ans+1ll*D2(i.first)*i.second)%P;
        for(auto i:B) ans=(ans+1ll*D2(i.first)*i.second)%P;
        return;
    }
    Solve(ls[p],L/2),Solve(rs[p],L/2);
}
int main(){
    n=rd(),k=rd();
    rep(i,1,n) { int l=rd(),r=rd(); Ins(rt,0,U-1,l,r); }
    m=Up(k),Solve(rt,U);
    // O(1)求得完全覆盖部分的答案
    int t=(2ll*D3(m/2)+1ll*D2(k-m/2+1)*m)%P;
    ans=(ans+1ll*cm*t)%P;
    printf("%d\n",ans);
}

「USACO 2021 US Open Platinum」Balanced Subsets

Posted on 2023-08-12 Edited on 2024-04-19 In CP 题解

「USACO 2021 US Open Platinum」Balanced Subsets

考虑题目给出的定义对应怎样的图形，显然是一个凸的封闭图形

不妨通过左右边线描述，从上到下

1.左边线先左移再右移

2.右边线先右移再左移

不妨直接令 \(dp_{i,l,r,f1,f2}\) 表示当前第 \(i\) 行，当前左右边线为 \(l,r(l\leq r,\forall j\in[l,r],a_{i,j}=G)\)

\(f1,f2\) 表示当前左右边线处于左移还是右移状态

以左边线为例，定义右移开始的时刻为第一个 \(l>l'\) 的时刻

容易得到转移，是一个前/后缀和的形式

那么对于 \([l,r]\) 两维分别做前缀和，然后 \(O(n^3)\) 转移即可

注意转移过程中要确保 \([l',r'],[l,r]\) 有交

~~以下是暴力二维前缀和+手艹9种转移的代码~~ QQ图片20210506190147.jpg

const int N=170,P=1e9+7;

int n;
int c[N];
char s[N];
int dp[N][N][2][2];
int F[N][N][2][2];
int S(int a,int b,int x1,int x2,int y1,int y2){
	x1--,y1--;
	return (0ll+F[x2][y2][a][b]-F[x1][y2][a][b]-F[x2][y1][a][b]+F[x1][y1][a][b])%P;
}

int main(){
	n=rd();
	int ans=0;
	rep(i,1,n) {
		scanf("%s",s+1);
		rep(j,1,n) c[j]=c[j-1]+(s[j]=='G');
		memset(F,0,sizeof F);
		rep(i,1,n) rep(j,1,n) rep(a,0,1) rep(b,0,1) {
			F[i][j][a][b]=(0ll+F[i-1][j][a][b]+F[i][j-1][a][b]-F[i-1][j-1][a][b]+dp[i][j][a][b])%P;
		}
		memset(dp,0,sizeof dp);
		rep(l,1,n) rep(r,l,n) if(c[r]-c[l-1]==r-l+1) {
			(dp[l][r][0][0]+=S(0,0,l,r,l,r))%=P;
			(dp[l][r][1][0]+=S(0,0,1,l-1,l,r))%=P;
			(dp[l][r][0][1]+=S(0,0,l,r,r+1,n))%=P;
			(dp[l][r][1][1]+=S(0,0,1,l-1,r+1,n))%=P;
	
			(dp[l][r][0][1]+=S(0,1,l,r,r,n))%=P;
			(dp[l][r][1][1]+=S(0,1,1,l-1,r,n))%=P;
			
			(dp[l][r][1][0]+=S(1,0,1,l,l,r))%=P;
			(dp[l][r][1][1]+=S(1,0,1,l,r+1,n))%=P;
			
			(dp[l][r][1][1]+=S(1,1,1,l,r,n))%=P;
		}
		rep(l,1,n) rep(r,l,n) if(c[r]-c[l-1]==r-l+1) dp[l][r][0][0]++;
		rep(l,1,n) rep(r,l,n) rep(a,0,1) rep(b,0,1) if(dp[l][r][a][b]) 
			ans=(ans+dp[l][r][a][b])%P;
	}
	Mod2(ans),printf("%d\n",ans);
}

「USACO 2020.12 Platinum」Sleeping Cows

Posted on 2023-08-12 Edited on 2024-04-19 In CP 题解

「USACO 2020.12 Platinum」Sleeping Cows

写容斥就输了。。

为每个牛棚考虑牛，从大到小，考虑每一个牛棚是否匹配

\(dp_{i,j,f}\) 表示后 \(i\) 个牛棚中有 \(j\) 个钦定要匹配但是还未匹配的牛棚， \(f=0/1\) 表示是否存在一个牛棚未选

每次移动 \(i\) ，会有一部分牛不能继续匹配

如果 \(f=1\) ，那么这部分牛必须被全部后面的 \(j\) 个牛棚中某一些匹配，否则不合法

如果 \(f=0\) ，这一部分可以随便匹配一定数量的牛

用组合数维护转移权值即可



const int N=3010,P=1e9+7;

int n;
int dp[N][N][2],I[N],J[N],C[N][N];
int A[N],B[N];
void Add(int &u,int x){
	u+=x,Mod1(u);
}
ll qpow(ll x,ll k=P-2) {
	ll res=1;
	for(;k;k>>=1,x=x*x%P) if(k&1) res=res*x%P;
	return res;
}

int main(){
	freopen("cow.in","r",stdin),freopen("cow.out","w",stdout);
	n=rd();
	rep(i,1,n) B[i]=rd();
	rep(i,1,n) A[i]=rd();
	sort(A+1,A+n+1),sort(B+1,B+n+1);
	int p=n;
	dp[n+1][0][0]=1;
	rep(i,0,n) rep(j,C[i][0]=1,i) C[i][j]=(C[i-1][j]+C[i-1][j-1])%P;
	rep(i,J[0]=1,n) J[i]=1ll*J[i-1]*i%P;
	I[n]=qpow(J[n]);
	drep(i,n,1) I[i-1]=1ll*I[i]*i%P;
	
	drep(i,n,0) {
		int c=0;
		while(p && B[p]>A[i]) c++,p--;
		rep(j,0,n-i) {
			rep(k,0,min(j,c)) {
				int t=1ll*dp[i+1][j][0]*C[c][k]%P*C[j][k]%P*J[k]%P;
				Add(dp[i][j-k][1],t),Add(dp[i][j-k+1][0],t);
				if(k==c) {
					t=1ll*dp[i+1][j][1]*C[j][k]%P*J[k]%P;
					Add(dp[i][j-k][1],t),Add(dp[i][j-k+1][1],t);
				}
			}
		}
	}
	int ans=(dp[0][0][0]+dp[0][0][1])%P;
	printf("%d\n",ans);
}

「USACO 2021 US Open Platinum」United Cows of Farmer John

Posted on 2023-08-12 Edited on 2024-04-19 In CP 题解

「USACO 2021 US Open Platinum」United Cows of Farmer John

考虑依次枚举右端点 \(i\) ，计算左边合法的方案数，设一个数 \(x\) 上次出现的位置为 \(lst_x\)

则 \(i\) 能够作为右端点的区间就是 \([lst_{a_i}+1,i-2]\)

考虑什么样的位置可以作为左端点，显然这个点在 \([1,i]\) 中是最后一次出现

我们将不妨这样的点权值设为 \(w_i=1\)

考虑一个点作为中间点贡献怎样的区间，同样的，这个点在 \([1,i]\) 中是最后一次出现

并且，能够贡献的区间 \(>\) 上一次出现的位置 \(lst_x\)

这个中间点能够匹配的左端点个数就是 \(\displaystyle \sum_{k=lst_{a_j}+1}^{j-1} w_k\)

现在我们要用数据结构动态修改某一个位置的 \(w_i\) ，增减 \([lst_{a_j}+1,j-1]\) 的区间，查询 \([lst_{a_i}+1,i-2]\)

不妨再为一个点增加点权 \(t_i\) ，此时我们要维护的操作

1.单点修改 \(w_i\)

2.区间修改 \(t_i\)

3.求 \(w_it_i\) 区间和

在线段树上每个节点维护 \(w_i\) 之和， \(w_it_i\) 之和，可以标记永久化 \(t_i\)

具体实现参考代码~~（实际写得很丑）~~

const int N=2e5+10,INF=1e9+10;


int n;
int lst[N],lst2[N],cnt;
ll s1[N<<2],s2[N<<2];
int t[N<<2];
// s1表示w之和，s2表示区间内部t[i]*w[i]之和，t[i]现在是永久化的标记
void Up(int p){
	s2[p]=s2[p<<1]+s2[p<<1|1];
	s1[p]=s1[p<<1]+s1[p<<1|1]+s2[p]*t[p];
}
void Upd(int p,int l,int r,int x){
	if(l==r) {
		s2[p]^=1,s1[p]=t[p]*s2[p];
		return;
	}
	int mid=(l+r)>>1;
	x<=mid?Upd(p<<1,l,mid,x):Upd(p<<1|1,mid+1,r,x);
	Up(p);
}

void Upd(int p,int l,int r,int ql,int qr,int x){
	if(ql>qr) return;
	if(ql<=l && r<=qr) {
		t[p]+=x,s1[p]+=x*s2[p];
		return;
	}
	int mid=(l+r)>>1;
	if(ql<=mid) Upd(p<<1,l,mid,ql,qr,x);
	if(qr>mid) Upd(p<<1|1,mid+1,r,ql,qr,x);
	Up(p);
}

struct Node{
	ll x,y;
	Node(ll x=0,ll y=0):x(x),y(y){  }
	Node operator + (const Node __) { return Node(x+__.x,y+__.y); }
};
Node Que(int p,int l,int r,int ql,int qr){
	if(ql>qr) return Node();
	if(ql<=l && r<=qr) return Node(s1[p],s2[p]);
	int mid=(l+r)>>1; Node res;
	if(ql<=mid) res=res+Que(p<<1,l,mid,ql,qr);
	if(qr>mid) res=res+Que(p<<1|1,mid+1,r,ql,qr);
	res.x+=res.y*t[p];
	return res;
}

int main(){
	n=rd();
	ll ans=0;
	rep(i,1,n) {
		int x=rd();
		if(lst[x]) {
			Upd(1,1,n,lst[x]),cnt--;
			Upd(1,1,n,lst2[x]+1,lst[x],-1);
		}
		Node t=Que(1,1,n,lst[x]+1,i-2);
		ans+=t.x;
		Upd(1,1,n,i),cnt++,Upd(1,1,n,lst[x]+1,i-1,1);
		lst2[x]=lst[x],lst[x]=i;
	}
	printf("%lld\n",ans);
}

「ROI 2019 Day1」运输 20/19

「ROI 2018 Day 2」无进位加法

\(O(N^2-N^3)\) ， 从高到低确定答案的每一个位

优化1：快速确定答案最高位的可能范围

优化2：快速维护 \(a_i\) 倒序

优化3

「USACO 2020 US Open Platinum」Exercise

Solution 1

Solution 2

最终优化

「ROI 2019 Day2」模式串查找 (口胡)

「ROI 2018 Day 1」量子隐形传态

如何存储距离

Lemma:

快速完成转移

「USACO 2020.12 Platinum」Spaceship

「USACO 2020.12 Platinum」Cowmistry

\(O(k\log k)\)

\(O(n\log k-n\log ^2 k)\)

关于实现

「USACO 2021 US Open Platinum」Balanced Subsets

「USACO 2020.12 Platinum」Sleeping Cows

「USACO 2021 US Open Platinum」United Cows of Farmer John

\(O(N^2-N^3)\) ，从高到低确定答案的每一个位