POJ-2065 SETI 高斯消元,扩展GCD
生活随笔
收集整理的這篇文章主要介紹了
POJ-2065 SETI 高斯消元,扩展GCD
小編覺得挺不錯的,現在分享給大家,幫大家做個參考.
該題題義是給定如下一個方程組:
F(1) = C1 (mod) P
F(2) = C2 (mod) P
F(3) = C3 (mod) P ...
其中F(1) = A(1,1)*x1 + A(1, 2)*x2 + A(1, 3)*x3...
其中A(i, j) = i ^ (j-1).
面對這樣一個方程,我們的做法就是先不管方程右端的(mod)P,因為F(i) 和 Ci 是同余的,那么在方程兩端乘以一個數是沒有影響的,代入某一個方程同樣是允許的。于是剩下的就是高斯消元。由于解是唯一的,因此解出來的系數矩陣一定是滿秩。接下來就可一用擴展gcd求解未知元了。
代碼如下:
#include <cstdlib> #include <cstring> #include <algorithm> #include <cstdio> using namespace std;typedef long long int Int64;char s[100];int length, MOD;Int64 x[100];int _pow(int a, int b) { // 快速冪int ret = 1;while (b) {if (b & 1) {ret *= a;ret %= MOD;}a *= a;a %= MOD;b >>= 1;}return ret; }struct Matrix {Int64 a[100][100];void init(){for (int i = 1; i <= length; ++i) {for (int j = 1; j <= length; ++j) {a[i][j] = _pow(i, j-1);}a[i][length+1] = s[i];}}void rswap(int x, int y, int s){ // 用于將后面的行交換到第一行上面來 Int64 t;for (int j = s; j <= length+1; ++j) {t = a[x][j];a[x][j] = a[y][j];a[y][j] = t;}}void multi(int x, Int64 M, int s){ // 同一行乘以某一數,當然這時為消元準備的for (int j = s; j <= length+1; ++j) {a[x][j] *= M;}}void relax(int x, int y, int s){ // 將某一方程整體迭代進另外一個表達式for (int j = s; j <= length + 1; ++j) {a[y][j] -= a[x][j];a[y][j] = (a[y][j] % MOD + MOD) % MOD;a[x][j] = (a[x][j] % MOD + MOD) % MOD;}}void print(){for (int i = 1; i <= length; ++i) {for (int j = 1; j <= length+1; ++j) {printf("%5I64d ", a[i][j]);}puts("");}puts("");} }M;Int64 GCD(Int64 a, Int64 b) {Int64 t;while (b) {t = b;b = a % b;a = t;}return a; }Int64 LCM(Int64 a, Int64 b) {return a / GCD(a, b) * b; }Int64 ext_gcd(Int64 a, Int64 b, Int64 &x, Int64 &y) {Int64 temp, ret;if (!b) {x = 1, y = 0;return a;}ret = ext_gcd(b, a % b, x, y);temp = x, x = y, y = temp - a/b*y;return ret; }void Gauss() {int i = 1, k;Int64 a, b, sum;memset(x, 0, sizeof (x));for (int j = 1; j <= length; ++j) {for (k = i; k <= length; ++k) {if (M.a[k][j]) break; }if (k > length) continue;if (k != i) {M.rswap(i, k, j);}for (k = i + 1; k <= length; ++k) {if (M.a[k][j]) {Int64 lcm = LCM(M.a[i][j], M.a[k][j]);a = lcm/M.a[i][j];b = lcm/M.a[k][j];if (a > 1) M.multi(i, a, j);if (b > 1) M.multi(k, b, j);M.relax(i, k, j);}}++i;}// M.print();// 一定會有唯一解,所以沒有進行無解判定for (k = i-1; k >= 1; --k) { // 當然這里可以暴力進行求解,畢竟給定的素數不是很大Int64 xx, yy, A, B, C = 0, L, G;for (int j = k+1; j <= length; ++j) {C += x[j] * M.a[k][j];}A = M.a[k][k], B = MOD;C = M.a[k][length+1] - C; G = ext_gcd(A, B, xx, yy);L = B / G;xx *= C / G; xx = (xx % L + L) % L;x[k] = xx;}for (int j = 1; j <= length; ++j) {printf(j == 1 ? "%I64d" : " %I64d", x[j]);}puts(""); }int main() {int N;scanf("%d", &N);while (N--) {scanf("%d %s", &MOD, s+1);length = strlen(s+1);for (int i = 1; i <= length; ++i) {if (s[i] != '*') {s[i] -= ('a' - 1);}else {s[i] = 0;}}M.init();Gauss();}return 0; }總結
以上是生活随笔為你收集整理的POJ-2065 SETI 高斯消元,扩展GCD的全部內容,希望文章能夠幫你解決所遇到的問題。
- 上一篇: POJ 3687 Labeling Ba
- 下一篇: 关于如何评价洗牌质量的猜想