如果题解区公式锅了，可以考虑直接去博客查看（（
Best Theorem 的板子题…
Best Theorem 本质上是用来统计有向图中欧拉回路数目定理。记 $t_x$ 表示以 $x$ 为根的树形图的数量，$G=\{V,E\}$， 那么有
$$\mathrm{ec}(G)= t_x \prod _{i\in V}(\deg_i-1)!$$
其中 $\deg_x$ 表示 $x$ 的度。如果图中有欧拉回路，那么存在 $\deg_{(in)x}=\deg_{(out)x}$，所以可以随便选一种来统计。
考虑如何统计出以某个点为根的生成树。发现其实从欧拉回路的数量上来证明，取哪个点都一样。所以就考虑把某个点删除之后的生成树数量，这就是以该点为根的生成树数量。这部分直接 Matrix Tree + 辗转相除 即可。
但是这题还没完。考虑起点为 $1$ ，所以以 $1$ 为开头的路，是可以循环同构的，而其他的点是不能循环同构的，也就是本质上，原来的best已经把从所有点出发的循环同构算过一遍了，唯独不算从根开始的全局同构。但现在要求算上，于是最后的结果需要乘上 $\deg_1$。
upd：
然后这个题因为叙述不严谨存在一些细节问题。首先用 Best 定理求解时原图必须存在欧拉回路，可以用出入度来判断；其次考虑如果要走完全部的边，不包含 1 的极大连通子图的大小至多为 1，否则必定有边走不到。其它细节稍注意一下就好了。
orz `JKlover`
```cpp
const int N = 110 ;
const int M = 500010 ;
const int P = 1000003 ;
typedef long long LL ;
int T ;
int n ;
int I[M] ;
int fa[M] ;
LL res ;
LL fac[M] ;
LL deg[M] ;
LL D[N][N] ;
LL K[N][N] ;
LL A[N][N] ;
LL Gauss(){
LL ans = 1 ;
for (int i = 1 ; i < n ; ++ i){
for (int j = i + 1 ; j < n ; ++ j){
while(K[j][i]){
LL t = K[i][i] / K[j][i] ;
for(int k = i ; k < n ; ++ k)
K[i][k] = (K[i][k] - t * K[j][k]) % P ;
swap(K[i], K[j]), ans = (-ans % P + P) % P ;
}
}
if (K[i][i]) (ans *= K[i][i]) %= P ;
}
return (ans % P + P) % P ;
}
void clear(){
memset(K, 0, sizeof(K)) ;
memset(A, 0, sizeof(A)) ;
memset(D, 0, sizeof(D)) ;
memset(I, 0, sizeof(I)) ;
for (int i = 1 ; i <= n ; ++ i) fa[i] = i ;
}
int find_f(int x){
return x == fa[x] ? x : fa[x] = find_f(fa[x]) ;
}
int main(){
cin >> T ; fac[0] = 1 ;
for (LL i = 1 ; i <= 500000 ; ++ i)
fac[i] = fac[i - 1] * i % P ;
while (T --){
scanf("%d", &n) ; clear() ;
for (int i = 1 ; i <= n ; ++ i){
int k, s ; cin >> s ; deg[i] = s ;
for (int j = 1 ; j <= s ; ++ j){
cin >> k ; A[i][k] ++, D[k][k] ++, ++ I[k] ;
if (find_f(i) != find_f(k)) fa[find_f(i)] = find_f(k) ;
}
}
for (int i = 1 ; i <= n ; ++ i)
for (int j = 1 ; j <= n ; ++ j)
K[i][j] = ((D[i][j] - A[i][j]) % P + P) % P ;
for (int i = 1 ; i <= n ; ++ i)
if (I[i] != deg[i]) goto lalala ;
for (int i = 1 ; i <= n ; ++ i)
if (find_f(i) != find_f(1) && deg[i]) goto lalala ;
res = deg[1] * Gauss() % P ;
for (int i = 1 ; i <= n ; ++ i)
if (deg[i]) (res *= fac[deg[i] - 1]) %= P ;
for (int i = 2 ; i <= n ; ++ i)
if (find_f(i) == find_f(1)) goto lelele ;
cout << fac[deg[1]] << endl ; continue ;
lalala : puts("0") ; continue ;
lelele : cout << res << endl ;
}
return 0 ;
}
```