Polyhedral Complex Derivation from Piecewise Trilinear Networks

Authors: jin-hwa  kim
Venue: NeurIPS_cc_2024_Conference
Published: 2024-04-30 06:55:38

Abstract: Recent advancements in visualizing deep neural networks provide insights into their structures and mesh extraction from Continuous Piecewise Affine (CPWA) functions. Meanwhile, developments in neural surface representation learning incorporate non-linear positional encoding, addressing issues like spectral bias; however, this poses challenges in applying mesh extraction techniques based on CPWA functions. Focusing on trilinear interpolating methods as positional encoding, we present theoretical insights and an analytical mesh extraction, showing the transformation of hypersurfaces to flat planes within the trilinear region under the eikonal constraint. Moreover, we introduce a method for approximating intersecting points among three hypersurfaces contributing to broader applications. We empirically validate correctness and parsimony through chamfer distance and efficiency, and angular distance, while examining the correlation between the eikonal loss and the planarity of the hypersurfaces.

Main paper: 
 introduction: 
 Recent advancements in visualizing deep neural networks [1] [2] [3] [4] significantly contribute to understanding their intricate structures. This progress provides valuable insights into the expressivity, robustness, training methodologies, and distinctive geometry of neural networks. By leveraging the inherent piecewise linearity in the Continuous Piecewise Affine (CPWA) functions, e.g., ReLU neural networks, each region of the input space is represented as a convex polyhedron. The assembly of these sets constructs a polyhedral complex that delineates the decision boundaries of neural networks [5] .Building upon these breakthroughs, we explore the exciting prospect of analytically extracting a mesh representation from neural implicit surface networks [6] [7] [8] . The extraction of a mesh not only offers a precise visualization and characterization of the network's geometry but also does so efficiently through the utilization of vertices and faces, in contrast to sampling-based methods [9, 10] . The sampling-based methods are limited by their black-box nature, which can obscure the underlying structure and lead to inefficiencies in capturing fine geometric details.These networks typically learn a signed distance function (SDF) by utilizing ReLU neural networks. However, recent approaches incorporate non-linear positional encoding techniques, such as trigonometric functions [11] or trilinear interpolations using hashing [12] or tensor factorization [13] , to mitigate issues like spectral bias [14] , ensuring fast convergence, and maintaining high-fidelity. Consequently, applying mesh extraction techniques based on CPWA functions becomes challenging.Focusing on the successful and widely-adopted trilinear interpolating methods, we present novel theoretical insights and a practical methodology for precise mesh extraction. Following a novel review The mesh is analytically extracted from piecewise trilinear networks, comprising both HashGrid [12] and ReLU neural networks, that have been trained to learn a signed distance function with the eikonal loss. We start with initial vertices and edges defined using the grid marks (1) , and its linear subdivision (2a); however, intersecting polygons (ref. Section 3.2) need bilinear subdivision (2b) and consequentially trilinear subdivision (2c). Note that the linear and bilinear subdivisions are specific instances of trilinear subdivisions with straightforward solutions.of the edge subdivision method [4] , an efficient mesh extraction technique for CPWA functions, within the framework of tropical geometry (Section 3), we demonstrate that ReLU neural networks incorporating the trilinear module are piecewise trilinear networks (Definition 4.1). Then, we establish that the hypersurface within the trilinear region becomes a plane under the eikonal [15] constraint (Theorem 4.5), a common practice during the training of SDF.Building upon this observation, we propose a method for approximating the determination of intersecting points among three hypersurfaces. Given the inherent curvature of these surfaces, achieving an exact solution is infeasible. In our approach, we utilize two hypersurfaces and a diagonal plane for a more feasible and precise approximation (Theorem 4.7). Additionally, we introduce a methodology for arranging vertices composed of faces to implicitly represent normal vectors. We experimentally confirm accuracy and simplicity using chamfer distance, efficiency, and angular distance, simultaneously exploring the relationship between the eikonal loss and the planarity of the hypersurface. In Figure 1 , we present the analytically generated skeleton and normal map produced by our method using piecewise trilinear networks composed of HashGrid and ReLU neural networks. For a more detailed view, please refer to Figure 10 in Appendix.Our contributions are summarized as follows:1. We present novel theoretical insights and a practical methodology for precise mesh extraction, from the piecewise trilinear networks using the HashGrid.2. We provide a theoretical analysis of the piecewise trilinear networks under the eikonal constraint revealing that, within a trilinear region, a hypersurface transforms into a plane.3. We validate its correctness and parsimony through chamfer distance and efficiency, and angular distance, showing the correlation between the eikonal loss and hypersurface planarity.Question: Do the main claims made in the abstract and introduction accurately reflect the paper's contributions and scope? Answer: [Yes] Justification: Our main claims made in the abstract and introduction are accurately reflected in theoretical and experimental expositions. Theorem 4.5 and Figure 5 provide the theoretical and experimental supports of the planarity constraint of the eikonal equation, Theorem 4.7 provides a plausible approximation for the curved edge subdivision, and Table 1 shows the consistent chamfer efficiency across multiple meshes from the Stanford 3D Scanning repository, while Figure 4 shows that the higher resolution of HashGrid is favorable to the planarity constraint, achieving lower chamfer distances. Guidelines:• The answer NA means that the abstract and introduction do not include the claims made in the paper. • The abstract and/or introduction should clearly state the claims made, including the contributions made in the paper and important assumptions and limitations. A No or NA answer to this question will not be perceived well by the reviewers. • The claims made should match theoretical and experimental results, and reflect how much the results can be expected to generalize to other settings. • It is fine to include aspirational goals as motivation as long as it is clear that these goals are not attained by the paper. 
