MIT-两个多项式相乘

问题描述

给定两个次数为 $n$ 的多项式：

$$\begin{array}{c}A(x)=a_0+a_{1}x+\cdots +a_n{x}^n\\ B(x)=b_0+b_{1}x+\cdots +b_n{x}^n\end{array}$$

求它们的积 $A(x)\cdot B(x)$。

---

算法实现

暴力算法

我们设 $A(x)={\sum }_{i=0}^n{a}_i{x}^i$ ，$B(x)={\sum }_{i=0}^n{b}_i{x}^i$ 和 $C(x)={\sum }_{k=0}^{2n}{c}_i{x}^i=A(x)B(x)$。其中 $c_k={\sum }_{i=0}^k{a}_i{b}_{k-i}$ ，$0\le k\le 2n$。

int* Brute_PolynomialMultiply(int a[], int b[], int n)
{int m = 2 * n;int c[m];for (int k = 0; k <= m; k ++ ){for (int i = 0; i <= k; i ++ )c[k] += a[i] * b[k - i];}return c;
}

时间复杂度 ：$O(n^2)$

---

分治算法

考虑到多项式乘法问题的结构性质，我们可以采用分治算法。该算法的核心思想是将大问题分解为多个较小的子问题，利用分治策略递归地解决这些子问题，从而减少计算量。

在多项式乘法中，我们可以将每个多项式分解为两部分，并递归地计算这些部分的乘积。具体来说，我们定义
$$\begin{array}{rcl}{A}_0(x)& =& a_0+a_{1}x+\cdots +a_{\lfloor \frac{n}{2}\rfloor -1}{x}^{\lfloor \frac{n}{2}\rfloor -1},\\ A_1(x)& =& a_{\lfloor \frac{n}{2}\rfloor }+a_{\lfloor \frac{n}{2}\rfloor +1}x+\cdots +a_n{x}^{n-\lfloor \frac{n}{2}\rfloor }.\end{array}$$

多项式$A(x)$ 可写成 $A(x)=A_0(x)+x^{n/2}{A}_1(x)$ 的形式 。

同样的，我们也可以像上述一样定义 $B_0(x)$ 和 $B_1(x)$，多项式 $B(x)$ 可写成下面的形式：

$$B(x)=B_0(x)+B_1(x)x^{\lfloor \frac{n}{2}\rfloor }$$

所以，我们可以通过以下方式计算它们的乘积：

$$A(x)B(x)=A_0(x)B_0(x)+A_0(x)B_1(x)x^{\lfloor \frac{n}{2}\rfloor }+A_1(x)B_0(x)x^{\lfloor \frac{n}{2}\rfloor }+A_1(x)B_1(x)x^{2\lfloor \frac{n}{2}\rfloor }$$

这样就将原来问题规模为 $n$ 的问题分割成 $4$ 个 $\frac{n}{2}$ 的子问题了，分别为：
$$\begin{array}{cc}U(x)=A_0(x)B_0(x),& V(x)=A_0(x)B_1(x),\\ W(x)=A_1(x)B_0(x),& Z(x)=A_1(x)B_1(x)\end{array}$$

// 多项式相加
vector<int> addPoly(const vector<int>& A, const vector<int>& B) {int n = max(A.size(), B.size());vector<int> result(n, 0);for (int i = 0; i < n; i++) {if (i < A.size()) result[i] += A[i];if (i < B.size()) result[i] += B[i];}return result;
}vector<int> DC_multiplyPoly(const vector<int>& A, const vector<int>& B) {int n = A.size();// 基本情况if (A.empty() || B.empty()) return {};if (n == 1) return {A[0] * B[0]};int mid = n / 2;// 分割多项式vector<int> A0(A.begin(), A.begin() + mid);  // [0, mid - 1)vector<int> A1(A.begin() + mid, A.end()); //vector<int> B0(B.begin(), B.begin() + mid);vector<int> B1(B.begin() + mid, B.end());// 递归计算vector<int> U = DC_multiplyPoly(A0, B0);vector<int> V = DC_multiplyPoly(A0, B1);vector<int> W = DC_multiplyPoly(A1, B0);vector<int> Z = DC_multiplyPoly(A1, B1);// 假设P0 = V(x) + W(x)vector<int> P0 = addPoly(V, W);// 合并结果vector<int> result(2 * n, 0);for (int i = 0; i < U.size(); i ++ ) result[i] += U[i];for (int i = 0; i < P0.size(); i ++ ) result[i + mid] += P0[i];for (int i = 0; i < Z.size(); i ++ ) result[i + 2 * mid] += Z[i];return result;
}

