Notes

Read these notes to help understand the method and how it is best implemented.

Technical

See the modules technical notes for:

  • Quick Access To All Demos : domain data requirements, sign flip when using the unified Poisson kernel, and the choice of rank for 2D and 3D.

  • Generator: unified vs standard kernel types, application notes on choosing \(\alpha\) and \(\beta\) when generating kernels, and warm starting.

  • Decomposition: 2D (SVD) vs 3D (Symmetric CP) decompositions, absorbing eigenvalues, the meaning and usage of rank for tensors in 3D, and what a safe rank means.

  • Mathlib: inverse vs forward Poisson equations, matrix-form vs vector-form linear solvers, Poisson filter computation and convolution.

  • Convergence: how adaptive truncation works.

  • IO Handler: how to load pre-generated Poisson filters database, csv and hlsli/c++ exports.

General

  • This work is based on MAC-Grids.

  • In the context of kernel and filter generation we use target Jacobi iteration and Poisson filters order interchangeably in the code. Order means the target Jacobi iteration of a filter/kernel.

  • See this section in Mathlib Technicals for the order of convolution when applying Poisson filters in 2D and 3D.

  • When solving an inverse Poisson equation (like Poisson pressure) if going with UNIFIED kernel instead of STANDARD you need the following sign flip on the right hand side \(b\) in \(Ax=b\) to make it work. This is the only downside of benefiting from a UNIFIED kernel.

    data_domain *= -1. if opts.kernel.kernel_type == com.PoissonKernelType.UNIFIED else 1.

  • The Jacobi convergence curve provides a lower bound on the Poisson filter possible solutions. This means if there is no numerical loss due to filter reduction or truncation, Poisson filters must exactly match the convergence behaviour of Jacobi within machine precision, which will be the best case for numerical quality, and worst case performance-wise. Any reduction or truncation results in degradation of the convergence behaviour in exchange for performance improvement.

  • Problem with SymCP Eigenvalues in tensor decomposition: contrary to 2D, eigenvalues are not necessarily sorted in terms of variance explained in 3D. When using grouped ranks, solutions might get worse before getting better when including higher ranks. See Safe Rank explanation in Decomposition Technicals for more details on how rank definition is different in 3D.

  • Remember there is no consensus on the definition of tensor eigenvalues. We use the definition that best suits our case, and remains consistent for both 2D and 3D filter computations. Such definition is only possible because the Poisson kernel is symmetric. See “All Real Eigenvalues of Symmetric Tensors”, 2014 and “A Spectral Theory for Tensors”, 2010 .

  • Eigenvalues from SymCP are real (no complex part).

Potential Applications and Future Work

  • Multi-grid smoother: multi-grid solvers have smoothing steps that are typically iterative. The performance of multi-grid solvers can be improved by using high-order Poisson filters in the smoothing step. The results will be interesting to observe as we might only need one or very a few multi-grid cycle iterations to converge.

  • Dirichlet boundary: Current method only supports Neumann boundaries. Addressing Dirichlet boundary treatment with Poisson filters will be highly desirable to expand the scope of our method to many application domains that need Dirichlet boundary treatment.

    For instance, one can consider using Poisson filters for Seamless Cloning. You can check out Fourier implementation of Poisson image editing, 2012 to see how the Dirichlet boundary can be transformed into a Neumann boundary problem when using spectral methods. This method is consistent with the spectral point of view that Poisson filters is based on.

  • Improving Mirror Marching algorithm: when dealing with Neumann boundaries, current method only marches parallel to spatial axes at no angle, which results in perfect solutions around walls but non-ideal artifacts around corners. A Scattering Marching algorithm that uses the object surface normals would likely remedy this artifact.

Suggested Reads/Codes