related work: 
 Mesh extraction in deep neural networks. A conventional black-box approach is to assess a given function at sampled points to derive the complex of the function, commonly through marching cubes [9] or marching tetrahedra [10, 16] . Despite the inclusion of mesh optimization techniques [17] in modern mesh packages, this exhaustive method often generates redundant complexities due to discretization, increasing the computational load by using unnecessary mesh elements (e.g., fragmented planars). Aiming for exact extraction, initial efforts focused on identifying linear regions [18] and employing Analytical Marching [2] to exactly reconstruct the zero-level isosurface. Contemporary approaches involve region subdivision, wherein the regions of the complex are progressively subdivided from neuron to neuron and layer to layer, allowing for the efficient calculation of the exponentially increasing linear regions [4, [19] [20] [21] . Nevertheless, many previous studies have explored CPWA functions consisting of fully connected layers and ReLU activation, primarily owing to the mathematical simplicity of hyperplanes defining linear regions. However, this restricts their applications in emerging domains of deep implicit surface learning that diverge from CPWA functions.Positional encoding. The spectral bias of a multi-layer perceptron (MLP) hinders effectively learning high frequencies, both theoretically and empirically [14] . To address this limitation, Fourier feature mapping, employing trigonometric projections similar to those used in the Transformer architecture [22] , proves successful in representing complex 3D objects and scenes. Furthermore, to address its slow convergence and enhance rendering quality, novel positional encoding methods using trilinear interpolations, e.g., HashGrid [12] or TensoRF [13] , are introduced. Both techniques generate positional features through trilinear interpolation among the eight nearest corner features derived from pre-defined 3D grids. These features are obtained by hashing from multi-resolution hash tables or factorized representations of a feature space, respectively. However, these approaches cause the function represented by neural networks to deviate from CPWA functions, as trilinear interpolation is not an affine transformation. Consequently, the complex regions are divided by intricate hypersurfaces.Eikonal equation and SDF. An eikonal 1 equation is a non-linear first-order partial differential equation governing the propagation of wavefronts. Notably, this finds application in the regularization of a signed distance function (SDF) in the context of neural surface modeling [6, 8, 23] . If Ω represents a subset of the space X equipped with the metric d, the SDF f (x) is defined as follows:EQUATIONwhere ∂Ω denotes the boundary of Ω. The metric d is defined (with a slight notational abuse) as:EQUATIONwhere d represents the shortest distance from x to ∂Ω in the Euclidean space of X . In the Euclidean space with a piecewise smooth boundary, the SDF exhibits differentiability almost everywhere, and its gradient adheres to the eikonal equation. Specifically, for points x on the boundary, ∇f(x) = N (x)where the Jacobian of f represents the outward normal vector, or equivalently, |∇f (x)| = 1.Question: Does the paper discuss the limitations of the work performed by the authors? Answer: [Yes] Justification: Please refer to Section 6. Guidelines:• The answer NA means that the paper has no limitation while the answer No means that the paper has limitations, but those are not discussed in the paper. • The authors are encouraged to create a separate "Limitations" section in their paper.• The paper should point out any strong assumptions and how robust the results are to violations of these assumptions (e.g., independence assumptions, noiseless settings, model well-specification, asymptotic approximations only holding locally). The authors should reflect on how these assumptions might be violated in practice and what the implications would be. • The authors should reflect on the scope of the claims made, e.g., if the approach was only tested on a few datasets or with a few runs. In general, empirical results often depend on implicit assumptions, which should be articulated. • The authors should reflect on the factors that influence the performance of the approach.For example, a facial recognition algorithm may perform poorly when image resolution is low or images are taken in low lighting. Or a speech-to-text system might not be used reliably to provide closed captions for online lectures because it fails to handle technical jargon. • The authors should discuss the computational efficiency of the proposed algorithms and how they scale with dataset size. • If applicable, the authors should discuss possible limitations of their approach to address problems of privacy and fairness. • While the authors might fear that complete honesty about limitations might be used by reviewers as grounds for rejection, a worse outcome might be that reviewers discover limitations that aren't acknowledged in the paper. The authors should use their best judgment and recognize that individual actions in favor of transparency play an important role in developing norms that preserve the integrity of the community. Reviewers will be specifically instructed to not penalize honesty concerning limitations.Guidelines:• The answer NA means that the paper does not include experiments.• The authors should answer "Yes" if the results are accompanied by error bars, confidence intervals, or statistical significance tests, at least for the experiments that support the main claims of the paper. • The factors of variability that the error bars are capturing should be clearly stated (for example, train/test split, initialization, random drawing of some parameter, or overall run with given experimental conditions). • The method for calculating the error bars should be explained (closed form formula, call to a library function, bootstrap, etc.) • The assumptions made should be given (e.g., Normally distributed errors).• It should be clear whether the error bar is the standard deviation or the standard error of the mean. • It is OK to report 1-sigma error bars, but one should state it. The authors should preferably report a 2-sigma error bar than state that they have a 96% CI, if the hypothesis of Normality of errors is not verified. • For asymmetric distributions, the authors should be careful not to show in tables or figures symmetric error bars that would yield results that are out of range (e.g. negative error rates). • If error bars are reported in tables or plots, The authors should explain in the text how they were calculated and reference the corresponding figures or tables in the text. 
preliminaries: 
 We present an overview of tropical geometry (Section 3.1) to provide a formal definition of the tropical hypersurface and the tropical algebra of neural networks. Then, we discuss the edge subdivision algorithm [4] for extracting the tropical hypersurface (Section 3.2). We supplement Appendix A.1 for less familiar readers. Later, we extend this idea to include considerations for trilinear interpolation in Section 4. You may choose to skip these sections if you are already acquainted with these concepts.Tropical 2 geometry [25, 26] describes polynomials and their geometric characteristics. As a skeletonized version of algebraic geometry, tropical geometry represents piecewise linear meshes using tropical semiring and tropical polynomial (ref. Appendix A.1). The set of points where a tropical polynomial f is non-differentiable, due to tropical sum ⊕ := max(•), is called tropical hypersurface:Definition 3.1 (Tropical hypersurface). The tropical hypersurface of a tropical polynomial f (x) = c 1 x α1 ⊕ • • • ⊕ c r x αr is defined as: V (f ) = {x ∈ R d : c i x αi = c j x αj = f (x) for some α i ̸ = α j }Although the neural networks ν is generally non-convex, it is known that the difference of two tropical signomials (allowing a real-valued α in a tropical polynomial) can represent ν [27] (ref.Appendix A.3). Given this observation, we proceed to examine the following proposition:Proposition 3.3 (Decision boundary of neural networks). The decision boundary for L-layer neural networks with ReLU activation ν (L) (x) = 0 is a subset of or equals with the region where two tropical monomials, or two arguments of max, equal (Definition 3.1), as follows:B ⊆ x : ν (i) j (x) = 0, i ∈ (1, . . . , L), j ∈ (1, . . . , H) if i ̸ = L else (1)where the subscript j denotes the j-th neuron, assuming the hidden size of neural networks is H, while the final L-th layer has a single output for rather discussions as an SDF (Definition 4.3).Remind that tropical operations satisfy the usual laws of arithmetic for recursive consideration.The piecewise hyperplanes corresponding to ν (i) j are the candidates shaping the decision boundary. It enables for-loop iterations over the neurons across layers, as in Algorithm 2, not visiting every exponentially growing linear region [18, 27] . The equality in Proposition 3.3 generally does not hold, but we can filter them by checking if ν (L) 1 (x) = 0. For more details, please refer to Berzins [4] . We note that our formal explanation through tropical geometry is unprecedented in the previous work [4] .Polyhedral3 complex derivation from a multitude of linear regions is inefficient and may be infeasible in some cases. Rather, Berzins [4] argue that tracking vertices and edges of the decision boundary sequentially considering the polynomial number of hyperplanes is particularly efficient, even enabling parallel tensor computations for edge subdivisions.Initially, we start with a unit cube of eight vertices and a dozen edges, denoted by V and E, respectively. For each piecewise linear, or folded in their term, hyperplane, we find the intersection points of each one of the edges and the hyperplane (if any) to add to V. The divided edges by the intersection are also added to E. Note that we also find new edges of polygons, the intersections of convex polyhedra representing corresponding linear regions [5] , and the folded hyperplane. Repeating these for all folded hyperplanes, before selecting all vertices x ⋆ where ν(x ⋆ ) = 0. Then, the valid edges connecting those vertices would represent the decision boundary.To elaborate with details, the order of subdividing is invariant as stated in Proposition 3.3; however, we must perform the subdivision with all hyperplanes in the previous layers in advance if we want to keep the current set of edges not crossing linear regions. One critical advantage of this principle is that we can find the intersection point by evaluating the output ratio of the two vertices of an edge, ν(x 0 ) and ν(x 1 ), which are proportional to the distances to the hyperplane by Thales's theorem [28] . The intersection point would simply be x0,1 = (1 -w)x 0 + wx 1 where w= |d 0 |/|d 0 -d 1 | and d k = ν (l) j (x k ), upon the fact that d 0 • d 1 < 0 if the hyperplane divides the edge.Updating edges needs two steps: dividing the edge into two new edges and adding new edges of intersectional polygons by a hyperplane for convex polyhedra representing linear regions. The latter is a tricky part that we need to find every pair of two vertices among { x•,• }, both within the same linear region and on the two common hyperplanes, obviously including the current hyperplane. They argue that the sign-vectors for the preactivations provide an efficient way to find them. Definition 3.4 (Sign-vectors). For the current hyperplane specified by the i * -th layer and its j * -th preactivation, we define the sign-vectors with a small positive constant ϵ:γ(x k ) ϕ(i,j) =      +1 if ν (i) j (x k ) > +ϵ 0 if |ν (i) j (x k )| ≤ +ϵ -1 if ν (i) j (x k ) < -ϵ (4)for all i and j where ϕ(i, j) ≤ ϕ(i * , j * ), and ϕ(i, j) ∈ {0} ∪ N is an ascending indexing function from lower layers, handling arbitrary hidden sizes of networks to vectorize a signed matrix.Notice that a hyperplane divides a space into two half-spaces, two linear regions. Collectively, if the sign-vectors for the i-th layer and its j-th preactivation of two vertices have the same values except for zeros, which are wild cards for matching, we can say that the two vertices are within the same linear region in the current step (i * , j * ). For the second, the sign-vector has at least three zeros specifying a point by at least three hyperplanes. If two zeros are at the same indices in the two sign-vectors, the two vertices are on the same two hyperplanes, forming an edge. Notice that our principle asserts edges should be within linear regions.Berzins [4] argue that the edge subdivision provides the optimal time and memory complexities of O(|V|), linear in the number of vertices, while it can be efficiently computed in parallel. According to Hanin and Rolnick [20] , the number of linear regions is known for o(N L /L!) where N is the total number of neurons, or ϕ(L, 1) in our notation. 
method: 
 4.1 Piecewise trilinear networks Definition 4.1 (Trilinear interpolation). Let a unit cube be in the D-dimensional space, where its corners are represented by H ∈ R 2 D ×F where F is a corner-feature dimension. Given weights w ∈ [0, 1] D , the trilinear interpolation function ψ is defined as:EQUATIONwhere the interpolating weight ω(i, w) ∈ R is defined using a left-aligned and zero-trailing binarization function B ∈ {0, 1} D (ref. Appendix B.1) as follows:EQUATIONwhich is the volume of the opposite subsection of the hypercube divided by the weight point w. For a general hypercube, scaling the weight for a unit hypercube gets the same result. Going forward, we assume D=3 for trilinear interpolation without explicitly stating it otherwise.Lemma 4.2 (Nested trilinear interpolation). Let a nested cube be inside of a unit cube, where its positions of the eight corners deviate -a j or +b j , such that a j , b j ≥ 0, from w j for each dimension j, using the notations of Definition 4.1. The eight corners of the nested cube are the trilinear interpolations of the eight corners of the unit cube. Then, the trilinear interpolation with the unit cube and the nested cube for w is identical.Proof. Notice that trilinear interpolation is linear for each dimension to prove it. The detailed proof can be found in Lemma D.1.Definition 4.3 (Piecewise trilinear networks). Let a positional encoding module be τ : R D → R F using the trilinear interpolation of spatially nearby learnable vectors on a three-dimensional grid, e.g., HashGrid [12] or TensoRF [13] , where τ is an instance of ψ in Definition 4.1. Usually, they transform the input coordinate x to w as a relative position inside a unit grid in a corresponding resolution.H is a learnable parameter specified in the corresponding method. Then, we define trilinear neural networks ν, by prepending τ to the neural networks ν from Definition 3.2. Here, the input and output dimensions of ν (L) are F and 1 (an SDF distance), respectively.EQUATIONNote that τ is gridwise trilinear. Using Lemma 4.2, ν is trilinear within the composite intersections of linear regions of ν and trilinear regions of τ . (We discuss how to access virtual cubic corner features where trilinear regions are not cubic in Section 5.) So, we regard ν as piecewise trilinear.In a trilinear space, τ projects a line to a curve except the lines on the grid. Notably, the diagonal line t = (t, t, t) ⊺ where t ∈ [0, 1] is projected to a cubic Bézier curve (ref. Proposition D.2). We aim to generalize for the curved edge subdivision for hypersurfaces. Lemma 4.4 (Curved edge of two hypersurfaces). In a piecewise trilinear region, let an edge (x 0 , x 7 ) be the intersection of two hypersurfaces νi (x) = νj (x) = 0, while x 0 and x 7 are on the two hypersurfaces. Then, the edge is defined as:{x : τ (x) = (1 -t)τ (x 0 ) + tτ (x 7) where x ∈ R D and t ∈ R}.Proof. The detailed proof using the piecewise linearity of ν is provided in Lemma D.4.Notice that it does not decrease the number of variables to find a line solution since the trilinear interpolation τ still forms hypersurfaces making complex cases. Yet, we theoretically demonstrate that the eikonal constraints on ν render them hyperplanes with linear solutions. Theorem 4.5 (Hypersurface and eikonal constraint). A hypersurface τ (x) = 0 intersects two points τ (x 0 ) = τ (x 7 ) = 0 while τ (x 1...6 ) ̸ = 0 for the remaining six points. These points form a cube, with x 0 and x 7 positioned on the diagonal of the cube. The hypersurface satisfies the eikonal constraint ∥∇τ (x)∥ 2 2 = 1 for all x ∈ [0, 1] 3 . Then, the hypersurface of τ (x) = 0 is a plane. Corollary 4.6 (Affine-transformed hypersurface and eikonal constraint). Let ν 0 (x) and ν 1 (x) be two affine transformations defined byW ⊺ i x + b i , i ∈ {0, 1}. A hypersurface ν 0 τ (x) = 0 passing two points ν 0 τ (x 0 ) = ν 0 τ (x 7 ) = 0 while ν 0 τ (x 1...6 ) ̸ = 0 satisfies the eikonal constraint ∥∇ν 1 τ (x) ∥ 2 2 = 1 for all x ∈ [0, 1] 3 .Then, the hypersurface of ν 0 τ (x) = 0 is a plane.The proofs using its partial derivatives can be found in Theorem D.5 and Corollary D.6.Since we aim for a practical polyhedral complex derivation that piecewise trilinear networks represent, we replace one of the hypersurfaces forming the curved edge with a diagonal plane in a piecewise trilinear region. We choose this option not only for its mathematical simplicity but also due to the characteristics of trilinear hypersurfaces, as illustrated in Figure 9 , Appendix G. Thus, the new vertices lie on at least two hypersurfaces, and the new edges exist on the same hypersurface, while the eikonal constraint minimizes any associated error from this approximation. This error is tolerated by the hyperparameter ϵ =1e-4 as specified in Definition 3.4 and Section 6. More discussions on this matter can be found in Appendix C.1. Theorem 4.7 (Intersection of two hypersurfaces and a diagonal plane). Let ν0 (x) = 0 and ν1 (x) = 0 be two hypersurfaces passing two points x 0 and x 7 such that ν0 (x 0 ) = ν1 (x 0 ) = 0 and ν0 (x 7 ) = ν1 (x 7 ) = 0, P i := ν0 (x i ) and Q i := ν1 (x i ), P α = P 0 ; P 1 ; P 4 ; P 5 , P β = P 2 ; P 3 ; P 6 ; P 7 , andX = [(1 -x) 2 ; x(1 -x); (1 -x)x; x 2 ]. Then, x ∈ [0, 1]of the intersection point of the two hypersurfaces and a diagonal plane of x = z is the solution of the following quartic equation:EQUATIONProof. Please refer to Theorem D.7 for the detailed proof, which uses linear algebra to rearrange two trilinear equations with two variables to get a solution.Note that finding roots of a polynomial is the eigenvalue decomposition of a companion matrix [29] , which can be parallelized (ref. Lemma D.8) . Given that the companion matrix remains small for a quartic equation in R 4×4 , the overall complexity remains O(|V|). For the detailed complexity analysis on the proposed method, please refer to Appendix E.Skeletonization selects the vertices such that: Faces are formed by edges lying on the same plane and sharing a common region. This is identified by evaluating the sign-vectors (ref. Definition 3.4 and Section 5). Note that faces are not inherently triangular; however, if needed, they can be triangulated via triangularization for further analysis.V ⋆ = {x | |ν(x)| ≤ ϵ, x ∈ V},Normals are conventionally described by the order of the vertices as each face has ambiguity with two opposite directions. Let x 0,1,2 be the vertices forming a triangular face. The normal is defined using cross-product as n = (x 1 -x 0 ) × (x 2 -x 0 ) where the normal is orthogonal to both vectors, with a direction given by the right-hand rule. Please refer to Appendix B.2 for the details. 
implementation: 
 We employ the HashGrid [12] for τ , while our discussion on its generalizability is provided in Appendix C.2. We illustrate our method in Algorithm 1, 2, and Figure 2 , inspired by Hanin and Rolnick [30] for its visualization of linear regions. We describe the implemental details of the algorithms in the following sections and the corresponding caption of Figure 6 .Initialization. While the original edge subdivision algorithm [4] starts with a unit cube of eight vertices and a dozen edges, we know that the unit cube is subdivided by orthogonal planes consisting of the grid. Let an input coordinate be x ∈ [0, 1] 3 , and there is a unit cube such that one of its corners is at the origin. Taking into account multi-resolution grids, we obtain the marks m i , such that the grid planes are x = m i , y = m i , and z = m i , following Algorithm 4 in Appendix. Notice that we carefully consider the offset s/2, which prevents the zero derivatives from aligning across all resolution levels (ref.Appendix A of Müller et al. [12] ). Leveraging the obtained marks, we derive the associated grid vertices and edges, which serve as the initial sets for V and E, respectively. Please refer to Figure 2a where its multi-resolution grids are visualized using randomized colors.Optimizing sign-vectors. The number of orthogonal planes may significantly surpass the number of neurons, making it impractical to allocate elements for the sign-vectors (Definition 3.4). To address this, leveraging the orthogonality of the grid planes, we store a plane index along with an indicator of whether the input is on the plane (γ = 0) or not (γ = 1).Piecewise trilinear region. In Theorem 4.7, we need the outputs of corners P and Q in a common piecewise trilinear region; however, some corners may be outside of the region. For this, we replace all ReLU activations using the mask m that: m(i) j = (ν (i) j (x 0 ) > ϵ) | (ν (i) j (x 7 ) > ϵ), where x 0 and x 7 are the vertices of an edge of interest, for the calculated cubic corners x 0...7 . This implies that ν exhibits linearity, except for instances where two vertices simultaneously deviate from linearity. For P and Q, we select the last two pairs of i and j such that γ ϕ(i,j) (x 0 ) = γ ϕ(i,j) (x 7 ) = 0, indicating two common hypersurfaces passing two points x 0 and x 7 . 
experiment: 
 Objective. For a given trilinear neural network with the eikonal constraint, our goal is to get a boundary mesh of the zero-set of the SDF. (Note that NeuS [6] showed how to convert the densityInput: Sets of vertices V and edges E (from Section 5), # of layer L, hidden size H, trilinear networksν for i = 0 to L -1 do for j = 0 to H -1 do V, E = CESub(V, E, ν, i, j) ▷ Algorithm 2 end for end for V, E = CESub(V, E, ν, L, 0) V ⋆ , E ⋆ = Skeletonize(V, E, ν, ϵ) ▷ Section 4.3 F ⋆ = ExtractFaces(V ⋆ , E ⋆ ) ▷ Section 4.3 return V ⋆ , E ⋆ , F ⋆ Algorithm 2 Curved Edge Subdivision (CESub) Input: Vertices V, Edges E, trilinear networks ν, layer index i, neuron index j d = ν(V; i, j), V new = ∅ for e in E do if d e0 • d e1 < 0 then P , Q = Corners(e) ▷ Section 5 x new = Intersect(P , Q, e) ▷ Theorem 4.7 V new ← V new + {x new } e ← [e 0 , Index(x new )] E ← E + {[e 1 , Index(x new )]} end if end for V ← V + V new E ← E+ FindPolygons(V new ) ▷ Section 3.2 return V, Enetworks of NeRFs to an SDF in the frameworks of volume rendering.) To evaluate this, we utilize marching cubes with excessive sampling to get a pseudo-ground truth mesh to remove Bayes error caused by underfitting. We explicitly specify the samplings in the results.Hyperparameters. We used the number of layers L of 3 and hidden size H of 16 for the networks, and ϵ of 1e-4 for the sign-vectors (ref. Definition 3.4) . The weight for the eikonal loss is 1e-2. For HashGrid, the resolution levels of 4, feature size of 2, base resolution N min of 2, and max resolution N max of 32 by default. N min and N max are doubled (x2) or quadrupled (x4) for Medium or Large settings in Figure 4 . We use the official Python package of tinycudnn4 for the HashGrid module.Chamfer distance and efficiency. The chamfer distance is a metric used to evaluate the similarity between two sets of sampled points from two meshes. Let S 0 and S 1 be the two sets of points. Specifically, the bidirectional chamfer distance (CD) is defined as follows:EQUATION)We randomly select 100K points on the faces through ray-marching from the origin, directing toward a randomly chosen point on a unit sphere.0K 60K 120K 180K 240K 300K 360K Table 1 and Figure 4 show the results for the Standford 3D Scanning repository [31] . Varying the number of samples in marching cubes [9] , we plot the baseline with respect to the number of generated vertices. Given that our method extracts vertices from the intersection points of hypersurfaces, it enables parsimonious representations using the lower number of vertices having better CD, summarized by chamfer efficiency (CE) 5 in Table 1 , and ours are located in the lower-left region compared with the baselines in Figure 4 . In particular, we observe a strong CE in a simple object, e.g., Drill Bit, which achieved a CE of 47.8 versus 19.2 (MC).Our method inherently holds an advantage over MC in capturing optimal vertices and faces, especially in cases where the target surfaces exhibit planar characteristics. When dealing with curved surfaces, the smaller grid size of the Large model in Figure 4 views a curved surface in a smaller interval, and it tends to approach the plane. From the plotting of Small to Table 1 : Chamfer distance (CD) and chamfer efficiency (CE) for the Stanford 3D Scanning repository [31] . The CD (×1e-6) is evaluated using the marching cubes (MC) with 256 3 grid samples as the reference ground truth for all measurements. Please refer to Table 2 and Table 3 for the full results, and Table 4 for standard deviation and time spent in Appendix (Ours took 0.97 ± 0.21 while MC also took within 1 sec for the bunny.) (11) . For the plot, we conduct experiments using the unit sphere.EQUATIONLarge, our CDs (colored crosses) are decreasing to zero, consistently with better CEs, confirming this speculation. The complete results of the Large models and its visualization can be found in Tables 6 to 8 in Appendix F.1 and the right side of Figure 10 in Appendix, respectively. We also compare ours with Marching Tetrahedra (MT) [10, 16] and Neural Dual Contour (NDC) [32] in Figure 3 , where our method achieved the lowest chamfer distance and the most efficiency method with respect to the number of vertices. For the rationale behind our selection, please refer to Appendix F.3.Angular distance. The angular distance measures how much the normal vectors deviate from the ground truths. Analogous to the chamfer distance, we sample 100K points on the surface and calculate the normal vectors as described in Section 4.3. The angular distance is defined as follows:EQUATIONIn Table 5 in Appendix, the angular distances for the Stanford bunny are shown. As we expected, our approach efficiently estimates the faces using a parsimonious number of vertices compared to the sampling-based method. We confirm this trend is consistent across other objects.Eikonal constraint for planarity. To assess the efficacy of the eikonal constraint, which enforces the planarity of hypersurfaces as outlined in Theorem 4.5, we quantify the flatness error associated with the planarity constraints presented in the proof of Theorem D.5. Let x 0 and x 7 be two vertices consisting of a diagonal edge, while x 1...6 are its remaining cubic corners. We obtain ν(x 1...6 ) in a piecewise trilinear region as described in Section 5. This evaluation is conducted as follows:EQUATIONIn Figure 5 , we empirically validate that the eikonal loss enforces planarity in the hypersurfaces within piecewise trilinear regions, as described in Theorem 4.5. The absence of the eikonal loss results in elevated planarity errors, leading to the construction of an inaccurate mesh. [32] ) suffer over-smoothing or inaccuracy of the actual surface in the Small networks within the nose, whereas our method reflects these details (see the boundaries of a nose) with consistent normals (in colors). MT efficiently demands SDF values by utilizing six tetrahedra within grids to get intermediate vertices. However, it is inefficient in terms of the number of vertices extracted. NDC faces a generalization issue for a zero-shot setting for the given networks, producing unsmooth surfaces.Visualization. We provide the qualitative analyses in Figure 6 comparing with MT [10, 16] and NDC [32] , along with MC. While the competitors exhibit over-smoothing or inaccuracies, particularly in the nose region, our method preserves surface details with consistent normals and outperforms in both accuracy and generalization. In addition, we present Figures 10 to 13 in Appendix G encompassing the comparison between the Small and Large models (Figure 10 ), visualizing the marching cubes varying the number of samples (Figure 11 ), detailed comparisons (Figure 12 ), and the gallery of the Stanford 3D Scanning meshes from the Large models (Figure 13 ).Limitations. Although we implement the algorithm using pre-defined parallel tensor operations in PyTorch [33] , achieving time and memory complexity of O(|V|) as detailed in Appendix E, our method inherently requires forward passing the intermediate vertices during edge subdivisions iterating over neurons, which depends on previous steps (Algorithm 1). For the Stanford bunny using the Small model, the intermediate counts for vertices and edges exceed 11K and 19K, respectively, although they are eventually reduced to 4.3K vertices for face extraction. This overhead is negligible for the Small models, taking 0.97 seconds; however, for the Large models, the process takes 4.15±0.61 seconds for 140.3±5.6K vertices with our optimized algorithm, compared to 1.31 and 11.88 seconds for marching cubes with 256 and 512 samples per axis, respectively. We note that parallelizing over multiple GPUs would be a reliable option to improve latency. For the discussions on satisfying the eikonal equation in practice and the generalization to other positional encodings using trilinear interpolation, please refer to Appendix C.1 and Appendix C.2, respectively. 
conclusions: 
 In conclusion, our exploration of novel techniques in visualizing trilinear neural networks unveils its intricate structures and extends mesh extraction applications beyond CPWA functions. Focused on trilinear interpolating methods, our study yields novel theoretical insights and a practical methodology for precise mesh extraction. Notably, Theorem 4.5 shows the transformation of hypersurfaces to planes within the trilinear region under the eikonal constraint. Empirical validation using chamfer distance and efficiency, and angular distance affirms the correctness and parsimony of our methodology. The correlation analysis between eikonal loss and hypersurface planarity confirms the theoretical findings. Our proposed method for approximating intersecting points among three hypersurfaces, coupled with a vertex arrangement methodology, broadens the applications of our techniques. This work establishes a foundation for future research and practical applications in the evolving landscape of deep neural network analysis and toward real-time rendering mesh extraction. A Tropical geometry and tropical algebra of neural networks 
experiments compute resources: 
 Question: For each experiment, does the paper provide sufficient information on the computer resources (type of compute workers, memory, time of execution) needed to reproduce the experiments?Answer: [Yes] Justification: Appendix F.1 provides the information on the computer resources.Guidelines:• The answer NA means that the paper does not include experiments.• The paper should indicate the type of compute workers CPU or GPU, internal cluster, or cloud provider, including relevant memory and storage. • The paper should provide the amount of compute required for each of the individual experimental runs as well as estimate the total compute. • The paper should disclose whether the full research project required more compute than the experiments reported in the paper (e.g., preliminary or failed experiments that didn't make it into the paper). 
code of ethics: 
 Question: Does the research conducted in the paper conform, in every respect, with the NeurIPS Code of Ethics https://neurips.cc/public/EthicsGuidelines? 
broader impacts: 
 Question: Does the paper discuss both potential positive societal impacts and negative societal impacts of the work performed? 
safeguards: 
 Question: Does the paper describe safeguards that have been put in place for responsible release of data or models that have a high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)? Answer: [NA] Justification: This work poses no such risks. Guidelines:• The answer NA means that the paper poses no such risks.• Released models that have a high risk for misuse or dual-use should be released with necessary safeguards to allow for controlled use of the model, for example by requiring that users adhere to usage guidelines or restrictions to access the model or implementing safety filters. • Datasets that have been scraped from the Internet could pose safety risks. The authors should describe how they avoided releasing unsafe images. • We recognize that providing effective safeguards is challenging, and many papers do not require this, but we encourage authors to take this into account and make a best faith effort. 12. Licenses for existing assets Question: Are the creators or original owners of assets (e.g., code, data, models), used in the paper, properly credited and are the license and terms of use explicitly mentioned and properly respected? Answer: [Yes] Justification: The Stanford 3D Scanning repository (MIT License) is used with a proper citation [31] , along with HashGrid [12] for our trilinear module. Guidelines:• The answer NA means that the paper does not use existing assets.• The authors should cite the original paper that produced the code package or dataset.• The authors should state which version of the asset is used and, if possible, include a URL.• Depending on the country in which research is conducted, IRB approval (or equivalent) may be required for any human subjects research. If you obtained IRB approval, you should clearly state this in the paper. • We recognize that the procedures for this may vary significantly between institutions and locations, and we expect authors to adhere to the NeurIPS Code of Ethics and the guidelines for their institution. • For initial submissions, do not include any information that would break anonymity (if applicable), such as the institution conducting the review. 
 
 Appendix: In this section, we present an overview of tropical geometry (Appendix A.1) to provide a formal definition of the tropical hypersurface (Appendix A.2). This hypersurface is significant as it serves as the decision boundary for ReLU neural networks optimizing an SDF, resulting in the formation of a polyhedral complex. For this, we delve into the examination of the tropical algebra of neural networks (Appendix A.3). Furthermore, we discuss the edge subdivision algorithm [4] for extracting the tropical hypersurface (Section 3.2), and we extend this discussion to include considerations for trilinear interpolation in Section 4. Given that Appendix A.1 to A.3 offer insights into tropical geometry, illustrating how each neuron contributes to forming the polyhedral complex, you may choose to skip these subsections if you are already acquainted with these concepts.Tropical geometry [25, 26] describes polynomials and their geometric characteristics. As a skeletonized version of algebraic geometry, tropical geometry represents piecewise linear meshes using tropical semiring and tropical polynomial.The tropical semiring is an algebraic structure that extends real numbers and -∞ with the two operations of tropical sum ⊕ and tropical product ⊙. Tropical sum and product are conveniently replaced with max(•) and +, respectively. Notice that -∞ is the tropical sum identity, zero is the tropical product identity, and these satisfy associativity, commutativity, and distributivity. For example, a classical polynomial x 3 + 2xy + y 4 would represent max(3x, x + y + 2, 4y). The classical polynomial is rewritten for a tropical polynomial by naively replacing + with ⊕ and • with ⊙, respectively. After this, we rewrite by the definitions of tropical operations as max(3x, x+y +2, 4y).Roughly speaking, this is the tropicalization of the polynomial, denoted tropical polynomial or simply f . Note that this example (ours with max notation) comes from the Wikipedia 6 .Additionally, we present a graphical diagram in Figure 7 illustrating how the tropical polynomial divides the 2D space for clearer understanding. Please see the caption for further details. The set of points where a tropical polynomial f is non-differentiable, mainly due to tropical sum, max(•), is called tropical hypersurface, denoted V (f ). Note that V is an analogy to the vanishing set of a polynomial.Definition A.1 (Tropical hypersurface, restated). The tropical hypersurface of a tropical polynomial f (x) = c 1 x α1 ⊕ • • • ⊕ c r x αr is defined as:V (f ) = {x ∈ R d : c i x αi = c j x αj = f (x) for some α i ̸ = α j }where the tropical monomial c i x αi with d-variate x is defined as , a i, d ) and tropical product ⊙ := +. In other words, this is the set of points where f (x) equals two or more monomials in f .c i ⊙ x ai,1 1 ⊙ • • • ⊙ x a i,d d with d-variate α i = (a i,1 , . . .Definition A.2 (ReLU neural networks, restated). Assuming input and output dimensions are consistent, L-layer neural networks with ReLU activation are a composition of functions defined as:EQUATIONwhere σ(x) = max(x, 0) andρ (i) (x) = W (i) x + b (i) .Although the neural networks ν is generally non-convex, it is known that the difference of two tropical signomials can represent ν [27] . A tropical signomial takes the same form as a tropical polynomial (ref. Definition 3.1), but it additionally allows the exponentials α to be real values rather than integers. With tropical signomials F , G, and H,EQUATIONEQUATIONwhere F (l) , G (l) , and H (l) are as follows:EQUATIONEQUATIONEQUATIONwhileEQUATION-where W+ and W-∈ R + . Proposition A.3 (Decision boundary of neural networks, restated). The decision boundary for Llayer neural networks with ReLU activation ν (L) (x) = 0 is where two tropical monomials, or two arguments of max, equal (Definition 3.1), as follows:B ⊆ x : ν (i) j (x) = H (i) j (x) -G (i) j (x) = 0, i ∈ (1, . . . , L), j ∈ (1, . . . , H) if i ̸ = L else (1)where the subscript j denotes the j-th neuron, assuming the hidden size of neural networks is H, while the L-th layer has a single output for rather discussions as an SDF. Remind that tropical operations satisfy the usual laws of arithmetic for recursive consideration.This implies that piecewise hyperplanes corresponding to the preactivations {ν (i) j } are the candidates shaping the decision boundary. Each linear region has the same pattern of zero masking by ReLU(•), effectively making it linear as an equation of a plane Ax + b = 0 for some A and b. In turn, if the masking pattern changes as an input is continually changed, it moves to another linear region. From a geometric point of view, this is why the number of linear regions grows polynomially with the hidden size H and exponentially with the number of layers L [18, 27] .B is zero-trailing binarization function, consisting of D = 3 bits. The zero-trailing means the number in bits is increasing left-aligned. For example, the number one is "100", while the numbers two and four are "010" and "001", respectively. The number three is "110". From the left, each bit indicates x, y, and z axis in our notation. This convention seems to be widely used in the CUDA implementation in C language. Notice that we use the j subscript of B j to select the axis. Naturally, i = 0 and i = 7 indicate the opposite two corner points.Algorithm 3 describes the algorithm in Section 4.3, which are sorting polygon vertices implicitly representing the direction of faces. If the viewing direction aligns opposite to the normal, the order of vertices follows a counter-clockwise arrangement. We determine the vertex order based on the normal derived from the Jacobian of the mean of vertices. We identify the counter-clockwise arrangement of vertices using dot product to compute relative angles θ ∈ [-π, π] and cross product for getting direction. Algorithm 4 describes the algorithm in Section 5, calculating the HashGrid marks from the hyperparameters of a HashGrid.Input: vertices x i , normal n, N = |{x i }| µ = i x i /N for i = 0 to N -1 do v i = (x i -µ)/∥x i -µ∥ 2 c i = ⟨v i , v 0 ⟩, d i = ⟨v i × v 0 , n⟩ θ i = c i (2 • 1 di≥0 -1) + 2 • 1 di<0 end for return ArgSort({-θ i }) Algorithm 4 HashGrid Marks Input: Resolution levels L, base resolution N min , max resolution N max m = {0, 1} b = log 2 (N max /N min )/(L -1) for l = 0 to L -1 do s = 1/ N min • exp 2 (l • b) -1 m ← m + Max Arange(-s/2, 1, s), 0 end for return Sort Unique(m)Theorem 4.7 handles where the eikonal equation does not make perfect hyperplanes. First, we observed the characteristics of trilinear hypersurfaces as illustrated in Figure 9 in Appendix. Among 64 edge scenarios within a trilinear region, a diagonal plane intersects with the hypersurface in nearly all cases. Moreover, we can take advantage of it since we can eliminate one variable by assuming the curved edge lies on the x = z plane. In this light, Theorem 4.7 used linear algebra to rearrange two trilinear equations to get Equations ( 71) and (74). As we mentioned in the manuscript, the new vertices lie on at least two hypersurfaces, and the new edges exist on the same hypersurface, while the eikonal constraint minimizes any associated error. This error is tolerated by the hyperparameter ϵ=1e-4 (Section 6) as specified in Definition 3.4.Figure 5 also suggests empirical evidence. Empirically, the flatness error, derived from the proof of Theorem D.5, is effectively controlled by the learning rate of the eikonal loss. In the experiments, we used 0.01 to give our results, while many other neural implicit surface learning methods were used to exploit a higher learning rate of 0.1 (e.g., Yariv et al. [7] ). Note that, in practice, the eikonal loss is only applied to near-surface using clamps (if the target position is obviously far from a surface, ignore it; it is wise since we are concerned with the surfaces, not wanting exact SDFs for all space) for effective learning, which is sufficient to our extraction method.Moreover, the smaller grid size of the Large model (ref. Section 6) in Figure 4 views a curved surface in a smaller interval, and it tends to regard it as a plane. From the plotting of Small to Large, our CD (colored crosses) is decreasing to zero, consistently with better CE, confirming this speculation.We can say that any encoding using trilinear interpolation (that) we defined can be applied in our method. Here, the definition is implicitly referred to Definition 4.1 (Trilinear interpolation), where the cubic corners are represented by H ∈ R D×F . Since we do not rely on a specific mapping function f : x → (H 0 , H 1 . . . H 7 ), we may apply our extraction method to any encoding using trilinear interpolation that we defined. Specifically, InstantNGP [12] used a hashing function to map a corner point to the arbitrary entry of hash tables, allowing gracious collisions with online adaptivity (refer to their Sec. 3), while TensoRF [13] can be interpreted that it used a mapping function to the summation of the multiplications of two vectors from the vector-matrix factorization. Therefore, if we have the full ranks for R 1 , R 2 , R 3 , in Eqn. 3 in Chen et al. [13] , among the VM factorization, which can be equivalently represented by the HashGrid with a sufficient number of hash table sizes.Then, the discussion should be directed toward their encoding characteristics, which impact our extraction method due to their online adaptivity (HashGrid) or low-rank factorization (TensoRF). When we train an SDF, the optimization to global minima is generally difficult to achieve. Instead, previous works on neural implicit surface learning methods resort to optimizing the SDF on nearsurface since we are interested in modeling the zero-set surfaces, allowing some errors for the points that are far from the surfaces. Practically, this is done by defining a loss function using clamps and ignoring far-surface points. For this reason, the online adaptivity of HashGrid or low-rank factorization of TensoRF can exploit this nature: focused learning of the surfaces adaptively allocating their learnable parameters to the near-surface, which means there might be enough room for the regularization by the eikonal equation. The investigation of the efficiency of representing the implicit neural surfaces using various positional encodings is seemingly out of the scope of this work since our main interest lies in the theoretical exposition of the mesh extraction from it and setting up a feasible bridge between theory and practice.Lemma D.1. (Nested trilinear interpolation, restated) Let a nested cube be inside of a unit cube, where its positions of the eight corners deviate -a j or +b j , such that a j , b j ≥ 0, from w j for each dimension j, using the notations of Definition 4.1. The eight corners of the nested cube are the trilinear interpolations of the eight corners of the unit cube. Then, the trilinear interpolation with the unit cube and the nested cube for w is identical.Proof. Let a left-aligned and zero-trailing binarization function be B(•) ∈ {0, 1} D and D = 3. The partial derivative of the trilinear interpolation with the nested cube Ĥ with respect to H 7 and H 0 would be:EQUATIONEQUATIONandEQUATIONEQUATIONrespectively. In a similar way, we can describe the other cases H 1 to H 6 of the nested trilinear weights since there are two cases for each dimension inside the product.Proposition D.2 (Cubic Bézier curve). Trilinear interpolation from Definition 4.1 with the weights on the diagonal line t = (t, t, t) ⊺ where t ∈ [0, 1] is a cubic Bézier curve.Proof. We rewrite the trilinear interpolation from Definition 4.1 with the weight t as follows:EQUATIONFour points P 0 := H 0 , P 1 := H1+H2+H4, P 2 := H3+H5+H6, and P 3 := H 7 define the cubic Bézier curve.As t is changed from zero to one, the curve starts at P 0 going toward P 1 and turns in the direction of P 2 before arriving at P 3 . Usually, the curve does not pass P 1 nor P 2 ; but these give the directional guide for the curve.Lemma D.3 (Hypersurfaces in a trilinear space). In a piecewise trilinear region where x lies, hypersurfaces of the piecewise trilinear networks are defined as follows, following the notation of Definition 4.1:EQUATION)where w ∈ [0, 1] D is the relative position in a grid cell and {x i } and {H i } are the corresponding corner positions and representations, respectively, in a trilinear space.Proof. In a piecewise trilinear region, we can use the linearity of ν as follows:EQUATIONwhich concludes the proof.Lemma D.4 (Curved edge of two hypersurfaces, restated). In a piecewise trilinear region, let an edge (x 0 , x 7 ) be the intersection of two hypersurfaces νi (x) = νj (x) = 0, while x 0 and x 7 are on the two hypersurfaces. Then, the edge is defined as:EQUATION)Proof. By Definition 4.3,EQUATIONwhile ν i is linear in a piecewise trilinear region. Since it similarly goes for ν j , the intersection line of two hyperplanes ν i and ν j is (1 -t)τ (x 0 ) + tτ (x 7 ) where t ∈ R. In other words, the edge is on the hypersurface ν i (1 -t)τ (x 0 ) + tτ (x 7 ) = 0 where t ∈ R. The edge is defined as follows:EQUATION)which concludes the proof.Notice that when τ (x) 1 + τ (x) 2 + τ (x) 4 = τ (x) 3 + τ (x) 5 + τ (x) 6 = 0, at least one connecting edge can be identified. This occurs along the diagonal line x = t • 1, where t ∈ [0, 1]. Along this line, τ (x) represents a cubic Bézier curve, as shown in Proposition D.2 with the control points of P 1 and P 2 set to zeros, which is the line.Theorem D.5 (Hypersurface and eikonal constraint, restated). A hypersurface τ (x) = 0 intersects two points τ (x 0 ) = τ (x 7 ) = 0 while τ (x 1...6 ) ̸ = 0 for the remaining six points. These points form a cube, with x 0 and x 7 positioned on the diagonal of the cube. The hypersurface satisfies the eikonal constraint ∥∇τ (x)∥ 2 2 = 1 for all x ∈ [0, 1] 3 . Then, the hypersurface of τ (x) = 0 is a plane.Proof. By the definitions of τ in Definition 4.3 and the eikonal constraint,EQUATIONwhereEQUATIONHowever, since τ (x 0 ) = τ (x 7 ) = 0 while τ (x 2 ) ̸ = 0 and τ (x 5 ) ̸ = 0, there is no solution τ (x 1...6 ) of ∂τ (x) ∂y = 0 for all x. Likewise, the same goes for ∂τ (x) ∂x and ∂τ (x) ∂z . Additionally, the three partial derivatives can be rewritten as follows:EQUATION)where a ̸ = 0, b ̸ = 0, and c ̸ = 0. By substituting B = -bC/c and A = -aC/c, we can show that -abC/c -abC/c = 0 ⇔ C = 0, which makes A = B = C = 0. Therefore, we rewrite A, B, and C as follows:EQUATIONEQUATIONEQUATIONrespectively. The following conditions would satisfy the eikonal constraint:-τ (x 2 ) -τ (x 4 ) + τ (x 6 ) = 0 (46) τ (x 1 ) -τ (x 3 ) -τ (x 5 ) = 0 (47) -τ (x 1 ) -τ (x 4 ) + τ (x 5 ) = 0 (48) +τ (x 2 ) -τ (x 3 ) -τ (x 6 ) = 0 (49) -τ (x 1 ) -τ (x 2 ) + τ (x 3 ) = 0 (50) +τ (x 4 ) -τ (x 5 ) -τ (x 6 ) = 0 (51)From Equation (48) and Equation (51),EQUATIONFrom Equation (47) and Equation (50),EQUATIONFrom Equation (47) and Equation (48),EQUATIONFrom Equation (50), Equation (51), and Equation (54),τ (x 1 ) + τ (x 2 ) + τ (x 4 ) = 0 (55) τ (x 3 ) + τ (x 5 ) + τ (x 6 ) = 0 (56)Using Equation (52), Equation (53), Equation (54), Equation (55), and Equation (56),τ (x) = x(1 -y)(1 -z) -(1 -x)yz τ (x 1 ) + (1 -x)y(1 -z) -x(1 -y)z τ (x 2 ) + (1 -x)(1 -y)z -xy(1 -z) τ (x 4 ) (57) = (x + α)τ (x 1 ) + (y + α)τ (x 2 ) + (z + α)τ (x 4 ) (58) = xτ (x 1 ) + yτ (x 2 ) + zτ (x 4 ) + α τ (x 1 ) + τ (x 2 ) + τ (x 4 ) (59) = xτ (x 1 ) + yτ (x 2 ) + zτ (x 4 ) = 0 (60) where α = 2xyz -xy -yz -zx.This is an equation of a plane where its normal vector is τ (x 1 ), τ (x 2 ), τ (x 4 ) . Again, to satisfy the eikonal constraint, τEQUATION) which concludes the proof.Corollary D.6 (Affine-transformed hypersurface and eikonal constraint, restated). Let ν 0 (x) and ν 1 (x) be two affine transformations defined byW ⊺ i x + b i , i ∈ {0, 1}. A hypersurface ν 0 τ (x) = 0 passing two points ν 0 τ (x 0 ) = ν 0 τ (x 7 ) = 0 while ν 0 τ (x 1...6 ) ̸ = 0 satisfies the eikonal constraint ∥∇ν 1 τ (x) ∥ 2 2 = 1 for all x ∈ [0, 1] 3 .Then, the hypersurface is a plane.Proof. Similarly to Theorem D.5, we investigate the Jacobian of the eikonal constraint to be zero as follows:EQUATION)Since ∂ν 1 /∂τ is a constant, using Lemma D.3 and replacing ν 1 and τ with ν 1 • ν -1 0 and ν 0 • τ , respectively, give us the same constaints for ν 0 τ (x 1...6 ) regardless of ν 1 .Theorem D.7 (Intersection of two hypersurfaces and a diagonal plane, restated). Let ν0 (x) = 0 and ν1 (x) = 0 be two hypersurfaces passing two points x 0 and x 7 such that ν0 (x 0 ) = ν1 (x 0 ) = 0 and ν0 (x 7 ) = ν1 (x 7 ) = 0, P i := ν0 (x i ) and Q i := ν1 (x i ), P α = P 0 ; P 1 ; P 4 ; P 5 , P β = P 2 ; P 3 ; P 6 ; P 7 , (65) andEQUATION)Then, x ∈ [0, 1] of the intersection point of the two hypersurfaces and a diagonal plane of x = z is the solution of the following quartic equation:EQUATION)Proof. We rearrange ν0 (x) using x = z as follows:ν0 (x) = (1 -y) • X ⊺ P 0 ; P 1 ; P 4 ; P 5 + y • X ⊺ P 2 ; P 3 ; P 6 ; P 7 (69) = (1 -y)X ⊺ P α + yX ⊺ P β = 0 (70)EQUATIONIn a similar way, we rewrite ν1 (x) as follows:EQUATIONNotice that this is a quartic equation of x. To find the coefficients of the general form of quadratic equation, we consider a mapping C from the interpolation coefficients τ i to polynomial coefficients c i :EQUATIONwhereτ 0 (1 -x) 2 + τ 1 x(1 -x) + τ 2 (1 -x)x + τ 3 x 2 = c 0 + c 1 x + c 2 x 2 .Then, the coefficients of the quartic equation of Equation ( 74) are as follows:EQUATIONwhere the anti-diagonal sum of the result is the quartic coefficients, i.e.,EQUATION)The solutions of the quartic equation such that x ∈ [0, 1] are the valid intersection points, and we choose one of them x if there are more than one solution. Note that we can compute a batch of quartic equations using the parallel computing of the eigenvalue decomposition with Lemma D.8.We determine y by Equation (71) using x, and z = x.We reproduce the method to find polynomial roots from companion matrix eigenvalues [29] . Lemma D.8 (Polynomial roots from companion matrix eigenvalues, reproduced). The companion matrix of the monic polynomial7 EQUATIONis defined as:EQUATION)The roots of the characteristic polynomial p(x) are the eigenvalues of C(p).We provide a detailed exposition of the complexity analysis in Section 3.2 to ensure the comprehensiveness of our work. Notice that our complexity analysis closely aligns with the approach outlined in Berzins [4] , as in their Appendix A.Algorithm 2 operates on the vertices and edges, V and E, respectively. Therefore, our initial focus is on the complexity associated with these inputs. Let |V|, |E|, and | Ê| ( Ê ⊂ E) be the number of vertices, edges, and splitting edges consisting of polygons (Section 3.2), at the iteration of i and j.1. During the initialization step (refer to Section 5), we determine |V| and |E| to construct a grid in R D . In the case of having M HashGrid marks, the quantities are defined as M D and D(M -1)M D-1 , respectively. In case of the large M , we consider vertex and edge pruning for efficient computations inspired by a pruning strategy proposed in Berzins [4] . If, for every future neuron, denoted as ϕ(i, j) > ϕ(i * , j * ), the two vertices of an edge exhibit identical signs, the edge can be safely pruned as it will not undergo splitting. We make the assumption of that M D < |V| to avoid an excessively dense grid, with a complexity of O(|V|).2. Evaluate ν(V; i, j) with a complexity of O(|V|), assuming negligible inference costs as constants. The intermediate results from the inference can be used for the sign-vectors (Definition 3.4).3. Identify splitting edges Ê by comparing two signs per edge, with a complexity of O(|E|).4. Obtain P and Q for Theorem 4.7 with a complexity of O(| Ê|) in a similar way. 5. Inference for Theorem 4.7 involves matrix-vector multiplications and the eigenvalue decomposition of a companion matrix for each edge, with a complexity of O(| Ê|). Considering the maximum size of 4 for a quartic equation, this operation is regarded as linear. 6. Find intersecting polygons by pairing two sign-vectors with two common zeros in a shared region of the new vertices resulting from splitting, with a complexity of O(| Ê|).In the study by Berzins [4] , a linear relationship between |V| and |E| for a mesh is empirically observed (see Figure 8a in Berzins [4] ). It is worth to note that the number of splitting edges can be effectively upper-bounded as | Ê| < |E|D/ϕ(i, j) using Theorem 5 in [20] . Consequently, all the steps outlined in Algorithm 2 demonstrate a linear dependency on the number of vertices, denoted as O(|V|).To obtain faces, we sort the vertices within each region based on the Jacobian with respect to the mean of the vertices (see Algorithm 3). Considering triangularization, where each face has three vertices, the complexity remains O(|V|), assuming the negligible cost of calculating the Jacobian, especially when performed in parallel.Table 2 and Table 3 are extended from the Stanford 3D Scanning dataset in Table 1 . Table 4 extended from the Stanford bunny of Table 1 , including standard deviations, and pseudo-ground truth results, obtained through marching cubes with 256 samples. Notice that we used the more excessive 512 samples for the Medium and Large models in Figure 4 . We found that 256 samples were not enough for those. Table 5 shows the angular distance results. Tables 2 to 5 use the Small model, while Table 8 uses the Large model (ref. Section 6) . We conducted all experiments with an NVIDIA V100 32GB.Table 2 : Chamfer distance (CD) and chamfer efficiency (CE) for the Part A of the Stanford 3D Scanning repository [31] using the Small model as the default setting (ref. Section 6) . The chamfer distance (×1e-6) is evaluated using the marching cubes (MC) with 256 3 grid samples as the reference ground truth for all measurements. The chamfer efficiency is defined as 100/( |V| × CD), providing a concise representation of the trade-off between the number of vertices and the chamfer distance.The results are averaged over three random seeds. Both benefit from QEM while our method efficiently retains better chamfer distances.Figure 8 : Chamfer distance for the Bunny with the Large model applying QEM (a mesh simplification). Both benefit from QEM while our method efficiently retains better chamfer distances.F.2 Mesh simplification using QEMIn Figure 8 , Quadric Error Metrics (QEM) [35] is a popular mesh simplification algorithm, it could potentially improve both MC and our method for efficiency. Overall, as shown, the mesh simplification favors our method. Notice that QEM entails a higher computational cost, at least O(n log n) where n is the number of vertices. Additionally, these methods can cause shape deformation or holes, sensitive to hyperparameters. Therefore, QEM may be applied with caution for further mesh simplification, acknowledging the potential for unpredictable side effects. We used MeshLab's Quadric Edge Collapse Decimation with the default settings to halve the number of vertices for each step.Marching Cubes (MC) remains the de facto sampling-based method. In contrast, our method is a white-box method grounded in the theoretical understanding of trilinear interpolation and polyhedral complex derivation. Here is our rationale for this evaluation:1. Analytical approaches for CPWA functions are limited by their inability to effectively manage spectral bias and their exponential computational cost, which struggles with exponentially growing linear regions (infeasible to run our benchmarks.) For example, the time complexity of Analytic Marching [2] is O((n/2) 2(L-1) |V P |n 4 L 2 ), which grows exponentially with the number of layers L and the width n of the networks, where |V P | represents the number of vertices per face (Section 4.1 from Lei and Jia [2] ). In contrast, our method has a linear time complexity with respect to the number of vertices, as discussed in Appendix E, allowing it to avoid visiting every linear region exponentially growing.2. We also exclude optimization-based methods as they typically rely on MC or its variants initialization in their pipelines (our method can also be integrated into them), such as MeshSDF [36] , Deep Marching Tetrahedra [10] , and FlexiCubes [37] . Additionally, these methods result in prolonged computational costs.We provided both rigorous quantitative Figure 3 and qualitative Figure 6 analyses on Marching Tetrahedra (MT) [9, 10] and Neural Dual Contour (NDC) [32] , which serve as alternative representative approaches. MT is a popular variant of MC that utilizes six tetrahedra within each grid to identify intermediate vertices. NDC is a data-driven approach (a pretrained model) that uses Dual Contour.In short, MT efficiently requires SDF values but is less efficient regarding the number of vertices extracted. (MT is better than MC in the same resolution, producing more vertices, though.) NDC faces a generalization issue for a zero-shot setting, producing unsmooth surfaces. We used the code from Deep Marching Tetrahedra [10] for MT and the official code of NDC. And PyMCubes for our modernized MC library.In Equation (11) and Figure 5 , we confirmed our hypothesis that the eikonal loss induces piecewise linearity in trilinear networks. To further support this, we also present the error analysis of the learned SDF. Although the Chamfer distance virtually evaluates the SDF values based on sampling, concerns may arise regarding noise in the pseudo-ground truth. For this reason, one might wish to verify whether the predicted SDFs converge toward near-zero.We randomly sampled 100K points on the surface of pretrained networks (the Stanford bunny Large model) and calculated the average of squares for predicted SDF (model's outputs) on the surface. As we expected, the predicted SDFs on the reconstructed surface converge to zero as the weight of the eikonal loss increases, and the number is close to zero ≈ 3e-8 as shown in Table 9 . Note that, for a Large model, fine-grained trilinear regions from the dense grid marks (from the more dense multi-resolution grids) allow relatively small errors of 97 × 1e-9 without the eikonal loss, although the eikonal loss of 1e-2 further decreases to 31 × 1e-9. Note that the stronger eikonal loss of 1e-1 proved ineffective due to over-regularization. In Table 10 , we explore the ablation study to determine whether model size impacts model fitting, which in turn affects the GT mesh (we used MC256 for the pseudo-ground truth of the Small models). With a smaller model size, the learned mesh becomes simpler due to underfitting. As our evaluation focuses on how accurately the extraction methods capture the underlying zero level-set decision boundary of the networks, the change in CD remains small. Notice that the model with the number of layers of 3 and the width of 64 has the longest runtime of 2.80, primarily because the number of intermediate vertices exceeds 80K. The time complexity is more sensitive to the number of vertices than to the number of layers or the width, which are only indirect factors. Note that the model capacity of the decoding networks is comparable to InstantNGP [12] , where the representational efficiency is driven by the HashGrid. In a face-to-face comparison, we examine marching cubes with 256 and 64 samplings, denoted as MC 256 and MC 64, respectively, and our method. Notably, MC 64 struggles to accurately represent pointy areas of meshes due to its coarse sampling. In contrast, our approach has an advantage in this aspect, as we analytically extract the mesh from the decision boundary apex points of networks. However, it is worth noting that our adaptive sampling introduces a disadvantage where some vertices are in close proximity; specifically, we observed that 10.1% of vertices are within 1e-3 proximity to others. We believe that optimizing the mesh can significantly improve the parsimony of the vertex count compared with marching cubes. Justification: We conform with the NeurIPS Code of Ethics.Guidelines:• The answer NA means that the authors have not reviewed the NeurIPS Code of Ethics.• If the authors answer No, they should explain the special circumstances that require a deviation from the Code of Ethics. • The authors should make sure to preserve anonymity (e.g., if there is a special consideration due to laws or regulations in their jurisdiction).Justification: There is no societal impact of the work performed. Guidelines:• The answer NA means that there is no societal impact of the work performed.• If the authors answer NA or No, they should explain why their work has no societal impact or why the paper does not address societal impact. • Examples of negative societal impacts include potential malicious or unintended uses (e.g., disinformation, generating fake profiles, surveillance), fairness considerations (e.g., deployment of technologies that could make decisions that unfairly impact specific groups), privacy considerations, and security considerations. • The conference expects that many papers will be foundational research and not tied to particular applications, let alone deployments. However, if there is a direct path to any negative applications, the authors should point it out. For example, it is legitimate to point out that an improvement in the quality of generative models could be used to generate deepfakes for disinformation. On the other hand, it is not needed to point out that a generic algorithm for optimizing neural networks could enable people to train models that generate Deepfakes faster. • The authors should consider possible harms that could arise when the technology is being used as intended and functioning correctly, harms that could arise when the technology is being used as intended but gives incorrect results, and harms following from (intentional or unintentional) misuse of the technology. • If there are negative societal impacts, the authors could also discuss possible mitigation strategies (e.g., gated release of models, providing defenses in addition to attacks, mechanisms for monitoring misuse, mechanisms to monitor how a system learns from feedback over time, improving the efficiency and accessibility of ML).In geometry, a polyhedron is a three-dimensional shape with polygonal faces, having straight edges and sharp vertices.https://github.com/NVlabs/tiny-cuda-nnCE is defined as 100/( |V| × CD), reflecting the trade-off between the number of vertices and the CD.https://en.wikipedia.org/wiki/Tropical_geometryThe nonzero coefficient of the highest degree is one. If you are dealing with the equation, you can make it by dividing the nonzero coefficient of the highest degree for both sides. 


