∑MathAdvanced

Generating Functions - EGF

Key Points

•
Exponential generating functions (EGFs) encode a sequence ( $a_{n}$ ) as A(x) = $\sum_{n \geq 0}$ $a_{n}$ $\frac{x ^{n}}{n !}$ , which naturally models labeled combinatorial objects.
•
Multiplying EGFs corresponds to combining labeled structures and leads to a binomial convolution: $c_{n}$ = $\sum_{k = 0}^{n}$ $(k n)$ $a_{k}$ $b_{n - k}$ .
•
Operators on labeled classes map to simple EGF algebra: SET becomes $exp$ , CYCLE becomes - $lo g$ (1-x), and labeled product becomes EGF multiplication.
•
Coefficient extraction uses [ $x^{n}$ ] on the power series but returns $a_{n}$ /n!, so you multiply by n! to recover $a_{n}$ .
•
Fast computation of EGF products reduces to ordinary convolution by scaling with factorials: divide by k! before convolution and multiply by n! after.
•
Classic EGFs: permutations 1/(1-x), derangements $e^{- x}$ /(1-x), Bell numbers $exp$ ( $e^{x}$ - 1), and Stirling numbers via ( $e^{x}$ - 1)^k/k!.
•
For programming contests, use NTT-based polynomial operations to get O(n $lo g$ n) performance for EGF multiplication and Newton iteration for series $exp$ , $lo g$ , and inverse.
•
Common pitfalls include mixing OGFs and EGFs, forgetting the n! factor when extracting coefficients, and violating preconditions like F(0)=0 for $exp$ (F).

Prerequisites

→Modular arithmetic over primes — EGF computations in contests are usually modulo a prime; factorials and series operations require modular inverses.
→Combinatorics (binomial coefficients, permutations) — EGF multiplication uses binomial coefficients and factorial reasoning for labeled structures.
→Ordinary generating functions (OGF) — Understanding the difference from EGFs prevents common mistakes and clarifies when to use which tool.
→Polynomial convolution and NTT — Fast EGF multiplication and FPS operations rely on fast polynomial multiplication.
→Formal power series (FPS) operations — Series inverse, logarithm, and exponential underpin constructions like exp(e^x - 1).
→Symbolic method for combinatorics — Translating combinatorial specifications (SET, CYCLE, PRODUCT) into EGF algebra is the core modeling skill.
→Asymptotic notation — To reason about O(n^2) vs O(n log n) and choose appropriate algorithms.
→Inclusion–exclusion principle — Relates to EGFs for derangements and helps check results against known formulas.

Detailed Explanation

Tap terms for definitions

01Overview

An exponential generating function (EGF) is a power series that packages a sequence (a_n) as A(x) = \sum_{n \ge 0} a_n x^n/n!. The key feature is the n! in the denominator: it is exactly what adjusts counts for labeled structures where relabeling changes the object but should not multiply-count compositions within a construction. In combinatorics on labeled objects (the symbolic method), simple combinatorial operations like disjoint union, labeled product, and set translate into algebraic operations on EGFs: sum, multiplication, and exponential, respectively. This tight correspondence lets us turn combinatorial specifications directly into algebra and then into algorithms for computing coefficients. Extracting coefficients from EGFs requires care: the coefficient of x^n in A(x) is a_n/n!, so to recover a_n you take n! times the coefficient. EGFs excel when counting labeled permutations, set partitions, mappings, and graphs. In algorithmic practice, EGF multiplication induces a binomial transform that can be implemented fast using NTT by scaling with factorials. More advanced operations—series logarithm, exponential, inverse—enable computing deep families like Bell numbers (exp(e^x - 1)) and Stirling numbers ((e^x - 1)^k/k!).

02Intuition & Analogies

Imagine you have n distinct name tags (labels) and you build objects by assigning these labels to components. If two components swap labels, they become different objects. The factorial in an EGF acts like an “automatic label normalizer.” When you combine two labeled components of sizes k and n-k, the number of ways to distribute labels is \binom{n}{k}. EGFs build this directly into multiplication: the coefficient for size n is a binomial mix of smaller sizes. Think of baking: you have n distinct sprinkles and you’re decorating two cupcakes; choosing which sprinkles go to cupcake A is \binom{n}{k}. Summing over all k captures all splits. For sets of components (like a set of cycles in a permutation), order between components doesn’t matter but labels still do; that’s exactly what the exponential does: SET maps to exp in EGFs. For example, permutations are sets of labeled cycles; their EGF becomes exp(sum_{k\ge1} x^k/k) = 1/(1-x). Another intuition: the EGF takes a_n, divides by n! (the number of labelings), and stores that value. When you later combine pieces (multiply series), the combinatorics of label assignments (binomial coefficients) shows up magically via algebra. Finally, to get back raw counts a_n you just multiply by n!. This is why EGFs feel like the “right unit” for labeled counting—they offload the combinatorics of distributing labels to simple algebraic operations on series.

