poj 2065SETI <高斯消元>

链接： http://poj.org/problem?id=2065

题意：

给你一个素数P(P<=30000)和一串长为n的字符串str[]。字母'*'代表0，字母a-z分别代表1-26，这n个字符所代表的数字分别代表f(1)、f(2)....f(n)。

定义: f (k) = ∑(0<=i<=n-1) a^iki (mod p) (1<=k<=n,0<=ai<P)

求a0、a1.....an-1

#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <algorithm>
using namespace std;
int P, N;
int d[100][100], x[100];
int Gcd( int a, int b )
{
    return b==0?a:Gcd( b, a%b );
} 
int Lcm( int a, int b )
{
    return a/Gcd( a, b )*b;
}
int P_M( int a, int b, int m )
{
    int ans=1, t=a%m;
    while( b ){
        if( b&1 )
            ans*=t;
        if( ans>m )ans%=m;
        b>>=1;
        t*=t;
        t%=m;
    }
    return ans; 
}
void Gauss(  )
{
    int i=1, j, p, k, t;
    for( j=1; j<=N; ++ j ){    
    
        for( p=i; p<=N; ++ p ){
            if( d[p][j] )break;
        }    
        
        if( p>N )continue;
        if( p!=i ){
        
            for( k=j;k<=N+1; ++ k ){
                t=d[p][j],d[p][j]=d[i][j], d[i][j]=t;
            }
        }
    
        for( p=i+1; p<=N; ++p ){
            if(d[p][j]){
                int f1=Lcm( d[p][j], d[i][j] )/d[p][j];
                int f2=Lcm( d[p][j], d[i][j] )/d[i][j];// 可以直接用f1=d[p][j],但是不能放在下面 一个for里面 
                for(k=j; k<=N+1; ++k ){    
                    d[p][k] =(d[p][k]*f1-d[i][k]*f2)%P;
                }
            }
        }
        ++i;
    }
    memset( x, 0, sizeof x );
    for( p=N; p>=1; --p ){
            int sum=0; 
            for(k=p+1; k<=N;++k ){
                sum = ((sum+x[k]*d[p][k])%P+P)%P;
            }
            for( k=0;k<=P; ++k ){
                if(( (d[p][p]*k)%P + P)%P== ((d[p][N+1]-sum)%P + P) % P){//同余方程求解 
                      x[p] = k;
                      break;
                }
            } 
    }
    for( int i=1; i<=N; ++ i ){
        printf( "%d ", x[i] );    
    }
    puts( "" );
}
int main( )
{
    int T;
    char s[100];
    scanf( "%d",&T );
    while( T-- ){
        scanf( "%d", &P );
        scanf( "%s", s);
        N=strlen( s );
        memset( d, 0, sizeof d );
        for( int i=0; i<N; ++ i ){
            if(s[i]=='*')d[i+1][N+1]=0;
            else d[i+1][N+1]=s[i]-'a'+1;
        }
        for( int i=1; i<=N; ++ i ){
            for( int j=1; j<=N; ++ j ){
                d[i][j]=P_M( i, j-1 , P);
            }
        }/*
        for( int i=1; i<=N; ++ i )
        {
            for(int j=1; j<=N+1; ++ j)    
            {
                printf( "%d ", d[i][j] );    
            }
            puts( "" );
        }*/
        Gauss( );    
        
    
    }
    return 0;
}