線性遞推數列

線性遞推

對于無限數列 $\{a_0,a_1,...\}$ 和有限非空數列 $\{r_0,r_1,...,r_{m?1}\}$ 。

若對于任意 $m?1\le n$ ，有 ${\sum }_{i=0}^{m?1}{a}_{n?i}{r}_i=0$ ，則稱數列 $r$ 為數列 $a$ 的線性遞歸式。

若 $r_0=1$，稱數列 $r$ 為數列 $a$ 的線性遞推式。

稱存在線性遞推式的無限數列為線性遞推數列。

$r_0=1$ 說明可以寫成 $a_n=?(a_{n?1}?r_1+a_{n?2}?r_2+...+a_{n?m+1}?r_{m?1})$。

對于有限數列 $\{a_0,a_1,...a_{n?1}\}$ 和有限非空數列 $\{r_0,r_1,...,r_{m?1}\}$ 。

若對于任意 $m?1\le p\le n$ ，有 ${\sum }_{i=0}^{m?1}{a}_{p?i}{r}_i=0$ ，則稱數列 $r$ 為數列 $a$ 的線性遞歸式。

若 $r_0=1$，稱數列 $r$ 為數列 $a$ 的線性遞推式。

稱這個線性遞推式的階數為它的長度 $?1$。（因為下標是從 $0$ 開始的）

稱數列 $a$ 階數最小的線性遞推式為數列 $a$ 的最短線性遞推式。

生成函數

對于有限數列 $\{a_0,...,a_{n?1}\}$，定義它的生成函數為多項式 $A(x)={\sum }_{i=0}^{n?1}{a}_i{x}^i$。

對于無限數列 $\{a_0,a_1,...\}$，定義它的生成函數為形式冪級數 $A(x)={\sum }_0^{\infty }{a}_i{x}^i$。

對于無限數列 $\{a_0,a_1,...\}$ 和有限非空數列 $\{r_0,r_1,...,r_{m?1}\}$ 。

設數列 $a$ 的生成函數為 $A$，數列 $r$ 的生成函數為 $R$。

數列 $r$ 為數列 $a$ 的線性遞歸式等價于存在次數不超過 $m?2$ 次的多項式 $C$，滿足 $AR+C=0$。

對于有限數列 $\{a_0,a_1,...,a_{n?1}\}$ 和有限非空數列 $\{r_0,r_1,...,r_{m?1}\}$ 。

設數列 $a$ 的生成函數為 $A$，數列 $r$ 的生成函數為 $R$。

數列 $r$ 為數列 $a$ 的線性遞歸式等價于存在次數不超過 $m?2$ 次的多項式 $C$，滿足 $AR+C\equiv 0(\mathrm{m}\mathrm{o}\mathrm{d}{x}^n)$。

Berlekamp-Massey

注意：上面的都是下標從 $0$ 開始，這里的算法流程下標是從 $1$ 開始的。

Berlekamp-Massey 算法是用來在 $O(n^2)$ 的時間里求一個數列的最短線性遞推式。

方法思想是增量法。

重申明確一下問題模型：給定 $n$ 個元素的數列 $a_1,...a_n$，求一個最短的數列 $r_1,...r_m$。要求滿足 ${?}_{m<i\le n}{a}_i={\sum }_{j=1}^m{a}_{i?j}{r}_j$。要求在 $O(n^2)$ 時間內解決。

注意1：這個遞推的下標可以看作是一種卷積形式，相加和為一個“定值”，這有利于你快速理解下面的一些下標變換與對應。

注意2：這個條件只需要 $>m$ 的下標做到，換言之如果被當成了 $r$ 中的一個“基”，是不需要滿足這個式子的，因為會有下標為負的情況。

注意3：$a_i=a[i]=a(i)$ 只是表達方式不同。

假設遞推式已經經過了 $t$ 次更新，第 $i$ 次更新的遞推式記為 $r(i)$（這是一個多項式），長度為 $m(i)$。初始時，定義 $r(0)=?,m(0)=0$。

一個一個考慮加入數列 $a$ 中元素，現在即將加入 $a_i$。假設現在遞推式長度為 $m=m(t)$。

