循环最长公共子序列/Cyclic Longest Common Subsequence

循环最长公共子序列

参考：https://arxiv.org/pdf/1208.0396.pdf

问题描述

我们需要求两个循环序列的最长公共子序列，对于循环序列来说，把任意一个元素作为起始元素都认为是同一个序列，例如，循环序列’algorithm’同时也可以写作’rithmalgo’或者’thmalgori’，对于正常的序列’algorithm’和’grammar’来说，最长公共子序列为’grm’，但是对于循环序列来说，它们的最长公共子序列是’grma’

预备知识

对于普通的非循环序列的最长公共子序列来说，记两个序列分别为$A,B$，$A$的元素个数为$m$，$B$的元素个数为$n$，例如$A$为$BBAA$，$B$为$ABABB$，可以将其转化为有向图，每个单元格可以看作图中的一个节点，向右或者向下可以连一条边，对于单元格$(i,j)$来说，如果$A[i]$和$B[j]$的元素相同，则从$(i-1,j-1)$向$(i,j)$斜着连一条边，如下图所示，记该有向图为$G_{A,B}$

• 寻找最长公共子序列相当于寻找从$(0,0)$到$(m,n)$的最短路径，因为最短路径相当于走最多的斜边，而斜边对应着相同的元素；当有多条最短路径时，任意一条都是符合条件的。

算法描述

我们把循环最长公共子序列问题用符号CLCS表示，普通的最长公共子序列用符号LCS表示，为了解决CLCS问题，我们的方法是，将其中一个字符串复制一遍，例如将$A$复制一遍后成为$AA$，则$AA$与$B$构成的有向图为$G_{AA,B}$（如下图所示），我们需要求的CLCS对应的最短路径是$G_{AA,B}中$从$(0,0)$到$(m,n)$，从$(1,0)$到$(m+1,n)$，…,，从$(m−1,0)$到$(2m−1,n)$中的其中一条；朴素的算法时间复杂度为$O(m^2n)$，因为我们需要寻找$m$条最短路径，即解决$m$次LCS问题。而由于所有的路径都是存在于同一张图$G_{AA,B}$上，因此我们可以复用一些最短路径计算的工作，从而减少计算量。

我们固定图$G$的方向，按照图示的方向，我们可以说某一条边在图中比另一条边低（下面的边比上面的边低），对于路径来说也是类似的，我们定义最低最短路径为：从起点到终点的所有最短路径中，最低的那条，即处于图中最下方的那条。

我们用$parent[i][j]$来表示最短路径中$(i,j)$节点的前一个节点，在同样是最短路径的情况下，我们优先选择$\leftarrow$，然后是$\nwarrow$，最后才选择$\uparrow$，这样就能获得一条以$(0,0)$为起始点的最低最短路径。

如果我们计算出了图$G_{AA,B}$的最低最短路径，首先，我们可以直接获取到从$(0,0)$到$(m,n)$的最短路径，但是在我们要获取从$(1,0)$到$(m+1,n)$的最短路径时，我们需要以$(1,0)$作为起点，求新的最低最短路径，如果之前的最低最短路径经过了$(1,0)$，则之前的最低最短路径仍然是新的最低最短路径，如果之前的最低最短路径没经过$(1,0)$，经过的点为$(1,k)$，则我们可以基于之前经过$(1,k)$的那条最低最短路径来更新新的最低最短路径。

假设最初的最低最短路径经过点$(1,k)$（$k$必定满足$parent[1][k]=\nwarrow$，且所有$j<k, parent[1][j]=\leftarrow$），我们记这条最低最短路径为$R$，我们记经过$(1,0)$的最低最短路径为$L$，我们记路径$P$上的斜边数量为$len(P)$，即公共子序列的长度，在这种情况下，$len(R)=len(L)+1$

我们记以$(i,j)$为终点的最短路径上的斜边数量为$len(i,j)$，即以$(i,j)$为终点的最长公共子序列的长度，容易得到：

$len(i,j-1)\le len(i,j)$，因为向右走斜边数量不会减少

$len(i,j)\le len(i,j-1) + 1$，因为多走了一列，且每一列的斜边最多只能选择一个，因此最多只能多一个斜边

而由于$R$是最低最短路径，且$len(R)=len(L)+1$，因此当$(i,j-1) \isin L$且$(i,j) \isin R$时，$len(i,j)=len(i,j-1)+1$，所以我们只需要将$R$路径上的$parent[i][j]$修改为$\leftarrow$即可得到新的最低最短路径；设原本的起点是$(root-1,0)$，可以得到伪代码如下，RE-ROOT部分的时间复杂度为$O(m+n)$

RE-ROOT操作之后，就能够很容易得到以任意$(i,0)$为起点的最低最短路径了，对应着以某个字符作为起始字符的最长公共子序列，整体算法的为代码如下，总体时间复杂度为$O(mn)$

例题

https://oj.uz/problem/view/IZhO13_rowords

代码实现

#include <bits/stdc++.h>
 
using namespace std;
 
