🔢 mathematics

Laplacian regularized eikonal equation with Soner boundary condition on polyhedral meshes

This paper proposes a cell-centered finite volume algorithm for solving a Laplacian regularized eikonal equation with Soner boundary conditions on polyhedral meshes, demonstrating second-order convergence and significant computational efficiency over time-dependent methods for large-scale or distant distance field calculations.

Original authors: Jooyoung Hahn, Karol Mikula, Peter Frolkovič

Published 2026-08-14

📖 7 min read🧠 Deep dive

CC BY 4.0

Original authors: Jooyoung Hahn, Karol Mikula, Peter Frolkovič

Original paper licensed under CC BY 4.0 (http://creativecommons.org/licenses/by/4.0/). ✨ This is an AI-generated explanation of the paper below. It is not written or endorsed by the authors. For technical accuracy, refer to the original paper. Read full disclaimer

Imagine you are standing in a vast, dark cave filled with jagged rocks, stalactites, and hidden chambers. You want to know exactly how far you are from the nearest wall or rock at every single point in the cave. This isn't just a game of "how far is it to the exit?" but a complex 3D map where every speck of dust needs a distance tag. In the world of science and engineering, this "distance map" is called a distance function. It's the secret sauce behind everything from designing safer cars that understand their surroundings, to simulating how fire spreads through a forest, to predicting how electrical signals race through a beating heart.

To create these maps, scientists use a mathematical rule called the Eikonal equation. Think of this equation as a set of instructions for a wave of light or sound spreading out from a source. The rule says, "The wave moves at a constant speed, and the distance it has traveled is simply the time it took." However, in the real world, things get messy. The cave walls might be weirdly shaped, or the source might be a tiny speck inside a huge room. If you try to solve this math problem on a computer using standard methods, the solution can get "stuck" or behave strangely near the walls, especially if the cave has sharp corners or weird shapes. This is where a special rule, known as the Soner boundary condition, comes in. It's like a traffic cop at the cave entrance, ensuring the wave doesn't try to sneak out of the cave in a way that breaks the laws of physics.

For a long time, the best way to solve this was to pretend the wave was moving forward in time, step by step, until it filled the whole cave. But if the cave is huge and the source is tiny, this "time-stepping" method is incredibly slow. It's like trying to fill a swimming pool by pouring in a single cup of water every second; you'd be waiting forever to see the far end get wet. This paper introduces a clever new trick to speed things up, turning a slow, step-by-step race into an instant, all-at-once calculation, even on the most complicated, blocky computer models of the world.

The Paper's Big Idea: A Smoother, Faster Way to Map the World

The authors of this paper, Jooyoung Hahn, Karol Mikula, and Peter Frolkovič, have developed a new numerical algorithm to solve the Eikonal equation on polyhedral meshes. If you imagine a 3D computer model as a giant Lego structure, a "polyhedral mesh" is just a fancy way of saying the structure is built from blocks that can have any number of sides, not just cubes. This is crucial because real-world objects (like car engines or human hearts) are rarely perfect cubes; they are complex shapes that need these irregular blocks to be modeled accurately.

The team's main innovation is solving a modified version of the Eikonal equation called the Laplacian regularized eikonal equation. Here is the magic trick: instead of letting the "distance wave" travel slowly over time, they add a "smoothing" ingredient (the Laplacian term) that acts like an infinite-speed messenger. This allows the distance information to instantly reach every corner of the domain, rather than waiting for a wave to physically travel there.

However, there's a catch. If you make the smoothing too strong, the map becomes blurry and inaccurate. If you make it too weak, the math becomes unstable and crashes. The authors figured out a "Goldilocks" strategy. They start with a strong smoothing effect to get a rough, stable map, and then gradually reduce the smoothing in a specific sequence. With each step, they use the previous result as a starting point for the next, finer calculation. This is like sculpting a statue: you first chip away the big chunks of stone with a heavy hammer (strong smoothing), and then you switch to a fine chisel (weak smoothing) to get the perfect details.

What They Found and Why It Matters

The researchers tested their method on a variety of scenarios, from simple spheres to complex, hollowed-out shapes with sharp corners. They ran these tests on four different levels of mesh detail, ranging from about 8,000 blocks to over 28 million blocks.

The Speed Boost:
The most exciting finding is the dramatic reduction in computational cost. When the region of interest is far away from the starting object, their new method is vastly faster than the traditional "time-stepping" approach. In one test case with a very fine mesh (over 8 million blocks), their algorithm was nearly 50 times faster than the older method to reach the same level of accuracy. In another case with 28 million blocks, the speedup was even more dramatic, reaching a ratio of nearly 69 times faster. This means that problems that used to take hours or days to solve could potentially be done in minutes.

The Accuracy:
The paper also checked how close their "smoothed" maps were to the true mathematical answer. For smooth shapes (like a perfect sphere), they found that their method achieves a second-order experimental order of convergence in the $L^1$ norm error. In plain English, this means that as they made the computer blocks smaller (increasing the mesh resolution), the error in their distance map dropped very quickly, proving the method is highly accurate for smooth problems. For shapes with sharp corners or singularities, the accuracy was slightly lower (closer to first-order), which is expected and consistent with previous research.

The "Soner" Safety Net:
A key part of their success was correctly applying the Soner boundary condition. Without this, the algorithm would try to calculate distances in directions that don't make physical sense, leading to errors. The authors showed that their method respects this condition perfectly, ensuring the distance map behaves correctly even at the boundaries of the domain.

The "How" Behind the Magic

The method relies on a technique called the cell-centered finite volume method. Imagine the 3D space divided into tiny cells (the polyhedral blocks). The algorithm calculates the average value of the distance function inside each cell and ensures that the "flow" of information across the walls of these cells is balanced.

To handle the tricky math of the non-linear equation, they used a linearization technique. They took a known, slightly imperfect solution and used it to guess the direction of the wave, turning a hard, non-linear problem into a series of easier, linear problems. They solved these linear problems iteratively, refining the guess each time.

Crucially, this method is designed for parallel computing. Because the algorithm only needs information from the immediate neighbors of a cell (a "1-ring" neighborhood), it can be easily split up across many computer processors. This makes it perfect for modern supercomputers that use domain decomposition to tackle massive problems.

The Bottom Line

This paper doesn't claim to have solved every possible distance-mapping problem in the universe. It explicitly notes that for very small regularization parameters (when the smoothing is almost gone), the math can become unstable, and finding the perfect parameter value is still an area for future study. However, for the specific goal of computing distance functions on complex, polyhedral meshes, the authors have demonstrated a robust, highly efficient, and accurate method.

By combining a vanishing viscosity approach (gradually removing the smoothing) with a smart boundary condition, they have created a tool that is significantly faster than current state-of-the-art methods for large-scale simulations. Whether it's helping engineers design better combustion engines or helping doctors understand heart rhythms, this new algorithm offers a way to map the invisible distances of our world with unprecedented speed and precision.

Technical Summary: Laplacian Regularized Eikonal Equation with Soner Boundary Condition on Polyhedral Meshes

Problem Statement
The paper addresses the numerical computation of a distance function from given objects (represented by curves, surfaces, or boundary parts $\Gamma$ ) within a three-dimensional computational domain $\Omega$ . The domain is discretized using polyhedral cells, a common requirement for industrial problems with complex boundary shapes. The goal is to find the viscosity solution of the Eikonal equation:
$|\nabla u(x)| = 1, \quad x \in \Omega \setminus \Gamma$
subject to $u(x) = 0$ on $\Gamma$ . A critical challenge identified is the necessity of the Soner boundary condition (or state constraint condition) on the domain boundary $\partial \Omega \setminus \Gamma$ :
$\nu(x) \cdot \nabla u(x) \geq 0$
where $\nu$ is the outward normal. Without this condition, particularly on non-convex domains, the numerical solution may fail to converge to the correct viscosity solution.

Existing methods face limitations in this context:

Fast Marching Method (FMM): While efficient, its heap data structure is difficult to maintain efficiently in parallel computations using domain decomposition.
Fast Sweeping Method (FSM): Effective for simple domains but relies on specific ordering that is not straightforward for general polyhedral meshes.
Time-Relaxed Formulations: Solving the time-dependent Eikonal equation (propagating the front over time) is robust but computationally expensive, especially when the region of interest is far from $\Gamma$ or when using fine meshes with many cells.

Methodology
The authors propose a numerical algorithm based on a cell-centered Finite Volume Method (FVM) to solve a Laplacian regularized Eikonal equation. Inspired by the vanishing viscosity method, the governing equation is modified to include a Laplacian term:
$-\epsilon \Delta u_\epsilon(x) + |\nabla u_\epsilon(x)| = 1, \quad x \in \Omega \setminus \Gamma$
with $u_\epsilon = 0$ on $\Gamma$ and the Soner boundary condition applied on $\partial \Omega \setminus \Gamma$ . Here, $\epsilon > 0$ is a regularization parameter.

Key methodological components include:

Linearization: The nonlinear term $|\nabla u_\epsilon|$ is linearized using a known approximation $u_{\epsilon'}$ (from a previous iteration with a larger $\epsilon$ ) to form a linear convection-diffusion equation. The convection velocity is defined as $v = \nabla u_{\epsilon'} / |\nabla u_{\epsilon'}|_\sigma$ .
Discretization on Polyhedral Meshes: The authors adapt the flux-balanced approximation for polyhedral cells. Faces are tessellated into triangles to handle non-planar faces. The discretization handles internal faces and boundary faces separately.
Soner Boundary Condition Implementation: The numerical scheme explicitly enforces the Soner condition by splitting boundary triangles into sets based on the sign of the normal flux. Terms corresponding to outflow on the boundary (where the condition might be violated) are set to zero, while inflow terms are handled via a modified inflow-based gradient.
Sequential Regularization Strategy: To ensure convergence to the viscosity solution, the algorithm solves a sequence of problems with decreasing regularization parameters $\epsilon_n$ . The parameter is chosen empirically as $\epsilon_n = (h_L)^{1/2n}$ , where $h_L$ is the characteristic mesh length. The solution from step $n-1$ serves as the initial guess for step $n$ .
Parallelization: The method utilizes a standard cell-centered FVM structure with a 1-ring face neighborhood for domain decomposition. This allows for straightforward parallel implementation without complex global data structures like heaps.

Key Contributions

Algorithm for Polyhedral Meshes: A robust cell-centered FVM scheme capable of solving the Eikonal equation with the Soner boundary condition on general polyhedral meshes, overcoming the geometric restrictions of methods requiring acute triangles or tetrahedra.
Handling Nonlinearity and Regularization: A specific strategy for linearizing the Eikonal equation and selecting a decreasing sequence of regularization parameters to balance numerical stability (avoiding oscillations from small $\epsilon$ ) and accuracy (avoiding excessive smoothing from large $\epsilon$ ).
Computational Efficiency: Demonstrating that the Laplacian regularized approach offers a significant reduction in computational cost compared to time-relaxed formulations, particularly for large meshes or distant regions of interest, due to the infinite propagation speed inherent in the elliptic regularization.

Numerical Results
The authors validated the method using ten test cases (EX1–EX10) involving spheres, planes, disks, and complex domains with varying $\Gamma$ configurations.

Convergence: For smooth solutions (e.g., distance from a sphere), the method achieves a second-order experimental order of convergence (EOC) in the $L^1$ norm. For non-smooth solutions (e.g., near corners or singularities), the EOC approaches 1, consistent with theoretical expectations.
Parameter Sensitivity: Experiments confirmed that the choice of $\epsilon$ is critical. A value that is too small leads to instability, while a value that is too large reduces accuracy. The sequence $\epsilon_n = (h_L)^{1/2n}$ with $n=5$ was found to yield optimal results.
Performance Comparison: When compared to a time-relaxed bidirectional Eikonal equation solver (implemented in the same codebase), the proposed method demonstrated dramatic speedups. The efficiency gain increased with the number of cells and the distance of the region of interest from $\Gamma$ . In some cases, the proposed method was nearly 70 times faster for the same error tolerance.

Significance and Claims
The paper claims that the proposed algorithm provides a practical and efficient solution for computing distance functions on complex, industrial-grade polyhedral meshes. Its primary significance lies in:

Parallel Scalability: The ability to implement the solver using standard domain decomposition with 1-ring neighborhoods, making it suitable for high-performance computing environments where memory consumption is a bottleneck.
Cost Reduction: The ability to dramatically reduce computational costs compared to time-marching methods, specifically when calculating distances in large domains or far from the source objects.
Robustness: The successful application of the Soner boundary condition within a finite volume framework to ensure the correct viscosity solution is obtained on non-convex domains.

The authors conclude that the method is ready for integration into standard cell-centered FVM codes and is particularly advantageous for applications requiring fine-scale 3D simulations, such as combustion modeling, turbulence modeling (wall distance), and cardiac electrophysiology.

The Paper's Big Idea: A Smoother, Faster Way to Map the World

What They Found and Why It Matters

The "How" Behind the Magic

The Bottom Line

Technical Summary: Laplacian Regularized Eikonal Equation with Soner Boundary Condition on Polyhedral Meshes

More like this