[论文解读] A Quadratic-Time Algorithm for General Multivariate Polynomial Interpolation

Michael Hecht, Bevan L. Cheeseman+2 more|arXiv (Cornell University)|Oct 30, 2017

Polynomial and algebraic computation参考文献 17被引用 3

一句话总结

本文提出 PIP-SOLVER，一种用于多元多项式插值的新型递归算法，无需求逆范德蒙德矩阵即可在 O(N²) 时间内计算系数，相比基于矩阵求逆的方法具有更优的数值稳定性和性能。该方法可高效求解大规模插值问题，保证精度并降低计算成本。

ABSTRACT

For $m,n \in \mathbb{N}$, $m\geq 1$ and a given function $f : \mathbb{R}^m\longrightarrow \mathbb{R}$ the polynomial interpolation problem (PIP) is to determine a \emph{generic node set} $P \subseteq \mathbb{R}^m$ and the coefficients of the uniquely defined polynomial $Q\in\mathbb{R}[x_1,\dots,x_m]$ in $m$ variables of degree $\mathrm{deg}(Q)\leq n \in \mathbb{N}$ that fits $f$ on $P$, i.e., $Q(p) = f(p)$, $\forall\, p \in P$. We here show that in general, i.e., for arbitrary $m,n \in \mathbb{N}$, $m \geq 1$, there exists an algorithm that determines $P$ and computes the $N(\mbox{m,n})=\#P$ coefficients of $Q$ in $\mathcal{O}\big(N(\mbox{m,n})^2\big)$ time using $\mathcal{O}\big(\mbox{m}N(\mbox{m,n})\big)$ storage, without inverting the occurring Vandermonde matrix. We provide such an algorithm, termed PIP-SOLVER, based on a recursive decomposition of the problem and prove its correctness. Since the present approach solves the PIP without matrix inversion, it is computationally more efficient and numerically more robust than previous approaches. We demonstrate this in numerical experiments and compare with previous approaches based on matrix inversion and linear systems solving.

研究动机与目标

• 开发一种无需依赖矩阵求逆的、数值稳定且计算高效的通用多元多项式插值（PIP）算法。
• 克服基于求解线性系统或求逆大型范德蒙德矩阵的现有方法所面临的高计算成本和数值不稳定性问题。
• 提供一种递归分解框架，实现在 O(N²) 时间和 O(mN) 存储内求解 PIP，其中 N = N(m,n) 为单项式的数量。
• 证明该算法即使在大规模问题下仍能保持高精度，且在运行时间和数值鲁棒性方面优于传统方法。
• 为未来扩展至其他基函数（如切比雪夫多项式或傅里叶基）以及分布式或并行计算环境奠定基础。

提出的方法

• 该算法通过沿超平面将多元插值问题递归分解为低维子问题，实现子解的独立计算。
• 构建一个通用节点集 P ⊆ ℝ^m，使得对所有 p ∈ P 满足插值条件 Q(p) = f(p)，且存在唯一总次数 ≤ n 的多项式 Q。
• 该方法完全避免了对 N×N 多元范德蒙德矩阵 V_{m,n} 的显式求逆，从而规避了经典方法中的计算瓶颈。
• 通过将插值投影到余维为 1 的超平面 H_ε 上，递归求解子问题，利用单项式基的结构和递归系数恢复机制。
• 通过避免病态矩阵运算，转而依赖递归且数值行为良好的计算，确保数值稳定性。
• 针对小规模子问题（N ≤ 5…8×10⁴）提出一种混合策略，结合递归方法与直接矩阵求逆，以在大规模场景下优化性能。

实验结果

研究问题

• RQ1 能否在不求逆范德蒙德矩阵的前提下，以 O(N²) 时间求解一般多元多项式插值问题？
• RQ2 沿超平面的插值问题递归分解如何实现数值稳定且高效的求解器？
• RQ3 如何最优地分解 PIP 以最小化条件数并最大化子问题中的数值精度？
• RQ4 该算法能否推广至任意基函数（如切比雪夫或傅里叶多项式）？
• RQ5 与基于矩阵求逆的求解器相比，该递归方法在问题规模增大时的性能和精度如何变化？

主要发现

• 所提出的 PIP-SOLVER 算法以 O(N²) 时间和 O(mN) 存储计算插值多项式的系数，其中 N = N(m,n) = (m+n choose m)。
• 该算法完全避免了矩阵求逆，从而突破了矩阵求逆的理论下限 O(N² log N)，在实践中优于 Coppersmith-Winograd 和基于 Strassen 的方法。
• 数值实验表明，PIP-SOLVER 在精度和运行时间方面显著优于基于线性系统求解或矩阵求逆的经典求解器。
• 改进不仅体现在预因子上，更体现在缩放行为上：该方法即使在极大规模问题下仍能保持恒定的数值精度。
• 由于计算成本降低和数值稳定性增强，该算法使此前无法求解的插值问题成为可能。
• 由于子问题相互独立，分布式内存并行实现是可行的，可进一步降低运行时间并提升内存可扩展性。