Decision: The paper proposes a new method for mesh extraction from piecewise trilinear networks, grounded in a theoretical grounding in trilinear interpolation and polyhedral complex derivation.  The method is shown to be more accurate and efficient in comparison to MC and other alternative approaches.

The paper received four reviews with unanimity in their (post discussion) ratings. The reviewers carefully engaged with the paper and the author response and participated in ensuing discussions. They found the paper to be well motivated, and a reasonable, theoretically grounded, extension of previous work on complex polyhedral complex extraction from ReLU networks to architectures with positional encodings implemented using trilinear interpolation. The reviewers particularly appreciated the theoretical grounding of the proposed approach, with interesting insights, and holding a lot of promise. 

The initial reviews were mixed around the borderline. The primary concerns are around (a) the organization of the paper “focusing on lower abstraction details”, lacking supporting figures and visualizations compromising readability, and making the paper difficult to understand; (b) for a mathematical paper, riddled with typos, imprecise statements, poor phrasing, and undefined terms (definition before usage); (b) with unconvincing experimental evaluation – both qualitatively and quantitatively; inadequate comparison with alternate approaches, and, (d) an inadequate limitations section. 

The authors addressed many of these concerns during the discussion resulting in two reviewers revising their ratings upward, resulting in a unanimous rating of ‘borderline accept’ from all reviewers. I recommend that the paper be accepted subject to the authors incorporating the clarifications and results from their rebuttal and responses, as well as other promised changes improving the presentation of the paper.

