NOI模拟赛 华为

在这里插入图片描述
考虑设 f [ x ] [ y ] f[x][y] f[x][y] 为从 ( x , y ) (x,y) (x,y) 出发走到 ( 0 , 0 ) (0,0) (0,0) 的期望步数。有 f [ x ] [ y ] = p ∗ f [ x + 1 ] [ y ] + q ∗ f [ x ] [ y + 1 ] + 1 f[x][y]=p*f[x+1][y]+q*f[x][y+1]+1 f[x][y]=pf[x+1][y]+qf[x][y+1]+1
推出 f [ x + 1 ] [ y ] = ( f [ x ] [ y ] − q ∗ f [ x ] [ y + 1 ] − 1 ) / p f[x+1][y]=(f[x][y]-q*f[x][y+1]-1)/p f[x+1][y]=(f[x][y]qf[x][y+1]1)/p 于是设第一行为主元,可以递推出第0行用这些主元表示的系数 g ( x , y ) g(x,y) g(x,y) 。设 G x ( t ) = ∑ i = 0 m − 1 g ( x , i ) t i G_x(t) = \sum_{i=0}^{m-1}g(x,i)t^i Gx(t)=i=0m1g(x,i)ti 。在循环卷积意义下有 G x ( t ) = G x − 1 ( t ) 1 − q t p G_x(t)=G_{x-1}(t) \frac{1-qt}{p} Gx(t)=Gx1(t)p1qt ,于是通过循环卷积可以快速求出第0行的系数,再利用 f [ x ] [ y ] = p ∗ f [ x + 1 ] [ y ] + q ∗ f [ x ] [ y + 1 ] + 1 f[x][y]=p*f[x+1][y]+q*f[x][y+1]+1 f[x][y]=pf[x+1][y]+qf[x][y+1]+1 f [ 0 ] [ 0 ] = 0 f[0][0]=0 f[0][0]=0 列出方程组,高斯消元即可算出主元,接着对于每个询问,同样做循环卷积算出系数,再带入算出的主元即可
复杂度 O ( m 3 + t m log ⁡ m log ⁡ n ) O(m^3+tm \log m \log n) O(m3+tmlogmlogn)

#include <bits/stdc++.h>
#define int long long
using std::array;
using std::cerr;
using std::cin;
using std::copy;
using std::cout;
using std::function;
using std::get;
using std::make_pair;
using std::make_tuple;
using std::map;
using std::mt19937;
using std::pair;
using std::reverse;
using std::set;
using std::sort;
using std::swap;
using std::tuple;
using std::uniform_int_distribution;
using std::unique;
using std::vector;
using ll = long long;
namespace qwq {
mt19937 eng;
void init(int Seed) { return eng.seed(Seed); }
int rnd(int l = 1, int r = 1000000000) { return uniform_int_distribution<int>(l, r)(eng); }
}  // namespace qwq
template <typename T>
inline void chkmin(T &x, T y) {
    if (x > y)
        x = y;
}
template <typename T>
inline void chkmax(T &x, T y) {
    if (x < y)
        x = y;
}
template <typename T>
inline T min(const T &x, const T &y) {
    return x < y ? x : y;
}
template <typename T>
inline T max(const T &x, const T &y) {
    return x > y ? x : y;
}
char buf[100000], *bufs, *buft;
#define gc() \
    ((bufs == buft && (buft = (bufs = buf) + fread(buf, 1, 100000, stdin))), bufs == buft ? -1 : *bufs++)

