A Matrix-Free Exact Newton Method

SIAM Journal on Scientific Computing, Volume 46, Issue 3, Page A1423-A1440, June 2024.
Abstract. A modification of Newton’s method for solving systems of [math] nonlinear equations is presented. The new matrix-free method is exact as opposed to a range of inexact Newton methods in the sense that both the Jacobians and the solutions to the linear Newton systems are computed without truncation. It relies on a given decomposition of a structurally dense invertible Jacobian of the residual into a product of [math] structurally sparse invertible elemental Jacobians according to the chain rule of differentiation. Inspired by the adjoint mode of algorithmic differentiation, explicit accumulation of the Jacobian of the residual is avoided. Prospective, generally applicable implementations of the new method can be based on similar ideas. Sparsity is exploited for the direct solution of the linear Newton systems. Optimal exploitation of sparsity yields various well-known computationally intractable combinatorial optimization problems in sparse linear algebra such as Bandwidth or Directed Elimination Ordering. The method is motivated in the context of a decomposition into elemental Jacobians with bandwidth [math] for [math]. In the likely scenario of [math], the computational cost of the standard Newton algorithm is dominated by the cost of accumulating the Jacobian of the residual. It can be estimated as [math], thus exceeding the cost of [math] for the direct solution of the linear Newton system. The new method reduces this cost to [math], yielding a potential improvement by a factor of [math]. Supporting run time measurements are presented for the tridiagonal case showing a reduction of the computational cost by [math]. Generalization yields the combinatorial Matrix-Free Exact Newton Step problem. We prove NP-completeness, and we present algorithmic components for building methods for the approximate solution. Potential applications of the matrix-free exact Newton method in machine learning of surrogates for computationally expensive nonlinear residuals are touched on briefly as part of various conclusions to be drawn.