UVa1457/LA4746 Decrypt Messages

题目链接

  本题是2009年icpc亚洲区域赛上海赛区的题目

题意

  假设从2000年1月1日00:00:00到现在经过了x秒，计算$x^{q}modp$，设答案为a。这里p严格大于x。已知p, q, a，求现在时刻x的所有可能值。
  提示：如果一个年份是4的倍数但不是100的倍数，或者这个年份是400的倍数，这个年份是闰年。闰年的二月有29天，其他年（平年）的二月只有28天。在本题中，如果年份除以10 的余数为5或者8，则这一年的最后还会有一个“闰秒”。比如，2005年12月31日23:59:59的下一秒是2005 年12月31日23:59:60，再下一秒才是2006 年1月1日00:00:00。

输入格式

  输入的第一行为数据组数T。每组数据包含一行，包含3 个整数p, q, a（2＜p≤1000000007，1＜q≤10，0≤a＜p，p 保证为素数）。

输出格式

  对于每组数据，输出所有可能的时间，按照时间顺序排列。如果无解，出“Transmission error”。

分析

  离散对数解高次模方程经典题目，详细思路《算法竞赛入门经典–训练指南》第444页有写，涉及原根、离散对数等知识，参见OI Wiki。

模方程$x^q\equiv a(\mathrm{mod}m)$。如果找到了p的一个原根m，只需设$x=m^y,a=m^z$，则方程变为$m^{qy}\equiv m^z(\mathrm{mod}m)$，即$qy\equiv z(\mathrm{mod}p-1)$。q是已知量，而z可以用大步小步算法得到（z是以m为底的a的离散对数，且解唯一），因此只需解这个模线性方程就可以得到y（注意，这一步可能多解），最后进行一次模取幂运算，得到x。

  注意，题面交代0≤a＜p，对a=0的情况无法换元$a=m^z$求离散对数，这时候x只有一个答案$x=0$，直接输出即可。由于模p可达$10^9$，取模前的乘法运算需要用long long承接一下。

测试数据生成

  给一份生成测试数据的python脚本：

# -*- coding: utf-8 -*-from random import randint, choiceT, Q, N = 1000000008, 10, 2000if __name__ == '__main__':f, primes = [False] * T, []for i in range(2, T):if not f[i]:print(i)for j in range(2*i, T, i):f[j] = Trueif i > 2:primes.append(i)with open("in.txt", "w") as f:f.write(f'{N+1}\n')f.write(f'{choice(primes)} {randint(2, Q)} 0\n')for _ in range(N):p, q = choice(primes), randint(2, Q)f.write(f'{p} {q} {randint(0, p-1)}\n')

  运行大约2～3分钟就得到测试输入，贴到uDebug就能得到正确输出。

AC 代码

#include <iostream>
#include <algorithm>
#include <iomanip>
#include <cmath>
#include <vector>
#include <map>
using namespace std;#define T 32
int c[T] = {0}, p, q, a, t, kase = 0;int gcd(int a, int b, int& x, int& y) {if (!b) {x = 1; y = 0; return a;} else {int g = gcd(b, a%b, y, x);y -= a/b*x;return g;}
}int inv(int a, int n) {int x, y;return gcd(a, n, x, y) == 1 ? (x + n) % n : -1;
}int pow_mod(long long a, int n) {int ans = 1;while (n) {if (n & 1) ans = ans*a % p;a = a*a % p; n >>= 1;}return ans;
}bool is_primitive_root(int m) {for (int i=2; i*i<p; ++i) if ((p-1) % i == 0 && (pow_mod(m, i) == 1 || pow_mod(m, (p-1)/i) == 1)) return false;return true;
}int bsgs(int a, int b) {int m = sqrt(p+.5), s = inv(pow_mod(a, m), p), e = 1; map<int, int> x;x[1] = 0;for (int i=1; i<m; ++i) {e = e*(long long)a % p;if (!x.count(e)) x[e] = i;}for (int i=0, j=(p+m-1)/m; i<j; ++i) {if (x.count(b)) return i*m + x[b];b = b*(long long)s %p;}return -1;
}int f(int x) {return (x%4 == 0 && (x%100 || x%400 == 0) ? 31622400 : 31536000) + (x%10 == 5 || x%10 == 8);
}void print(int x) {int y = T-1, m = 12, d, hh, mm, ss;for (int i=1; i<T; ++i) if (c[i] >= x) {y = c[i] == x ? i : i-1;break;}if ((x -= c[y]) == 0) {cout << y+2000 << ".01.01 00:00:00" << endl;} else {bool f = y%4 == 0 && (y%100 || y%400 == 0);for (int i=1, t=0; i<13; ++i) {int s = (i==2 ? 28+f : (i==4 || i==6 || i==9 || i==11 ? 30 : 31)) * 86400;if (i==12 || t+s > x) {m = i; x -= t;break;}t += s;}d = min((x + 86399) / 86400, 31); x -= 86400*(d-1);hh = min(x / 3600, 23); x -= 3600*hh;mm = min(x / 60, 59); ss = x - 60*mm;cout << y+2000 << '.' << setw(2) << m << '.' << setw(2) << d << ' '<< setw(2) << hh << ':' << setw(2) << mm << ':' << setw(2) << ss << endl;}
}void solve() {cin >> p >> q >> a;cout << "Case #" << ++kase << ':' << endl;if (a == 0) {cout << "2000.01.01 00:00:00" << endl;return;}int m, x, y;for (m=2; m<p; ++m) if (is_primitive_root(m)) break;int b = bsgs(m, a), g = gcd(q, p-1, x, y);if (b % g) {cout << "Transmission error" << endl;return;}vector<int> ans; x = (x*(long long)b/g % (p-1) + p-1) % (p-1); y = (p-1) / g;for (int i=0; i<g; ++i) ans.push_back(pow_mod(m, x+i*y));sort(ans.begin(), ans.end());g = unique(ans.begin(), ans.end()) - ans.begin();for (int i=0; i<g; ++i) print(ans[i]);
}int main() {for (int i=1; i<T; ++i) c[i] = c[i-1] + f(i-1);cin >> t; cout << setfill('0');while (t--) solve();return 0;
}