Reviews: review: paper_summary: paper_summary: The work addresses the problem of exact iso-surface extraction given an implicit surface represented by a neural network. The idea behind this line of work is to closely analyze the architecture of the neural network to efficiently extract the true (or at least an approximate of) learned iso-surface. Previous work titled 'Polyhedral Complex Extraction from ReLU Networks using Edge Subdivision' develops the analysis for ReLU networks. However, many newer works use local latent codes to achieve better reconstruction quality. These latent codes are usually obtained via trilinear interpolation from a grid of learnable latent codes. The previous method does not trivially apply to iso-surface extraction in this case. The current work continues the analysis for ReLU networks that now have trilinear positional encoding as input. The paper first shows that under the eikonal constraint, trilinear regions decode into flat planes. Then, it derives a polynomial defining the intersection of two hypersurfaces. With this, the paper suggests an algorithm for zero-level set extraction given the trained neural network and the positional encoding.
strength_weakness: strength_weakness: Here are the main concerns I found in the paper in order of importance (with the first being the most important):

* I found the paper very hard to understand. I am quite familiar with neural implicit representations (though not with Polyhedral Complex Extraction), but it still took me at least three full passes before I started understanding the general idea of the paper. The paper is self-contained and formal, which is great, but I would consider simplifying the language for the sake of understandability.

