POJ 2065 SETI 高斯消元解线性同余方程

题意：

给出mod的大小，以及一个不大于70长度的字符串。每个字符代表一个数字，且为矩阵的增广列。系数矩阵如下

1^0 * a₀ + 1^1 * a₁ + ... + 1^(n-1) * a_n-1　　=　　f(1)

2^0 * a₀ + 2^1 * a₁ + ... + 2^(n-1) * an-1　   =　　f(2)

........

n^0 * a₀ + n^1 * a₁ + ... + n^(n-1) * a_n-1 =　　f(n)

快速幂取模下系数矩阵

#include <cstdio>
#include <cstring>
#include <algorithm>
using namespace std;

const int MAXN = 100;
int a[MAXN][MAXN], x[MAXN];
int MOD;

void debug(int n, int m)
{
    for(int i=0; i<n; i++)
    {
        for(int j=0; j<m; j++)
            printf("%d ", a[i][j]);
        printf("  %d
", a[i][m]);
    }
    puts("**************************************************************");
}

int powermod(int a,int b)
{
    int ans=1;
    a%=MOD;
    while(b)
    {
        if(b&1) ans=ans*a%MOD;
        a=a*a%MOD;
        b=b>>1;
    }
    return ans;
}

int gcd(int a,int b)               //递归算法
{
    return b ? gcd(b, a%b) : a;
}

int lcm(int a, int b)
{
    return a*b/gcd(a,b);
}

int Guass(int equ,int var)
{
//    debug(equ, var);
    int row,col;
    row=col=0;
    while(row<equ && col<var)
    {
        //列非零主
        int r=row;
        for(int i=row; i<equ; i++)
            if(a[i][col]!=0)
            {
                r=i;
                break;
            }
        if(r!=row)
        {
            for(int j=col; j<var+1; j++)
                swap(a[row][j],a[r][j]);
        }
        if(a[row][col]==0)//说明有自由变元
        {
            col++;
            continue;
        }
        //消元
        for(int i=row+1; i<equ; i++)
        {
            if(a[i][col]==0) continue;
            int l = lcm(a[row][col],a[i][col]);
            int ta = l/a[row][col];
            int tb = l/a[i][col];
            for(int j=col; j<var+1; j++)
                a[i][j] = ((tb*a[i][j] - ta*a[row][j]) % MOD + MOD) %MOD;
        }
//        debug(equ, var);
        row++;
        col++;
    }
//    for(int i=row; i<equ; i++)
//        if(a[i][var]!=0) return -1;
//    if(row < var) return 1;
    for(int i=row-1; i>=0; i--)
    {
        int tmp = a[i][var];
        for(int j=i+1; j<var; j++)
            tmp = ((tmp - x[j]*a[i][j])%MOD + MOD)%MOD;
        while(tmp%a[i][i]) tmp += MOD;
        x[i] = tmp/a[i][i]%MOD;
    }
    return 0;
}

char s[100];

int main()
{
//    freopen("in.txt", "r", stdin);
    int t;
    scanf("%d", &t);
    while(t--)
    {
        scanf("%d%s", &MOD, s);
        int n = strlen(s);
        memset(a, 0, sizeof(a));
        for(int i=0; i<n; i++)
            for(int j=0; j<n; j++)
                a[i][j] = powermod(i+1,j);
        for(int i=0; i<n; i++)
            a[i][n] = s[i]=='*'? 0:s[i]-'a'+1;
        Guass(n, n);
        for(int i=0; i<n; i++)
        {
            if(x[i]<0) x[i] += MOD;
            printf("%d%c", x[i], i==n-1? '
':' ');
        }
    }
    return 0;
}