时间复杂度：$O(n^2)$

假设 $n=2^h$ 且 $T(1)=O(1)$，$T(n)$ 是执行一次两个多项式均为 $n$ 次的DC_multiplyPoly所需的时间。

• vector<int> U = DC_multiplyPoly(A0, B0)：$T(\frac{n}{2})$
• vector<int> V = DC_multiplyPoly(A0, B1)：$T(\frac{n}{2})$
• vector<int> W = DC_multiplyPoly(A1, B0)：$T(\frac{n}{2})$
• vector<int> Z = DC_multiplyPoly(A1, B1)：$T(\frac{n}{2})$
• vector<int> P0 = addPoly(V, W)：$O(n)$
• …
• for (int i = 0; i < Z.size(); i ++ ) result[i + 2 * mid] += Z[i]：$O(n)$
• vector<int> A0(A.begin(), A.begin() + mid)：$O(1)$
• …
• vector<int> B1(B.begin() + mid, B.end())：$O(1)$

$$\begin{array}{rl}{T}_n& =4T(\frac{n}{2})+cn+O(n)\\ & =4\left[4T\left(\frac{n}{2^2}\right)+c\frac{n}{2}+O(\frac{n}{2})\right]+cn+O(n)\\ & =4^{2}T\left(\frac{n}{2^2}\right)+(1+2)cn+O(n)\\ & =4^2\left[4T\left(\frac{n}{2^3}\right)+c\frac{n}{2^2}+O(\frac{n}{2^2})\right]+(1+2)cn+O(n)\\ & =4^{3}T\left(\frac{n}{2^3}\right)+(1+2+2^2)cn+O(n)\\ & =...\\ & =4^{i}T\left(\frac{n}{2^i}\right)+\sum _{j=0}^{i-1}{2}^{j}cn+O(n)\\ & =...\\ & =4^{h}T\left(\frac{n}{2^h}\right)+\sum _{j=0}^{h-1}{2}^{j}cn+O(n)\\ & =(2^h{)}^{2}T(1)+\frac{1\cdot (1-2^h)}{1-2}cn+O(n)\\ & =n^{2}T(1)+cn(n-1)+O(n)\\ & =O(n^2)+c{n}^2-cn+O(n)\\ & =O(n^2)\end{array}$$

---

Karatsuba（优化的分治算法）

由上述的基础分治算法可以知道 $A(x)B(x)$ 多项式乘积的形式为：

$$A(x)B(x)=A_0(x)B_0(x)+[A_0(x)B_1(x)+A_1(x)B_0(x)]x^{\frac{n}{2}}+A_1(x)B_1(x)x^{2\lfloor \frac{n}{2}\rfloor }$$

我们在基础的分治算法的基础上，重新定义：

$$\begin{array}{rl}Y(x)& =[A_0(x)+A_1(x)]\times [B_0(x)+B_1(x)]\\ U(x)& =A_0(x)B_0(x)\\ Z(x)& =A_1(x)B_1(x)\end{array}$$

故有

$$A_0(x)B_1(x)+A_1(x)B_0(x)=Y(x)-U(x)-Z(x)$$

即最后 $A(x)B(x)$ 的多项式为：

$$A(x)B(x)=U(x)+[Y(x)-U(x)-Z(x)]x^{\lfloor \frac{n}{2}\rfloor }+Z(x)\times x^{2\lfloor \frac{n}{2}\rfloor }$$

我们用了一个小的数学技巧，将原本需要分割成 4 个 $\frac{n}{2}$ 的子问题，变成了 3 个 $\frac{n}{2}$ 个子问题。

