[HNOI2008]Cards

Time Limit: 10
Sec  Memory Limit:
162 MB

Description
小春现在很清闲,面对书桌上的N张牌,他决定给每张染色,目前小春只有3种颜色:红色,蓝色,绿色.他询问Sun有多少种染色方案,Sun很快就给出了答案.进一步,小春要求染出Sr张红色,Sb张蓝色,Sg张绝色.他又询问有多少种方案,Sun想了一下,又给出了正确答案.

最后小春发明了M种不同的洗牌法,这里他又问Sun有多少种不同的染色方案.两种染色方法相同当且仅当其中一种可以通过任意的洗牌法(即可以使用多种洗牌法,而每种方法可以使用多次)洗成另一种.Sun发现这个问题有点难度,决定交给你,答案可能很大,只要求出答案除以P的余数(P为质数).

第一行输入五个整数:Sr,Sb,Sg,M,P(M<=60,M+1N=Sr+Sb+Sg.接下来M行,每行描述一种洗牌法,每行有N个用空格隔开的整数X1X2...Xn.

恰为1到N的一个排列,表示使用这种洗牌法,第I位变成原来的Xi位的牌.

输入数据保证任意多次洗牌都可以用M种洗牌法中的一种代替,且对每种洗牌法,都存在一种洗牌法使得能回到

原状态.

20%的数据保证Sr+Sb+Sg<=20

100%的数据保证Max(Sr,Sb,Sg)<=20

Input
第一行输入五个整数:Sr,Sb,Sg,M,P(M<=60,M+1接下来M行,每行描述一种洗牌法,每行有N个用空格隔开的整数X1X2...Xn.

恰为1到N的一个排列,表示使用这种洗牌法,第I位变成原来的Xi位的牌.

输入数据保证任意多次洗牌都可以用M种洗牌法中的一种代替,且对每种洗牌法

都存在一种洗牌法使得能回到原状态.

20%的数据保证Sr+Sb+Sg<=20

100%的数据保证Max(Sr,Sb,Sg)<=20

Output
不同染法除以P的余数

Sample
Input
1 1 1 2 7

2 3 1

3 1 2
Sample
Output
2

Hint
有2 种本质上不同的染色法RGB
和RBG，使用洗牌法231 一次可得GBR 和BGR，使用洗牌法312 一次

可得BRG 和GRB。

 做这题，要用到了群论。扩展欧几里得求乘法逆元。。。补得我惨惨的。。至今没弄透彻。。。
首先要用到一个定理：
一个群的状态总数是是其每个置换群置换的不变量加起来除以置换群数目。
这里用到DP。对于每一个置换，找到与原序列构成的循环节。
例如：
1 2 3 4 5 6 7 8 9
2 1 4 6 7 9 3 8 5
其中1 2构成一个循环节，3 4 5 6 7
9构成一个循环节，8单独成循环节。若要使置换的群不变，只要使循环节内的所有卡片的颜色相同即可。

这样，以每一个循环节为阶段，记录每个环内数的个数s[i]，开三维数组做DP求出置换所得的不变量。
最后用扩展欧几里得求乘法逆元。
状态转移方程：

f[i,j,k]:=f[i-s[kk],j,k]+f[i,j-s[kk],k]+f[i,j,k-s[kk]]; 
(0<=i<=sr
0<=j<=sb
0<=k<=sg)
 
AC
CODE
 

program

hy_1004;

var

f,ff:

array

[

0..20

,

0..20

,

0..20

]

of

longint

;

a:

array

[

1..61

,

1..60

]

of

longint

;

s:

Array

[

1..60

]

of

longint

;

p:

array

[

1..60

]

of

boolean

;

ans,size,n,m,s1,s2,s3,mo,x,y:

longint

;

//============================================================================

procedure

init;

var

i,j:

longint

;

begin

readln(s1,s2,s3,m,mo);

n:=s1+s2+s3;
inc(m);

for

i:=

1

to

m-

1

do

for

j:=

1

to

n

do

read(a[i,j]);

for

i:=

1

to

n

do

a[m,i]:=i;

end

;

//============================================================================

function

work(x:

longint

):

longint

;

var

i,j,k,kk,y:

longint

;

begin

for

i:=

1

to

n

do

p[i]:=

false

;

for

i:=

0

to

s1

do

for

j:=

0

to

s2

do

for

k:=

0

to

s3

do

f[i,j,k]:=

0

;

size:=

0

;

for

i:=

1

to

n

do

if

not

(p[i])

then

begin

inc(size);
s[size]:=

1

;
p[i]:=

true

; y:=i;

while

not

(p[a[x,y]])

do

begin

p[a[x,y]]:=

true

;

inc(s[size]);
y:=a[x,y];

end

;

end

;

f[

0

,

0

,

0

]:=

1

;

for

kk:=

1

to

size

do

begin

for

i:=s1

downto

0

do

for

j:=s2

downto

0

do

for

k:=s3

downto

0

do

begin

     
4000
     

if

i+s[kk]<=s1

then

f[i+s[kk],j,k]:=(f[i+s[kk],j,k]+f[i,j,k])

mod

mo;

if

j+s[kk]<=s2

then

f[i,j+s[kk],k]:=(f[i,j+s[kk],k]+f[i,j,k])

mod

mo;

if

k+s[kk]<=s3

then

f[i,j,k+s[kk]]:=(f[i,j,k+s[kk]]+f[i,j,k])

mod

mo;

end

;

end

;
work:=f[s1,s2,s3];

end

;

//============================================================================

procedure

DP;

var

i:

longint

;

begin

ans:=

0

;

for

i:=

1

to

m

do

ans:=(ans+work(i))

mod

mo;

end

;

//============================================================================

function

extended_gcd(a,b:

longint

;

var

x,y:

longint

):

longint

;

var

t:

longint

;

begin

if

b=

0

then

begin

extended_gcd:=a;

x:=

1

;
y:=

0

;

end

else

begin

extended_gcd:=extended_gcd(b,a

mod

b,x,y);

t:=x;
x:=y; y:=t-(a

div

b)*y;

end

;

end

;

//============================================================================

procedure

consult;

var

tmp,t1,t2:

longint

;

begin

tmp:=extended_gcd(m,mo,x,y);

t1:=mo

div

tmp; t2:=m

div

tmp;

while

(x<=

0

)

or

(y>

0

)

do

begin

inc(x,t1);
dec(y,t2);

end

;

writeln

((ans*x)

mod

mo);

end

;

//============================================================================

begin

init;

DP;

consult; 

end

.