hessian-vector products

In machine learning, Hessian and other curvature matrices suffer from high computational complexity, lack of unified representation, and poor scalability, severely hindering applications such as second-order optimization, uncertainty quantification, and model compression. To address this, we propose a unified abstraction paradigm based on linear operators—specifically, matrix-vector multiplication interfaces—and systematically develop the first scalable, user-friendly curvature matrix representation framework. Building upon this paradigm, we introduce *curvlinops*, an open-source PyTorch library that supports automatic differentiation, incorporates structural priors (e.g., sparsity, Kronecker factorization), and maintains compatibility with major deep learning frameworks. Experiments demonstrate that our approach drastically simplifies implementation across curvature-driven tasks, enabling efficient, memory-bounded, and scalable computation even on models with tens of billions of parameters. This advances the practical adoption and democratization of curvature-based methods in large-scale deep learning.

Accurate estimation of the Hessian spectrum for large-scale foundation models has long been hindered by prohibitive computational costs, with existing approaches relying on small models or strong structural approximations. This work proposes a sharded local finite difference method compatible with Fully Sharded Data Parallel (FSDP), combined with stochastic Lanczos quadrature, enabling the first efficient and accurate computation of the true Hessian spectrum on hundred-billion-parameter language models. The method incurs only constant-level training overhead and reveals that widely used block-diagonal approximations severely break down even in moderately sized LLMs. Furthermore, it establishes a scalable paradigm for Hessian spectral analysis, supporting numerical stability studies under fp32/bf16 precision and end-to-end scaling law modeling.

Second-order optimization is hindered by the prohibitive computational cost of Hessian evaluation. This paper proposes an asynchronous, client-separated second-order optimization framework that, for the first time, decouples gradient and curvature (Hessian) computations across distinct clients and updates them asynchronously—eliminating manual hyperparameter tuning required by conventional lazy-Hessian approaches. Integrating cubic regularization with an inexact Hessian update mechanism, the framework effectively alleviates computational latency and communication bottlenecks in high-dimensional settings. We establish rigorous theoretical guarantees of strong convergence. Empirical evaluations on both synthetic and real-world datasets demonstrate consistent superiority over standard and lazy-Hessian baselines, achieving up to √τ real-time speedup while preserving the rapid convergence rate of second-order methods and enhancing practical training efficiency.

For large-scale optimization problems lacking second-order derivative information—such as large-scale data fitting and curve modeling—this paper proposes a parameterized compact Hessian approximation method that replaces the dense Hessian with a low-rank structure, significantly improving scalability. The method introduces the first unified parameterized framework that encompasses and generalizes classical compact quasi-Newton formulas, including BFGS and SR1, and incorporates a flexible vector selection mechanism to better adapt to heterogeneous problem structures. It integrates limited-memory techniques, low-rank matrix decomposition, and quasi-Newton update theory, with task-specific designs for eigenvalue computation, tensor decomposition, and nonlinear regression. Experimental results demonstrate that the proposed method maintains convergence accuracy while substantially reducing memory consumption and computational complexity compared to standard dense or conventional compact Hessian approximations.

Multi-shape Laplacian eigenbases are inherently incompatible due to spectral misalignment across non-isometric shapes, hindering consistent cross-shape functional representation and transfer. Method: This paper proposes a correspondence-guided coupled quasi-harmonic basis construction method. By approximately jointly diagonalizing Laplacian operators and incorporating geometrically consistent correspondence priors—such as stable region indicator functions—it constructs shape-shared, frequency-aligned quasi-harmonic bases. Contribution/Results: The method overcomes the limitation of conventional per-shape independent spectral bases and establishes, for the first time, a geometrically consistent multi-shape frequency-domain representation. It significantly improves accuracy and consistency in cross-shape function transfer across applications including shape editing, pose transfer, non-rigid matching, and shape similarity assessment. The framework provides a scalable theoretical foundation and practical toolset for multi-shape spectral analysis and collaborative processing.

Computing pointwise mean curvature on high-dimensional data manifolds typically incurs a prohibitive $O(m^4)$ complexity due to the explicit construction of the shape operator, hindering scalability. This work proposes an efficient estimation method that circumvents the explicit formation of the Hessian matrix by leveraging an exact algebraic identity. The approach combines covariance eigendecomposition, the cyclic property of the trace operator, and truncated singular value decomposition, while analytically approximating the expectation of nullspace eigenvectors under the Haar measure. This strategy substantially reduces computational complexity, achieving 50–300× speedup on real-world datasets with negligible loss in accuracy, thereby rendering mean curvature a practical geometric feature for machine learning applications.

This work addresses the instability and lack of convergence guarantees in training neural networks on Riemannian manifolds, which arise from neglecting manifold curvature. To this end, the authors propose a curvature-aware training framework that replaces Euclidean derivative supervision with a covariant Sobolev loss, aligns gradients via parallel transport, and incorporates a Laplace–Beltrami regularizer to enhance stability. For the first time, the framework explicitly integrates manifold curvature and the Jacobian of parallel transport into Sobolev-based training, enabling the derivation of geometry-dependent smoothness constants. Building on this, the authors establish curvature-aware linear and quadratic convergence theories. Empirical validation on Lie groups such as SO(3) and SE(3), as well as tasks involving surface imaging and physics-informed learning, demonstrates stable training dynamics and confirms the predicted theoretical convergence rates.

This work aims to unify and extend the applicability of Newton-type optimization algorithms by incorporating curvature information into gradient updates in a more general manner. To this end, the authors propose the Generalized Quadratic Gradient (GQG) framework, which abstracts the common structural properties of existing methods and expresses the update rule as the integration of the gradient with any positive-definite curvature matrix satisfying the stationarity condition of a local quadratic model. This framework transcends the conventional reliance on specific Hessian approximations—such as diagonal matrices or BFGS—and establishes a universal optimization paradigm applicable to arbitrary positive-definite curvature matrices. Grounded in a generalized analysis of local quadratic models and quasi-Newton theory, this study provides a rigorous theoretical foundation and methodological guidance for designing more flexible and efficient curvature-aware optimization algorithms.

This work addresses optimization problems defined over products of simplices, such as low-rank learning of discrete multivariate probability distributions and function data registration based on the Square-Root Velocity Function (SRVF) representation. To tackle the inherent constraints, the authors propose a smooth reparameterization that is strictly convex element-wise, transforming the constrained problem into an unconstrained optimization over a Riemannian manifold. The resulting problem is solved via Riemannian gradient descent (RGD). Theoretical analysis shows that this reparameterization maps second-order KKT points on the manifold to weak second-order KKT points of the original problem, ensuring theoretical soundness while enhancing computational efficiency. Experiments demonstrate that RGD significantly outperforms projected gradient descent (PGD), achieving more accurate shape-preserving registration in functional data and efficiently solving probability tensor decomposition tasks.

This work investigates the generalization ability of nonlinear least-squares models with ridge regularization at local minima attained after training. Leveraging average algorithmic stability analysis, it characterizes the local geometry of parameter space through the empirical Jacobian Gram matrix and a residual-curvature term. For the first time, an effective dimension is defined with respect to the trained model rather than the initialization point, and combined with curvature and covering complexity to derive a generalization error bound that depends on the geometric structure of learned features rather than the number of parameters. The framework explicitly links this bound to intrinsic data manifold dimensions or the number of stable activation regions in ReLU networks. Experiments confirm Jacobian contraction during training, the tightness of residual-curvature linearization, and a strong alignment between the proposed bound and observed generalization gaps.