• Bluestein's Algorithm


    用下降幂搞搞事情

    [ frac{m(m-1) + k(k+1) -(m-k)(m-k-1)}{2} \ =frac{m^2 - m + k^2 + k - m^2 - k^2 + 2mk - k + m}{2} \ = mk \ y_{m} = sum_{k=0}^n A_{k}omega^{mk}\=omega^{frac{m(m-1)}{2}}sum_{k=0}^n A_komega^{frac{k(k+1)}{2}}omega^{frac{-(m-k)(m-k-1)}{2}} ]

    后面是一个卷积形式,可以快速实现求 (f(x),f(x^1), f(x^2) .. f(x^n)) 以及循环卷积意义下的多项式快速幂。

    例题 UOJ500 任意基DFT

    观察发现

    [q_m = q_0x^m+sum_{i=0}^{m-1} xy^i =(frac{y}{x-1}+q_0)x^m-frac{y}{x-1}\ ]

    ​ 设 (c = frac{y}{x-1}) ,那么有

    [f(q_m) = sum_{i=0}^n a_i((c+q_0)x^m-c)^i \ =sum_{i=0}^n a_isum_{j=0}^i{ichoose j}(c+q_0)^jx^{mj}(-c)^{i-j} \ =sum_{i=0}^n(c+q_0)^ix^{mi}sum_{j=i}^n{i choose j}(-c)^{i-j}a_i ]

    (b_j = (c+q_0)^i sum_{j=i}^n{i choose j}(-c)^{i-j}a_i) ,这是一个卷积形式,可以用一遍 ( ext{FFT}) 求出所有的 (b) ,那么

    [f(q_m) = sum_{i=0}^nb_i x^{mi} ]

    再跑一遍 bluestein's 对 (1-m) 求解即可,注意用线性递推求幂次来优化常数。

    code

    /*program by mangoyang*/
    #pragma GCC optimize("Ofast", "inline")
    #include<bits/stdc++.h>
    #define inf (0x7f7f7f7f)
    #define Max(a, b) ((a) > (b) ? (a) : (b))
    #define Min(a, b) ((a) < (b) ? (a) : (b))
    typedef long long ll;
    using namespace std;
    template <class T>
    inline void read(T &x){
        int ch = 0, f = 0; x = 0;
        for(; !isdigit(ch); ch = getchar()) if(ch == '-') f = 1;
        for(; isdigit(ch); ch = getchar()) x = x * 10 + ch - 48;
        if(f) x = -x;
    }
    const int N = 1 << 21, mod = 998244353, G = 3;
    int a[N], b[N], A[N], B[N], C[N], d[N], e[N], js[N], inv[N], n, q, q0, qx, qy;
    inline void up(int &x, int y){
    	x = x + y >= mod ? x + y - mod : x + y;
    }
    inline int Pow(int a, int b){
    	int ans = 1;
    	for(; b; b >>= 1, a = 1ll * a * a % mod)
    		if(b & 1) ans = 1ll * ans * a % mod;
    	return ans;
    }
    namespace poly{
    	int rev[N], len, lg;
    	inline void init(int lenth){
    		for(len = 1, lg = 0; len <= lenth; len <<= 1, lg++);
    		for(int i = 0; i < len; i++)
    			rev[i] = (rev[i>>1] >> 1) | ((i & 1) << (lg - 1));
    	}
    	inline void dft(int *a, int sgn){
    		for(int i = 0; i < len; i++)
    			if(i < rev[i]) swap(a[i], a[rev[i]]);
    		for(int k = 2; k <= len; k <<= 1){
    			int w = Pow(G, (mod - 1) / k);
    			if(sgn == -1) w = Pow(w, mod - 2);
    			for(int i = 0, now; i < len; i += k){
    				now = 1;
    				for(int j = i; j < i + (k >> 1); j++){
    					int x = a[j], y = 1ll * a[j+(k>>1)] * now % mod;
    					a[j] = x + y >= mod ? x + y - mod : x + y;
    					a[j+(k>>1)] = x - y < 0 ? x - y + mod : x - y;
    					now = 1ll * now * w % mod;
    				}
    			}
    		}
    		if(sgn == -1){
    			int x = Pow(len, mod - 2);
    			for(int i = 0; i < len; i++) a[i] = 1ll * a[i] * x % mod;
    		}
    	}
    }
    int main(){
    	read(n), read(q), js[0] = 1, inv[0] = 1;
    	for(int i = 1; i <= n; i++){
    		js[i] = 1ll * js[i-1] * i % mod;
    		inv[i] = Pow(js[i], mod - 2);
    	}
    	for(int i = 0; i <= n; i++) read(a[i]);
    	read(q0), read(qx), read(qy);
    	int c = 1ll * qy * Pow(qx - 1, mod - 2) % mod;
    	for(int i = 0, ss = 1; i <= n; i++){
    		A[n-i] = 1ll * js[i] * a[i] % mod;
    		B[i] = 1ll * inv[i] * ss % mod;
    		ss = 1ll * ss * (mod - c) % mod;
    	}
    	poly::init(n + n + 1);
    	poly::dft(A, 1), poly::dft(B, 1);
    	for(int i = 0; i < poly::len; i++)
    		C[i] = 1ll * A[i] * B[i] % mod;
    	poly::dft(C, -1);
    	for(int i = 0; i <= n; i++){
    		b[i] = 1ll * C[n-i] * inv[i] % mod;
    		b[i] = 1ll * b[i] * Pow(c + q0, i) % mod;
    	}
    	memset(A, 0, sizeof(A));
    	memset(B, 0, sizeof(B));
    	memset(C, 0, sizeof(C));
    	d[0] = 1, e[0] = 1;
    	int ix = Pow(qx, mod - 2);
    	for(int i = 1, ss = 1, tt = 1; i <= n + q + 1; i++){
    		ss = 1ll * qx * ss % mod;
    		tt = 1ll * ix * tt % mod;
    		d[i] = 1ll * d[i-1] * ss % mod;
    		e[i] = 1ll * e[i-1] * tt % mod;
    	}
    	for(int i = 0; i <= n; i++)
    		A[i] = 1ll * b[i] * d[i] % mod;
    	for(int i = -n + 1; i <= q; i++)
    		B[i+n] = i <= 0 ? e[abs(i)] : e[i-1];
    	poly::init(n + n + q);
    	poly::dft(A, 1), poly::dft(B, 1);
    	for(int i = 0; i < poly::len; i++)
    		C[i] = 1ll * A[i] * B[i] % mod;
    	poly::dft(C, -1);
    	int ans = 0;
    	for(int i = 1; i <= q; i++)
    		ans ^= 1ll * C[i+n] * d[i-1] % mod;
    	cout << ans << endl;
    	return 0;
    }
    
  • 相关阅读:
    牛客网编程练习(2018校招真题编程题汇总)------字符串价值
    牛客网编程练习(2018校招真题编程题汇总)------排序
    牛客网编程练习(2018校招真题编程题汇总)------回文素数
    牛客网编程练习(2018校招真题编程题汇总)------判断题
    牛客网编程练习(2018校招真题编程题汇总)------删除重复字符
    mysql5.7出现大量too many connections及too many open files错误,且配置最大连接数未生效
    commons-lang3之StringUtils
    commons-lang3 事件机制 <EventListenerSupport>
    springboot文件上传下载简单使用
    Redis5.0.4复制
  • 原文地址:https://www.cnblogs.com/mangoyang/p/12989428.html
Copyright © 2020-2023  润新知