1. 程式人生 > 其它 >Berlekamp-Massey

Berlekamp-Massey

\(BM\)演算法(\(Berlekamp-Massey\))

給定一個長度為\(n\)的序列,在\(O(n^2)\)的時間內求出序列的最短遞推式(前提是要能通過這\(n\)項求出至少一個遞推式)

\(\forall i>m,a_i=\sum_{j=1}^mf_ja_{i-j}\)

考慮增量法,設當前遞推式為\(f_{cnt}\),check到了\(a_i\)

\(dlt_i=a_i-\sum_{j=1}^{m_{cnt}}f_{cnt,j}a_{i-j}\)

  • \(dlt_i=0\),說明\(f_{cnt}\)合法
  • 否則考慮在\(f_{cnt}\)的基礎上,加一些東西,構造出一個新的\(f\)

\(fail_i\)表示第\(i\)個遞推式失效的位置,對於之前求出的一個遞推式\(f_{prv}\)有(設\(x=fail_{prv}\)):

\(dlt_{x}=a_x-\sum_{j=1}^{m_{prv}}f_{prv,j}a_{x-j}\)

把下標代換一下可得到(設\(i-x=k\)

\(dlt_x=a_{i-k}-\sum_{j=1}^{m_{prv}}f_{prv,j}a_{i-k-j}\)

於是可以求出一個\(f'\)使得\(\sum_{j=1}^{m'}f'_ja_{i-j}=dlt_x\)

但是我們想加上的東西是\(dlt_i\),所以係數乘一下就行

於是就構造出來\(f_{cnt+1}=f_{cnt}+\frac{dlt_i}{dlt_x}f'\)

對應一下下標關係就是加上一個\((0\cdots0,k,-kf_{prv,1},-kf_{prv,2}\cdots)\),前面是\(i-fail_{prv}-1\)\(0\)

想求出最短的遞推式,就取\(i-fail_{prv}+len(prv)\)最小的就行了

for(register int i=1; i<=n; i++){
        // for(register int j=0; j<dp[cnt].size(); j++) printf("%d ",dp[cnt][j]);puts("");
        if(!cnt){
            if(a[i]){
                dlt[i] = a[i];
                fail[cnt++] = i;
                dp[cnt].resize(i,0);
            }
            continue;
        }
        int &sum = dlt[i], num = dp[cnt].size();
        sum = a[i]; fail[cnt] = i;
        for(register int j=0; j<num; j++) sum = dec(sum,1ll*a[i-j-1]*dp[cnt][j]%MOD);
        if(!sum) continue;
        int prv = cnt-1, mnlen = i-fail[prv]+(int)dp[prv].size();
        for(register int j=0; j<cnt; j++)
            if(i-fail[j]+(int)dp[j].size()<mnlen)
                mnlen = i-fail[j]+(int)dp[j].size(), prv = j;
        int k = 1ll*dlt[i]*ksm(dlt[fail[prv]])%MOD;
        dp[cnt+1] = dp[cnt]; cnt++;
        while(dp[cnt].size()<mnlen) dp[cnt].push_back(0);
        iinc(dp[cnt][i-fail[prv]-1],k);
        num = dp[prv].size();
        for(register int j=0; j<num; j++) ddec(dp[cnt][i-fail[prv]+j],1ll*k*dp[prv][j]%MOD);
    }