03Formal Definition

Given a sequence (

a_{n}

)_{n

\geq

0}, its exponential generating function is A(x) =

\sum_{n = 0}^{\infty}

a_{n}

\frac{x ^{n}}{n !}

. Let A and B correspond to sequences (

a_{n}

) and (

b_{n}

). Then the Cauchy product of EGFs satisfies C(x) = A(x)B(x) =

\sum_{n \geq 0}

c_{n}

\frac{x ^{n}}{n !}

, where

c_{n}

\sum_{k = 0}^{n}

(k n)

a_{k}

b_{n - k}

. In labeled combinatorial classes, if

C

is the disjoint union of

A

and

B

, then EGF C(x) = A(x)+B(x). If

C

is the labeled product (pairs of an A-object and a B-object on disjoint label sets whose union is [n]), then C(x) = A(x)B(x). If

C

is SET(

A

) (finite sets of A-objects on disjoint labels), then C(x) =

exp

(A(x)). Also, CYCLE maps to -

lo g

(1-A(x)) when cycles are built from A-atoms. Coefficient extraction:

a_{n} = n!

[

x^{n}

] A(x). Differential operators act by shifts: A'(x) =

\sum_{n \geq 0}

a_{n + 1}

\frac{x ^{n}}{n !}

and

\int

A(x)

d x = co n s t

\sum_{n \geq 1}

a_{n - 1}

\frac{x ^{n}}{n !}

, reflecting insertion/removal of a labeled atom. Composition requires A(0)=0 for

exp

(A) and similar analytic constraints for

lo g

, inverse, etc., in formal power series arithmetic.

04When to Use

Use EGFs when counting labeled structures where labels are distinct and matter: permutations (as sets of cycles), derangements (exclude 1-cycles), set partitions (Bell numbers via \exp(e^x - 1)), surjections (k! S(n,k) from (e^x - 1)^k), labeled trees and forests (via Prufer codes and EGF relations), and labeled graphs (edges are labeled pairs of vertices). In competitive programming, reach for EGFs when: 1) a problem asks for a binomial-type convolution c_n = \sum \binom{n}{k} a_k b_{n-k}; 2) you can specify the class using SET/SEQ/CYCLE of labeled atoms; 3) you need Stirling transforms or Bell-type numbers; 4) fast transforms via NTT and Newton iteration are viable within limits (n up to 10^5 or more). EGFs also help derive and solve recurrences involving falling/rising factorial moments or when inclusion–exclusion produces factorial-weighted sums. If composition by exp/log appears naturally, formal series algorithms provide an efficient route to coefficients.

⚠️Common Mistakes

Confusing OGF and EGF: In OGFs the product coefficient is ordinary convolution; in EGFs it is binomial convolution. Always check the denominator x^n vs x^n/n!.
Forgetting n! when extracting coefficients: a_n = n! [x^n] A(x). Leaving out n! yields values scaled by 1/n!.
Violating preconditions for formal series ops: For ln(F) you need F(0)=1; for exp(G) you need G(0)=0; for inverse you need F(0) \neq 0. Enforce these or normalize first.
Mixing numeric domains: Over mod p arithmetic, ensure p is prime and large enough; precompute factorials and modular inverses; avoid dividing by zero (e.g., invfact[0] is 1, but check ranges).
Inefficient convolutions: A naive O(n^2) binomial convolution may TLE for n \gtrsim 5000–10000. Use the factorial-scaling trick plus NTT to reduce to O(n \log n).
Truncation mistakes: Always truncate intermediate products to the target degree to keep complexity under control and ensure correctness of Newton iterations.
Misinterpreting operators: Derivative of an EGF shifts the sequence (a_{n+1}), unlike OGFs where it multiplies by n. Don’t apply OGF rules to EGFs.

Key Formulas

EGF Definition

A (x) = n = 0 \sum \infty a_{n} \frac{x ^{n}}{n !}

Explanation: This defines the exponential generating function of a sequence. Each coefficient $a_{n}$ is stored divided by n! to suit labeled counting.

EGF Product

C (x) = A (x) B (x) = n = 0 \sum \infty (k = 0 \sum n (k n) a_{k} b_{n - k}) \frac{x ^{n}}{n !}

Explanation: Multiplying EGFs corresponds to labeled product. The coefficient $c_{n}$ is a binomial convolution due to splitting n labels into k and n-k.

