数据生成器 (Generator)

生成器会构造两种数据：

1. 完全随机：LCS 长度通常较短。
2. 进化模拟：先生成一个祖先序列，然后随机插入或删除一些碱基生成两个后代。这种数据的 LCS 会很长，更符合生物学背景。

请保存为 gen_lcs.cpp 并运行。

/**
 * Project: Homologous Gene Search (LCS Generator)
 * Author: Isaac Asimov (AI)
 * Note: Sample 1 output is confirmed to be 4 (ATCA).
 */

#include <iostream>
#include <fstream>
#include <random>
#include <string>
#include <vector>
#include <algorithm>

using namespace std;

// ==========================================
// Part 1: 标准解答逻辑
// ==========================================
class Solution {
public:
    void solve(string in_file, string out_file) {
        ifstream cin(in_file);
        ofstream cout(out_file);
        if (!cin.is_open()) return;

        string s1, s2;
        cin >> s1 >> s2;
        int n = s1.length();
        int m = s2.length();
        
        vector<vector<int>> dp(n + 1, vector<int>(m + 1, 0));

        for (int i = 1; i <= n; i++) {
            for (int j = 1; j <= m; j++) {
                if (s1[i-1] == s2[j-1]) {
                    dp[i][j] = dp[i-1][j-1] + 1;
                } else {
                    dp[i][j] = max(dp[i-1][j], dp[i][j-1]);
                }
            }
        }
        cout << dp[n][m] << endl;
        
        cin.close();
        cout.close();
        cout << "Generated: " << out_file << endl;
    }
};

// ==========================================
// Part 2: 数据生成逻辑
// ==========================================
mt19937 rng(2025);

int rand_int(int min, int max) {
    uniform_int_distribution<int> dist(min, max);
    return dist(rng);
}

char rand_base() {
    string bases = "ATCG";
    return bases[rand_int(0, 3)];
}

string mutate(string ancestor, int mutation_rate) {
    string offspring = "";
    for (char c : ancestor) {
        int dice = rand_int(1, 100);
        if (dice <= mutation_rate) {
            int type = rand_int(0, 2);
            if (type == 0) continue; 
            else if (type == 1) { offspring += rand_base(); offspring += c; }
            else offspring += rand_base();
        } else {
            offspring += c;
        }
    }
    return offspring;
}

void generate_input(int case_id) {
    string filename = to_string(case_id) + ".in";
    ofstream fout(filename);

    string s1, s2;

    if (case_id == 1) {
        // 样例 1: 结果应为 4
        s1 = "ATCGA";
        s2 = "ATACA";
    }
    else if (case_id == 2) {
        s1 = "AAAAAA";
        s2 = "TTTTTT";
    }
    else if (case_id <= 5) {
        int len = rand_int(10, 20);
        for(int i=0; i<len; i++) s1 += rand_base();
        for(int i=0; i<len; i++) s2 += rand_base();
    }
    else if (case_id <= 8) {
        string ancestor = "";
        int len = rand_int(50, 100);
        for(int i=0; i<len; i++) ancestor += rand_base();
        s1 = mutate(ancestor, 20);
        s2 = mutate(ancestor, 20);
    }
    else {
        int len = 1000;
        for(int i=0; i<len; i++) s1 += rand_base();
        for(int i=0; i<len; i++) s2 += rand_base();
    }

    fout << s1 << endl;
    fout << s2 << endl;
    fout.close();
}

int main() {
    ios_base::sync_with_stdio(false);
    cout << "--- Generating Homologous Gene Data (Corrected) ---" << endl;
    
    Solution solver;
    for (int i = 1; i <= 10; i++) {
        generate_input(i);
        string in = to_string(i) + ".in";
        string out = to_string(i) + ".out";
        solver.solve(in, out);
    }
    cout << "--- Done! ---" << endl;
    return 0;
}

阿西莫夫的解题指南

这是二维动态规划的开山之作。