設 $\Delta (i)=a_i?{\sum }_{j=1}^{m}a[i?j]?r(i)[j]$。

• $\Delta (i)=0$，證明目前遞推式仍然滿足條件，合法不用修改，繼續考慮加入下一個。

• $\Delta (i)=0$。設 $\text{fail}(t)=i$，表示 $r(t)$ 遞推式第一次失效的位置在 $i$。

  • $t=0$。

    意味著 $a_i$ 前面的數都為 $0$。那么不管遞推數列怎么配求和都是 $0$。顯然只能將 $i$ 這個下標也被遞推式包含。

    所以新的遞推式直接就是 $r(t+1)=\{\underset{i}{0,0,...,0}\}$。

    這樣就不需要 $a_i$ 滿足上面求和的式子了，它已經作為遞推式的一個“基”了。

  • $t=0$。

    考慮構造一個遞推式 $R$，滿足 ${?}_{|R|<k<i}{\sum }_{j=1}^{|R|}{a}_{k?j}{R}_j=0$，${\sum }_{j=1}^{|R|}{a}_{i?j}{R}_j=\Delta (i)$。

    尋找之前某次失效的遞推式，$0\le p<t$，顯然這個遞推式的失效位置為 $\text{fail}(p)$。

    同時，設 $\omega =\frac{\Delta (i)}{\Delta (\text{fail}(p))}$，則 $R=\{\underset{i?\text{fail}(p)?1}{0,0,\dots ,0},\omega ,?\omega ?r(p)[1],?\omega ?r(p)[2],...,?\omega ?r(p)[|r(p)|]\}$。

    令 $r(t+1)=r(t)+R$ 即可。（加法遵循多項式加法原則，即對位系數相加）

    ---

    中場休息——請確保上面的下標之間關系理清楚了，再繼續往下看正確性

    ---

    為什么這樣是對的？理論解釋不清楚，不妨眼見為實，我們直接帶進去算。

    • 對于 $a_i$。

      $R[1,i?\text{fail}(p)?1]$ 都為 $0$，對應 $a[i?1,\text{fail}(p)+1]$ 相乘和為 $0$。

      真正乘起來有值的是從 $R[i?\text{fail}(p)]?a[\text{fail}(p)]$ 項開始的。

      并且后面對應相乘，$r(p)[1]?a[\text{fail}(p)?1],r(p)[2]?a[\text{fail(p)-2}]......$

      等等！！！

      這相乘求和不就是第 $p$ 個遞推式失效的位置對應的表達式嗎？轉化一下即 $=a_{\text{fail(p)}}?\Delta (fail(p))$。

      $?\omega ?a_{\text{fail(p)}}?\omega ?(a_{\text{fail(p)}}?\Delta (\text{fail}(p)\left)\right)=\omega ?\Delta (\text{fail}(p))=\Delta (i)$。

      發現 $R$ 把 $r(t)$ 失效的差距部分恰好補上了。

      由此可知 $r(t+1)$ 對 $a_i$ 成立。

    • 對于 $m\le k<i$。其實推導與上面差不多。

      $R[1,i?\text{fail}(p)?1]$ 都為 $0$，對應 $a[k?1,k?i+\text{fail}(p)+1]$ 相乘和為 $0$。

      真正乘起來有值的是從 $R[i?\text{fail}(p)]?a[k?i+\text{fail}(p)]$ 項開始的。

      并且后面對應相乘，$r(p)[1]?a[k?i+\text{fail}(p)?1],r(p)[2]?a[k?i+\text{fail(p)-2}]......$

      轉化一下即 $=a_{k?i+\text{fail}(p)}?\Delta (k?i+\text{fail}(p))$。

      因為 $k<i$，所以 $k?i+\text{fail}(p)<\text{fail}(p)$。

      還記得 $\text{fail}(p)$ 的定義嗎？——是 $r(p)$ 第一個失效的位置。換言之，在此位置之前 $r(p)$ 都是成立的。

      所以 $\Delta (k?i+\text{fail}(p))=0$。

      $?\omega ?a_{k?i+\text{fail(p)}}?\omega ?(a_{k?i+\text{fail(p)}}?\Delta (k?i+\text{fail}(p)\left)\right)=\omega ?\Delta (k?i+\text{fail}(p))=0$。

    綜上，我們構造的這個 $R$ 本質上只起到了對 $a_i$ 補充 $\Delta (i)$ 的效果，對于其余 $k$ 的貢獻都是 $0$。

    這是利用了 $r(p)$ 在 $1～\text{fail}(p)?1$ 都滿足關系，而在 $\text{fail}(p)$ 相差 $\Delta (p)$ 的性質。

    ---

    此外我們還希望遞推式的長度越短越好，也就是說 $\max (m(t),i?\text{fail}(p)+m(p))$ 最短。

    貪心地只需要動態維護最短的 $i?\text{fail}(p)+m(p)$，每次算出 $r(t+1)$ 時都與之前的 $r(p)$ 比一下誰更短即可。