// 多项式相减 A - B
vector<int> subPoly(const vector<int>& A, const vector<int>& B) {int n = max(A.size(), B.size());vector<int> result(n, 0);for (int i = 0; i < n; i++) {if (i < A.size()) result[i] += A[i];if (i < B.size()) result[i] -= B[i];}return result;
}vector<int> Karatsuba_multiplyPoly(const vector<int>& A, const vector<int>& B) {int n = A.size();// 基本情况if (A.empty() || B.empty()) return {};if (n == 1) return {A[0] * B[0]};int mid = n / 2;// 分割多项式vector<int> A0(A.begin(), A.begin() + mid);  // [0, mid - 1)vector<int> A1(A.begin() + mid, A.end()); //vector<int> B0(B.begin(), B.begin() + mid);vector<int> B1(B.begin() + mid, B.end());// 递归计算 修改点在这里！！！！vector<int> U = Karatsuba_multiplyPoly(A0, B0);vector<int> Z = Karatsuba_multiplyPoly(A1, B1);vector<int> Y = Karatsuba_multiplyPoly(addPoly(A0, A1), addPoly(B0, B1));// 假设P0 = Y(x) - U(x) - Z(x)P0 = subPoly(subPoly(Y, U), Z);// 合并结果vector<int> result(2 * n, 0);for (int i = 0; i < U.size(); i ++ ) result[i] += U[i];for (int i = 0; i < P0.size(); i ++ ) result[i + mid] += P0[i];for (int i = 0; i < Z.size(); i ++ ) result[i + 2 * mid] += Z[i];return result;
}

时间复杂度：$O(n^{{\log }_{2}3})\approx O(n^{1.585})$

假设 $n=2^h$ 且 $T(1)=O(1)$，$T(n)$ 是执行一次两个多项式均为 $n$ 次的Karatsuba_multiplyPoly所需的时间。

• vector<int> U = DC_multiplyPoly(A0, B0)：$T(\frac{n}{2})$
• vector<int> Z = DC_multiplyPoly(A1, B1)：$T(\frac{n}{2})$
• vector<int> Y = Karatsuba_multiplyPoly(addPoly(A0, A1), addPoly(B0, B1))：$T(\frac{n}{2})$
• addPoly(V, W)：$O(n)$
• subPoly(V, W)：$O(n)$
• …
• for (int i = 0; i < Z.size(); i ++ ) result[i + 2 * mid] += Z[i]：$O(n)$
• vector<int> A0(A.begin(), A.begin() + mid)：$O(1)$
• …
• vector<int> B1(B.begin() + mid, B.end())：$O(1)$

$$\begin{array}{rl}{T}_n& =3T(\frac{n}{2})+cn+O(n)\\ & =3\left[3T\left(\frac{n}{2^2}\right)+c\frac{n}{2}+O(\frac{n}{2})\right]+cn+O(n)\\ & =3^{2}T\left(\frac{n}{2^2}\right)+(1+\frac{3}{2})cn+O(n)\\ & =3^2\left[3T\left(\frac{n}{2^3}\right)+c\frac{n}{2^2}+O(\frac{n}{2^2})\right]+(1+\frac{3}{2})cn+O(n)\\ & =3^{3}T\left(\frac{n}{2^3}\right)+(1+\frac{3}{2}+[\frac{3}{2}{]}^2)cn+O(n)\\ & =...\\ & =3^{h}T\left(\frac{n}{2^h}\right)+\sum _{j=0}^{h-1}[\frac{3}{2}{]}^{j}cn+O(n)\\ & =3^{h}T(1)+\frac{c}{2}(3^h-n)+O(n)\\ & =(2^{{\log }_{2}3}{)}^{h}T(1)+\frac{c}{2}((2^{{\log }_{2}3}{)}^h-n)+O(n)\\ & =2^{h{\log }_{2}3}T(1)+\frac{c}{2}(2^{h{\log }_{2}3}-n)+O(n)\\ & =(2^h{)}^{{\log }_{2}3}T(1)+\frac{c}{2}((2^h{)}^{{\log }_{2}3}-n)+O(n)\\ & =O(n^{{\log }_{2}3})+\frac{c}{2}(n^{{\log }_{2}3}-n)+O(n)\\ & =O(n^{{\log }_{2}3})\end{array}$$