#define O(x) cerr << #x << " : " << x << '\n'
const double Pi = acos(-1);
const int MAXN = 262144, MOD = 998244353, inv2 = (MOD + 1) / 2, I32_INF = 0x3f3f3f3f;
const long long I64_INF = 0x3f3f3f3f3f3f3f3f;
auto Ksm = [](int x, int y) -> int {
    if (y < 0) {
        y %= MOD - 1;
        y += MOD - 1;
    }
    int ret = 1;
    for (; y; y /= 2, x = (long long)x * x % MOD)
        if (y & 1)
            ret = (long long)ret * x % MOD;
    return ret;
};
auto Mod = [](int x) -> int {
    if (x >= MOD)
        return x - MOD;
    else if (x < 0)
        return x + MOD;
    else
        return x;
};
template <const int N_num, const int M_num>
struct Graph {
    int H[N_num];
    struct Edge {
        int to, lac;
    } e[M_num];
    inline void add_edge(int x, int y) {
        e[*H] = { y, H[x] };
        H[x] = (*H)++;
    }
    inline void init() {
        memset(H, -1, sizeof H);
        *H = 0;
    }
};
#define go(x, y) for (int i = x.H[y], v; (v = x.e[i].to) && ~i; i = x.e[i].lac)
inline int ls(int k) { return k << 1; }
inline int rs(int k) { return k << 1 | 1; }
using ull = unsigned long long;

