Trong lĩnh vực khoa học máy tính và toán học ứng dụng, các thuật toán đa thức đóng một vai trò quan trọng. Bài viết này sẽ giới thiệu hai kỹ thuật mạnh mẽ để thực hiện phép nhân đa thức một cách hiệu quả: Biến đổi Fourier nhanh (FFT) và Biến đổi số học nhanh (NTT).
- Biến đổi Fourier nhanh (FFT)
1.1 Biểu diễn đa thức và nhân đa thức
Một đa thức có thể được biểu diễn theo hai cách chính:
Biểu diễn theo hệ số Một đa thức bậc $n$ thường được viết dưới dạng tổng các lũy thừa của $x$ với các hệ số tương ứng: $P(x) = \sum_{i=0}^n a_i x^i$. Trong đó, $a_i$ là các hệ số của đa thức $P(x)$. Để nhân hai đa thức $A(x)$ và $B(x)$ tạo ra $C(x) = A(x) \cdot B(x)$, hệ số $c_k$ của $C(x)$ được tính bằng phép chập: [ c_k = \sum_{i=0}^k a_i b_{k-i} ] Phép tính này yêu cầu độ phức tạp thời gian là $O(n^2)$, với $n$ là bậc của đa thức.
Biểu diễn theo điểm-giá trị Một đa thức bậc $n$ được xác định duy nhất bởi $n+1$ cặp điểm-giá trị $(x_0, y_0), (x_1, y_1), \dots, (x_n, y_n)$, trong đó $y_j = P(x_j)$. Ưu điểm lớn của dạng biểu diễn này là phép nhân đa thức trở nên cực kỳ đơn giản. Nếu ta có $A(x)$ và $B(x)$ dưới dạng điểm-giá trị (với cùng các điểm $x_i$), thì $C(x_i) = A(x_i) \cdot B(x_i)$. Để nhân hai đa thức bậc $n$, ta cần $2n+1$ điểm (vì bậc của đa thức kết quả có thể lên tới $2n$). Sau khi có các giá trị điểm của $A(x)$ và $B(x)$, phép nhân chỉ mất $O(n)$ thời gian. Như vậy, thách thức chính là làm thế nào để chuyển đổi giữa hai dạng biểu diễn này một cách hiệu quả?
- Thuật toán chuyển đổi từ biểu diễn hệ số sang biểu diễn điểm-giá trị được gọi là Biến đổi Fourier rời rạc (DFT).
- Thuật toán ngược lại, từ biểu diễn điểm-giá trị sang biểu diễn hệ số, được gọi là Biến đổi Fourier rời rạc nghịch đảo (IDFT). Quy trình nhân hai đa thức bằng FFT có thể hình dung như sau:
- Hệ số sang điểm-giá trị: Áp dụng DFT cho $A(x)$ và $B(x)$ để chuyển chúng sang dạng điểm-giá trị.
- Nhân điểm-giá trị: Nhân các giá trị tương ứng của $A(x)$ và $B(x)$ để thu được các điểm-giá trị của $C(x)$.
- Điểm-giá trị sang hệ số: Áp dụng IDFT cho các điểm-giá trị của $C(x)$ để chuyển nó trở lại dạng hệ số.
1.2 Căn đơn vị và tính chất
FFT sử dụng các số phức đặc biệt gọi là căn đơn vị (roots of unity).
- Số phức: Một số phức có dạng $z = a + bi$, trong đó $a, b$ là số thực và $i^2 = -1$.
- Căn đơn vị bậc $N$: Là các nghiệm của phương trình $z^N = 1$. Trong mặt phẳng phức, các căn đơn vị bậc $N$ nằm trên đường tròn đơn vị (bán kính 1) và chia đường tròn thành $N$ phần bằng nhau. Căn đơn vị thứ $k$ của bậc $N$ được ký hiệu là $\omega_N^k$ và có thể biểu diễn dưới dạng Euler: $\omega_N^k = e^{i \frac{2\pi k}{N}} = \cos\left(\frac{2\pi k}{N}\right) + i \sin\left(\frac{2\pi k}{N}\right)$.
Các căn đơn vị có một số tính chất quan trọng:
- Tính chất nhân: $\omega_N^j \cdot \omega_N^i = \omega_N^{i+j}$
- Tính chất giảm bậc: $\omega_{2N}^{2k} = \omega_N^k$. (Điều này có thể thấy rõ từ công thức Euler: $e^{i \frac{2\pi (2k)}{2N}} = e^{i \frac{2\pi k}{N}}$).
- Tính chất chu kỳ: $\omega_N^{k+N} = \omega_N^k$. (Tương đương với việc quay thêm một vòng tròn đầy đủ).
- Tính chất đối xứng: $\omega_N^{k + N/2} = -\omega_N^k$. (Tương đương với việc quay nửa vòng tròn).
1.3 Thuật toán FFT (Chia để trị)
Để tính $P(\omega_N^k)$ cho tất cả $k = 0, \dots, N-1$ một cách hiệu quả, FFT sử dụng chiến lược chia để trị. Giả sử ta có đa thức $P(x) = \sum_{i=0}^{N-1} a_i x^i$, trong đó $N$ là lũy thừa của 2. Ta có thể chia $P(x)$ thành hai đa thức con, một chứa các hệ số ở chỉ số chẵn và một chứa các hệ số ở chỉ số lẻ: [ P_{\text{chẵn}}(x) = a_0 + a_2x + a_4x^2 + \dots ] [ P_{\text{lẻ}}(x) = a_1 + a_3x + a_5x^2 + \dots ] Từ đó, $P(x)$ có thể được viết lại như sau: [ P(x) = P_{\text{chẵn}}(x^2) + x \cdot P_{\text{lẻ}}(x^2) ] Bây giờ, ta sẽ thay các căn đơn vị vào biểu thức này:
- Trường hợp 1: $k < N/2$ $P(\omega_N^k) = P_{\text{chẵn}}((\omega_N^k)^2) + \omega_N^k \cdot P_{\text{lẻ}}((\omega_N^k)^2)$ Sử dụng tính chất $\omega_N^{2k} = \omega_{N/2}^k$, ta có: $P(\omega_N^k) = P_{\text{chẵn}}(\omega_{N/2}^k) + \omega_N^k \cdot P_{\text{lẻ}}(\omega_{N/2}^k)$
- Trường hợp 2: $k \ge N/2$ (đặt $k' = k - N/2$, vậy $k' < N/2$) $P(\omega_N^{k'+N/2}) = P_{\text{chẵn}}((\omega_N^{k'+N/2})^2) + \omega_N^{k'+N/2} \cdot P_{\text{lẻ}}((\omega_N^{k'+N/2})^2)$ Sử dụng tính chất $\omega_N^{2(k'+N/2)} = \omega_N^{2k'+N} = \omega_N^{2k'} = \omega_{N/2}^{k'}$ và $\omega_N^{k'+N/2} = -\omega_N^{k'}$, ta có: $P(\omega_N^{k'+N/2}) = P_{\text{chẵn}}(\omega_{N/2}^{k'}) - \omega_N^{k'} \cdot P_{\text{lẻ}}(\omega_{N/2}^{k'})$ Như vậy, $P(\omega_N^k)$ và $P(\omega_N^{k+N/2})$ có thể được tính từ các giá trị của $P_{\text{chẵn}}(\omega_{N/2}^k)$ và $P_{\text{lẻ}}(\omega_{N/2}^k)$ chỉ trong $O(1)$ thời gian sau khi các phép gọi đệ quy đã hoàn tất. Điều này dẫn đến độ phức tạp $O(N \log N)$ cho DFT. Đối với IDFT, ta chỉ cần thay thế $\omega_N^k$ bằng $\omega_N^{-k}$ (tức là $\cos(-x) + i \sin(-x)$) và chia kết quả cuối cùng cho $N$.
1.4 Triển khai FFT
Phiên bản đệ quy của FFT có thể gây tốn kém do chi phí gọi hàm và sao chép mảng. Một phương pháp triển khai hiệu quả hơn là FFT dạng lặp, sử dụng kỹ thuật "đảo bit" (bit-reversal permutation) để sắp xếp lại các phần tử đầu vào. Sau đó, các phép tính được thực hiện từ dưới lên (từ các nhóm nhỏ nhất đến nhóm lớn nhất).
#include <iostream>
#include <vector>
#include <complex> // Thư viện cho số phức
#include <cmath> // Thư viện cho các hàm toán học như cos, sin, acos
#include <algorithm> // Thư viện cho hàm swap
#include <iomanip> // Để điều chỉnh định dạng đầu ra số thực
// Sử dụng namespace std để tránh viết std:: nhiều lần
using namespace std;
// Hằng số PI
const double PI = acos(-1.0);
// Hàm đọc số nguyên nhanh. Có thể dùng cin/scanf thông thường nếu không cần tối ưu I/O.
inline int docSoNguyen() {
int giaTri = 0, dau = 1;
char kyTu = getchar();
while (kyTu < '0' || kyTu > '9') {
if (kyTu == '-') dau = -1;
kyTu = getchar();
}
while (kyTu >= '0' && kyTu <= '9') {
giaTri = (giaTri << 3) + (giaTri << 1) + kyTu - '0'; // Tương đương giaTri * 10 + kyTu - '0'
kyTu = getchar();
}
return giaTri * dau;
}
// Định nghĩa kiểu dữ liệu số phức
using SoPhuc = complex<double>;
// Mảng dùng để lưu thứ tự đảo bit của các chỉ số
vector<int> mangDaoBit;
int soBitDaoBit; // Số bit cần thiết để biểu diễn chỉ số lớn nhất
// Hàm tính và lưu trữ mảng đảo bit
void tinhMangDaoBit(int N) {
mangDaoBit.resize(N);
soBitDaoBit = 0;
while ((1 << soBitDaoBit) < N) soBitDaoBit++; // Tìm số bit nhỏ nhất b sao cho 2^b >= N
// Tính thứ tự đảo bit cho từng chỉ số
for (int i = 0; i < N; ++i) {
mangDaoBit[i] = 0;
for (int j = 0; j < soBitDaoBit; ++j) {
// Nếu bit thứ j của i là 1, đặt bit thứ (soBitDaoBit - 1 - j) của mangDaoBit[i] là 1
if ((i >> j) & 1) {
mangDaoBit[i] |= (1 << (soBitDaoBit - 1 - j));
}
}
}
}
// Hàm thực hiện Biến đổi Fourier nhanh (FFT)
// N: Kích thước của mảng đầu vào (phải là lũy thừa của 2)
// mangHeSo: Vector chứa các hệ số của đa thức (hoặc các giá trị điểm)
// cheDoBienDoi: 1 cho DFT (Biến đổi Fourier xuôi), -1 cho IDFT (Biến đổi Fourier ngược)
void thucHienFFT(int N, vector<SoPhuc>& mangHeSo, int cheDoBienDoi) {
// Sắp xếp lại các phần tử của mảng theo thứ tự đảo bit
for (int i = 0; i < N; ++i) {
if (i < mangDaoBit[i]) { // Chỉ swap nếu chỉ số hiện tại nhỏ hơn chỉ số đảo bit của nó để tránh swap 2 lần
swap(mangHeSo[i], mangHeSo[mangDaoBit[i]]);
}
}
// Vòng lặp chính của FFT: kết hợp các nhóm nhỏ hơn thành các nhóm lớn hơn
// 'doDaiNhom' đại diện cho kích thước hiện tại của các nhóm đang được kết hợp
for (int doDaiNhom = 1; doDaiNhom < N; doDaiNhom <<= 1) { // doDaiNhom tăng theo lũy thừa của 2: 1, 2, 4, ..., N/2
// canDonViGoc là căn đơn vị cơ sở cho cấp độ kết hợp hiện tại
// Nó là e^(i * 2*PI / (2*doDaiNhom)) hoặc e^(-i * 2*PI / (2*doDaiNhom)) tùy thuộc vào cheDoBienDoi
SoPhuc canDonViGoc(cos(PI / doDaiNhom), sin(PI / doDaiNhom) * cheDoBienDoi);
// Duyệt qua các khối (block) trong mảng
// Mỗi khối có kích thước 2 * doDaiNhom
for (int i = 0; i < N; i += 2 * doDaiNhom) {
SoPhuc canDonViHienTai(1, 0); // canDonViHienTai bắt đầu từ 1 và nhân với canDonViGoc trong mỗi bước
// Duyệt qua các cặp phần tử trong khối hiện tại
// doDaiNhom là nửa kích thước của khối, tương ứng với số phép kết hợp trong mỗi khối
for (int j = 0; j < doDaiNhom; ++j) {
// 'x' là phần tử của nửa đầu khối, 'y' là phần tử của nửa sau khối
SoPhuc x = mangHeSo[i + j];
SoPhuc y = canDonViHienTai * mangHeSo[i + j + doDaiNhom];
// Áp dụng phép biến đổi "bướm" (butterfly operation)
mangHeSo[i + j] = x + y;
mangHeSo[i + j + doDaiNhom] = x - y;
canDonViHienTai *= canDonViGoc; // Cập nhật căn đơn vị cho cặp tiếp theo
}
}
}
}
int main() {
// Tối ưu hóa I/O, đặc biệt hữu ích trong các cuộc thi lập trình
ios_base::sync_with_stdio(false);
cin.tie(NULL);
int bacDaThucA, bacDaThucB; // Bậc của hai đa thức A và B
bacDaThucA = docSoNguyen();
bacDaThucB = docSoNguyen();
// Khởi tạo vector số phức cho các hệ số của đa thức
// Kích thước ban đầu là bac + 1 vì đa thức bậc bac có bac + 1 hệ số
vector<SoPhuc> daThucA(bacDaThucA + 1);
vector<SoPhuc> daThucB(bacDaThucB + 1);
// Đọc các hệ số của đa thức A và lưu vào phần thực của số phức
for (int i = 0; i <= bacDaThucA; ++i) daThucA[i].real(docSoNguyen());
// Đọc các hệ số của đa thức B và lưu vào phần thực của số phức
for (int i = 0; i <= bacDaThucB; ++i) daThucB[i].real(docSoNguyen());
// Xác định kích thước N cho FFT. N phải là lũy thừa của 2 và lớn hơn hoặc bằng (bacA + bacB + 1).
// Bậc của đa thức kết quả là bacA + bacB, nên cần ít nhất bacA + bacB + 1 điểm.
int kichThuocFFT = 1;
while (kichThuocFFT <= bacDaThucA + bacDaThucB) {
kichThuocFFT <<= 1; // Nhân 2 cho đến khi đủ lớn
}
// Thay đổi kích thước vector để phù hợp với kichThuocFFT, các phần tử mới sẽ được khởi tạo mặc định là 0
daThucA.resize(kichThuocFFT);
daThucB.resize(kichThuocFFT);
// Tính mảng đảo bit một lần
tinhMangDaoBit(kichThuocFFT);
// Thực hiện DFT (Biến đổi Fourier xuôi) cho cả hai đa thức
thucHienFFT(kichThuocFFT, daThucA, 1);
thucHienFFT(kichThuocFFT, daThucB, 1);
// Nhân các giá trị điểm-giá trị tương ứng của hai đa thức
// Đây là bước nhân đa thức trong miền tần số (point-value form)
for (int i = 0; i < kichThuocFFT; ++i) {
daThucA[i] = daThucA[i] * daThucB[i];
}
// Thực hiện IDFT (Biến đổi Fourier ngược) để chuyển kết quả về dạng hệ số
thucHienFFT(kichThuocFFT, daThucA, -1);
// In các hệ số của đa thức kết quả
// Các giá trị sau IDFT cần được chia cho kichThuocFFT và làm tròn về số nguyên gần nhất
for (int i = 0; i <= bacDaThucA + bacDaThucB; ++i) {
// Sử dụng floor(x + 0.5) hoặc round(x) để làm tròn số thực.
// std::fixed và std::setprecision là để xử lý lỗi làm tròn nhỏ của số thực.
cout << (int)(floor(daThucA[i].real() / kichThuocFFT + 0.5)) << " ";
}
cout << endl; // In một ký tự xuống dòng ở cuối
return 0; // Kết thúc chương trình thành công
}
- Biến đổi số học nhanh (NTT)
2.1 Vấn đề về độ chính xác và Giới thiệu NTT
FFT sử dụng các số phức và các phép tính dấu phẩy động, điều này có thể dẫn đến các lỗi làm tròn nhỏ, tích lũy và ảnh hưởng đến độ chính xác của kết quả, đặc biệt với các đa thức có hệ số lớn. Để khắc phục hạn chế này, Biến đổi số học nhanh (NTT) ra đời. NTT thực hiện các phép biến đổi tương tự FFT nhưng hoàn toàn trong số học modulo. Thay vì sử dụng các căn đơn vị phức, NTT sử dụng các căn nguyên thủy (primitive roots) modulo một số nguyên tố $p$. Điều này đảm bảo rằng tất cả các phép tính đều là số nguyên và chính xác.
2.2 Căn nguyên thủy và tính chất
- Căn nguyên thủy modulo $p$: Một số nguyên $g$ được gọi là căn nguyên thủy của số nguyên tố $p$ nếu các lũy thừa $g^1, g^2, \dots, g^{p-1}$ tạo thành một hoán vị của $1, 2, \dots, p-1$ khi lấy modulo $p$. (Nói cách khác, $g$ là một phần tử sinh của nhóm nhân $(\mathbb{Z}/p\mathbb{Z})^*$).
- Thay thế căn đơn vị: Trong NTT, chúng ta định nghĩa "căn đơn vị" $\omega_N$ tương tự như trong FFT, nhưng trong môi trường modulo. Cụ thể, $\omega_N = g^{\frac{p-1}{N}} \pmod p$.
Để NTT hoạt động, các "căn đơn vị" này phải có các tính chất tương tự như trong FFT. Điều này đúng với các căn nguyên thủy:
- $\omega_N^N \equiv (g^{\frac{p-1}{N}})^N \equiv g^{p-1} \equiv 1 \pmod p$ (theo Định lý Fermat nhỏ).
- $\omega_{2N}^{2k} \equiv (g^{\frac{p-1}{2N}})^{2k} \equiv (g^{\frac{p-1}{N}})^k \equiv \omega_N^k \pmod p$.
- $\omega_N^{k+N} \equiv (g^{\frac{p-1}{N}})^{k+N} \equiv (g^{\frac{p-1}{N}})^k \cdot (g^{\frac{p-1}{N}})^N \equiv \omega_N^k \cdot g^{p-1} \equiv \omega_N^k \cdot 1 \equiv \omega_N^k \pmod p$.
- $\omega_N^{k + N/2} \equiv (g^{\frac{p-1}{N}})^{k+N/2} \equiv (g^{\frac{p-1}{N}})^k \cdot (g^{\frac{p-1}{N}})^{N/2} \equiv \omega_N^k \cdot g^{\frac{p-1}{2}} \pmod p$. Từ định lý Fermat nhỏ, $g^{p-1} \equiv 1 \pmod p$, suy ra $g^{\frac{p-1}{2}} \equiv \pm 1 \pmod p$. Vì $g$ là căn nguyên thủy, $g^{\frac{p-1}{2}}$ không thể là $1 \pmod p$ (nếu không thì bậc của $g$ sẽ nhỏ hơn $p-1$), vậy $g^{\frac{p-1}{2}} \equiv -1 \pmod p$. Do đó, $\omega_N^{k+N/2} \equiv -\omega_N^k \pmod p$. Để $\omega_N$ được định nghĩa và có đủ $N$ giá trị khác nhau, $N$ phải là ước của $p-1$. Trong các ứng dụng FFT/NTT, $N$ thường là lũy thừa của 2. Do đó, ta cần chọn một số nguyên tố $p$ sao cho $p-1$ có một thừa số là lũy thừa lớn của 2.
Các modulo NTT phổ biến:
- $p = 998244353$. Đây là một số nguyên tố có dạng $119 \cdot 2^{23} + 1$. Căn nguyên thủy của nó là $g=3$.
- $p = 1004535809$. Đây là một số nguyên tố có dạng $479 \cdot 2^{21} + 1$. Căn nguyên thủy của nó là $g=3$.
2.3 Triển khai NTT
Cấu trúc của NTT rất giống FFT, chỉ khác ở chỗ các phép cộng, trừ, nhân được thực hiện modulo $p$, và căn đơn vị được thay bằng lũy thừa của căn nguyên thủy. Phép chia cho $N$ trong IDFT được thay bằng phép nhân với nghịch đảo modulo của $N$ (tức là $N^{p-2} \pmod p$ theo định lý Fermat nhỏ).
#include <iostream>
#include <vector>
#include <algorithm> // Thư viện cho hàm swap
// Sử dụng namespace std để tránh viết std:: nhiều lần
using namespace std;
// Hằng số modulo và căn nguyên thủy cho NTT
// MOD = 998244353 là một số nguyên tố có dạng C * 2^K + 1, với K=23
const int MOD = 998244353;
// CAN_NGUYEN_THUY = 3 là một căn nguyên thủy của MOD
const int CAN_NGUYEN_THUY = 3;
// Hàm đọc số nguyên nhanh. Có thể dùng cin/scanf thông thường.
inline int docSoNguyen() {
int giaTri = 0, dau = 1;
char kyTu = getchar();
while (kyTu < '0' || kyTu > '9') {
if (kyTu == '-') dau = -1;
kyTu = getchar();
}
while (kyTu >= '0' && kyTu <= '9') {
giaTri = (giaTri << 3) + (giaTri << 1) + kyTu - '0';
kyTu = getchar();
}
return giaTri * dau;
}
// Hàm lũy thừa theo modulo (base^exp % MOD)
long long luyThuaModulo(long long coSo, long long soMu) {
long long ketQua = 1;
coSo %= MOD; // Đảm bảo cơ số nằm trong khoảng [0, MOD-1]
while (soMu > 0) {
if (soMu % 2 == 1) ketQua = (ketQua * coSo) % MOD;
coSo = (coSo * coSo) % MOD;
soMu /= 2;
}
return ketQua;
}
// Hàm nghịch đảo modulo (tính x^(MOD-2) % MOD)
// Dùng để thực hiện phép chia A/B bằng cách nhân A * B^(MOD-2) % MOD
long long nghichDaoModulo(long long n) {
return luyThuaModulo(n, MOD - 2);
}
// Mảng dùng để lưu thứ tự đảo bit của các chỉ số
vector<int> mangDaoBitNTT;
int soBitDaoBitNTT; // Số bit cần thiết để biểu diễn chỉ số lớn nhất
// Hàm tính và lưu trữ mảng đảo bit
void tinhMangDaoBitNTT(int N) {
mangDaoBitNTT.resize(N);
soBitDaoBitNTT = 0;
while ((1 << soBitDaoBitNTT) < N) soBitDaoBitNTT++;
for (int i = 0; i < N; ++i) {
mangDaoBitNTT[i] = 0;
for (int j = 0; j < soBitDaoBitNTT; ++j) {
if ((i >> j) & 1) {
mangDaoBitNTT[i] |= (1 << (soBitDaoBitNTT - 1 - j));
}
}
}
}
// Hàm thực hiện Biến đổi số học nhanh (NTT)
// N: Kích thước của mảng đầu vào (phải là lũy thừa của 2)
// mangHeSo: Vector chứa các hệ số của đa thức
// cheDoBienDoi: 1 cho DFT (Biến đổi xuôi), 0 cho IDFT (Biến đổi ngược)
void thucHienNTT(int N, vector<long long>& mangHeSo, int cheDoBienDoi) {
// Sắp xếp lại các phần tử của mảng theo thứ tự đảo bit
for (int i = 0; i < N; ++i) {
if (i < mangDaoBitNTT[i]) {
swap(mangHeSo[i], mangHeSo[mangDaoBitNTT[i]]);
}
}
// Lấy nghịch đảo của căn nguyên thủy nếu đang thực hiện IDFT
long long canNguyenThuyDao = nghichDaoModulo(CAN_NGUYEN_THUY);
// Vòng lặp chính của NTT
for (int doDaiNhom = 1; doDaiNhom < N; doDaiNhom <<= 1) {
// canDonViGoc là căn nguyên thủy cho cấp độ kết hợp hiện tại
// Nếu cheDoBienDoi = 1, dùng CAN_NGUYEN_THUY. Nếu cheDoBienDoi = 0, dùng canNguyenThuyDao
long long canDonViGoc = luyThuaModulo(cheDoBienDoi ? CAN_NGUYEN_THUY : canNguyenThuyDao, (MOD - 1) / (doDaiNhom * 2));
// Duyệt qua các khối trong mảng
for (int i = 0; i < N; i += 2 * doDaiNhom) {
long long canDonViHienTai = 1; // canDonViHienTai bắt đầu từ 1
// Duyệt qua các cặp phần tử trong khối hiện tại
for (int j = 0; j < doDaiNhom; ++j) {
long long x = mangHeSo[i + j];
long long y = (canDonViHienTai * mangHeSo[i + j + doDaiNhom]) % MOD;
// Áp dụng phép biến đổi "bướm" (butterfly operation)
mangHeSo[i + j] = (x + y) % MOD;
mangHeSo[i + j + doDaiNhom] = (x - y + MOD) % MOD; // Đảm bảo kết quả không âm
canDonViHienTai = (canDonViHienTai * canDonViGoc) % MOD; // Cập nhật căn đơn vị
}
}
}
// Nếu là IDFT, chia tất cả các hệ số cho N (bằng cách nhân với nghịch đảo modulo của N)
if (cheDoBienDoi == 0) {
long long invN = nghichDaoModulo(N);
for (int i = 0; i < N; ++i) {
mangHeSo[i] = (mangHeSo[i] * invN) % MOD;
}
}
}
int main() {
// Tối ưu hóa I/O
ios_base::sync_with_stdio(false);
cin.tie(NULL);
int bacDaThucA, bacDaThucB; // Bậc của hai đa thức A và B
bacDaThucA = docSoNguyen();
bacDaThucB = docSoNguyen();
vector<long long> heSoDaThucA(bacDaThucA + 1);
vector<long long> heSoDaThucB(bacDaThucB + 1);
// Đọc các hệ số của đa thức A
for (int i = 0; i <= bacDaThucA; ++i) heSoDaThucA[i] = docSoNguyen();
// Đọc các hệ số của đa thức B
for (int i = 0; i <= bacDaThucB; ++i) heSoDaThucB[i] = docSoNguyen();
// Xác định kích thước N cho NTT. N phải là lũy thừa của 2 và đủ lớn.
int kichThuocNTT = 1;
while (kichThuocNTT <= bacDaThucA + bacDaThucB) {
kichThuocNTT <<= 1;
}
// Thay đổi kích thước vector để phù hợp với kichThuocNTT
heSoDaThucA.resize(kichThuocNTT);
heSoDaThucB.resize(kichThuocNTT);
// Tính mảng đảo bit một lần
tinhMangDaoBitNTT(kichThuocNTT);
// Thực hiện NTT (DFT) cho cả hai đa thức
thucHienNTT(kichThuocNTT, heSoDaThucA, 1);
thucHienNTT(kichThuocNTT, heSoDaThucB, 1);
// Nhân các giá trị điểm-giá trị (modulo)
for (int i = 0; i < kichThuocNTT; ++i) {
heSoDaThucA[i] = (heSoDaThucA[i] * heSoDaThucB[i]) % MOD;
}
// Thực hiện NTT ngược (IDFT) để lấy lại hệ số
thucHienNTT(kichThuocNTT, heSoDaThucA, 0);
// In kết quả
for (int i = 0; i <= bacDaThucA + bacDaThucB; ++i) {
cout << heSoDaThucA[i] << " ";
}
cout << endl;
return 0;
}