zoukankan      html  css  js  c++  java
  • uoj50【UR#3】链式反应

    • 题解:

      • 令$a(x)$为破坏死光的$EFG$,$f(x)$为方案的$EGF$:$f(x) = x + int   frac{1}{2} f^2(x) a(x)   dt$;
      • 注意到$f(0)=0$,所以考虑如何解:$f'(x) =  frac{1}{2} a(x) f(x)^2 + 1$
      • 设$g(f) = 1 + frac{1}{2}af^2$,即求解$f' = g(f)$;
      • 主要思想是牛顿迭代,假设已经求得$f  equiv f_{0}   pmod {x^{n}} $:
      • 泰勒展开得:$f' equiv g(f_{0}) + g'(f_{0})(f - f_{0}) pmod {x^{2n}} $
      • $f' - g'(f_{0})f equiv g(f_{0}) - g'(f_{0})f_{0} $
      • 进一步令:$r equiv e ^ {int   - g'(f_{0}) dt}$,(注意有$r' = -g'(f_{0}) r$)
      • 两边同时乘以$r$ 得到:
      • $f'r - fr' equiv  (g(f_{0}) - g'(f_{0})f_{0}) r$
      • $(fr)' equiv (g(f_{0}) - g'(f_{0})f_{0}) r$
      • $ f equiv frac{  int   (1 - frac{1}{2}af_{0}^{2})r dt }{r} pmod {x^{2n}}$
      • 就可以求了;
      1 #include<bits/stdc++.h>
      2 #define mod 998244353
      3 #define Run(i,l,r) for(int i=l;i<=r;++i)
      4 #define Don(i,l,r) for(int i=l;i>=r;--i)
      5 const int N=1<<19;
      6 using namespace std;
      7 int n,a[N],fac[N],ans[N];
      8 char s[N],ps[1000000],*pp=ps;
      9 int pw(int x,int y){
     10     if(y<0)y+=mod-1;
     11     int re=1;
     12     for(;y;y>>=1,x=1ll*x*x%mod)if(y&1)re=1ll*re*x%mod;
     13     return re;
     14 }
     15 void push(char x){
     16     if(pp==ps+1000000)fwrite(ps,1,1000000,stdout),pp=ps;
     17     *pp++=x;
     18 }
     19 void write(int x){
     20     static int top,sta[20];
     21     if(!x){push('0'),push('
    ');return;}
     22     while(x)sta[++top]=x%10,x/=10;
     23     while(top)push(sta[top--]^'0');
     24     push('
    ');
     25 }
     26 void flush(){fwrite(ps,1,pp-ps,stdout);}
     27 namespace poly{
     28     int g=3,iv[N],rev[N],L;
     29     void init(int l){for(int i=1;i<=l;++i)iv[i]=pw(i,mod-2);}
     30     void cls(int*A,int l,int r){for(int i=l;i<r;++i)A[i]=0;}
     31     void cpy(int*A,int*B,int l){for(int i=0;i<l;++i)A[i]=B[i];}
     32     void der(int*A,int l){Run(i,0,l-2)A[i]=1ll*(i+1)*A[i+1]%mod;A[l-1]=0;}
     33     void dif(int*A,int l){Don(i,l-1,1)A[i]=1ll*iv[i]*A[i-1]%mod;A[0]=0;}
     34     void ntt(int*A,int l,int f){
     35         for(L=0;(1<<L)<l;++L);
     36         for(int i=0;i<l;++i){
     37             rev[i]=(rev[i>>1]>>1)|((i&1)<<(L-1));
     38             if(i<rev[i])swap(A[i],A[rev[i]]);
     39         }
     40         for(int i=1;i<l;i<<=1){
     41             int wn=pw(g,f*(mod-1)/(i<<1));
     42             for(int j=0;j<l;j+=i<<1){
     43                 int w=1;
     44                 for(int k=0;k<i;++k,w=1ll*w*wn%mod){
     45                     int x=A[j+k],y=1ll*w*A[j+k+i]%mod;
     46                     A[j+k]=(x+y)%mod,A[j+k+i]=(x-y+mod)%mod;
     47                 }
     48             }
     49         }
     50         if(!~f){for(int i=0;i<l;++i)A[i]=1ll*A[i]*iv[l]%mod;}
     51     }
     52     void inv(int*A,int*B,int l){
     53         static int t[N];
     54         if(l==1){B[0]=pw(A[0],mod-2);return;}
     55         int len=l<<1;
     56         inv(A,B,l>>1);cls(B,l,len);
     57         cpy(t,A,l);cls(t,l,len);
     58         ntt(t,len,1);ntt(B,len,1);
     59         for(int i=0;i<len;++i)B[i]=1ll*B[i]*(2-1ll*t[i]*B[i]%mod+mod)%mod;
     60         ntt(B,len,-1);cls(B,l,len);
     61     }
     62     void ln(int*A,int*B,int l){
     63         static int t[N];
     64         int len=l<<1;
     65         inv(A,B,l);
     66         cpy(t,A,l);cls(t,l,len);
     67         der(t,l);
     68         ntt(B,len,1);ntt(t,len,1);
     69         for(int i=0;i<len;++i)B[i]=1ll*B[i]*t[i]%mod;
     70         ntt(B,len,-1);cls(B,l,len);
     71         dif(B,l);
     72     }
     73     void exp(int*A,int*B,int l){
     74         static int t[N];
     75         if(l==1){B[0]=1;return;}
     76         int len=l<<1;
     77         exp(A,B,l>>1);cls(B,l,len);
     78         ln(B,t,l);
     79         for(int i=0;i<l;++i)t[i]=(A[i]-t[i]+mod)%mod;
     80         t[0]++;
     81         ntt(B,len,1);ntt(t,len,1);
     82         for(int i=0;i<len;++i)B[i]=1ll*B[i]*t[i]%mod;
     83         ntt(B,len,-1);cls(B,l,len);
     84     }
     85     void solve(int*A,int l){
     86         static int t[N],r[N];
     87         if(l==1){A[0]=0;return;}
     88         int len=l<<1;
     89         solve(A,l>>1);cls(A,l,len);
     90         cpy(t,a,l);cls(t,l,len);
     91         ntt(A,len,1);ntt(t,len,1);
     92         for(int i=0;i<len;++i){
     93             int tmp=A[i];
     94             A[i]=(mod-1ll*iv[2]*t[i]%mod*A[i]%mod*A[i]%mod)%mod;
     95             t[i]=(mod-1ll*t[i]*tmp%mod)%mod;
     96         }
     97         ntt(A,len,-1);cls(A,l,len);A[0]++;
     98         ntt(t,len,-1);cls(t,l,len);dif(t,l);
     99         exp(t,r,l);inv(r,t,l);
    100         ntt(A,len,1);ntt(r,len,1);
    101         for(int i=0;i<len;++i)A[i]=1ll*A[i]*r[i]%mod;
    102         ntt(A,len,-1);cls(A,l,len);
    103         dif(A,l);
    104         ntt(A,len,1);ntt(t,len,1);
    105         for(int i=0;i<len;++i)A[i]=1ll*A[i]*t[i]%mod;
    106         ntt(A,len,-1);cls(A,l,len);
    107     }
    108 }
    109 int main(){
    110 //    freopen("uoj50.in","r",stdin);
    111 //    freopen("uoj50.out","w",stdout);
    112     scanf("%d%s",&n,s+1);
    113     for(int i=fac[0]=1;i<=n;++i){
    114         fac[i]=1ll*fac[i-1]*i%mod;
    115         if(s[i]-'0')a[i-1]=pw(fac[i-1],mod-2);
    116     }
    117     int len=1;for(;len<=n;len<<=1);
    118     poly::init(len<<1);
    119     poly::solve(ans,len);
    120     for(int i=1;i<=n;++i)write(1ll*ans[i]*fac[i]%mod);
    121     flush(); 
    122     return 0;
    123 }
    View Code
  • 相关阅读:
    Win搭建JAVA环境
    Python JSON存储数据
    XML
    模块5
    模块4
    模块3
    模块2
    模块
    开放封闭原则
    函数续
  • 原文地址:https://www.cnblogs.com/Paul-Guderian/p/10487214.html
Copyright © 2011-2022 走看看