Adaptive Multigrid
This page documents developments on the design and implementation of an adaptive multigrid solver that allows for the efficient solution of Poissonlike systems of equations, discretized over an quad/octtree.
This material is based upon work supported by the National Science Foundation under Grant No. 1422325:
 Title: III: CGV: Small: Designing an Adaptive Method for Solving Large Linear Systems of Equations in Two and ThreeDimensional Space
 PI: Michael Kazhdan
 Duration: 09/2014  08/2017
 Students: Fabian Prada
Goals and Challenges
In numerous fields, ranging from image processing to scientific simulation, solving a large Poisson system is an essential step in integrating locally defined properties into a consistent global solution. Though the linear system is global (with properties at one point in space affecting the solution at points far away), in many cases the details of the solution are only required at particular locations. As a result, computational effort expended on estimating the details of the solution away from these regions of interest is wasted. The proposed research seeks to design, implement, and make publicly available a new solver for addressing this problem  providing a way to adapt traditional hierarchical solvers so that computation is only focused where it is needed.
The motivation for this approach stems from the fact that in solving a Poisson equation, longdistance effects are dominated by lower frequencies, so highresolution components of the solution need only be computed in regions of interest. Using experiences from prior work in surface reconstruction and large image processing, the PI will explore a modification of traditional multigrid, providing a hierarchical solver even when the function spaces do not nest. To do this the multigrid solver will be modified in two ways. First, to compensate for the fact that the coarser solution cannot be upsampled to the finer resolutions, the solution of the linear system will be represented as a linear combination of functions at all levels of the hierarchy, not just the finest one. Second, instead of using traditional restriction and prolongation operators to transition between successive levels of the hierarchy, the solver will proceed by iterating through the levels of the hierarchy (finetocoarse and then coarsetofine as in a typical V or Wcycle), relaxing the system within each level, and adjusting the constraints at the other levels, accounting for the components of the solution that have already been met. (For example, in the prolongation phase, the coarser solution will be incorporated into the finer levels by adjusting the righthandside, rather than updating the current estimate of the solution.) The research will explore both the general design of such a system, as well as domainspecific implementations that facilitate its integration in targeted scientific and entertainment disciplines. Results will be made publicly available in the form C/C++ source code for constructing adapted quadtrees/octrees, formulating the system constraints, and evaluating the computed solution, as well as documentation and applicationspecific implementations.
Direct Results
Initial research will focus on extending the cascadic solver used in Poisson Surface Reconstruction. This includes:
 [Version 7.0] Extending the support for function representation/evaluation using the secondorder Bspline basis. (This has been used to provide color intepolation support in the context of surface reconstruction.)


Geometry
 Geometry + Color

Notes:
 For performing color interpolation, we "splat" the color values associated with a sample using a kernel whose width is inversely proportional to the sampling density and whose weight is proportional. This provides a color field that adapts to the sampling density and is defined without requiring an explicit solve.
 [Version 8.0] Extending the function representation to BSplines of different degrees. (This supports a parameterizable tradeoff between the accuracy of higher order elements and efficiency/sparsity of linear systems derived from lower order elements.)





Degree: 1
 Degree: 2
 Degree: 3
 Degree: 4

Time: 3.7 s
 Time: 4.9 s
 Time: 9.6 s
 Time: 18.0 s

Memory: 75 MB
 Memory: 141 MB
 Memory: 374 MB
 Memory: 1016 MB

Triangles: 149 K
 Triangles: 148 K
 Triangles: 148 K
 Triangles: 148 K

Notes:
 The slightly choppier look of the degree1 reconstruction may be due to the fact that we cannot use a secondorder fit when the implicit function is represented using piecewise trilinear elements.
 Although the system is defined using degreed elements, the normal field is still represented using secondorder elements so that similar constraints are presented to all four systems, and results can be compared more fairly.
 [Version 9.0] Extending to higher order systems, supporting more general integral and interpolation constraints, and supporting additional boundary conditions. (This enables the implementation of the SSD reconstruction algorithm which uses a biLaplacian regularizer using our system, providing a more efficient a less memory intensive solution.)




Poisson Reconstruction
 SSD Reconstruction (Calakli and Taubin)
 SSD Reconstruction (Ours)

Time: 110 s
 Time: 354 s
 Time: 92 s

Memory: 2120 MB
 Memory: 6011 MB
 Memory: 2225 MB

Triangles: 4.3 M
 Triangles: 4.0 M
 Triangles: 4.3 M

Notes:
 Our SSD implementation turns out to be slightly faster than our implementation of Poisson Reconstruction. This is because the constraints (excepting interpolation) are trivial and do not need to be computed explicitly.
 Initial evaluation indicates that our SSD reconstruction better reproduces some of the detail in the input. This is likely due to the manner in which default parameters are set.
 [Version 10.0] Extension of the cascadic solver to a fully multigrid solver that supports constraints expressed both in terms of the finite element basis and pointwise values. (This provides a solver with multigrid performance and linear time complexity.)
Depth
 VCycle: 1
 VCycle: 2
 VCycle: 3
 VCycle: 4
 VCycle: 5
 VCycle: 6
 VCycle: 7

0
 8.3e17
 8.3e17
 8.3e17
 8.3e17
 8.3e17
 8.3e17
 8.3e17

1
 5.8e16
 4.2e17
 4.2e17
 4.2e17
 4.2e17
 4.2e17
 4.2e17

2
 7.1e7
 8.7e11
 1.1e14
 7.3e16
 7.4e16
 7.5e16
 7.4e16

3
 6.6e4
 1.2e5
 2.9e7
 7.5e9
 2.0e10
 5.5e12
 1.5e13

4
 3.7e3
 8.2e4
 2.4e4
 7.3e5
 2.2e5
 6.6e6
 2.0e6

5
 1.6e2
 4.2e3
 1.7e3
 9.9e4
 6.7e4
 4.8e4
 3.5e4

6
 8.5e3
 1.7e3
 3.5e4
 7.1e5
 1.9e5
 1.1e5
 7.8e6

7
 6.7e3
 1.1e3
 1.8e4
 3.1e5
 5.2e6
 1.0e6
 3.8e7

8
 4.3e4
 1.9e5
 3.9e6
 7.0e7
 1.2e7
 1.9e8
 2.7e9

Residuals at the different levels of the multigrid solver as a function of the number of vcycles. Results were obtained using the Stanford Bunny model as input.

Notes:
 Our full multigrid solver currently only supports the unscreened version of the Poisson Reconstruction algorithm. This is because that algorithm assigns depthdependent interpolation weights to achieve more efficient convergence. Though that works well in practice, it violates the underlying assumption of multigrid that there is a single (continuous) energy that is being minimized at all levels of the multiresolution hierarchy.
Indirect Results
In addition to applications in surface reconstruction, we have also begun exploring the use of hierarchical solvers within the broader domain of geometry processing. This has lead to important new results in solving the opticalflow problem over surfaces. This was an important step of the processing that allowed for the computation of motion graphs for unstructued texture meshes.
Contact Info
Any opinions, findings, and conclusions or recommendations expressed in this material are those of the author(s) and do not necessarily reflect the views of the National Science Foundation
Last updated 08/23/2016