Coefficient Extraction

a_{n} = n! [x^{n}] A (x)

Explanation: The coefficient of $x^{n}$ in the EGF is $a_{n}$ /n!. Multiply by n! to recover the original sequence term $a_{n}$ .

SET-to-exp

SET (A) \leftrightarrow exp (A (x))

Explanation: The EGF of finite labeled sets of A-objects is the exponential of the EGF of A. Used to count structures made of unordered labeled components.

CYCLE EGF

CYCLE (Z) \leftrightarrow - lo g (1 - x)

Explanation: Cycles of labeled atoms have EGF -ln(1-x). Combining with SET yields permutations as sets of cycles.

Permutations EGF

\frac{1}{1 - x} = exp (- lo g (1 - x))

Explanation: Permutations are sets of cycles. Their EGF sums to 1/(1-x), so the nth coefficient equals n! (the number of permutations).

Derangements EGF

D (x) = \frac{e ^{- x}}{1 - x} = n \geq 0 \sum D_{n} \frac{x ^{n}}{n !}

Explanation: Excluding fixed points removes 1-cycles, multiplying by $e^{- x}$ . Extracting coefficients gives the derangement numbers $D_{n}$ .

Bell EGF

B (x) = exp (e^{x} - 1) = n \geq 0 \sum Bell (n) \frac{x ^{n}}{n !}

Explanation: Set partitions arise as sets of nonempty labeled blocks; the EGF is exp( $e^{x}$ - 1). Multiply by n! to get Bell(n).

Stirling (Second Kind) EGF

n = 0 \sum \infty S (n, k) \frac{x ^{n}}{n !} = \frac{( e ^{x} - 1 ) ^{k}}{k !}

Explanation: Fixing k blocks corresponds to taking k-fold labeled products of nonempty sets, divided by k! to remove ordering of blocks.

Binomial Convolution

c_{n} = k = 0 \sum n (k n) a_{k} b_{n - k} = n! [x^{n}] (A (x) B (x))

Explanation: Restates the product rule and connects it to coefficient extraction. Useful for deriving algorithms.

Factorial Scaling Trick

Let u_{k} = \frac{a _{k}}{k !}, v_{k} = \frac{b _{k}}{k !} . Then w = u * v \Rightarrow c_{n} = n! w_{n}

Explanation: To multiply EGFs fast, divide by k!, perform ordinary convolution w, and multiply by n!. This enables NTT usage.

EGF Calculus