* The work is positioned as a better way of iso-surface extraction, yet there is no visualization of extracted meshes in the main paper. It is very hard to say if the mesh extracted with the proposed method is actually better than the ones extracted with Marching Cubes. In my view, MC64 seems to be actually closer to MC256 than the mesh extracted with 'ours' in Figure 7. A much more detailed visual comparison of extracted meshes seems necessary in the main paper given the claims.

* The method is shown to work only with relatively small neural networks: the experiments are performed on a model with three layers and a hidden size of 16, while DeepSDF, for instance, operates on eight layers, each with a width of 512. It seems that it can be used only as a source of understanding the behavior of these architectures, not as a practical method of iso-surface extraction.

* Comparison with only Marching Cubes seems very limited; it would be great to see a comparison with 'Analytic Marching' as well.

I will consider adjusting the final rating if weakness #2 and question #1 are addressed in the rebuttal
questions: questions: * I am quite confused about the following (probably due to my limited understanding of the math in the paper): A standard ReLU network with a constant latent code as input (i.e., DeepSDF) can be seen as a trivial trilinear positional encoding (i.e., the entire domain is trilinear). Given that such networks are usually trained to satisfy the eikonal constraint, why aren’t they limited to always produce flat plane reconstructions, given Corollary 4.6?

* It seems that a natural measure of how accurately the algorithm extracts the learned iso-surface would be:

        1. Sample N points on the surface of the extracted mesh
        2. Predict SDF on each of them
        3. Calculate sum of squared SDF values