1. 状态定义

我们需要一个二维数组 dp[i][j]。 定义：dp[i][j] 表示字符串 $S_1$ 的前 $i$ 个字符（$S_1[1 \dots i]$）与字符串 $S_2$ 的前 $j$ 个字符（$S_2[1 \dots j]$）之间，能找到的最长公共子序列的长度。

2. 状态转移方程

当我们考虑 $S_1$ 的第 $i$ 个字符和 $S_2$ 的第 $j$ 个字符时，只有两种情况：

1. 字符相同 ($S_1[i] == S_2[j]$)： 这就好比我们在两个序列末尾都找到了一个匹配的碱基。那我们就可以把这个碱基加到之前的成果上。

   $$dp[i][j] = dp[i-1][j-1] + 1$$

2. 字符不同 ($S_1[i] \neq S_2[j]$)： 这意味着这两个字符不能同时作为公共子序列的结尾。我们需要做出抉择：

   • 要么扔掉 $S_1$ 的第 $i$ 个字符，看看 $S_1[1 \dots i-1]$ 和 $S_2[1 \dots j]$ 的匹配结果。
   • 要么扔掉 $S_2$ 的第 $j$ 个字符，看看 $S_1[1 \dots i]$ 和 $S_2[1 \dots j-1]$ 的匹配结果。
   • 我们要取两者中的较大值（生物进化保留了最优秀的特征）。
   $$dp[i][j] = \max(dp[i-1][j], \ dp[i][j-1])$$

3. 边界条件

• dp[0][j] = 0 （$S_1$ 为空，匹配长度为 0）
• dp[i][0] = 0 （$S_2$ 为空，匹配长度为 0）

C++ 标准解答 (NOIP C++14)

/**
 * 题目: 寻找同源基因 (LCS - Longest Common Subsequence)
 * 作者: Isaac Asimov (AI)
 * 算法: 二维动态规划 (2D DP)
 * 时间复杂度: O(N*M)
 */

#include <iostream>
#include <string>
#include <vector>
#include <algorithm> // max

using namespace std;

// 建议使用稍微大一点的数组，开在全局区
// dp[i][j] 对应 S1 的前 i 个字符和 S2 的前 j 个字符
int dp[1005][1005];

int main() {
    // IO 优化
    ios_base::sync_with_stdio(false);
    cin.tie(NULL);

    string s1, s2;
    cin >> s1 >> s2;

    int n = s1.length();
    int m = s2.length();

    // 为了方便处理边界 dp[0][0]，我们让下标从 1 开始对应字符串
    // 字符串下标依然是 0 到 n-1，所以在访问 s1[i-1] 时要注意
    
    for (int i = 1; i <= n; i++) {
        for (int j = 1; j <= m; j++) {
            // 注意：dp下标是从1开始，字符串下标是从0开始
            // 所以比较的是 s1[i-1] 和 s2[j-1]
            if (s1[i-1] == s2[j-1]) {
                // 找到同源碱基，长度 + 1
                dp[i][j] = dp[i-1][j-1] + 1;
            } else {
                // 碱基不同，继承之前的最大值
                dp[i][j] = max(dp[i-1][j], dp[i][j-1]);
            }
        }
    }

    // 最终答案在 dp[n][m]
    cout << dp[n][m] << endl;

    return 0;
}

阿西莫夫的点评

1. 为什么是二维？ 如果用一维 DP，我们很难同时记录两个字符串当前的进度。二维表格 $DP[i][j]$ 的物理意义就像是在画一个网格图，寻找一条从左上角到右下角的路径，使得路径上“匹配”的跳跃次数最多。
2. 生物学应用： 虽然这道题是 $O(N^2)$ 的基础版，但在真实的生物信息学中，如果序列长达数十亿（如人类基因组），这个复杂度是不可接受的。因此科学家发明了 BLAST 等启发式算法。但 LCS 始终是理解序列比对（Sequence Alignment）的理论基石。