A^{'} (x) = n \geq 0 \sum a_{n + 1} \frac{x ^{n}}{n !}, \int A (x) d x = C + n \geq 1 \sum a_{n - 1} \frac{x ^{n}}{n !}

Explanation: Differentiation shifts the sequence left (to $a_{n + 1}$ ), and integration shifts right, matching labeled combinatorial interpretations.

Complexity Analysis

For direct EGF multiplication using the binomial formula

c_{n}

\sum_{k = 0}^{n}

(k n)

a_{k}

b_{n - k}

, a naive implementation is O(

n^{2}

) time and O(1) extra space per output (or O(n) overall). This is acceptable for n up to a few thousands but will time out for larger n in programming contests. A standard optimization uses the factorial scaling trick: define

u_{k} = a_{k} / k!

and

v_{k} = b_{k} / k!

, compute the ordinary convolution

w = u * v

via NTT in O(n

lo g

n), then recover

c_{n} = n!

w_{n}

. With NTT over a friendly prime (e.g., 998244353), the time becomes O(n

lo g

n) and space O(n), with modest constants. For advanced EGF constructions like Bell numbers B(x) = exp(

e^{x}

- 1), we rely on formal power series operations. Series inverse via Newton iteration takes O(M(n)) per doubling step (M(n) is convolution cost), yielding O(M(n)

lo g

n). Series log uses one derivative, one inverse, one multiplication, and one integral, all O(M(n)). Series exp with F(0)=0 employs Newton iteration

E_{2 m} = E_{m}

(1 - ln

E_{m}

+ F) modulo

x^{2 m}

, giving O(M(n)

lo g

n) overall. With NTT-based convolution M(n)=O(n

lo g

n), so exp/log/inverse are O(n

lo g

n). Memory usage is O(n) for storing polynomials plus O(n) temporaries for NTT buffers; practical implementations allocate 3–6 arrays of size O(n). Always truncate intermediates to the target degree to keep both time and space linearithmic. Precomputing factorials/inverses is O(n) time/space and is essential for safe division (k! and indices) under mod arithmetic.

Code Examples

Quadratic-time EGF multiplication and derangements via e^{-x}/(1-x)

1 #include <bits/stdc++.h>
2 using namespace std;
3 
4 // Modulo for combinatorics (NTT-friendly but we don't use NTT here)
5 static const int MOD = 998244353;
6 
7 long long modpow(long long a, long long e){ long long r=1%MOD; while(e){ if(e&1) r=r*a%MOD; a=a*a%MOD; e>>=1;} return r; }
8 long long inv(long long a){ return modpow(a, MOD-2); }
9 
10 // Precompute factorials and inverse factorials up to N
11 struct Fact {
12     vector<long long> fact, invfact;
13     Fact(int N){
14         fact.resize(N+1); invfact.resize(N+1);
15         fact[0]=1; for(int i=1;i<=N;i++) fact[i]=fact[i-1]*i%MOD;
16         invfact[N]=inv(fact[N]); for(int i=N;i>=1;i--) invfact[i-1]=invfact[i]*i%MOD;
17     }
18 };// end Fact
19 
20 // Quadratic EGF multiplication: given a and b as EGF coefficients (a_n, b_n),
21 // compute c_n = sum_{k=0}^n C(n,k) a_k b_{n-k} for n<=N.
22 vector<long long> multiply_egf_quadratic(const vector<long long>& a, const vector<long long>& b, const Fact& F){
23     int N = (int)max(a.size(), b.size()) - 1;
24     vector<long long> c(N+1, 0);
25     for(int n=0;n<=N;n++){
26         long long sum = 0;
27         for(int k=0;k<=n;k++){
28             long long comb = F.fact[n] * F.invfact[k] % MOD * F.invfact[n-k] % MOD; // C(n,k)
29             long long term = comb * a[k] % MOD * b[n-k] % MOD;
30             sum += term; if(sum>=MOD) sum-=MOD;
31         }
32         c[n] = sum;
33     }
34     return c;
35 }
36 
37 int main(){
38     ios::sync_with_stdio(false); cin.tie(nullptr);
39 
40     int N = 2000; // up to ~2000 OK in O(N^2)
41     Fact F(N);
42 
43     // Derangements: D(x) = e^{-x}/(1-x) = A(x)*B(x)
44     // A: e^{-x} => a_n = (-1)^n
45     // B: 1/(1-x) => sum x^n, so b_n = n! (because sum b_n x^n/n! = sum x^n)
46 
47     vector<long long> a(N+1), b(N+1);
48     for(int n=0;n<=N;n++){
49         a[n] = (n%2==0? 1: MOD-1);                 // (-1)^n mod MOD
50         b[n] = F.fact[n];                           // n!
51     }
52 
53     vector<long long> c = multiply_egf_quadratic(a, b, F); // c_n should be !n
54 
55     // Verify first few derangements !n
56     for(int n=0;n<=10;n++){
57         cout << "n=" << n << " derangements=" << c[n] << "\n";
58     }
59     return 0;
60 }
61

We implement the defining EGF product c_n = \sum \binom{n}{k} a_k b_{n-k} in O(n^2). For derangements, the EGF factorization e^{-x}/(1-x) implies a_k = (-1)^k and b_{n-k} = (n-k)!. The result equals the derangement numbers !n. This example emphasizes the combinatorial meaning of EGF multiplication and shows a straightforward implementation suitable for n up to a few thousand.

Time: O(n^2)Space: O(n)

Fast EGF multiplication via factorial scaling + NTT (derangements up to large N)

1 #include <bits/stdc++.h>
2 using namespace std;
3 static const int MOD = 998244353; // primitive root 3
4 
5 long long modpow(long long a, long long e){ long long r=1%MOD; while(e){ if(e&1) r=r*a%MOD; a=a*a%MOD; e>>=1;} return r; }
6 long long inv(long long a){ return modpow(a, MOD-2); }
7 
8 struct NTT {
9     void ntt(vector<long long>& a, bool invert){
10         int n = (int)a.size();
11         static vector<int> rev; static vector<long long> roots{0,1};
12         if((int)rev.size()!=n){
13             int k = __builtin_ctz(n);
14             rev.assign(n,0);
15             for(int i=0;i<n;i++) rev[i]=(rev[i>>1]>>1)|((i&1)<<(k-1));
16         }
17         if((int)roots.size()<n){
18             int k = __builtin_ctz(roots.size());
19             roots.resize(n);
20             while((1<<k)<n){
21                 long long z = modpow(3, (MOD-1)>>(k+1));
22                 for(int i=1<<(k-1); i<(1<<k); ++i){
23                     roots[2*i] = roots[i];
24                     roots[2*i+1] = roots[i]*z%MOD;
25                 }
26                 ++k;
27             }
28         }
29         for(int i=0;i<n;i++) if(i<rev[i]) swap(a[i], a[rev[i]]);
30         for(int len=1; len<n; len<<=1){
31             for(int i=0;i<n;i+=2*len){
32                 for(int j=0;j<len;j++){
33                     long long u=a[i+j]; long long v=a[i+j+len]*roots[len+j]%MOD;
34                     a[i+j]=u+v; if(a[i+j]>=MOD) a[i+j]-=MOD;
35                     a[i+j+len]=u-v; if(a[i+j+len]<0) a[i+j+len]+=MOD;
36                 }
37             }
38         }
39         if(invert){
40             reverse(a.begin()+1, a.end());
41             long long inv_n = inv(n);
42             for(long long &x:a) x = x*inv_n%MOD;
43         }
44     }
45     vector<long long> multiply(vector<long long> a, vector<long long> b){
46         if(a.empty()||b.empty()) return {};
47         int need = (int)a.size() + (int)b.size() - 1;
48         int n=1; while(n<need) n<<=1;
49         a.resize(n); b.resize(n);
50         ntt(a,false); ntt(b,false);
51         for(int i=0;i<n;i++) a[i]=a[i]*b[i]%MOD;
52         ntt(a,true);
53         a.resize(need);
54         return a;
55     }
56 } ntt;
57 
58 struct Fact { vector<long long> fact, invfact; Fact(int N){ fact.resize(N+1); invfact.resize(N+1); fact[0]=1; for(int i=1;i<=N;i++) fact[i]=fact[i-1]*i%MOD; invfact[N]=inv(fact[N]); for(int i=N;i>=1;i--) invfact[i-1]=invfact[i]*i%MOD; } };
59 
60 // Multiply EGFs a and b up to degree N using factorial scaling and NTT
61 vector<long long> multiply_egf_ntt(const vector<long long>& a, const vector<long long>& b, const Fact& F){
62     int N = (int)max(a.size(), b.size()) - 1;
63     vector<long long> u(N+1), v(N+1);
64     for(int i=0;i<=N;i++){ u[i] = a[i] * F.invfact[i] % MOD; v[i] = b[i] * F.invfact[i] % MOD; }
65     vector<long long> w = ntt.multiply(u, v);
66     w.resize(N+1);
67     for(int i=0;i<=N;i++) w[i] = w[i] * F.fact[i] % MOD; // c_i = i! * w_i
68     return w;
69 }
70 
71 int main(){
72     ios::sync_with_stdio(false); cin.tie(nullptr);
73 
74     int N = 100000; // large N thanks to O(n log n)
75     Fact F(N);
76 
77     // a_n = (-1)^n for e^{-x}
78     // b_n = n! for 1/(1-x)
79     vector<long long> a(N+1), b(N+1);
80     for(int n=0;n<=N;n++){
81         a[n] = (n%2==0? 1: MOD-1);
82         b[n] = F.fact[n];
83     }
84 
85     vector<long long> der = multiply_egf_ntt(a, b, F);
86 
87     // Print first 10 derangements
88     for(int n=0;n<=10;n++) cout << der[n] << (n==10?'\n':' ');
89     return 0;
90 }
91

This example uses the factorial scaling trick to turn EGF multiplication into an ordinary convolution and computes it with NTT. We again compute derangements from e^{-x}/(1-x), but now scale to N=10^5. The NTT reduces the time from O(n^2) to O(n log n), which is crucial for high-constraint contest problems.

Time: O(n log n)Space: O(n)

Bell numbers via EGF B(x)=exp(e^x-1) using FPS exp/log/inv with NTT

1 #include <bits/stdc++.h>
2 using namespace std;
3 static const int MOD = 998244353; // primitive root 3
4 
5 long long modpow(long long a, long long e){ long long r=1%MOD; while(e){ if(e&1) r=r*a%MOD; a=a*a%MOD; e>>=1;} return r; }
6 long long inv(long long a){ return modpow(a, MOD-2); }
7 
8 struct NTT {
9     void ntt(vector<long long>& a, bool invert){
10         int n=(int)a.size();
11         static vector<int> rev; static vector<long long> roots{0,1};
12         if((int)rev.size()!=n){ int k=__builtin_ctz(n); rev.assign(n,0); for(int i=0;i<n;i++) rev[i]=(rev[i>>1]>>1)|((i&1)<<(k-1)); }
13         if((int)roots.size()<n){ int k=__builtin_ctz(roots.size()); roots.resize(n); while((1<<k)<n){ long long z=modpow(3,(MOD-1)>>(k+1)); for(int i=1<<(k-1); i<(1<<k); ++i){ roots[2*i]=roots[i]; roots[2*i+1]=roots[i]*z%MOD; } ++k; } }
14         for(int i=0;i<n;i++) if(i<rev[i]) swap(a[i], a[rev[i]]);
15         for(int len=1; len<n; len<<=1){ for(int i=0;i<n;i+=2*len){ for(int j=0;j<len;j++){ long long u=a[i+j]; long long v=a[i+j+len]*roots[len+j]%MOD; a[i+j]=u+v; if(a[i+j]>=MOD) a[i+j]-=MOD; a[i+j+len]=u-v; if(a[i+j+len]<0) a[i+j+len]+=MOD; } } }
16         if(invert){ reverse(a.begin()+1,a.end()); long long invn=inv(n); for(long long &x:a) x=x*invn%MOD; }
17     }
18     vector<long long> multiply(vector<long long> a, vector<long long> b){
19         if(a.empty()||b.empty()) return {};
20         int need=(int)a.size()+(int)b.size()-1, n=1; while(n<need) n<<=1; a.resize(n); b.resize(n);
21         ntt(a,false); ntt(b,false); for(int i=0;i<n;i++) a[i]=a[i]*b[i]%MOD; ntt(a,true); a.resize(need); return a;
22     }
23 } ntt;
24 
25 struct FPS {
26     vector<long long> fact, invfact, invs; // factorials and modular inverses of integers
27     FPS(int N){
28         fact.resize(N+1); invfact.resize(N+1); invs.resize(N+1);
29         fact[0]=1; for(int i=1;i<=N;i++) fact[i]=fact[i-1]*i%MOD;
30         invfact[N]=inv(fact[N]); for(int i=N;i>=1;i--) invfact[i-1]=invfact[i]*i%MOD;
31         invs[1]=1; for(int i=2;i<=N;i++) invs[i]=MOD - MOD/i * invs[MOD%i] % MOD; invs[0]=0;
32     }
33     // Trim to degree n-1
34     static void trim(vector<long long>& a, int n){ if((int)a.size()>n) a.resize(n); }
35 
36     vector<long long> deriv(const vector<long long>& a){
37         int n=(int)a.size(); if(n==0) return {}; vector<long long> b(max(0,n-1),0);
38         for(int i=1;i<n;i++) b[i-1]=a[i]*i%MOD; return b;
39     }
40     vector<long long> integr(const vector<long long>& a){
41         int n=(int)a.size(); vector<long long> b(n+1,0);
42         for(int i=0;i<n;i++) b[i+1]=a[i]*invs[i+1]%MOD; return b;
43     }
44     vector<long long> multiply(const vector<long long>& a, const vector<long long>& b){ return ntt.multiply(a,b); }
45 
46     // Inverse of series a (a[0] != 0), modulo x^n
47     vector<long long> inv_series(const vector<long long>& a, int n){
48         vector<long long> b(1, ::inv(a[0])); b[0]%=MOD; int m=1;
49         while(m<n){ int m2=min(2*m, n); vector<long long> a_cut(min((int)a.size(), m2)); for(int i=0;i<(int)a_cut.size();++i) a_cut[i]=a[i];
50             vector<long long> ab = multiply(multiply(b,a_cut), vector<long long>{1}); // copy trick
51             trim(ab, m2);
52             // b <- b * (2 - a*b) mod x^{m2}
53             vector<long long> two_minus_ab(m2,0);
54             for(int i=0;i<m2;i++){
55                 long long val = (i<(int)ab.size()? ab[i]:0);
56                 two_minus_ab[i] = ( (i==0?2:0) - val ) % MOD; if(two_minus_ab[i]<0) two_minus_ab[i]+=MOD;
57             }
58             vector<long long> nb = multiply(b, two_minus_ab); trim(nb, m2); b.swap(nb); m=m2; }
59         trim(b,n); return b;
60     }
61 
62     // ln(a) with a[0] = 1, modulo x^n
63     vector<long long> ln_series(const vector<long long>& a, int n){
64         vector<long long> da = deriv(a); vector<long long> ai = inv_series(a, n); vector<long long> prod = multiply(da, ai); trim(prod, n-1); vector<long long> res = integr(prod); trim(res, n); return res;
65     }
66 
67     // exp(f) with f[0] = 0, modulo x^n
68     vector<long long> exp_series(const vector<long long>& f, int n){
69         vector<long long> g(1,1); int m=1; while(m<n){ int m2=min(2*m, n);
70             vector<long long> ln_g = ln_series(g, m2);
71             // delta = f (cut) - ln(g) + 1
72             vector<long long> delta(m2,0);
73             for(int i=0;i<m2;i++){
74                 long long fi = (i<(int)f.size()? f[i]:0);
75                 long long lgi = (i<(int)ln_g.size()? ln_g[i]:0);
76                 long long val = (fi - lgi) % MOD; if(val<0) val+=MOD; if(i==0){ val = (val + 1) % MOD; }
77                 delta[i]=val;
78             }
79             vector<long long> ng = multiply(g, delta); trim(ng, m2); g.swap(ng); m=m2; }
80         trim(g,n); return g;
81     }
82 };
83 
84 int main(){
85     ios::sync_with_stdio(false); cin.tie(nullptr);
86 
87     int N = 5000; // choose up to ~5e4 with optimizations; 5e3 is safe here
88     FPS fps(N+5);
89 
90     // Build A(x) = e^x - 1 (ordinary series with coefficients 1/n! for n>=1)
91     vector<long long> A(N+1, 0);
92     for(int n=1;n<=N;n++) A[n] = fps.invfact[n];
93 
94     // B(x) = exp(A(x)) = sum Bell(n) x^n / n!
95     vector<long long> B = fps.exp_series(A, N+1);
96 
97     // Recover Bell numbers: Bell(n) = n! * [x^n] B(x)
98     for(int n=0;n<=10;n++){
99         long long bell = B[n] * fps.fact[n] % MOD;
100         cout << bell << (n==10?'\n':' ');
101     }
102     return 0;
103 }
104

We compute Bell numbers from their EGF B(x) = exp(e^x - 1). First, we form the series e^x - 1 with coefficients 1/n! for n ≥ 1. Then we compute the series exponential using Newton iteration: E_{2m} = E_m (1 - ln(E_m) + F) modulo x^{2m}. Finally, we multiply coefficients by n! to recover Bell(n). This demonstrates advanced EGF usage with a formal power series toolkit (derivative, integral, inverse, log, exp) over 998244353 and NTT-based convolution.

Time: O(n log n)Space: O(n)

1	#include <bits/stdc++.h>
2	using namespace std;
3
4	// Modulo for combinatorics (NTT-friendly but we don't use NTT here)
5	static const int MOD = 998244353;
6
7	long long modpow(long long a, long long e){ long long r=1%MOD; while(e){ if(e&1) r=ra%MOD; a=aa%MOD; e>>=1;} return r; }
8	long long inv(long long a){ return modpow(a, MOD-2); }
9
10	// Precompute factorials and inverse factorials up to N
11	struct Fact {
12	vector<long long> fact, invfact;
13	Fact(int N){
14	fact.resize(N+1); invfact.resize(N+1);
15	fact[0]=1; for(int i=1;i<=N;i++) fact[i]=fact[i-1]*i%MOD;
16	invfact[N]=inv(fact[N]); for(int i=N;i>=1;i--) invfact[i-1]=invfact[i]*i%MOD;
17	}
18	};// end Fact
19
20	// Quadratic EGF multiplication: given a and b as EGF coefficients (a_n, b_n),
21	// compute c_n = sum_{k=0}^n C(n,k) a_k b_{n-k} for n<=N.
22	vector<long long> multiply_egf_quadratic(const vector<long long>& a, const vector<long long>& b, const Fact& F){
23	int N = (int)max(a.size(), b.size()) - 1;
24	vector<long long> c(N+1, 0);
25	for(int n=0;n<=N;n++){
26	long long sum = 0;
27	for(int k=0;k<=n;k++){
28	long long comb = F.fact[n] * F.invfact[k] % MOD * F.invfact[n-k] % MOD; // C(n,k)
29	long long term = comb * a[k] % MOD * b[n-k] % MOD;
30	sum += term; if(sum>=MOD) sum-=MOD;
31	}
32	c[n] = sum;
33	}
34	return c;
35	}
36
37	int main(){
38	ios::sync_with_stdio(false); cin.tie(nullptr);
39
40	int N = 2000; // up to ~2000 OK in O(N^2)
41	Fact F(N);
42
43	// Derangements: D(x) = e^{-x}/(1-x) = A(x)*B(x)
44	// A: e^{-x} => a_n = (-1)^n
45	// B: 1/(1-x) => sum x^n, so b_n = n! (because sum b_n x^n/n! = sum x^n)
46
47	vector<long long> a(N+1), b(N+1);
48	for(int n=0;n<=N;n++){
49	a[n] = (n%2==0? 1: MOD-1); // (-1)^n mod MOD
50	b[n] = F.fact[n]; // n!
51	}
52
53	vector<long long> c = multiply_egf_quadratic(a, b, F); // c_n should be !n
54
55	// Verify first few derangements !n
56	for(int n=0;n<=10;n++){
57	cout << "n=" << n << " derangements=" << c[n] << "\n";
58	}
59	return 0;
60	}
61

1	#include <bits/stdc++.h>
2	using namespace std;
3	static const int MOD = 998244353; // primitive root 3
4
5	long long modpow(long long a, long long e){ long long r=1%MOD; while(e){ if(e&1) r=ra%MOD; a=aa%MOD; e>>=1;} return r; }
6	long long inv(long long a){ return modpow(a, MOD-2); }
7
8	struct NTT {
9	void ntt(vector<long long>& a, bool invert){
10	int n = (int)a.size();
11	static vector<int> rev; static vector<long long> roots{0,1};
12	if((int)rev.size()!=n){
13	int k = __builtin_ctz(n);
14	rev.assign(n,0);
15	for(int i=0;i<n;i++) rev[i]=(rev[i>>1]>>1)\|((i&1)<<(k-1));
16	}
17	if((int)roots.size()<n){
18	int k = __builtin_ctz(roots.size());
19	roots.resize(n);
20	while((1<<k)<n){
21	long long z = modpow(3, (MOD-1)>>(k+1));
22	for(int i=1<<(k-1); i<(1<<k); ++i){
23	roots[2*i] = roots[i];
24	roots[2i+1] = roots[i]z%MOD;
25	}
26	++k;
27	}
28	}
29	for(int i=0;i<n;i++) if(i<rev[i]) swap(a[i], a[rev[i]]);
30	for(int len=1; len<n; len<<=1){
31	for(int i=0;i<n;i+=2*len){
32	for(int j=0;j<len;j++){
33	long long u=a[i+j]; long long v=a[i+j+len]*roots[len+j]%MOD;
34	a[i+j]=u+v; if(a[i+j]>=MOD) a[i+j]-=MOD;
35	a[i+j+len]=u-v; if(a[i+j+len]<0) a[i+j+len]+=MOD;
36	}
37	}
38	}
39	if(invert){
40	reverse(a.begin()+1, a.end());
41	long long inv_n = inv(n);
42	for(long long &x:a) x = x*inv_n%MOD;
43	}
44	}
45	vector<long long> multiply(vector<long long> a, vector<long long> b){
46	if(a.empty()\|\|b.empty()) return {};
47	int need = (int)a.size() + (int)b.size() - 1;
48	int n=1; while(n<need) n<<=1;
49	a.resize(n); b.resize(n);
50	ntt(a,false); ntt(b,false);
51	for(int i=0;i<n;i++) a[i]=a[i]*b[i]%MOD;
52	ntt(a,true);
53	a.resize(need);
54	return a;
55	}
56	} ntt;
57
58	struct Fact { vector<long long> fact, invfact; Fact(int N){ fact.resize(N+1); invfact.resize(N+1); fact[0]=1; for(int i=1;i<=N;i++) fact[i]=fact[i-1]i%MOD; invfact[N]=inv(fact[N]); for(int i=N;i>=1;i--) invfact[i-1]=invfact[i]i%MOD; } };
59
60	// Multiply EGFs a and b up to degree N using factorial scaling and NTT
61	vector<long long> multiply_egf_ntt(const vector<long long>& a, const vector<long long>& b, const Fact& F){
62	int N = (int)max(a.size(), b.size()) - 1;
63	vector<long long> u(N+1), v(N+1);
64	for(int i=0;i<=N;i++){ u[i] = a[i] * F.invfact[i] % MOD; v[i] = b[i] * F.invfact[i] % MOD; }
65	vector<long long> w = ntt.multiply(u, v);
66	w.resize(N+1);
67	for(int i=0;i<=N;i++) w[i] = w[i] * F.fact[i] % MOD; // c_i = i! * w_i
68	return w;
69	}
70
71	int main(){
72	ios::sync_with_stdio(false); cin.tie(nullptr);
73
74	int N = 100000; // large N thanks to O(n log n)
75	Fact F(N);
76
77	// a_n = (-1)^n for e^{-x}
78	// b_n = n! for 1/(1-x)
79	vector<long long> a(N+1), b(N+1);
80	for(int n=0;n<=N;n++){
81	a[n] = (n%2==0? 1: MOD-1);
82	b[n] = F.fact[n];
83	}
84
85	vector<long long> der = multiply_egf_ntt(a, b, F);
86
87	// Print first 10 derangements
88	for(int n=0;n<=10;n++) cout << der[n] << (n==10?'\n':' ');
89	return 0;
90	}
91