template <typename T1, typename T2> istream &operator>>(istream &is, pair<T1, T2> &pa) { is >> pa.first >> pa.second; return is; }
template <typename T> istream &operator>>(istream &is, vector<T> &vec) { for (auto &v : vec) is >> v; return is; }
template <typename T1, typename T2> ostream &operator<<(ostream &os, const pair<T1, T2> &pa) { os << "(" << pa.first << "," << pa.second << ")"; return os; }
template <typename T> ostream &operator<<(ostream &os, const vector<T> &vec) { os << "["; for (auto v : vec) os << v << ","; os << "]"; return os; }
template <typename T> ostream &operator<<(ostream &os, const deque<T> &vec) { os << "deq["; for (auto v : vec) os << v << ","; os << "]"; return os; }
template <typename T> ostream &operator<<(ostream &os, const set<T> &vec) { os << "{"; for (auto v : vec) os << v << ","; os << "}"; return os; }
template <typename T> ostream &operator<<(ostream &os, const multiset<T> &vec) { os << "{"; for (auto v : vec) os << v << ","; os << "}"; return os; }
template <typename T> ostream &operator<<(ostream &os, const unordered_set<T> &vec) { os << "{"; for (auto v : vec) os << v << ","; os << "}"; return os; }
template <typename T> ostream &operator<<(ostream &os, const unordered_multiset<T> &vec) { os << "{"; for (auto v : vec) os << v << ","; os << "}"; return os; }
template <typename TK, typename TV> ostream &operator<<(ostream &os, const unordered_map<TK, TV> &mp) { os << "{"; for (auto v : mp) os << v.first << "=>" << v.second << ","; os << "}"; return os; }
template <typename TK, typename TV> ostream &operator<<(ostream &os, const map<TK, TV> &mp) { os << "{"; for (auto v : mp) os << v.first << "=>" << v.second << ","; os << "}"; return os; }
template <typename T> void resize_array(vector<T> &vec, int len) { vec.resize(len); }
template <typename T, typename... Args> void resize_array(vector<T> &vec, int len, Args... args) { vec.resize(len); for (auto &v : vec) resize_array(v, args...); }
template <typename T1, typename T2> pair<T1, T2> operator+(const pair<T1, T2> &l, const pair<T1, T2> &r) { return make_pair(l.first + r.first, l.second + r.second); }
template <typename T1, typename T2> pair<T1, T2> operator-(const pair<T1, T2> &l, const pair<T1, T2> &r) { return make_pair(l.first - r.first, l.second - r.second); }
long long gcd(long long a, long long b) { return b ? gcd(b, a % b) : a; }
mt19937 mrand(random_device{}());
int rnd(int x) { return mrand() % x; }

#define rep(i, a, n) for (int i = a; i < (n); i++)
#define per(i, a, n) for (int i = (n)-1; i >= a; i--)
#define pb push_back
#define mp make_pair
#define all(x) (x).begin(), (x).end()
#define fi first
#define se second
#define sz(x) ((int)(x).size())
typedef vector<int> vi;
typedef long long ll;
typedef unsigned int uint;
typedef unsigned long long ull;
typedef pair<int, int> pii;
typedef double db;
#if DEBUG
#define dbg(x) cerr << #x << " = " << (x) << " (L" << __LINE__ << ") " << __FILE__ << endl;
#else
#define dbg(x)
#endif

class Solution {
public:
    void Solve() {
        string s,t;
        while(cin>>s>>t) {
            s += s;
            // dbg(s);
            int n = sz(s), m=sz(t);
 
            auto clcs = [&] () {
                vector<vi> dp(n+1, vi(m+1));
                vector<vi> parent(n+1, vi(m+1));
                rep(i,1,n+1) parent[i][0] = 2;
 
                rep(i,0,n) {
                    rep(j,0,m) {
                        if (s[i] == t[j]) dp[i+1][j+1] = dp[i][j] + 1;
                        else dp[i+1][j+1] = max(dp[i][j+1], dp[i+1][j]);
 
                        if (dp[i+1][j+1] == dp[i+1][j]) parent[i+1][j+1] = 0;
                        else if (s[i] == t[j] && dp[i+1][j+1] == dp[i][j] + 1) parent[i+1][j+1] = 1;
                        else parent[i+1][j+1] = 2;
                    }
                }
 
                auto lcs = [&] (int root) {
                    int x = root + n/2;
                    int y = m;
                    int cnt = 0;
                    while (x > root && y > 0) {
                        if (parent[x][y] == 1) {
                            cnt++;
                            x--; y--;
                            // dbg(mp(x,y));
                            // dbg(cnt);
                        } else if (parent[x][y] == 0) y--;
                        else x--;
                    }
                    return cnt;
                };
 
                auto reroot = [&] (int root) {
                    int x = root;
                    int y = 0;
 
                    while (y <= m && parent[x][y] != 1) y++;
                    if (y > m) return;
 
                    parent[x][y] = 0;
                    while (x < n && y < m) {
                        if (parent[x+1][y] == 2) {
                            x++;
                            parent[x][y] = 0;
                        } else if (parent[x+1][y+1] == 1) {
                            x++; y++;
                            parent[x][y] = 0;
                        } else {
                            y++;
                        }
                    }
 
                    while (x < n && parent[x+1][y] == 2) {
                        x++;
                        parent[x][y] = 0;
                    }
                };
 
                int ans = 0;
                ans = max(ans, lcs(0));
                rep(root, 1, n/2) {
                    reroot(root);
                    ans = max(ans, lcs(root));
                }
                // dbg(ans);
                return ans;
            };
 
            int ans = clcs();
            reverse(all(s));
            ans = max(ans, clcs());
            cout << ans << endl;
        }
    }
private:
};

void set_io(const string &name = "") {
    ios::sync_with_stdio(false);
    cin.tie(nullptr);
#if FILE_IO
    if (!name.empty()) {
        freopen((name + ".in").c_str(), "r", stdin);
        freopen((name + ".out").c_str(), "w", stdout);
    }
#endif
}

int main() {
    set_io("input");
    Solution().Solve();
 
    return 0;
}