it should be 0 if we extracted the learned surface precisely. Would be interesting to see this metric. Can we confirm Corollary 4.6 with it?
limitations: limitations: The limitations are well discussed. I would also add a discussion on the depth and width of the neural networks for which the method is practical to use.
paper_summary: paper_summary: The paper presents a theoretical and practical framework to extract meshes from piecewise trilinerar networks using hash grid representation. They further show that hypersurfaces within the trilinear regions become planes under the eikonal constraint. A method for approximating the intersection between 3 hypersurfaces is proposed (approximated as 2 hypersurfaces and a diagonal). Evaluation is performed against marching cubes resolutions using chamfer distance, chamfer efficiency and angular distance.
strength_weakness: strength_weakness: - While there are some improvements over marching cubes, the improvements are not quite substantial, and come at the cost of a higher runtimes
- No comparisons to any alternate approaches over marching cubes are provided -- e.g. marching cubes + QEM decimation etc.
- No rigorous qualititative comparisons are provided, which could be important since chamfer does not always capture the shape quality
questions: questions: Would be nice if the authors can add some qualitative comparisons and also against more than just marching cubes. Further current evaluation is only performed on a handful of meshes, would be nice if an evaluation on a bigger dataset can be performed
limitations: limitations: Yes
paper_summary: paper_summary: The paper extends the idea of extracting the polyhedral complex from ReLU networks to extracting a more general cell complex from ReLU networks with a positional encoding, specifically, trilinear encodings, such as HashGrid and TensoRF, commonly used in neural field representations of geometries or scenes.
It first theoretically extends the concept of piecewise linear to piecewise trilinear networks. It then proposes a modified *curved* edge subdivision algorithm to extract a polyhedral mesh approximating the neural implicit surface. The key assumption governing the approximation error is the eikonal loss, which ensures that the isosurfaces of the trilinear functions are close to planar.
The algorithm is validated on several shapes qualitatively and quantitatively, comparing to marching cubes in terms of the mesh size, accuracy, and time. It demonstrates sparser meshes at an increased computational cost.
strength_weakness: strength_weakness: Despite the method being able to handle volumetric complex of the input space, the the work only considers the "skeletonized" complex, i.e., the boundary mesh. At least briefly describing the former would make for a more rounded story. Instead, the skeletonized method is framed as a potential alternative to marching cubes, which unfortunately is not convincing from the results. The only demonstrated advantage is in terms of the mesh size. This argument would be more convincing given a comparison to marching cubes with some basic mesh post-processing. The work would also benefit from comparing the same procedure to piecewise linear networks and the underlying edge subdivision method in terms of training time, extraction time, mesh size, and accuracy.

