[NOI2017]泳池

Description
  有一个 1001×n 1001 × n 的的网格,最下面一行连着海滩,每个格子有 q q 的概率是安全的,1q的概率是危险的。
  定义一个矩形是合法的当且仅当:

  - 这个矩形中每个格子都是安全的
  - 必须紧贴网格的下边界

  问你最大的合法子矩形大小为 k k 的概率是多少。

  n109,k1000

Solution
  “快喊救护车来,这里有一棵白菜被题目虐死了!!!”
  这题目做了一天。。(我好菜菜啊

  首先拆开询问,等价于是用小于等于 k k 的答案减去小于等于k1的答案。
  我们来一个 根本想不出的 DP D P ,令 f[i][j] f [ i ] [ j ] 表示考虑一共有 j j 列,现在从最靠近海滩的那一行(称为第1行)到第 i i 行的情况不考虑,默认为一定安全。然后我们在第i+1行放至少一个危险格子,保证当前最大合法子矩形合法的概率。同时注意保证 ijk i ∗ j ≤ k ,即就算下面的 i i 行全是安全的仍然满足条件。
  那么有:
  

f[i][j]=t=0j1u>ivi(f[u][t]qt(ui))(1q)(f[v][jt1]q(jt1)(vi))

  意思是在第 i+1 i + 1 行枚举一个危险格子,它的坐标是 (i+1,t+1) ( i + 1 , t + 1 ) ,它的左边和右边我们分别枚举,并中间有多少空的就应该乘 q q 的那么多次方。
  (为什么一个是等于一个是大于等于呢?因为若存在几个连续高度为一的危险格子,那么就会把方案算重,这么写就保证只在连续的那一段最后一个高度为一的地方计入答案)
  前缀和优化,我们记F[i][j]=uif[u][j],枚举 i i j的时间为 k k ,所以我们复杂度可以做到O(k2)

  定义 h[i] h [ i ] 为有 i i 列的答案是多少,于是有
  

h[i]=f[1][i]qi+j=1if[1][j1]qj1(1q)h[ij]ik

  

h[i]=j=1k+1f[1][j1]qj1(1q)h[ij]i>k h [ i ] = ∑ j = 1 k + 1 f [ 1 ] [ j − 1 ] ∗ q j − 1 ∗ ( 1 − q ) ∗ h [ i − j ] − − − i > k

  那么求解出 h[n] h [ n ] 即为答案。


  然而并没有介素,还过不了呐筒子们啊!矩阵乘法优化直接炸飞起~
  
  来来来,先看一下这个戳我戳我:),学习一下常系数齐次线性递推优化矩阵快速幂,找了一上午资料最后看见这个感觉整棵菜都好起来了。然后这个会写那个链接所补充的模板题就好。
  我们可以先推出 ik i ≤ k h[i] h [ i ] 们,之后那个式子, h[ij] h [ i − j ] 的前面都是常数对吧,直接就是个其次线性方程了,就用那个方法做就行。暴力多项式取模是 O(k2) O ( k 2 ) 的,总复杂度 O(k2logn) O ( k 2 l o g n ) 可以通过此题了。
  
  好像还有 FFT F F T 优化多项式取模?!留坑。

Source

//2018-4-28
//miaowey
//
#include <bits/stdc++.h>
using namespace std;

#define Set(a, v) memset(a, v, sizeof(a))
#define For(i, a, b) for(int i = (a); i <= (int)(b); ++i)
#define Forr(i, a, b) for(int i = (a); i >= (int)(b); --i)

#define N (1000 + 5)
const int P = 998244353;

inline int Mod(int a){
    if(a < 0) return a + P;
    return a >= P? a - P: a;
}
inline int Mul(int a, int b){
    return (long long)a * b % P;
}

inline int Pow(int a, int anum){
    int ret = 1;
    while(anum){
        if(anum & 1) ret = Mul(ret, a);
        a = Mul(a, a); anum >>= 1;
    }
    return ret;
}

int n, x, y, p, q, powq[N], f[N][N], g[N], h[N], b[N], c[N], tmp[N << 1];

void ExpMul(int *A, int *B, int k){
    For(i, 0, k) For(j, 0, k) tmp[i + j] = Mod(tmp[i + j] + Mul(A[i], B[j]));

    Forr(i, k << 1, k){
        For(j, 0, k - 1) tmp[i + j - k] = Mod(tmp[i + j - k] + Mul(tmp[i], g[k - j]));
        tmp[i] = 0;
    }

    For(i, 0, k - 1) A[i] = tmp[i], tmp[i] = 0;
}

void ExpPow(int *A, int anum, int *ret, int k){
    while(anum){
        if(anum & 1) ExpMul(ret, A, k);
        ExpMul(A, A, k); anum >>= 1;
    }
}

int Solve(int k){
    if(!k) return Pow(p, n);

    Set(f, 0); f[k + 1][0] = 1;

    Forr(i, k, 1){
        int m = min(n, k / i);

        f[i][0] = 1;
        For(j, 0, m) g[j] = Mul(Mul(f[i + 1][j], powq[j]), p);

        For(j, 1, m){
            For(t, 0, j - 1) f[i][j] = Mod(f[i][j] + Mul(g[t], f[i][j - t - 1]));
            f[i][j] = Mod(f[i][j] + Mul(f[i + 1][j], powq[j]));
        }
    }

    ++k;
    For(i, 1, k) g[i] = Mul(Mul(f[1][i - 1], powq[i - 1]), p);

    h[0] = 1;
    For(i, 1, k - 1){
        h[i] = Mul(powq[i], f[1][i]);
        For(j, 1, i) h[i] = Mod(h[i] + Mul(h[i - j], g[j]));
    }

    if(n < k) return h[n];

    int ret = 0;

    Set(b, 0); Set(c, 0); b[0] = c[1] = 1;
    ExpPow(c, n, b, k);

    For(i, 0, k - 1) ret = Mod(ret + Mul(b[i], h[i]));

    return ret;
}

int main(){
#ifndef ONLINE_JUDGE
    freopen("pool.in", "r", stdin);
    freopen("pool.out", "w", stdout);
#endif

    int k;

    scanf("%d%d%d%d", &n, &k, &x, &y);  
    q = Mul(x, Pow(y, P - 2)); p = Mod(1 - q);

    powq[0] = 1;
    For(i, 1, 1000) powq[i] = Mul(powq[i - 1], q);

    printf("%d\n", Mod(Solve(k) - Solve(k - 1)));

    return 0;
}
  • 0
    点赞
  • 1
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值