int fac[MAXN], ifac[MAXN];
namespace POLY {
int SZ, R[MAXN], W[MAXN], INV[MAXN + 1];
void POLYINIT() {
    INV[1] = 1;
    for (int i = 2; i <= MAXN; ++i) INV[i] = (long long)(MOD - MOD / i) * INV[MOD % i] % MOD;
}
void INIT(int len) {
    if (SZ == len)
        return;
    SZ = len;
    for (int i = 1; i < len; ++i) R[i] = (R[i >> 1] >> 1) | (i & 1 ? (len >> 1) : 0);
    int wn = Ksm(3, (MOD - 1) / len);
    W[len >> 1] = 1;
    for (int i = (len >> 1) + 1; i < len; ++i) W[i] = (long long)W[i - 1] * wn % MOD;
    for (int i = (len >> 1) - 1; i > 0; --i) W[i] = W[i << 1];
}
unsigned long long c[MAXN];
void Ntt(vector<int> &F, int limit, int type) {
    copy(F.begin(), F.begin() + limit, c);
    for (int i = 1; i < limit; ++i)
        if (i < R[i])
            swap(c[i], c[R[i]]);
    for (int o = 2, j = 1; o <= limit; o <<= 1, j <<= 1) {
        for (int i = 0; i < limit; i += o) {
            for (int k = 0; k < j; ++k) {
                unsigned long long OI = c[i + j + k] * W[k + j] % MOD;
                c[i + j + k] = c[i + k] + MOD - OI;
                c[i + k] += OI;
            }
        }
    }
    if (type == -1) {
        reverse(c + 1, c + limit);
        int inv = INV[limit];
        for (int i = 0; i < limit; ++i) c[i] = c[i] % MOD * inv % MOD;
    }
    for (int i = 0; i < limit; ++i) F[i] = c[i] % MOD;
}
int w;
typedef std::pair<int, int> complex;
complex operator+(complex &x, complex &y) {
    return std::make_pair((x.first + y.first) % MOD, (x.second + y.second) % MOD);
}
complex operator-(complex &x, complex &y) {
    return std::make_pair((x.first - y.first + MOD) % MOD, (x.second - y.second + MOD) % MOD);
}
complex operator*(complex &x, complex &y) {
    return std::make_pair(((long long)x.first * y.first + (long long)w * x.second % MOD * y.second) % MOD,
                          ((long long)x.first * y.second + (long long)y.first * x.second) % MOD);
}
complex ksm(complex x, int y) {
    complex res(1, 0);
    for (; y; y /= 2, x = x * x)
        if (y & 1)
            res = res * x;
    return res;
}
bool check(int x) { return Ksm(x, (MOD - 1) / 2) == 1; }
int Sqrt(int x) {
    if (!x)
        return x;
    else {
        int a = qwq::rnd() % MOD;
        while (check(((long long)a * a + MOD - x) % MOD)) a = qwq::rnd() % MOD;
        w = ((long long)a * a + MOD - x) % MOD;
        complex b(a, 1);
        int ans1(ksm(b, (MOD + 1) / 2).first);
        return min(ans1, MOD - ans1);
    }
}

int M, N;

struct Poly {
    vector<int> v;
    int &operator[](const int &pos) { return v[pos]; }
    int len() { return v.size(); }
    void set(int l) { return v.resize(l); }
    void adjust() {
        while (v.size() > 1 && !v.back()) v.pop_back();
    }
    void rev() { reverse(v.begin(), v.end()); }
    void Ntt(int L, int type) {
        int limit = 1 << L;
        INIT(limit);
        set(limit);
        POLY::Ntt(v, limit, type);
    }
    void Squ() {
        int L = ceil(log2(len())) + 1, limit = 1 << L;
        Ntt(L, 1);
        for (int i = 0; i < limit; ++i) v[i] = (long long)v[i] * v[i] % MOD;
        Ntt(L, -1);
        adjust();
    }
    void operator+=(Poly &x) {
        if (len() < x.len())
            set(x.len());
        for (int i = 0; i < x.len(); ++i) v[i] = Mod(v[i] + x[i]);
        adjust();
    }
    void operator-=(Poly &x) {
        if (len() < x.len())
            set(x.len());
        for (int i = 0; i < x.len(); ++i) v[i] = Mod(v[i] - x[i]);
        adjust();
    }
    Poly operator*(Poly &x) {
        Poly ret, tmp0 = *this, tmp1 = x;
        int L = ceil(log2(tmp0.len() + tmp1.len() - 1)), n = 1 << L;
        tmp0.Ntt(L, 1);
        tmp1.Ntt(L, 1);
        ret.set(n);
        for (int i = 0; i < n; ++i) ret[i] = (long long)tmp0[i] * tmp1[i] % MOD;
        ret.Ntt(L, -1);
        for (int i = 0; i < M; ++i) ret[i] = (ret[i] + ret[i + M]) % MOD;
        ret.set(M);
        return ret;
    }
    Poly operator-(Poly &x) {
        Poly ret;
        ret.set(max(len(), x.len()));
        for (int i = 0; i < len(); ++i) ret[i] = v[i];
        for (int i = 0; i < x.len(); ++i) ret[i] = Mod(ret[i] - x[i]);
        return ret;
    }
    Poly operator+(Poly &x) {
        Poly ret;
        ret.set(max(len(), x.len()));
        for (int i = 0; i < len(); ++i) ret[i] = v[i];
        for (int i = 0; i < x.len(); ++i) ret[i] = Mod(ret[i] + x[i]);
        return ret;
    }
    void operator*=(Poly &x) {
        Poly tmp = x;
        int L = ceil(log2(len() + x.len() - 1)), n = 1 << L;
        Ntt(L, 1);
        tmp.Ntt(L, 1);
        for (int i = 0; i < n; ++i) v[i] = (long long)v[i] * tmp[i] % MOD;
        Ntt(L, -1);
        for (int i = 0; i < M; ++i) v[i] = (v[i] + v[i + M]) % MOD;
        set(M);
    }
    Poly GetInv(int deg = -1) {
        if (deg == 1)
            return { { Ksm(v[0], MOD - 2) } };
        Poly ret = GetInv((deg + 1) / 2), tmp;
        int L = ceil(log2(deg)) + 1, n = 1 << L, mx = min(len(), deg);
        tmp.set(deg);
        for (int i = 0; i < mx; ++i) tmp[i] = v[i];
        tmp.Ntt(L, 1);
        ret.Ntt(L, 1);
        for (int i = 0; i < n; ++i) ret[i] = (2 - (long long)tmp[i] * ret[i] % MOD + MOD) * ret[i] % MOD;
        ret.Ntt(L, -1);
        ret.set(deg);
        return ret;
    }
    pair<Poly, Poly> operator%(Poly &x) {
        if (x.len() > len())
            return { { { 0 } }, *this };
        Poly tmp0 = *this, tmp1 = x;
        tmp0.rev();
        tmp1.rev();
        tmp1 = tmp1.GetInv(len() - x.len() + 1);
        tmp0 *= tmp1;
        tmp0.set(len() - x.len() + 1);
        tmp0.rev();
        tmp1 = tmp0 * x;
        Poly ret = *this - tmp1;
        ret.set(x.len() - 1);
        return { tmp0, ret };
    }
    Poly Dif(int deg = -1) {
        Poly tmp;
        tmp.set(max(len() - 1, 1ll));
        for (int i = 0; i < len() - 1; ++i) tmp[i] = v[i + 1] * (i + 1LL) % MOD;
        if (~deg)
            tmp.set(deg);
        return tmp;
    }
    Poly GetSqrt(int deg = -1) {
        if (deg == 1)
            return { { POLY::Sqrt(v[0]) } };
        Poly ret = GetSqrt((deg + 1) / 2), tmp0 = ret.GetInv(deg), tmp1;
        int L = ceil(log2(deg)) + 1, mx = min(len(), deg);
        tmp1.set(deg);
        for (int i = 0; i < mx; ++i) tmp1[i] = v[i];
        tmp0 *= tmp1;
        tmp0.set(deg);
        ret += tmp0;
        for (auto &i : ret.v) i = (long long)i * inv2 % MOD;
        return ret;
    }
    Poly Int(int deg = -1) {
        Poly tmp;
        tmp.set(len() + 1);
        for (int i = 1; i < tmp.len(); ++i) tmp[i] = (long long)v[i - 1] * INV[i] % MOD;
        if (~deg)
            tmp.set(deg);
        return tmp;
    }
    Poly GetLn(int deg = -1) {
        Poly tmp0 = Dif(deg), tmp1 = GetInv(deg);
        tmp0 *= tmp1;
        tmp0.set(deg);
        return tmp0.Int(deg);
    }
    Poly GetExp(int deg = -1) {
        if (deg == 1)
            return { { 1 } };
        Poly tmp0 = GetExp((deg + 1) / 2), tmp1 = tmp0.GetLn(deg);
        for (int i = 0; i < deg; ++i) tmp1[i] = Mod(v[i] - tmp1[i]);
        tmp1[0] = Mod(tmp1[0] + 1);
        tmp0 *= tmp1;
        tmp0.set(deg);
        return tmp0;
    }
    Poly GetKsm(int deg, int K) {
        Poly tmp0 = GetLn(deg);
        for (auto &i : tmp0.v) i = (long long)i * K % MOD;
        return tmp0.GetExp(deg);
    }
};
vector<int> getmulpointvalue(Poly &a, vector<int> &x) {
    static Poly tmp[MAXN * 4];
    function<void(int, int, int)> get = [&](int u, int l, int r) -> void {
        if (l == r) {
            tmp[u] = { { Mod(-x[l]), 1 } };
            return;
        }
        int mid = (l + r) / 2;
        get(ls(u), l, mid);
        get(rs(u), mid + 1, r);
        tmp[u] = tmp[ls(u)] * tmp[rs(u)];
    };
    get(1, 0, x.size() - 1);
    vector<int> ret(x.size());
    function<void(int, int, Poly, int)> solve = [&](int l, int r, Poly f, int u) -> void {
        if (l == r) {
            ret[l] = f[0];
            return;
        }
        int mid = (l + r) / 2;
        solve(l, mid, (f % tmp[ls(u)]).second, ls(u));
        solve(mid + 1, r, (f % tmp[rs(u)]).second, rs(u));
    };
    solve(0, x.size() - 1, (a % tmp[1]).second, 1);
    return ret;
}
vector<int> Interpolation(vector<pair<int, int>> &x) {
    static Poly tmp[MAXN * 4];
    vector<int> tmpx;
    for (auto &i : x) tmpx.push_back(i.first);
    function<void(int, int, int)> get = [&](int u, int l, int r) -> void {
        if (l == r) {
            tmp[u] = { { Mod(-tmpx[l]), 1 } };
            return;
        }
        int mid = (l + r) / 2;
        get(ls(u), l, mid);
        get(rs(u), mid + 1, r);
        tmp[u] = tmp[ls(u)] * tmp[rs(u)];
    };
    get(1, 0, tmpx.size() - 1);
    Poly ret = tmp[1];
    ret = ret.Dif();
    vector<int> tmpy = getmulpointvalue(ret, tmpx);
    for (int i = 0; i < x.size(); ++i) tmpy[i] = (long long)Ksm(tmpy[i], MOD - 2) * x[i].second % MOD;
    function<Poly(int, int, int)> solve = [&](int u, int l, int r) -> Poly {
        if (l == r)
            return { { tmpy[l] } };
        int mid = (l + r) / 2;
        Poly tmp0 = solve(ls(u), l, mid) * tmp[rs(u)];
        Poly tmp1 = solve(rs(u), mid + 1, r) * tmp[ls(u)];
        return tmp0 + tmp1;
    };
    return solve(1, 0, x.size() - 1).v;
}
// Poly t[MAXN * 4];
}  // namespace POLY
using namespace POLY;