Similarly, while the method is evaluated for different positional encodings, the neural network size remains fixed (3x16) making it hard to judge how the method scales in practice.

The topic and the proposed method are nuanced and technical, so it does not help that there are many imprecise statements, including in theorems, and poor phrasing, which render the manuscript hard to follow (see questions and suggestions below). The authors also overlook the chance to help the reader with more visual aids.

Overall, I think the work has great potential which is not fully realized and publication-ready due to the unrefined presentation and somewhat uninspiring application.
questions: questions: - In Theorem 4.5, is it true that $\tau$ is a trilinear function? While $\tau$ is defined elsewhere in 4.3, this assumption is missing in the theorem itself which is confusing. You should also clarify in the manuscript how trilinearity is used in the proof, e.g. in (41) - (43).
- WIth that, could you clarify Theorem 4.5? My best understanding is that it states that the only trilinear function which satisfies the eikonal equation is a plane, but a counterexample to this is easy to conceive.
- In Line 221, do you perhaps mean $ |\tilde{\nu}(\mathbf{x})| \leq \epsilon$? Otherwise, you select the edges of the sub-level set, instead of the level set.
- Why do you use the $256^3$ marching cubes mesh as a ground-truth instead of using the input mesh on whose SDF you train?
- Can you comment on the extraction of the volumetric mesh? Is my understanding correct that this is what you have before the skeletonization step? Does this and your method more broadly offer a way to explain the success of trilinear positional encodings in terms of the underlying complex?
limitations: limitations: The limitations paragraph focuses solely on the runtime of the method. The section would benefit from reiterating and discussing the reliance on the eikonal constraint, which is known to be hard to satisfy globally. I would also urge the authors to consider discussing the limitations of the experiments as described in the weaknesses.


**Errors and suggestions**

- The motivation of discussing the normals in 4.3 is unclear. 
- All the left-hand terms in 40-43 should be squared. The rhs in (40) is the square of the gradient norm, not the gradient norm.
- The work would benefit from a discuss on the size of the NN and the encoding for practical applications. To my surprise, [12] uses similarly sized networks to your 3x16 with decent performance.

- Line 19: "polyhedral" -> 'polyhedron"
- Line 50: "This provides a theoretical exposition of the eikonal constraint" is unclear
- Line 59: you should give some examples of where marching cubes are used in a learning pipeline and explain what you mean by redundancy (probably flat regions which are described by many triangles?)
- Equation 2: I would not recommend overloading the notation for $d$
- Definition 3.2: what do you mean by "without loss of generality" when defining the NN?
- Line 109: "subscription" -> "subscript"
- Line 110: "for rather discussions as an SDF" is entirely unclear
- Lines 112-114 are very unclear
- Line 127: "subdivising" -> "subdividing"
- Definition 4.1: what is $\psi$? You should also define/explain $F$ before using it. Why is the definition given for general $D$ if this is a trilinear interpolation?
- Line 203: "also the" -> "also due to the" ?
- Line 248: ", similarly" is not needed
- Line 259: Can you be more precise with the objective? The goal is to get a boundary mesh of the zero-set of the SDF. 
- Line 262: reformulate to use a verb
- Line 276: CE is not defined in the main text.
- Line 294: "plenary" -> "planarity"?
- Line 560: do not write "by the way"
- Overall, please check the text with a spell and/or style checker
paper_summary: paper_summary: The paper proposes a theoretical link between the tropical algebraic interpretation of ReLU-based neural networks and surfaces represented by signed distance functions trained on those networks. Based on that theoretical interpretation, a surface extraction algorithm is proposed.
strength_weakness: strength_weakness: (1) I do not like the current presentation of the paper. I understand that it is theoretical, but its focus on lower abstraction details from the very beginning makes it a niche paper in NeurIPS. For example, the Experiments section focus on the details of the empirical settings, moving the actual results to an Appendix. Another example is that there is little effort in providing an overview of the method. Overall this is the core design of the paper: focus on the definitions and details before properly motivating them. In its current form, this paper is more aligned with an applied mathematics journal. This is the core reason of my current rating and fixing that problem is not possible without rewriting the paper. However, I think the theoretical framework is sound so I am not recommending rejection. All other following weaknesses derive from this presentation problem.

(2) This paper needs graphical intuitions to assist unexperienced readers. For example, showing the plot of a simple tropical polynomial, with the hypersurface points highlighted, would give a good intuition of the rationale behind Definition 3.1 (i.e. the linking points between each piece-wise linear region must be where two monomials are equal). If space for the images is a problem, I would advise to include them in the supplementary material, with references in the main text. This more didactical approach would increase the public that could appreciate the paper.

(3) Another example is the subdivision scheme. Showing an schematic image of the subdivision would help understanding.

(4) I also think a better overview image is needed. Figure 1 does not make a good overview because its abstraction level is too high. A better image should refer to core symbols in the definitions. Not having a proper overview makes understanding the motivations behind some definitions more difficult. For example, $\psi$ and $H$ are defined in Definition 4.1, and I could not find the motive for them in the main text. They are only used at Proposition D.2, at the Appendix. This kind of graphical intuition is even more important in a geometrical problem, which is the case of this paper.

(5) The extracted surfaces have noise and discretization problems, probably because of the grid-based approach and interpolation. Grid-based methods are known to introduce discretization problems (Neural Geometric Level of Detail - NGLOD, Instant Neural Geometric Primitives - Instant-NGP).

(6) The evaluation is lacking. The paper only compares against classic marching cubes. Being a classic algorithm, several improvements were introduced by more recent papers over the years.
questions: questions: (1) Typo in line 426: $x^3 + 2xy + y^4$ should be $x^3 \oplus 2xy \oplus y^4$. This is important because an unexperienced reader may get stuck in there.

(2) I believe the notation $v_j^{(i)}(x)$ in Proposition 3.3 is a little bit ambiguous. It is not immediately obvious that the $j$-th component is from the result of the layer computation $v^{(i)}(x)$ or if it is the $j$-th weight in layer $(i)$. In the second case $v_j^{(i)}(x)$ would have a weird interpretation, but I can see readers being stuck in there. Changing the notation probably would have consequences, but at least saying that the component $j$ is from the resulting vector $v^{(i)}(x)$ in the commentary that follows the definition would avoid that ambiguity.

(3) Also in Proposition 3.3 I would add an intuition commentary afterwards saying that $v_j^{(i)}(x)$ is a tropical polynomial, thus $v_j^{(i)}(x) = 0$ is an hyperplane composed of the points where two monomials are equal, intuitively connecting Proposition 3.3 with Definition 3.1.
limitations: limitations: Yes, section 7 describes limitations.

Comments: title: Thanks and commitment to revisions
Follow-Up on Rebuttal Feedback
comment: As the author-reviewer discussion ends, we sincerely thank our reviewers for their thorough assessments and detailed, constructive feedback, which have been invaluable in improving our work. We are committed to revising our manuscript and appendix according to our discussions and will prioritize this.

Thank you once again for your dedication, time, and effort. 

We truly appreciate the value of your contributions.


On behalf of all authors
I thank the authors for the clarification. I appreciate the additional table reporting the runtimes for different model sizes. I would encourage you to include this in the appendix.

Given that some of my concerns have been addressed and given the overall commitment of the authors to improve the work, I will increase my rating to borderline accept.

However, my main concern remains with potentially publishing misleading information. I urge the authors to take time to consider the subdivision process and redo Figure C for the final version. Currently, Figure C *is* wrong - in the generic case, the final decision boundary does not overlap with intermediate neuron decision boundaries. Each neuron in a ReLU network subdivides some existing pieces. The composition with the positional encoding changes the shape of these pieces but not the fact that these are still subdivided in a generic manner. 
Please consult existing works describing this, e.g. 
- Figure 4 in Humayun et al. SplineCam: Exact Visualization and Characterization of Deep Network Geometry and Decision Boundaries. CVPR 2023
- Figure 3 in Jordan et al. Provable Certificates for Adversarial Examples: Fitting a Ball in the Union of Polytopes. NeurIPS 2019.

You can also verify this yourself using the approach I described previously. 


Similarly, you cannot simply state that _Small_ performs worse so you are just not going to report the results.
Dear authors, thank you for you rebuttal. I appreciate the implementation improvement which greatly improves the time taken for meshing by the approach. However, I still think an evaluation an evaluation on a bigger set of meshes can strengthen the statistical significance of the quantitative evaluations, since the current evaluation is performed on very limited set of meshes. Further, visually I still find the benefits of the approach not too substantial, e.g. in Fig. F, it is mentioned that the improved normals would help lighting, but the reader has no idea how significant this improvement is going to be. We see that QEM helps both proposed and MC, I still wonder (a) how "MC + QEM" vs "Ours + QEM" looks qualitatively for smaller and larger meshes, and
(b) how "MC + QEM" compares against the proposed approach both qualitatively and quantitatively, to see how helpful the proposed approach is.  

I will keep my initial rating of borderline accept
Dear Reviewer 1zQb, as we approach the conclusion of the reviewer-author discussion period, do you have any remaining questions, or have your concerns been fully addressed? Your confirmation and reconsideration of your rating would be greatly appreciated.
**Figure C** 

We would like to clarify a potential misunderstanding. While Figure C is a schematic figure, we believe it accurately illustrates the characteristics of our method. The final decision boundary (0 level-set, shown in red) *may exclude* some intermediate neuron decision boundaries, as demonstrated. In Definition 3.1 (p.3), we introduced the concept of the **tropical** hypersurface, corresponding to the nonlinear boundaries formed by all neurons (including an SDF). We then referred to Proposition 3.3 (Decision boundary of neural networks), which stated that the decision boundary is a tropical hypersurface; however, the reverse is not true. Not all points on the tropical hypersurface constitute the decision boundary at the zero level-set. We formally denoted this relationship as $\mathcal{B} \subseteq ...$, rather than as $\mathcal{B} = ...$ in Proposition 3.3. And, in Lines 112-113, we stated that "the preactivations, $\nu_j^{(i)}$ are the *candidates* shaping the decision boundary." For more rigorous explanations on this matter, please refer to Proposition 6.1. of Zhang et al. (ICML 2018) "Tropical Geometry of Deep Neural Networks."

Notice that there is one more detail. The bottom line of the bunny is from a grid line, not from neuron's folded hypersurface. Geometrically, the grids of the HashGrid is non-linear boundaries which efficiently divided the Euclidean space into a lattice structure. To achieve this via only neurons, we need more neurons, which may explain why HashGrid (Muller et al., 2022) does not need many neurons or layers in decoding networks. 

Your constructive feedback has highlighted areas that may lead to misunderstandings, and we are fully committed to improving the clarity of our presentation.

Thank you for suggesting the figure from Hanin & Rolnick (2019). We will consider visualizing our concepts using this reference. Sadly, we cannot update the rebuttal PDF in the reviewer-author discussion period.


