1. 從一道模板題說起多項式對數(shù)函數(shù)ln到底是什么如果你在洛谷、Codeforces或者任何一個算法競賽社區(qū)混跡過一段時間大概率會刷到過“P4725 【模板】多項式對數(shù)函數(shù)多項式 ln”這道題。它就像算法競賽選手在多項式領(lǐng)域的一個“成人禮”標志著從只會做加減乘除的“小學生”進階到開始觸及多項式更深刻運算的“中學生”。但很多人在第一次接觸時都會懵多項式還能取對數(shù)這玩意兒有什么用難道是把1 2x 3x^2丟進計算器按ln鍵嗎顯然不是。這里的“多項式對數(shù)函數(shù)”是一個形式冪級數(shù)上的形式運算。它解決的核心問題是給定一個常數(shù)項為1的多項式或形式冪級數(shù)A(x)求另一個多項式B(x)使得在形式冪級數(shù)的意義下exp(B(x)) A(x)。這里 exp 是指數(shù)函數(shù)。換句話說B(x) 就是 A(x) 的“形式對數(shù)”。這個運算在組合數(shù)學、生成函數(shù)、多項式算法中有著極其重要的地位。比如當你用生成函數(shù)刻畫一個組合結(jié)構(gòu)時對其取 ln 往往對應(yīng)著將連通分量拆解出來在多項式牛頓迭代求解中l(wèi)n 和 exp 是一對關(guān)鍵的基礎(chǔ)算子。這道模板題之所以經(jīng)典是因為它完美地將多項式求導、積分、求逆、乘法這幾個基礎(chǔ)操作串聯(lián)了起來形成了一個完整的算法鏈條。網(wǎng)上能找到的題解和代碼很多但大多只給出了“怎么做”的步驟和代碼對于“為什么這么做”、“每一步背后的數(shù)學原理是什么”、“實現(xiàn)時有哪些一踩就炸的坑”卻語焉不詳。我這篇文章就想結(jié)合我多次實現(xiàn)和調(diào)試的經(jīng)驗把這些隱藏在水面下的東西徹底講透。我們不止要會套模板更要理解這個模板的每一顆螺絲釘是怎么擰上去的。2. 核心公式推導為什么求ln變成了求導、求逆和積分幾乎所有教程都會直接甩給你這個公式 若 A(x) 1 a_1 x a_2 x^2 ...且 A(0)1則ln(A(x)) ∫ [A(x) / A(x)] dx這個公式是整套算法的基石。我們來一步步拆解它看它到底是怎么來的。2.1 從形式微分的定義出發(fā)首先我們得認同對形式冪級數(shù)也可以定義“導數(shù)”。這很直觀對于多項式 A(x) ∑_{i0}^{n} a_i x^i其形式導數(shù) A(x) 就是 ∑_{i1}^{n} i * a_i x^{i-1}。就是把每一項的指數(shù)拿下來當系數(shù)然后指數(shù)減一?，F(xiàn)在考慮我們想求的 B(x) ln(A(x))。這里 ln 是一個形式運算。我們對這個等式兩邊同時關(guān)于 x 求形式導數(shù)利用鏈式法則左邊B(x) 右邊d/dx [ln(A(x))] A(x) / A(x) 這里直接類比了實數(shù)域上 ln(f(x)) 的導數(shù)為 f(x)/f(x)于是我們得到了一個關(guān)鍵等式B(x) A(x) / A(x)。2.2 從微分到積分得到了 B(x) 的導數(shù)那么 B(x) 本身自然就是對其積分B(x) ∫ B(x) dx ∫ [A(x) / A(x)] dx注意這里積分會有一個積分常數(shù) C。因為是不定積分。那么 C 是多少我們需要利用初始條件A(0) 1。我們希望 B(x) 也是一個形式冪級數(shù)并且通常定義 ln(1) 0。所以 B(0) ln(A(0)) ln(1) 0。另一方面我們對 ∫ [A(x) / A(x)] dx 求出的結(jié)果其常數(shù)項就是積分產(chǎn)生的常數(shù) C。為了讓 B(0) 0我們必須令 C 0。所以最終公式里我們直接寫為定積分形式從 0 積到 x或者理解為取不定積分后忽略常數(shù)項因為常數(shù)項為0。在實現(xiàn)時我們做不定積分然后手動將結(jié)果的常數(shù)項設(shè)為0即可。2.3 公式的可行性分析這個公式將 ln 運算轉(zhuǎn)化為了三個我們已知能做的操作求導 (A(x))O(n) 復雜度極其簡單。求逆 (1 / A(x))這里需要計算 A(x) 的乘法逆元。這需要用到多項式求逆算法通常使用牛頓迭代法復雜度 O(n log n)。積分 (∫ ... dx)O(n) 復雜度是求導的逆過程同樣簡單。所以整個多項式 ln 的算法復雜度就卡在了多項式求逆這一步為 O(n log n)。這也就是為什么多項式求逆是多項式全家桶里更基礎(chǔ)的一個模板。注意這個公式成立有一個絕對的前提A(x) 的常數(shù)項必須為 1。為什么 從數(shù)學上看ln(A(x)) 要想展開成形式冪級數(shù)必須在 x0 處有定義且 ln(A(0)) 需要是一個有限值我們?nèi)?。如果 A(0)0ln(0) 無定義如果 A(0) 是其他非1常數(shù) c那么 ln(A(x)) ln(c) ln(1 (A(x)-c)/c)。這里 ln(c) 是一個實數(shù)常數(shù)但我們的多項式是在某個模數(shù)如998244353的有限域上運算的ln(c) 在這個域里可能沒有定義除非 c 是模數(shù)的原根相關(guān)。為了保證運算純粹在模意義下進行且結(jié)果是一個多項式常數(shù)項為0最方便且通用的約定就是要求 A(0)1。這樣 ln(1)0一切都很干凈。3. 算法步驟拆解與零基礎(chǔ)實現(xiàn)指南理解了公式我們來把算法步驟徹底細化。假設(shè)我們有多項式 A(x)其次數(shù)為 n-1通常我們處理長度為 n 的數(shù)組下標 0 到 n-1 對應(yīng)次數(shù) 0 到 n-1 的系數(shù)且滿足 A[0] 1。3.1 第一步計算 A(x) 的導數(shù) A(x)這一步是熱身。設(shè) A(x) a0 a1x a2x^2 ... a_{n-1}x^{n-1}。 那么 A(x) a1 2a2x 3a3*x^2 ... (n-1)*a_{n-1}*x^{n-2}。在代碼中這就是一個簡單的循環(huán)// 假設(shè)系數(shù)存儲在數(shù)組 a 中長度為 n for (int i 1; i n; i) { da[i-1] 1LL * a[i] * i % mod; // da 存儲導數(shù)系數(shù) } // da 的有效長度變?yōu)?n-1注意邊界導數(shù)的次數(shù)比原多項式低一次。3.2 第二步計算 A(x) 的乘法逆元 B(x) 1 / A(x)這是整個算法的核心和性能瓶頸。我們需要求一個多項式 B(x)使得 A(x) * B(x) ≡ 1 (mod x^n)。這里mod x^n的意思是我們只關(guān)心乘積的前 n 項0 到 n-1 次更高次的項可以忽略。求逆通常使用牛頓迭代法。其思想是假設(shè)我們已經(jīng)求出了在模 x^{ceil(m/2)} 意義下的逆元 B_0(x)如何快速得到在模 x^m 意義下的逆元 B(x)推導過程略涉及泰勒展開結(jié)論是迭代公式為B(x) ≡ B_0(x) * (2 - A(x) * B_0(x)) (mod x^m)實際操作時我們采用遞歸或迭代倍增的方式初始條件當 n1 時A(x) 只有一個常數(shù)項 a0。由前提 a0 1所以在模 x^1 意義下其逆元就是 1。假設(shè)我們已經(jīng)求出在模 x^{ceil(n/2)} 意義下的逆元 B_0(x)。目標計算模 x^n 意義下的逆元 B(x)。根據(jù)公式我們需要計算T(x) A(x) * B_0(x) (mod x^n)// 注意這里模數(shù)要提升到 x^nT(x) (2 - T(x)) (mod x^n)// 對 T(x) 的每一項做 2 - t_i 運算B(x) T(x) * B_0(x) (mod x^n)這個過程需要多項式乘法。利用 NTT快速數(shù)論變換可以將乘法優(yōu)化到 O(n log n)。由于這部分是獨立模板代碼較長。其關(guān)鍵點在于每次迭代時對于 A(x) 我們只需要前 m 項當前目標長度對于 B_0(x) 我們知道它在模 x^{m/2} 下是精確的。計算A(x) * B_0(x)時結(jié)果長度會增長但我們只取前 m 項。然后進行2 - T(x)的系數(shù)運算最后再乘一次 B_0(x) 并取前 m 項。3.3 第三步計算 C(x) A(x) * B(x)現(xiàn)在我們有了導數(shù)da長度為 n-1和逆元b長度為 n。將它們相乘。注意da的次數(shù)是 n-2b的次數(shù)是 n-1它們的乘積次數(shù)最高為 (n-2)(n-1)2n-3。但我們最終只需要前 n-1 項因為下一步積分后我們要得到 n 項結(jié)果。所以我們可以只計算到長度至少為 n-1 的卷積。設(shè)dc da * b。我們?nèi)c的前 n-1 項。注意dc[0]對應(yīng)的是A(x)*B(x)的常數(shù)項。3.4 第四步對 C(x) 積分得到最終結(jié)果積分是導數(shù)的逆運算。如果C(x) c0 c1*x c2*x^2 ...那么它的積分∫ C(x) dx C c0*x (c1/2)*x^2 (c2/3)*x^3 ...其中 C 是積分常數(shù)。我們已經(jīng)知道結(jié)果的常數(shù)項必須為 0。所以我們計算 對于 i 從 0 到 n-2res[i1] dc[i] * inv(i1) % mod其中inv(i1)是 i1 在模 mod 下的乘法逆元需要預處理。 而res[0] 0。這樣得到的res就是ln(A(x))的前 n 項系數(shù)。3.5 完整流程圖示與復雜度分析輸入: A(x), 滿足 A[0]1, 次數(shù)界 n 輸出: B(x) ln(A(x)) mod x^n 1. 求導: DA(x) derivative(A(x)) // O(n) 2. 求逆: IA(x) inverse(A(x), n) // O(n log n) 使用牛頓迭代NTT 3. 乘法: C(x) DA(x) * IA(x) mod x^{n-1} // O(n log n) NTT乘法結(jié)果取前n-1項 4. 積分: B(x) integral(C(x)) // O(n) 常數(shù)項設(shè)為0 5. 返回 B(x)總時間復雜度由兩次 O(n log n) 的操作主導即求逆和乘法。空間上需要一些臨時數(shù)組進行變換和計算。4. 實戰(zhàn)代碼剖析從模塊構(gòu)建到邊界處理光說不練假把式。下面我結(jié)合一個典型的基于 NTT模數(shù) 998244353原根為 3的實現(xiàn)來逐塊解析代碼并指出那些容易寫錯、調(diào)試到崩潰的細節(jié)。4.1 基礎(chǔ)工具函數(shù)快速冪與逆元const int mod 998244353, g 3; // 原根 int qpow(int a, int b) { int res 1; while (b) { if (b 1) res 1LL * res * a % mod; a 1LL * a * a % mod; b 1; } return res; }qpow是標準快速冪。inv函數(shù)通常直接調(diào)用qpow(a, mod-2)但頻繁調(diào)用時建議預處理 1~n 的逆元。4.2 核心NTT 與多項式乘法這是所有多項式操作的基礎(chǔ)。代碼較長但結(jié)構(gòu)固定。關(guān)鍵點在于rev數(shù)組的蝴蝶變換以及三層循環(huán)的迭代實現(xiàn)。這里我強調(diào)幾個易錯點長度與界限NTT 要求長度是 2 的冪。每次進行多項式乘法前必須計算lim 1, bit 0; while (lim n m) lim 1, bit;。然后初始化 rev 數(shù)組。結(jié)果清零對于長度為lim的數(shù)組一定要確保lim范圍內(nèi)的數(shù)據(jù)是有效的或者在使用前清空。特別是多次調(diào)用時舊數(shù)據(jù)可能殘留。逆變換后的縮放NTT 逆變換后每一項需要乘以lim的逆元。int invlim qpow(lim, mod-2); for (int i0; ilim; i) a[i]1LL*a[i]*invlim%mod;4.3 多項式求逆的實現(xiàn)細節(jié)這是最難寫對的部分。我給出一個相對清晰的迭代版本框架void poly_inv(int *a, int *b, int n) { // 計算 b(x)使得 a(x)*b(x) ≡ 1 (mod x^n) static int tmp[N]; // 臨時數(shù)組需要足夠大如4倍n b[0] qpow(a[0], mod-2); // 初始條件常數(shù)項逆元 for (int len 2; (len 1) n; len 1) { // 當前目標是求出模 x^len 下的逆元 int lim len 1; // 乘法需要長度 // 將 a 的前 len 項拷貝到 tmp并做 NTT for (int i 0; i len; i) tmp[i] i n ? a[i] : 0; for (int i len; i lim; i) tmp[i] b[i] 0; ntt(tmp, lim, 1); ntt(b, lim, 1); // 根據(jù)公式 B B0 * (2 - A * B0) 計算 for (int i 0; i lim; i) { b[i] 1LL * b[i] * (2 - 1LL * tmp[i] * b[i] % mod mod) % mod; } ntt(b, lim, -1); // 重要b 中 len 之后的項是無意義的必須清零防止影響下一輪 for (int i len; i lim; i) b[i] 0; } // 最后確保 b 只有前 n 項有效后面清零如果傳入的n不是2的冪 for (int i n; i lim; i) b[i] 0; }踩坑實錄1清零清零清零這是多項式題最經(jīng)典的錯誤。在牛頓迭代的每一輪結(jié)束后b數(shù)組在len之后的系數(shù)必須手動設(shè)為0。因為 NTT 逆變換后這些位置可能留有上一輪或計算過程中的垃圾值。下一輪循環(huán)時我們會把整個b數(shù)組包括后面的垃圾值做 NTT這些垃圾值會污染整個頻域?qū)е陆Y(jié)果完全錯誤。這個 bug 非常隱蔽因為小數(shù)據(jù)時可能因為長度不夠碰不到垃圾內(nèi)存而僥幸正確大數(shù)據(jù)一定掛。4.4 多項式 ln 的完整實現(xiàn)集成了求導、求逆、乘法和積分。void poly_derivative(int *a, int *da, int n) { for (int i 1; i n; i) da[i-1] 1LL * a[i] * i % mod; da[n-1] 0; // 導數(shù)長度減一最后一位可置0 } void poly_integral(int *a, int *ia, int n) { ia[0] 0; // 常數(shù)項為0 // 預處理1~n的逆元 inv[i] for (int i 1; i n; i) ia[i] 1LL * a[i-1] * inv[i] % mod; } void poly_ln(int *a, int *res, int n) { // 前提檢查a[0] 必須為 1 assert(a[0] 1); static int da[N], ia[N], tmp[N]; // 1. 求導 poly_derivative(a, da, n); // da 長度為 n-1 // 2. 求逆 poly_inv(a, ia, n); // ia 是 A(x) 的逆長度為 n // 3. 乘法da * ia int lim 1, bit 0; while (lim (n-1) n) lim 1, bit; // 結(jié)果需要前 n-1 項 for (int i 0; i lim; i) tmp[i] (i n-1) ? da[i] : 0; for (int i 0; i lim; i) { // 注意ia 只有前 n 項有效后面在 poly_inv 中已清零 // 但這里為了安全可以只拷貝前 n 項后面置0 if (i n) ia[i] ia[i]; else if (i lim) ia[i] 0; // 確保 ia 在 lim 長度內(nèi)有效 } // 這里需要調(diào)用一個標準的 NTT 乘法函數(shù)輸入 tmp 和 ia結(jié)果存回 tmp ntt_mul(tmp, ia, lim); // 假設(shè)這個函數(shù)處理了 NTT 變換、點乘、逆變換和縮放 // 現(xiàn)在 tmp 的前 n-1 項是 da * ia 的結(jié)果 // 4. 積分 poly_integral(tmp, res, n); // 積分后長度變?yōu)?n }踩坑實錄2長度對齊與數(shù)組越界在調(diào)用poly_inv時我們傳入的長度是n它會計算模x^n的逆。但注意poly_inv內(nèi)部可能按2的冪分配內(nèi)存。我們傳給poly_ln的數(shù)組a其有效長度就是n。但在求導后da有效長度是n-1。進行乘法da * ia時da長度n-1ia長度n卷積結(jié)果長度至少需要(n-1)(n-1)2n-2才能保證前n-1項精確。我們設(shè)置的lim必須大于等于這個值。同時要確保傳入 NTT 乘法的數(shù)組在lim長度內(nèi)都有定義要么是有效系數(shù)要么是0。任何未初始化的值都會導致錯誤。5. 調(diào)試技巧與常見問題排查就算你完全理解了算法第一遍代碼也幾乎不可能一次 AC。以下是我總結(jié)的排查清單5.1 結(jié)果完全不對/隨機數(shù)檢查 NTT 的正確性這是根源。寫一個簡單的測試比如計算(12x) * (13x)看結(jié)果是不是1 5x 6x^2。確保正變換、點乘、逆變換、縮放每一步都正確。檢查數(shù)組清零如 4.3 節(jié)所述在poly_inv的每一輪迭代后必須清零b數(shù)組len之后的部分。在poly_ln中調(diào)用 NTT 前確保tmp和ia在lim范圍內(nèi)的數(shù)據(jù)是干凈的。檢查長度計算lim是否足夠大while (lim n m) lim 1中的n和m是否正確對于da * ian n-1da的長度m nia的長度所以lim需要至少2n-1的下一個2的冪。5.2 結(jié)果前幾項對后面錯檢查求逆的邊界poly_inv函數(shù)是否保證了結(jié)果嚴格只有前n項有效在倍增過程中我們計算的是模x^len的逆但len可能大于n。函數(shù)最后需要把n之后的系數(shù)清零。檢查積分用的逆元表inv[i]是否預處理正確inv[i]是i在模mod下的逆元通常用線性遞推inv[i] mod - 1LL * (mod/i) * inv[mod%i] % mod來求。確保inv[1] 1。5.3 常數(shù)項不為0檢查輸入確認輸入多項式a[0]是否真的為 1。題目可能不保證需要自己先判斷或處理。檢查積分函數(shù)poly_integral是否將res[0]設(shè)為了 05.4 性能問題TLENTT 的蝴蝶變換rev數(shù)組是否預處理每次乘法都重新計算會超時。不必要的拷貝在poly_inv和poly_ln中盡量減少大數(shù)組的memcpy操作。使用指針和就地計算。乘法優(yōu)化對于da * ia我們只需要前n-1項?？梢允褂谩鞍朐诰€卷積”的思路進行優(yōu)化但模板題通常不需要標準的 NTT 乘法即可通過。5.5 一個實用的調(diào)試方法對拍寫一個暴力版本的poly_ln用于小數(shù)據(jù)范圍比如 n 10的驗證。 暴力版本可以模擬形式冪級數(shù)的運算先預處理逆元然后根據(jù)定義ln(A(x)) ∑_{k1} (-1)^{k-1} * (A(x)-1)^k / k。因為 A(x)-1 的常數(shù)項為0所以這個級數(shù)在模 x^n 意義下是有限的只需要算到 kn-1 即可。 用這個暴力程序去驗證你的 NTT 優(yōu)化版本在小數(shù)據(jù)n5,6,7...下的結(jié)果是否一致。這是定位問題最有效的方式。6. 從模板到應(yīng)用ln 在生成函數(shù)中的意義搞懂了實現(xiàn)我們再來聊聊它到底有什么用這樣下次遇到問題你才能想到用它。6.1 組合意義的連接集合與連通分量這是 ln 最經(jīng)典的應(yīng)用。假設(shè)我們有一個組合類 A比如所有的圖其指數(shù)生成函數(shù)EGF為 A(x)。那么A(x)的 expB(x) exp(A(x))通常代表了由 A 中的“連通”對象任意組合而成的“所有”對象。例如A(x) 是連通圖的 EGF那么 B(x) 就是所有圖的 EGF。反過來A(x)的 lnC(x) ln(B(x))就代表了從“所有”對象中提取出“連通”分量。例如已知所有圖的 EGF B(x)那么 ln(B(x)) 就是連通圖的 EGF。在很多計數(shù)問題中我們更容易求出所有方案的生成函數(shù)而想要得到連通方案的生成函數(shù)就需要對其取 ln。6.2 多項式牛頓迭代中的角色牛頓迭代是求解多項式方程F(G(x)) 0的強大工具。例如求exp、求sqrt開根、求復合逆函數(shù)等。 在推導這些迭代式時ln和exp經(jīng)常作為一對互逆的運算出現(xiàn)。例如求G(x) exp(F(x))可以轉(zhuǎn)化為方程ln(G(x)) - F(x) 0然后應(yīng)用牛頓迭代。此時poly_ln就成了迭代過程中必須調(diào)用的子程序。6.3 形式微分與形式積分的工具ln的公式本身完美結(jié)合了求導、求逆和積分。這使得它成為學習多項式形式運算的一個優(yōu)秀案例。掌握了它你就掌握了處理形式冪級數(shù)的一整套基本工具鏈。7. 總結(jié)與擴展思考實現(xiàn)一個poly_ln就像搭樂高。你需要準備好“求導”、“求逆”其內(nèi)部又需要“NTT乘法”、“積分”這幾個基礎(chǔ)模塊然后按照公式∫ (A / A) dx把它們正確地拼接起來。其中多項式求逆是最復雜、最容易出錯的一環(huán)務(wù)必理解其牛頓迭代的倍增思想并牢記迭代后清零的紀律。在競賽中poly_ln很少單獨出題它往往是更大問題的一塊拼圖。比如你需要先對某個生成函數(shù)取 ln進行一些操作再取 exp。因此將它寫對、寫熟封裝成一個可靠的函數(shù)是進軍更高級多項式算法如指數(shù)函數(shù)、三角函數(shù)、快速冪、復合逆的必經(jīng)之路。最后關(guān)于常數(shù)項不為1的情況理論上可以通過提取公因式解決若 A(0) c ≠ 0則 ln(A(x)) ln(c) ln(A(x)/c)。但 ln(c) 在模意義下需要離散對數(shù)來求解這超出了普通多項式模板的范圍。所以模板題和常見應(yīng)用都默認常數(shù)項為1。寫多項式代碼是對耐心和細心的雙重考驗。一個符號的錯誤、一次忘記的清零都可能導致調(diào)試數(shù)小時。但一旦你徹底征服了它那種對復雜算法了如指掌、對每一行代碼都充滿自信的感覺是無與倫比的。希望這篇超詳細的拆解能幫你少走些彎路真正把這塊硬骨頭啃下來。