非常神仙的状压DP+容斥原理。
首先，给出一个状压方程：$f_S$表示点集为$S$的情况下，整个点集构成强连通图的方案数。
这个DP方程还是比较容易想到的，但是没有办法正常转移，考虑通过容斥原理进行转移。
对于一个点集，它无法构成强连通分量的方案，就是我们选择一个出度为$0$的强连通分量，这个强连通分量并不包含整体的方案，就是无法构成的方案数，也就是缩点后的图是一个至少两个节点的DAG。
那么，我们可以钦定一个点集$j,jsubset S$作为出度为$0$的强连通分量，那么可以得到，这样其他的不是从这个点集连向其他点集的边是随意选取的，也就是$2^p$种方案。
但是由于我们不知道$S-j$的部分中有没有出现出度为$0$的强连通分量，如果出现，那么就重复计算了。
所以考虑容斥，$t_i$表示$i$点集构成点数大于1的DAG的方案数，$t_S=sumlimits_{jsubseteq S,j e 0}{(-1)^{|j|-1} imes 2^{way_{i-j,j}}} imes t_{i-j}$，其中$way_{j,i-j}$表示$i-j$点集向$j$点集链接的边，$|j|$表示选择的强连通分量集的强连通分量数量。
上面的式子应该会很好理解，就不多解释了，那么根据上面的式子，我们提出一个建设性的想法，$g_i$表示$i$的点集中，构成奇数个互不连通强连通分量-构成偶数个互不连通强连通分量的方案数，也就是说在状态中直接包含容斥系数，那么可以给出方程：$g_S=f_S-sumlimits_{jsubseteq S,uin j}g_{S-j} imes f_j$
这个式子的正确性关键在于$g_S$构成的强连通分量是互不连通(因为出度全部为$0$)的，所以这样转移显然是正确的。
那么接下来考虑$f_S$如何更新。
同样，我们回到上面的式子，$t_S$可以表达为：$t_S=sumlimits_{jsubset i,j e 0}g_j imes 2^{sum_S-w_j}$，$w_j$表示$j$连向$S$的边数。
然后，$f_S=2^{sum_S}-t_S$
然后我们就解决这道题啦！撒花撒花✿✿ヽ(°▽°)ノ✿
建议搭配下面文档食用
https://files-cdn.cnblogs.com/files/Winniechen/BZOJ3812.pdf
#include <cstdio>
#include <algorithm>
#include <cmath>
#include <cstring>
#include <cstdlib>
#include <queue>
#include <iostream>
#include <bitset>
using namespace std;
#define N (1<<15)+20
#define ll long long
#define mod 1000000007
int f[N],g[N],cnt[N],in[N],out[N],n,m,b[N],w[N],sum[N];
void dfs(int i,int j)
{
if(i&(j-1))dfs(i,i&(j-1));
w[j]=w[j-(j&-j)]+cnt[in[j&-j]&i];
}
int main()
{
scanf("%d%d",&n,&m);b[0]=1;
for(int i=1,x,y;i<=m;i++)
{
scanf("%d%d",&x,&y);x--,y--;
in[1<<x]|=1<<y,out[1<<y]|=1<<x;
b[i]=(b[i-1]<<1)%mod;
}
for(int S=1;S<1<<n;S++)
{
int x=S&-S,s=S^x;cnt[S]=cnt[s]+1;sum[S]=sum[s]+cnt[in[x]&S]+cnt[out[x]&S];
dfs(S,S);f[S]=b[sum[S]];
for(int j=s;j;j=s&(j-1))g[S]=(g[S]-(ll)f[S^j]*g[j])%mod;
for(int j=S;j;j=S&(j-1))f[S]=(f[S]-(ll)g[j]*b[sum[S]-w[j]])%mod;
g[S]=(g[S]+f[S])%mod;
}
printf("%d
",(f[(1<<n)-1]+mod)%mod);
}