In many applications, one seeks to approximate integration against a positive measure of interest by a positive discrete measure: a numerical quadrature rule with positive weights. One common desired discretization property is moment preservation over a finite dimensional function space, e.g., bounded-degree polynomials. Carath\'{e}odory's theorem asserts that if there is any finitely supported quadrature rule with more nodes than the dimension of the given function space, one can form a smaller (and hence more efficient) positive, nested, quadrature rule that preserves the moments of the original rule. We describe an efficient streaming procedure for Carath\'{e}odory-Steinitz pruning, a numerical procedure that implements Carath\'{e}odory's theorem for this measure compression. The new algorithm makes use of Givens rotations and on-demand storage of arrays to successfully prune very large rules whose storage complexity only depends on the dimension of the function space. This approach improves on a naive implementation of Carath\'{e}odory-Steinitz pruning whose runtime and storage complexity are quadratic and linear, respectively, in the size of the original measure. We additionally prove mathematical stability properties of our method with respect to a set of admissible, total-variation perturbations of the original measure. Our method is compared to two alternate approaches with larger storage requirements: non-negative least squares and linear programming, and we demonstrate comparable runtimes, with improved stability and storage robustness. Finally, we demonstrate practical usage of this algorithm to generate quadrature for discontinous Galerkin finite element simulations on cut-cell meshes.
翻译:在许多应用中,人们希望通过一个正离散测度——即具有正权重的数值积分法则——来逼近对某个感兴趣的正测度的积分。一种常见的理想离散化性质是在有限维函数空间(例如有界次数的多项式空间)上保持矩不变。Carathéodory定理断言:如果存在一个节点数多于给定函数空间维数的有限支撑求积法则,那么我们可以构造一个更小(因而更高效)的正的、嵌套的求积法则,使其保持原法则的矩。我们描述了一种用于Carathéodory-Steinitz剪枝的高效流式处理过程,该数值过程为实现测度压缩的Carathéodory定理提供了具体实现。新算法利用Givens旋转和按需存储数组的策略,能够成功地对存储复杂度仅依赖于函数空间维数的超大规模法则进行剪枝。相比Carathéodory-Steinitz剪枝的朴素实现(其时间复杂度和存储复杂度分别与原测度规模呈二次和线性关系),本方法有所改进。我们还从数学上证明了该方法对于原测度在总变差意义下的一类容许扰动具有稳定性。我们将本方法与两种存储需求更大的替代方法——非负最小二乘和线性规划——进行了比较,结果显示在运行时间相当的情况下,本方法具有更好的稳定性和存储鲁棒性。最后,我们展示了该算法在为切割网格上的间断Galerkin有限元模拟生成求积法则时的实际应用。