Mete Yurtoglu Google, 747 6th Street South, Kirkland, WA 98033 e-mail: myurtoglu@me.com Molly Carton Mechanical Engineering, University of Washington, Seattle, WA 98195-2600 e-mail: mcarton@uw.edu Duane Storti Mechanical Engineering, University of Washington, P.O. Box: 352600, Seattle, WA 98195-2600 e-mail: storti@uw.edu Treat All Integrals as Volume Integrals: A Unified, Parallel, Grid-Based Method for Evaluation of Volume, Surface, and Path Integrals on Implicitly Defined Domains We present a unified method for numerical evaluation of volume, surface, and path inte- grals of smooth, bounded functions on implicitly defined bounded domains. The method avoids both the stochastic nature (and slow convergence) of Monte Carlo methods and problem-specific domain decompositions required by most traditional numerical integra- tion techniques. Our approach operates on a uniform grid over an axis-aligned box con- taining the region of interest, so we refer to it as a grid-based method. All grid-based integrals are computed as a sum of contributions from a stencil computation on the grid points. Each class of integrals (path, surface, or volume) involves a different stencil for- mulation, but grid-based integrals of a given class can be evaluated by applying the same stencil on the same set of grid points; only the data on the grid points changes. When functions are defined over the continuous domain so that grid refinement is possible, grid-based integration is supported by a convergence proof based on wavelet analysis. Given the foundation of function values on a uniform grid, grid-based integration meth- ods apply directly to data produced by volumetric imaging (including computed tomogra- phy and magnetic resonance), direct numerical simulation of fluid flow, or any other method that produces data corresponding to values of a function sampled on a regular grid. Every step of a grid-based integral computation (including evaluating a function on a grid, application of stencils on a grid, and reduction of the contributions from the grid points to a single sum) is well suited for parallelization. We present results from a paral- lelized CUDA implementation of grid-based integrals that faithfully reproduces the out- put of a serial implementation but with significant reductions in computing time. We also present example grid-based integral results to quantify convergence rates associated with grid refinement and dependence of the convergence rate on the specific choice of difference stencil (corresponding to a particular genus of Daubechies wavelet). [DOI: 10.1115/1.4039639] 1 Introduction Evaluation of a definite surface integral is typically defined as the sum of contributions from a finely divided polygonal approxi- mation to the domain of integration [1], and the statement extends to cover line and volume integrals if polygon is generalized to n- dimensional polytope. Methods for evaluating integrals based on such a definition require a finely divided set of polytopes custom- ized to closely approximate the particular domain of integration. In cases where the domain is specified in terms of polytopes, the polytopes can be finely subdivided and such an approach is per- fectly reasonable. In contrast, this paper focuses on the case when the domain is defined by an implicit function (or by an approxima- tion based on a sampling of function values) so that a finely divided domain approximation is not available without further computation. In the field of engineering design, evaluation of definite inte- grals is often associated with determining mass or surface proper- ties of a solid model, and there is a considerable literature related to this task dating back to Sarraga [2] and Lee and Requicha [3,4] who nicely set the problem context. In part I of their work, Lee and Requicha [3] noted that the existing literature on computa- tional multiple integration focuses on problems with a compli- cated integrand but a simple domain, and characterized mass property computation as the converse problem with relatively sim- ple integrands on complicated domains. They also identified key issues: (1) the modeling scheme used to describe the solid domi- nates the design of integration algorithms and (2) dominant errors are typically associated with errors in representing the domain. Their treatment starts with consideration of solids represented in terms of their polyhedral boundaries, and states a method that combines the Divergence Theorem with integration over the para- metrized polygonal facets. While stating that integrals on regions with curved boundaries can be approximated by obtaining polygo- nal approximations of the curved surfaces or by applying Green’s theorem on the trimmed patch boundaries in the parameter plane, they are careful to note that that error estimates are not available for such methods. They go on to discuss Monte Carlo methods that employ computations at a large number N of sample points to provide reliable but slow ðOð ffiffiffiffiffiffiffi ðNÞ p ÞÞ convergence and, in part II [4], methods based on conversion to cellular approximations. Manuscript received September 8, 2017; final manuscript received March 9, 2018; published online April 26, 2018. Assoc. Editor: Yong Chen. Journal of Computing and Information Science in Engineering JUNE 2018, Vol. 18 / 021013-1 Copyright V C 2018 by ASME