[考试反思]数学专题测试4：深度
高斯消元专题测试
联赛后第一次爆零？？？我也不知道发生了什么
感觉和放假的关系不大啊。。。
虽说题比平时难而且数据范围出了一点锅，但是爆零有点过了？？
一眼看T3，MatrixTree啊肯定的，然后想到之前做的那个不等式的那道题直觉是容斥。
结果又想到了那个专题里幻想乡的那道题，于是更加笃定是容斥，结果思路就阻塞了。
20分的部分分好像很难写，就往后放了放，结果最后并没有时间写。
然后看起来比较简单的是T2，直到看到数据范围之前我都认为应该不会太难。
基础式子都会写，然后组合数忘记特判不合法的0了，于是把10分送掉了，可能也是因为用不习惯键盘打的时候着急了？
剩下的部分分完全看不出来它想让我干什么。。。只是因为数据范围少给了一个条件
考试的时候差点就在群里问这题是不是无解了。。。数据范围毫无破绽，如果s=1e18，m=1e9，其余为0，那么这就显然不可做了。
我还以为是我too naive于是就没有问。。。谁知道它就真的出锅了。
然而skyh猜测数据范围拿到了60分。。。猜测数据范围。。。学到了，也许吧。。。
其实最开始想到了这档部分分的打法，挺显然的，但是数据范围里没说，就没写呗。。。
看了眼T1，想起之前做的某道猜拳的题，期望dp呗。这道题要的是相对胜率更大呗，那就这么写呗
写了一半发现不是很对，虽说感觉能骗一点分但是能拿多少我也不知道。还多测，乱写可能就爆零了。
于是求稳，转向30分的部分分，大型分类讨论。
结果一直写到考试结束，没写完，交了，只有第一档部分分写了。
然后测试点又锅了没有第一档的部分分。
所以，就莫名其妙的爆零了？
嗯，就这么爆零了。
考试结束后没几分钟就涨了20分，，为何如此弱智。
T1：猜拳游戏
大意：猜拳，一回合n场，其中A决策固定（出的概率已知），赢的场多者赢得回合。A领先m2回合或B领先m1回合时取胜。求B最优决策时胜率。
$n \le 1000,m1,m2 \le 100$
比较明显的是，对于不同回合的相同局面采取的决策一定相同（局面是指已进行场数及本回合分差）。
平局的回合没有影响，但是对决策转移很碍事，考虑如何去掉。
平局了就再来一局呗，所以平局的那部分概率可以按比例分给双方。如果原来胜率为p败率是q那么在实际赢得一局的概率$\frac{p}{p+q}$
化简一下$\frac{p}{q+p}=\frac{1}{\frac{q}{p}+1}$。所以我们现在的目的就是最大化$\frac{q}{p}$
然后是一个比较经典的01分数规划思路，二分答案观察是否合法，所以现在设二分值为mid。
那么你可以认为，赢一回合积1分，输一回合扣mid分，平局不变，是否存在最优决策使总积分为正。
这个直接当成普通的期望dp倒着做回去就行了。
设$dp[i][j]$表示目前在第$i$场，已经净胜$j$场（可负）的最优积分。转移就是在石头剪刀布三种决策里取最大积分。
初状态就是$dp[n+1][x>0]=1,dp[n+1][0]=0,dp[n+1][x<0]=-mid$。答案就是$dp[1][0]$
然后我们求出了最大单回合胜率，问题在于净胜场达到$m1/m2$了。
画出转移图，点数不多，直接高斯消元解就可以。
居然卡记搜的常数弄成T40我还以为我写错了，差评
#include<bits/stdc++.h>
using namespace std;
int n,m1,m2;double p[4][1001],K,P,A[202][202],dp[1002][2006];
unordered_map<int,double>M[1111];
#define f(a,x) p[a][r]*dp[r+1][v+x]
double sch(){
for(int i=0;i<=n;++i)dp[n+1][i]=-K;for(int i=n+2;i<=n+n+2;++i)dp[n+1][i]=1;dp[n+1][n+1]=0;
for(int r=n;~r;--r)for(int v=n+n+1;v;--v)dp[r][v]=max(max(f(1,0)+f(2,-1)+f(3,1),f(1,1)+f(2,0)+f(3,-1)),f(1,-1)+f(2,1)+f(3,0));
return dp[1][n+1];
}
void Gauss(int n){
for(int i=0;i<=n;++i){
int np=i;
for(int j=i+1;j<=n;++j)if(fabs(A[j][i])>fabs(A[i][np]))np=j;swap(A[np],A[i]);
if(fabs(A[i][i])<1e-6)continue;
for(int j=n+1;j>=i;--j)A[i][j]/=A[i][i];
for(int j=0;j<=n;++j)if(j!=i)for(int k=n+1;k>=i;--k)A[j][k]-=A[j][i]*A[i][k];
}
}
int main(){
cin>>n>>m1>>m2;
if(n+m1+m2==0)return 0;
for(int i=1;i<=n;++i)for(int j=1;j<=3;++j)cin>>p[j][i],p[j][i]/=100;
double l=0,r=1e6;
while(K=(l+r)/2,r-l>1e-6){
if(sch()>0)l=K;else r=K;
for(int i=1;i<=n;++i)M[i].clear();
}
P=1-1/(1+K);
memset(A,0,sizeof A);
A[m2][m1+m2+1]=1;
for(int i=0;i<m1+m2-1;++i)A[i][i+1]=P-1;
for(int i=2;i<=m1+m2;++i)A[i][i-1]=-P;
for(int i=0;i<=m1+m2;++i)A[i][i]+=1;
Gauss(m1+m2);
printf("%.5lf\n",A[m1+m2][m1+m2+1]);
return main();
}T2：inequ
不会。
T3：生成树
大意：无向图边有3种色，蓝边不超过b绿边不超过g的生成树数量。$n \le 40,m \le 1000$
这除了矩阵树还能是啥？问题在于怎么处理特殊限制。
然而这次并不是容斥，思维僵化了。
怎么在矩阵树里对不同种类的边加以区分？加个边权如何？（怎么想到的？？？）
在矩阵树里的边权可以通过变元矩阵树解决啊，但是这会改变生成树的数量。
为了方便我们将没有限制的红边作为基准，它的边权为1，蓝边为x绿边为y。
那么对于一个用上了a条蓝边b条绿边的生成树，它所贡献的方案树就是$x^ay^b$
而对于一种特定的边权方案我们用矩阵树定理所求出的总方案数就是$\sum\limits_{a=0}^{n-1} \sum\limits_{b=0}^{n-1} r_{a,b}x^ay^b $
其中$r_{a,b}$就是用了啊条蓝边b条绿边的生成树数量。我们要求的答案就是$\sum\limits_{i =0}^{b} \sum\limits_{j=0}^{g} r_{i,j}$
我们现在能对于给定的$x,y$用矩阵树定理求出和式的值，能计算出$x^ay^b$这个数是几，我们要求的是$r$这个系数。
而$x,y$是我们自己随便定的，所以只要我们定若干组$x,y$我们就能列出许多等式，用它们解出所有$r$的值。
其实这个过程就是插值的过程，然而二维的插值并不会，所以就草率的高斯消元将就一下就可以了。时间复杂度$O(n^6)$
#include<bits/stdc++.h>
using namespace std;
#define mod 1000000007
int A[44][44],B[1666][1666],n,m,a[4][44][44],lb,lc,re[44][44],cnt,X[1666],Y[1666],ans;
int qp(int b,int t,int a=1){for(;t;t>>=1,b=1ll*b*b%mod)if(t&1)a=1ll*a*b%mod;return a;}
int Matrix_Tree(){
int a=1;
for(int i=1;i<n;++i){
int iv=qp(A[i][i],mod-2);a=1ll*a*A[i][i]%mod;
for(int j=1;j<n;++j)A[i][j]=1ll*A[i][j]*iv%mod;
for(int j=i+1;j<n;++j)for(int k=n-1;k>=i;--k)A[j][k]=(A[j][k]-1ll*A[j][i]*A[i][k]%mod+mod)%mod;
}
return a;
}
void Gauss(int n){
for(int i=1;i<=n;++i){
if(!B[i][i])for(int j=i+1;j<=n;++j)if(B[j][i]){swap(B[i],B[j]);break;}
int iv=qp(B[i][i],mod-2);
for(int j=1;j<=n+1;++j)B[i][j]=1ll*B[i][j]*iv%mod;
for(int j=1;j<=n;++j)if(i!=j)for(int k=n+1;k>=i;--k)B[j][k]=(B[j][k]-1ll*B[j][i]*B[i][k]%mod+mod)%mod;
}
}
int main(){
cin>>n>>m>>lb>>lc;
for(int i=1,x,y,k;i<=m;++i)cin>>x>>y>>k,a[k][x][x]++,a[k][y][y]++,a[k][x][y]--,a[k][y][x]--;
for(int x=0;x<n;++x)for(int y=0;x+y<n;++y)re[x][y]=++cnt,X[cnt]=x,Y[cnt]=y;
for(int x=0;x<n;++x)for(int y=0;x+y<n;++y){
for(int p=1;p<=cnt;++p)B[re[x][y]][p]=1ll*qp(x,X[p])*qp(y,Y[p])%mod;
for(int p=1;p<=n;++p)for(int q=1;q<=n;++q)A[p][q]=(a[1][p][q]+1ll*a[2][p][q]*x+1ll*a[3][p][q]*y)%mod;
B[re[x][y]][cnt+1]=Matrix_Tree();
}Gauss(cnt);
for(int x=0;x<=lb;++x)for(int y=0;y<=lc;++y)ans=(ans+B[re[x][y]][cnt+1])%mod;
cout<<ans<<endl;
}