**MC + QED**

Unfortunately, the *Small* models have coarser grid intervals in the HashGrid, which puts them at a disadvantage compared to the *Large* models in terms of our analytical approximation method. (Denser grid intervals result in smaller errors when approximating curved hypersurfaces.) Additionally, using the *Small* models is less representative of real-world scenarios than the *Large* models.


**Theorem 5**

Yes, that appears to be an unintended typo. We will correct it.

**Model Size**

Referring to `3x16` could be misleading. Our model has three layers, structured as 3 (coordinates) -> 8 (hash grid) -> 16 -> 16 -> 1.

One limitation of this ablation study is that model size impacts fitting, affecting the GT mesh (MC256 for the pseudo-GT). However, the major model capacity comes from the hash tables (explained above in the initial author feedback). Therefore, the changes in CD are small. For completeness, we provide the following tables:

For the Bunny Small models,

| Layers | Width | Vertices    | CD (x 1e-6) | Time (sec) | Note       |
|--------|-------|-------------|-------------|------------|------------|
| 3      | 4     | 4547        | 666         | 0.12       |            |
| 3      | 8     | 4661        | 750         | 0.19       |            |
| 3      | 16    | 4628        | 724         | 0.65       | Ours       |
| 3      | 32    | 4567        | 754         | 0.57       |            |
| 3      | 64    | 5588        | 759         | 2.86       |            |
| 2      | 4     | 4565        | 738         | 0.11       |            |
| 2      | 8     | 4863        | 768         | 0.13       |            |
| 2      | 16    | 4566        | 734         | 0.12       |            |
| 2      | 32    | 4577        | 754         | 0.22       |            |
| 2      | 64    | 4611        | 739         | 0.50       | InstantNGP |

Notice that the model with layer=3 and width=64 took 2.86, because of the large number of intermediate (not final) vertices. The time complexity is more sensitive to the number of intermediate vertices than to the number of layers or the width.

W.r.t. the number of layers:

| Layers | Width | Vertices    | CD (x 1e-6) | Time (sec) | Note       |
|--------|-------|-------------|-------------|------------|------------|
| 2      | 16    | 4566        | 734         | 0.12       |            |
| 3      | 16    | 4628        | 724         | 0.65       | Ours       |
| 4      | 16    | 4569        | 759         | 0.51       |            |
| 5      | 16    | 4579        | 749         | 0.59       |            |

The number of layers has a limited impact on runtime, especially for 3 to 5.

**Comparison to CPWA**

I appreciate your insight into how it could enhance the overall work. However, we have found that the experiments require a major code revision due to integration issues with HashGrid. Although, we will consider this in our revision.
I thank the authors for their detailed response. While some of my minor concerns are addressed, the presentation and clarity remain my main concern.


**Figure C**
I think this schematic is very misleading. The final decision boundary (0 level-set, red) is not composed of some intermediate neuron decision boundaries in the way you show. I would strongly encourage you to visualize some intermediate steps from your method accurately depicting a realistic subdivision.

If this is not possible in the available time due to how your code is, a simple way to visualize the partition is to pass the coordinates of the pixels through the NN and plot the pre-activations of the ReLU neurons (e.g. Figures 2 or 7 in Hanin, B. et al, Deep ReLU networks have surprisingly few activation patterns, NeurIPS 2019.)


**MC + QED**
I appreciate the additional MC+QED experiment, but I slightly question the setup - this is reported for the _Large_ model, which is mainly discussed in the Appendix (Table 6). There, your method has better CD than any MC, so of course QEM will perform better on your mesh. I would appreciate if you could report MC+QEM for the _Small_ (Table 1) which is used in the main paper and where finer MC meshes (128 or 196) have better CD but worse vertex counts.


**Theorem 5**
Could you please clarify Line 559 in the proof? Do you perhaps mean for _all_ $\mathbf{x}$?


**Model size**
My criticism is not that the size itself it small, but rather that the method is assessed just for a _single size_. E.g. Fig. 10 of instantNGP [12] reports multiple widths and depths of the MLP, so it would be reassuring to see how the quality/runtime behaves for models of size other than 3x16.


**Comparison to CPWA**
I understand the motivation of the spectral bias and positional encoding. However, because you are treating this topic from the perspective of surface meshes, contrasting the meshes obtained with and without an encoding is in my view a worthwhile consideration. Describing and explaining the qualitative and quantitative differences between these meshes could add to a unique interpretation of positional encodings and strengthen the unique aspect of your work. The whole point of increasing the spectral bias is to represent higher frequency components - are these present in the mesh relative to the plain ReLU NN? Training a Standford bunny with ReLU certainly isn't hard (the method you reference in [4] and the code therein already provides this). I do not expect the authors to perform the experiments and the analysis in the short available time and this is not critical for my current assessment, but as a discussion point, I do maintain that this could add to your unique work overall.


**Q4** 
Thank you, this clarifies my concern on the ground-truth.
Dear Reviewer FmNY, 

As we approach the end of the reviewer-author discussion period, we want to ensure that our feedback has sufficiently addressed your major concerns. We recognize the importance of this opportunity to address the reviewers' raised concerns, and we are committed to providing any necessary clarifications. Please let us know if you have further questions or need additional clarification.

Best regards,

On behalf of all authors
Thank you for your thoughtful and detailed feedback. We truly appreciate the time you invested in reviewing our manuscript.
Thank you.
I don't have any more questions and will adjust my rating. I think the paper is interesting and should be accepted, given the promised changes in the final version.
We understand that, for the extracted mesh using MC64 or MC256, we can measure the average of predicted SDFs and how much they deviate from the learned surface of networks.

In the same procedure, the corresponding MC numbers are as follows:

| MC Sample | Average of Squared SDFs |
|-|-|
| 64 | 2308 x 1e-9 |
| 256 | 66 x 1e-9 |

We confirmed that the numbers are higher than the mesh using our method (31 x 1e-9).
Thank you a lot for such a fast response.

Do you happen to also have these numbers for MC64 and MC256?
Regarding **Fig. A** (closer looks at bunny's nose), we understand you are referring to the small facets along the left-side of the nose line in the highlighted region. We assume these narrow facets are embedded in the facets of MC256, while MC64 fails to capture the left and right edge of the nose and represents this with more facets cutting these two sides of edges (*not aligned* with nose lines and rather assign more facets in here).

Regarding the quantitative evaluation, we are happy to report the predicted SDF on the surface as you suggested. We randomly sampled 100K points on the surface of pretrained networks (a bunny Large model) and calculated the *average* of squares for predicted SDF (model's outputs) on the surface. As we expected, the predicted SDFs on the reconstructed surface converge to zero as the weight of the eikonal loss increases, and the number is close to zero $\approx 3e-8$. Note that, for a large model, fine-grained trilinear regions from the dense grid marks (from the more dense multi-resolution grids) allow relatively small errors of $97 \times 1e-9$ without the eikonal loss, although the eikonal loss of $1e-2$ further decreases to $31 \times 1e-9$.

| Eikonal Loss Weight   |   Average of Squared SDFs |
|-----------------------|----------------------- |
| 0 | 97 x 1e-9 |
| 1e-4 | 91 x 1e-9 |
| 1e-3 | 68 x 1e-9 |
| 1e-2 | 31 x 1e-9 |

Plotting the predicted SDF values on the reconstructed surfaces would be beneficial. We could visualize this in a manner similar to **Fig. F**, which illustrated the normal errors on the surfaces. We anticipate that regions with high curvature will be particularly noteworthy. Notice that, sadly, we cannot update the rebuttal PDF during the author-reviewer discussion period. 

Apologies for overlooking your insights, and we hope this response addresses your concern.
Dear authors, thank you a lot for addressing my questions, it clarified some of my concerns.

Q1: My misunderstanding was coming from treating the "trilinear region" in line 50 as the trilinear region of the PE, not a trilenar region of PE + neural network. Your response clarified it. Thank you a lot.

W2 + Q2: I am still a bit confused about the metrics and visual comparisons in the paper. Even in the provided zoom-in to the nose of the Stanford bunny, it is still hard to see if the proposed method is closer to MC256. It is also hard to see the oversampling of the vertices mentioned in the rebuttal: the proposed method has only 20% fewer vertices, and the highlighted region contains the small facets along the nose line that are not present in MC64.

Regarding the quantitative evaluation, I agree that in the limit of MC resolution currently used metric should converge to the one proposed in Q2. However, MC256^3 has only 4 subdivisions of each edge from MC64^3, which is still quite low resolution and introduces additional noise into the metric. Wouldn't it be more direct and transparent to directly report a sum of squares for predicted SDF on the surface? In my understanding, since the proposed method (almost) exactly finds the intersections of the hyperplanes, the predicted SDF on the entire reconstructed surface should be zero (or with a very small error from the approximations made when searching for the intersection). This is not the case for Marching Cubes which has a fixed grid over the domain.

Plotting the predicted SDF value over the reconstructed surfaces should also give much better insights into where MC makes mistakes vs the proposed method and how severe they are.
Dear Reviewer MSAr,

Thank you for your thorough review and the constructive feedback you provided.

We would like to kindly inquire whether our responses to Weakness #2 and Question #1 have addressed your concerns sufficiently, as you mentioned you might consider adjusting the final rating based on our rebuttal. Given the limited time for the author-reviewer discussion period, we would also appreciate it if you could let us know if there are any additional questions or clarifications needed.

Thank you again for your time and effort.

Best regards,

On behalf of all authors
Thank you for your prompt response. 

We humbly ask that you please inform us if you have any remaining major issues that could affect the rating so we can address them during the author-reviewer discussion period. In particular, we have included Figures B and C in the rebuttal PDF to address your main concerns.
*Below are our comments on your questions, demonstrating our commitment to clarifying our paper.*
 
---

**The motivation of discussing the normals in 4.3 is unclear.**

- To define a face, we use a list of vertex indices. To determine the normal vector of the face (considering the ambiguity of two opposite directions), we implicitly rely on the order of the vertex indices as described in the paper.

**Line 50: "This provides a theoretical exposition of the eikonal constraint" is unclear**

- This indicates Thm 4.5 and Coro. 4.6 which stated "that within the trilinear region, the hypersurface transforms into a plane" under the eikonal constraint. 

**Line 59: you should give some examples of where marching cubes are used in a learning pipeline and explain what you mean by redundancy (probably flat regions which are described by many triangles?)**

- MeshSDF (Remelli et al., 2020) will be the example. Yes, your answer is right.

**Line 110: "for rather discussions as an SDF" is entirely unclear**

- An SDF has a single scalar output as a definition. In the manuscript, we later discuss the networks as an SDF in Sec. 4, specifically Def. 4.3. (Line 175). We will make sure to clarify this as you asked.

**Lines 112-114 are very unclear**

- We hope **Fig. C (c) in the rebuttal PDF** helps your understanding. Our method performs a single pass over neurons to find tropical hypersurface (non-linear boundaries; refer to **Fig. B in the rebuttal PDF**).

**Definition 4.1: what is $\psi$? You should also define/explain $F$ before using it.**

- In Def. 4.1, $\psi$ is a trilinear interpolation function given weight $\mathrm{w}$ and features $\mathrm{H}$. This is implemented using the HahsGrid as $\tau$ in Def. 4.3. In Lines 172-173, we explained that "they transform the input coordinate $\mathrm{x}$ to $\mathrm{w}$ with a multi-resolution scale and find a relative position within a grid cell." The features $\mathrm{H}$ are the hash table entries derived from the hashing function, where the input to the hash function is a given grid vertex and the output is the index in the hash table. We will provide a clear explanation in the revision.

- $F$ is the output dimension of $\tau$, positional encoding module. (Line 170)

**Line 259: Can you be more precise with the objective? The goal is to get a boundary mesh of the zero-set of the SDF.**

- By the definition of SDF, the zero-set of the SDF naturally indicates surfaces. But we respect your concern and will revise it.
Thank you for considering the review and the efforts to answer the comments. I have no further questions.