Poly f, g;

Poly ksm(Poly x, int y) {
    Poly res;
    res.set(M);
    res[0] = 1;
    while (y) {
        // cout<<res[0]<<' '<<x[0]<<'\n';
        if (y & 1)
            res *= x;

        // cout<<res[0]<<'\n';
        x *= x;  //
        y >>= 1;
    }
    return res;
}

ll val[MAXN], mat[2100][2100];

void solve() {
    for (int i = 1; i <= M; i++) {
        for (int j = i; j <= M; j++)
            if (mat[j][i]) {
                if (j != i)
                    swap(mat[i], mat[j]);
                break;
            }
        ll inv = Ksm(mat[i][i], MOD - 2);
        for (int j = i + 1; j <= M; j++) {
            ll mul = inv * (MOD - mat[j][i]) % MOD;
            (mat[j][0] += mat[i][0] * mul) %= MOD;
            for (int k = i; k <= M; k++) (mat[j][k] += mat[i][k] * mul) %= MOD;
        }
    }
    for (int i = 1; i <= M; i++) mat[i][0] = MOD - mat[i][0];
    for (int i = M; i; i--) {
        ll v = mat[i][0] * Ksm(mat[i][i], MOD - 2) % MOD;
        val[i - 1] = v;
        for (int j = i - 1; j > 0; j--) (mat[j][0] += MOD - v * mat[j][i] % MOD) %= MOD;
    }
}