一共遞推 $O(n)$ 次，每次修改 $O(n)$ 次，時間復雜度為 $O(n^2)$。

代碼實現中，$r$ 用的是 $\text{vector}$，下標從 $0$ 開始，所以有些下標會與上面的推導略有差異。且不保證一定正確！！！

//std參考code #include <bits/stdc++.h> using namespace std; #define mod 1000000007 #define int long long #define maxn 1005 int n, cnt; int a[maxn], fail[maxn], delta[maxn]; vector < int > r[maxn];int qkpow( int x, int y ) {int ans = 1;while( y ) {if( y & 1 ) ans = ans * x % mod;x = x * x % mod;y >>= 1;}return ans; }signed main() {scanf( "%lld", &n );for( int i = 1;i <= n;i ++ ) scanf( "%lld", &a[i] );for( int i = 1;i <= n;i ++ )if( ! cnt ) {if( a[i] ) {fail[cnt ++] = i;delta[i] = a[i];r[cnt].resize( i, 0 );}continue;}else {fail[cnt] = i;delta[i] = a[i];for( int j = 0;j < r[cnt].size();j ++ ) ( delta[i] -= a[i - j - 1] * r[cnt][j] ) %= mod;if( ! delta[i] ) continue;// int len = 0x7f7f7f7f, p;int len = i - fail[cnt - 1] + r[cnt - 1].size(), p = cnt - 1;for( int j = 0;j < cnt;j ++ )if( i - fail[j] + r[j].size() < len )len = i - fail[j] + r[j].size(), p = j;int omega = delta[i] * qkpow( delta[fail[p]], mod - 2 ) % mod;r[cnt + 1] = r[cnt];cnt ++;while( r[cnt].size() < len ) r[cnt].push_back( 0 );( r[cnt][i - fail[p] - 1] += omega ) %= mod;for( int j = 0;j < r[p].size();j ++ )( r[cnt][i - fail[p] + j] -= omega * r[p][j] ) %= mod; }printf( "%lld\n", r[cnt].size() );for( int i : r[cnt] ) printf( "%lld ", ( i + mod ) % mod );return 0; } /* Input 1 7 1 2 4 9 20 40 90 Output 1 4 0 0 10 0Input 2 18 2 4 8 16 32 64 128 256 512 2 4 8 16 32 64 128 256 512 Output 2 0 0 0 0 0 0 0 0 1 */

BM算法和線性遞推的應用

當我們確信這個題是 矩陣加速 / 線性遞推 的題且矩陣很難構造。

那么可以 暴力遞推/打表 出超過兩倍矩陣階數的長度，然后 BM+線性遞推 得到這個序列的任意項。

不知道矩陣階數就 能找多長,找多長 了。

某些特殊情況下也可以代替 DFA的最小化 。

打表題也可以試著找一下規律，因為只能找線性遞推的規律，所以大部分情況要靈活處理，比如說先差分再找規律之類的。

創作挑戰賽新人創作獎勵來咯，堅持創作打卡瓜分現金大獎

總結

以上是生活随笔為你收集整理的[线性代数学习笔记] 线性递推数列及 Berlekamp-Massey 算法的详细推导过程的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。