任意值域最长公共子序列问题

先上例题

这道题的值域是(26)，但是我们将要讲的方法是不限制值域的

暴力

先考虑一发朴素的(O(nm))算法

(f[i][j]=max(f[i-1][j],f[i][j-1]))

(f[i][j]=max(f[i][j],f[i-1][j-1]+1)(a[j]==b[i]))

注意这里定义是(a[j]==b[i])

我们发现(f[i][j])相较于(f[i-1][j-1])最多大(1)

那么我们也许可以构造一个(01)矩阵(M)，使得任意(f[i][j]=sumlimits_{k=1}^{n}M[i][k])

构造M矩阵

假设(a)串为：(abbabab)，(b)串为:(aabab)

先将(a)串翻转：(bababba)

我们把每个(b)中出现过的每个字母在(a)中的位置用二进制串表示出来

('a':0101001)

('b':1010110)

我们考虑用增量法（也许是叫增量法？）来构造(M)矩阵

先将空矩阵列出

我们称呼矩阵的第(i)行为(r[i])（指右边标的那个数字(i)）

假设我们已经有了(r[1])和(r[2])

下面来介绍以下怎么通过(r[2])得到(r[3])

首先我们将(r[2])分组，具体操作为每个(1)作为一组的开头（第一个数字也作为一组的开头）

(0 0 0 | 1 0 0 | 1)

然后我们将(b[i])对应的位置的字母的出现位置写出，按照同样的分组

(1 0 1 | 0 1 1 | 0)

对于每一组，我们取最右边的(1)作为(1)，其他位置都变成(0)

(r[3]=0 0 1| 0 0 1|1)

类推可以得到整个(M)

原理：

我们再看刚才的转移过程

(r[2]=0 0 0 | 1 0 0 | 1)

(1 0 1 | 0 1 1 | 0)

最右边的(1)表示从右往左看到最右边的竖线为止，已经完成了(LCS)长度为(1)的匹配

那么我们如果想根据这个长度为(1)的匹配去得到当前状态下长度为(2)的匹配，显然只要找到这个最右的(1)左边第一个可以成为(1)的数字，也就是下一组最右边的(1)

所以这样做是正确的

位运算加速递推

我们发现这样做复杂度还是(O(nm))，没卵用

但是这样做我们把运算搬到了二进制上，可以使用位运算大幅度减少常数

还是用(r[2])得到(r[3])的过程演示：先得到(b[2])位置的字母对应的位置二进制串或上(r[2])

设(X=1 0 1 1 1 1 1)

把(r[2])向左移一位，并将最低位设置为(1)

得到$Y=(0 0 1 0 0 1 1)

令(Z=X-Y)得到（数字意义上的减法）

(Z=100 110 0)

再令(Z=X xor Z)

得到(Z=001 001 1)

然后(r[3]=X and Z=001 001 1)

你问我为什么这样做我也不知道

总结一下(r[i]=X and ((X-(r[i-1]<<1|1))xor X))

如果我们用(bitset)的话，就可以做到优秀的(O(frac{nm}{w}))的复杂度啦

这样就可以通过(2e4)到(3e4)的数据

压位bitset优化

不过我们还可以上压位(bitset)

据说(bitset)常数巨大，我们可以考虑用(unsigned long long)模拟(bitset)操作，可以把(64)位(bitset)压成一位，这样不仅支持了加减法，还可以(O(frac{n}{64}))实现位运算

就完全可以跑过(5e4)的数据了

#include<bits/stdc++.h>
using namespace std;
namespace red{
#define int long long
#define ls(p) (p<<1)
#define rs(p) (p<<1|1)
#define mid ((l+r)>>1)
#define y1 qwq
#define id(x,y) (((x)-1)*m+(y))
	inline int read()
	{
		int x=0;char ch,f=1;
		for(ch=getchar();(ch<'0'||ch>'9')&&ch!='-';ch=getchar());
		if(ch=='-') f=0,ch=getchar();
		while(ch>='0'&&ch<='9'){x=(x<<1)+(x<<3)+ch-'0';ch=getchar();}
		return f?x:-x;
	} 
	const int N=50010,M=1100,inf=0x3f3f3f3f;
	int n,m;
	char  a[N],b[N];
	struct Bitset
	{
		unsigned long long a[M];
	}ans,x,y,A[N];
	Bitset operator & (Bitset a,Bitset b)
	{
		for(int i=0;i<M;++i) a.a[i]&=b.a[i];
		return a;
	}
	Bitset operator ^ (Bitset a,Bitset b)
	{
		for(int i=0;i<M;++i) a.a[i]^=b.a[i];
		return a;
	}
	Bitset operator | (Bitset a,Bitset b)
	{
		for(int i=0;i<M;++i) a.a[i]|=b.a[i];
		return a;
	}
	void operator-=(Bitset &a, Bitset b)
	{
    	    unsigned long long last = 0;
    	    for (int i = 0; i < M; i++)
       		    if(a.a[i] >= b.a[i] + last)
           		    a.a[i]-=b.a[i]+last,last = 0;
       		    else
    	  		    a.a[i]-=b.a[i]+last,last = 1;
	}
	void operator<<=(Bitset &a, Bitset b)  // a=b*2+1
	{
    	    unsigned long long last = 1;
    	    for (int i = 0; i < M; i++)
		    {
       		    unsigned long long x = b.a[i] >> 63;
       		    a.a[i] = (b.a[i] << 1 | last);
        	    last = x;
    	    }
	}
	inline void insert(Bitset &a,int x)
	{
		a.a[x>>6]|=1ull<<(x&63);
	}
	inline int count(Bitset A)
	{
		int ans=0;
		for(int i=0;i<M;++i)
		{
			int sum=0;
			while(A.a[i])
			{
				if(A.a[i]&1) ++sum;
				A.a[i]>>=1;
			}
			ans+=sum;
		}
		return ans;
	}
	inline void main()
	{
		scanf("%s%s",a+1,b+1);
		n=strlen(a+1),m=strlen(b+1);
		for(int i=1;i<=n;++i) insert(A[a[i]-'a'],i);
		for(int i=1;i<=m;++i)
		{
			y<<=ans;
			ans=ans|A[b[i]-'a'];
			x=ans;
			x-=y;
			ans=ans&(ans^x);
		}
		printf("%lld
",count(ans));
	}
}
signed main()
{
	red::main();
	return 0;
}

参考

唐文斌国家集训队论文

大佬题解