inline int read() {
    int x = 0, f = 1;
    char ch = getchar();
    while (ch < '0' || ch > '9') {
        if (ch == '-')
            f = -1;
        ch = getchar();
    }
    while (ch >= '0' && ch <= '9') {
        x = x * 10 + ch - '0';
        ch = getchar();
    }
    return x * f;
}

signed main() {
    freopen("huawei.in", "r", stdin);
    freopen("huawei.out", "w", stdout);
    POLYINIT();
    N = read(), M = read();
    int p = read(), q = read();
    int s = Ksm(p + q, MOD - 2);
    (p *= s) %= MOD, (q *= s) %= MOD;
    int ip = Ksm(p, MOD - 2), iq = Ksm(q, MOD - 2);
    f.set(M);
    f[0] = ip, f[1] = (MOD - q) * ip % MOD;
    f = ksm(f, N - 1);
    // cout<<f[0]<<" "<<f[1]<<'\n';
    for (int i = 1; i <= M; i++) {
        mat[i][0] = (MOD - ip) * (N - 1) % MOD;
        for (int j = 1; j <= M; j++) mat[i][j] = f[(i - j + M) % M];
        if (i > 1) {
            mat[i][0] += MOD - 1 - (MOD - ip) * (N - 1) % MOD * q % MOD;
            mat[i][0] %= MOD;
            for (int j = 1; j <= M; j++) (mat[i][j] += MOD - f[(i - j + M - 1) % M] * q % MOD) %= MOD;
            (mat[i][i] += MOD - p) %= MOD;
        }
    }
    /*	for(int i=1;i<=M;i++){
                    for(int j=0;j<=M;j++)
                            printf("%lld ",mat[i][j]);
                    puts("");
            }*/
    solve();
    int t = read();
    while (t--) {
        int x = read(), y = read();
        ll ans = (MOD - ip) * (N - 1 - x) % MOD;
        f.set(0);
        f.set(M);
        f[0] = ip, f[1] = (MOD - q) * ip % MOD;
        f = ksm(f, N - 1 - x);
        for (int i = 0; i < M; i++) (ans += f[(y - i + M) % M] * val[i]) %= MOD;
        printf("%lld\n", (ans % MOD + MOD) % MOD);
    }
    return 0;
}

  • 0
    点赞
  • 0
    收藏
    觉得还不错? 一键收藏
  • 0
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

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

抵扣说明:

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

余额充值