We generalize the Brezzi-Rappaz-Raviart approximation theorem, which allows to obtain existence and a priori error estimates for approximations of solutions to some nonlinear partial differential equations. Our contribution lies in the fact that we typically allow for nonlinearities having merely Lipschitz regularity, while previous results required some form of differentiability. This is achieved by making use of the theory of metrically regular mappings, developed in the context of variational analysis. We apply this generalization to derive some quasi-optimal error estimates for finite element approximations to solutions of viscous Hamilton-Jacobi equations and second order mean field games systems.
Bifurcation characterizes the qualitative changes in parameterized dynamical systems and is one of the major topics in the field. In this work, we study combinatorial bifurcations within the framework of combinatorial dynamical systems—a young but already well-established theory. We introduce the Conley–Morse persistence barcode, a compact algebraic descriptor of combinatorial bifurcations. This barcode captures structural changes in a dynamical system at the level of Morse decompositions and provides a characterization of the nature of observed transitions in terms of the Conley index. The construction of the Conley–Morse persistence barcode builds upon ideas from topological persistence. Specifically, we consider a persistence module obtained from the Conley index of invariant sets indexed over a poset. Using gentle algebras, we prove that this module decomposes into simple intervals (bars) and compute them by adapting the zigzag persistence algorithm to our purpose.
We present a general framework, treating Lipschitz domains in Riemannian manifolds, that provides conditions guaranteeing the existence of norming sets and generalized local polynomial reproductions—a powerful tool used in the analysis of various mesh-free methods and a mesh-free method in its own right. As a key application, we prove the existence of smooth local polynomial reproductions on compact subsets of algebraic manifolds in ℝ^n with Lipschitz boundary. These results are then applied to derive new findings on the existence, stability, regularity, locality, and approximation properties of shape functions for a coordinate-free moving least squares approximation method on algebraic manifolds, which operates directly on point clouds without requiring tangent plane approximations. There are two appendices: the first derives high order Markov inequalities for polynomials on algebraic manifolds and the second gives instructions for calculating the dimension of the space of degree m polynomials restricted to a real algebraic variety.
In 1946, S. Ulam invented Monte Carlo method, which has since become the standard numerical technique for making statistical predictions for long-term behaviour of dynamical systems. We show that this, or in fact any other numerical approach can fail for the simplest non-linear discrete dynamical systems given by the logistic maps f_a(x)=ax(1-x) of the unit interval. We show that there exist computable real parameters a∈ (0,4) for which almost every orbit of f_a has the same asymptotical statistical distribution in [0,1], but this limiting distribution is not Turing computable.
Coupling arguments are a central tool for bounding the deviation between two stochastic processes, but traditionally have been limited to Wasserstein metrics. In this paper, we apply the shifted composition rule—an information-theoretic principle introduced in our earlier work [3]—in order to adapt coupling arguments to the Kullback–Leibler (KL) divergence. Our framework combines the strengths of two previously disparate approaches: local error analysis and Girsanov’s theorem. Akin to the former, it yields tight bounds by incorporating the so-called weak error, and is user-friendly in that it only requires easily verified local assumptions; and akin to the latter, it yields KL divergence guarantees and applies beyond Wasserstein contractivity. We apply this framework to the problem of sampling from a target distribution π . Here, the two stochastic processes are the Langevin diffusion and an algorithmic discretization thereof. Our framework provides a unified analysis when π is assumed to be strongly log-concave (SLC), weakly log-concave (WLC), or to satisfy a log-Sobolev inequality (LSI). Among other results, this yields KL guarantees for the randomized midpoint discretization of the Langevin diffusion. Notably, our result: (1) yields the optimal O(√(d)/ε ) rate in the SLC and LSI settings; (2) is the first result to hold beyond the 2-Wasserstein metric in the SLC setting; and (3) is the first result to hold in any metric in the WLC and LSI settings.
We construct fully-discrete schemes for the Benjamin–Ono, Calogero–Sutherland DNLS, and cubic Szegő equations on the torus, which are exact in time with spectral accuracy in space. We prove spectral convergence for the first two equations, of order K^-s+1 in L^2 norm for initial data in H^s(𝕋) , s>1 , with an error constant depending linearly on the final time instead of exponentially. These schemes are based on explicit formulas, which have recently emerged in the theory of nonlinear integrable equations. Numerical simulations show the strength of the newly designed methods both at short and long time scales, thanks to the remarkable fact that the computational cost of the method is independent of the final time. These schemes open doors for the understanding of the long-time dynamics of integrable equations.
In [3], the authors used the Legendre transform to give a tractable method for studying Topological Data Analysis (TDA) in terms of sums of Gaussian kernels. In this paper, we prove a variant for sums of cosine similarity-based kernel functions, which requires considering the more general “c-transform” from optimal transport theory [18]. We then apply these methods to a point cloud arising from a recent breakthrough study, which exhibits a toroidal structure in the brain activity of rats [12]. A key part of this application is that the transport map and transformed density function arising from the theorem replace certain delicate preprocessing steps related to density-based denoising and subsampling.
Transformers are deep neural network architectures that underpin the recent successes of large language models. Unlike more classical architectures that can be viewed as point-to-point maps, a Transformer acts as a measure-to-measure map implemented as specific interacting particle system on the unit sphere: the input is the empirical measure of tokens in a prompt and its evolution is governed by the continuity equation. In fact, Transformers are not limited to empirical measures and can in principle process any input measure. As the nature of data processed by Transformers is expanding rapidly, it is important to investigate their expressive power as maps from an arbitrary measure to another arbitrary measure. To that end, we provide an explicit choice of parameters that allows a single Transformer to match N arbitrary input measures to N arbitrary target measures, under the minimal assumption that every pair of input-target measures can be matched by some transport map.
We develop a parabolic inf-sup theory for a modified TraceFEM semi-discretization in space of the heat equation posed on a stationary hypersurface embedded in ℝ^n . We consider the normal derivative volume stabilization and add an L^2 -type stabilization to the time derivative. We assume that the representation of and the integration over the surface are exact, however, all our results are independent of how the surface cuts the bulk mesh. For any mesh for which the method is well-defined, we establish necessary and sufficient conditions for inf-sup stability of the proposed TraceFEM in terms of H^1 -stability of a stabilized L^2 -projection and of an inverse inequality constant that accounts for the lack of conformity of TraceFEM. Furthermore, we prove that the latter two quantities are bounded uniformly for a sequence of shape-regular and quasi-uniform bulk meshes. We derive several consequences of uniform discrete inf-sup stability, namely uniform well-posedness, discrete maximal parabolic regularity, parabolic quasi-best approximation, convergence to minimal regularity solutions, and optimal order-regularity energy and L^2 L^2 error estimates. We show that the additional stabilization of the time derivative restores optimal conditioning of time-discrete TraceFEM typical of fitted discretizations.
Gaussian random fields play an important role in many areas of science and engineering. In practice, they are often simulated by sampling from a high-dimensional multivariate normal distribution, which arises from the discretisation of a suitable precision operator. Existing methods such as Cholesky factorization and Gibbs sampling become prohibitively expensive on fine meshes due to their high computational cost. In this work, we revisit the Multigrid Monte Carlo (MGMC) algorithm developed by Goodman Sokal (Physical Review D 40.6, 1989) in the quantum physics context. While the authors of this paper conclude that MGMC does not overcome critical slowing down in simulations of field theories near phase transitions, we demonstrate here that it has the potential to significantly accelerate sampling in spatial statistics. The class of Gaussian Random Fields we consider includes those with Matérn covariance, but is more general in that it also allows for non-stationary covariance functions. To show that MGMC can overcome the limitation of existing methods, we establish a grid-size-independent convergence theory based on the link between linear solvers and samplers for multivariate normal distributions, drawing on standard multigrid convergence arguments. We then apply this theory to linear Bayesian inverse problems. This application is achieved by extending the standard multigrid theory to operators with a low-rank perturbation. Moreover, we develop a novel bespoke random smoother which takes care of the low-rank updates that arise in constructing posterior moments. In particular, we prove that Multigrid Monte Carlo is algorithmically optimal in the limit of the grid-size going to zero. Numerical results support our theory, demonstrating that Multigrid Monte Carlo can be significantly more efficient than alternative methods when applied in a Bayesian setting.
In this paper we introduce Crouzeix-Raviart elements of general polynomial order k and spatial dimension d≥ 2 for simplicial finite element meshes. We give explicit representations of the non-conforming basis functions and prove that the conforming companion space, i.e., the conforming finite element space of polynomial order k is contained in the Crouzeix-Raviart space. We prove a direct sum decomposition of the Crouzeix-Raviart space into (a subspace of) the conforming companion space and the span of the non-conforming basis functions. Degrees of freedom are introduced which are bidual to the basis functions and give rise to the definition of a local approximation/interpolation operator. In two dimensions or for k=1 , these degrees of feedom can be split into simplex and ( d-1) dimensional facet integrals in such a way that, in a basis representation of Crouzeix-Raviart functions, all coefficients which correspond to basis functions related to lower-dimensional faces in the mesh are determined by these facet integrals. It will also be shown that such a set of degrees of freedom does not exist in higher space dimension and k>1 .
We numerically investigate the nearly self-similar blowup of the generalized axisymmetric Navier–Stokes equations. First, we rigorously derive the axisymmetric Navier–Stokes equations with swirl in both odd and even dimensions, marking the first such derivation for dimensions greater than three. Building on this, we generalize the equations to arbitrary positive real-valued dimensions, preserving many known properties of the 3D axisymmetric Navier–Stokes equations. To address scaling instability, we dynamically vary the space dimension to balance advection scaling along the r and z directions. A major contribution of this work is the development of a novel two-scale dynamic rescaling formulation, leveraging the dimension as an additional degree of freedom. This approach enables us to demonstrate a one-scale self-similar blowup with solution-dependent viscosity. Notably, the self-similar profile satisfies the axisymmetric Navier–Stokes equations with constant viscosity. We observe that the effective dimension is approximately 3.188 and appears to converge toward 3 as background viscosity diminishes. Furthermore, we introduce a rescaled Navier–Stokes model derived by dynamically rescaling the axial velocity in 3D. This model retains essential properties of 3D Navier–Stokes. Our numerical study shows that this rescaled Navier–Stokes model with two constant viscosity coefficients exhibits a nearly self-similar blowup with maximum vorticity growth on the order of O(10^30) .
This paper considers the problem of understanding the behavior of a general class of accelerated gradient methods on smooth nonconvex functions. Motivated by some recent works that have proposed effective algorithms, based on Polyak’s heavy ball method and the Nesterov accelerated gradient method, to achieve convergence to a local minimum of nonconvex functions, this work proposes a broad class of Nesterov-type accelerated methods and puts forth a rigorous study of these methods encompassing the escape from saddle points and convergence to local minima through both an asymptotic and a non-asymptotic analysis. In the asymptotic regime, this paper answers an open question of whether Nesterov’s accelerated gradient method (NAG) with variable momentum parameter avoids strict saddle points almost surely. This work also develops two metrics of asymptotic rates of convergence and divergence, and evaluates these two metrics for several popular standard accelerated methods such as the NAG and Nesterov’s accelerated gradient with constant momentum (NCM) near strict saddle points. In the non-asymptotic regime, this work provides an analysis that leads to the “linear” exit time estimates from strict saddle neighborhoods for trajectories of these accelerated methods as well the necessary conditions for the existence of such trajectories. Finally, this work studies a sub-class of accelerated methods that can converge in convex neighborhoods of nonconvex functions with a near optimal rate to a local minimum and at the same time this sub-class offers superior saddle-escape behavior compared to that of NAG.
Classical results in asymptotic statistics show that the Fisher information matrix controls the difficulty of estimating a statistical model from observed data. In this work, we introduce a companion measure of robustness of an estimation problem: the radius of statistical efficiency (RSE) is the size of the smallest perturbation to the problem data that renders the Fisher information matrix singular. We compute RSE up to numerical constants for a variety of testbed problems, including principal component analysis, generalized linear models, phase retrieval, bilinear sensing, and matrix completion. Interestingly, we observe a precise reciprocal relationship between RSE and the intrinsic complexity/sensitivity of the problem instance, paralleling the classical Eckart-Young theorem in numerical analysis. To establish our results, we develop theory for spectral functions of measures that extends well-known results from matrix analysis and eigenvalue optimization-a contribution that may be of interest beyond our immediate findings.
We show that the cohomology of the Regge complex in three dimensions is isomorphic to $\mathcal{H}^{{\scriptscriptstyle \bullet}}_{dR}(\Omega)\otimes\mathcal{RM}$, the infinitesimal-rigid-body-motion-valued de~Rham cohomology. Based on an observation that the twisted de~Rham complex extends the elasticity (Riemannian deformation) complex to the linearized version of coframes, connection 1-forms, curvature and Cartan's torsion, we construct a discrete version of linearized Riemann-Cartan geometry on any triangulation and determine its cohomology.
This paper concerns the long-standing question of representing (totally) anti-symmetric functions in high dimensions. We propose a new ansatz based on the composition of an odd function with a fixed set of anti-symmetric basis functions. We prove that this ansatz can exactly represent every anti-symmetric and continuous function and the number of basis functions has efficient scaling with respect to dimension (number of particles). The singular locus of the anti-symmetric basis functions is precisely identified.
We introduce a unified framework of symmetric resonance based schemes which preserve central symmetries of the underlying PDE. We extend the resonance decorated trees approach introduced in arXiv:2005.01649 to a richer framework by exploring novel ways of iterating Duhamel's formula, capturing the dominant parts while interpolating the lower parts of the resonances in a symmetric manner. This gives a general class of new numerical schemes with more degrees of freedom than the original scheme from arXiv:2005.01649. To encapsulate the central structures we develop new forest formulae that contain the previous class of schemes and derive conditions on their coefficients in order to obtain symmetric schemes. These forest formulae echo the one used in Quantum Field Theory for renormalising Feynman diagrams and the one used for the renormalisation of singular SPDEs via the theory of Regularity Structures. These new algebraic tools not only provide a nice parametrisation of the previous resonance based integrators but also allow us to find new symmetric schemes with remarkable structure preservation properties even at very low regularity.
We introduce and analyze a numerical approximation of the porous medium equation with fractional potential pressure introduced by Caffarelli and V\'azquez: \[ \partial_t u = \nabla \cdot (u^{m-1}\nabla (-\Delta)^{-\sigma}u) \qquad \text{for} \qquad m\geq2 \quad \text{and} \quad \sigma\in(0,1). \] Our scheme is for one space dimension and positive solutions $u$. It consists of solving numerically the equation satisfied by $v(x,t)=\int_{-\infty}^xu(x,t)dx$, the quasilinear non-divergence form equation \[ \partial_t v= -|\partial_x v|^{m-1} (- \Delta)^{s} v \qquad \text{where} \qquad s=1-\sigma, \] and then computing $u=v_x$ by numerical differentiation. Using upwinding ideas in a novel way, we construct a new and simple, monotone and $L^\infty$-stable, approximation for the $v$-equation, and show local uniform convergence to the unique discontinuous viscosity solution. Using ideas from probability theory, we then prove that the approximation of $u$ converges weakly-$*$, or more precisely, up to normalization, in $C(0,T; P(\mathbb{R}))$ where $P(\mathbb{R})$ is the space of probability measures under the Rubinstein-Kantorovich metric.The analysis include also fundamental solutions where the initial data for $u$ is a Dirac mass. Numerical tests are included to confirm the results. Our scheme seems to be the first numerical scheme for this type of problems.
I introduce the concept of a persistence diagram (PD) bundle, which is the space of PDs for a fibered filtration function (a set $\{f_p: \mathcal{K}^p \to \mathbb{R}\}_{p \in B}$ of filtrations that is parameterized by a topological space $B$). Special cases include vineyards, the persistent homology transform, and fibered barcodes for multiparameter persistence modules. I prove that if $B$ is a smooth compact manifold, then for a generic fibered filtration function, $B$ is stratified such that within each stratum $Y \subseteq B$, there is a single PD "template" (a list of "birth" and "death" simplices) that can be used to obtain the PD for the filtration $f_p$ for any $p \in Y$. If $B$ is compact, then there are finitely many strata, so the PD bundle for a generic fibered filtration on $B$ is determined by the persistent homology at finitely many points in $B$. I also show that not every local section can be extended to a global section (a continuous map $s$ from $B$ to the total space $E$ of PDs such that $s(p) \in \textrm{PD}(f_p)$ for all $p \in B$). Consequently, a PD bundle is not necessarily the union of "vines" $\gamma: B \to E$; this is unlike a vineyard. When there is a stratification as described above, I construct a cellular sheaf that stores sufficient data to construct sections and determine whether a given local section can be extended to a global section.
Generalised hardness of approximation (GHA) is the phenomenon that one can easily compute an ϵ -approximation to a solution of a computational problem for ϵ> ϵ _1 > 0 , but for ϵ < ϵ _1 (the approximation threshold) it suddenly becomes hard, for example, non-computable or intractable (non-polynomial time). In this paper we demonstrate the phenomenon that GHA happens when using AI techniques for solving inverse problems, namely training neural networks (NNs) to optimally perform on the training data. In particular, for any non-zero underdetermined linear inverse problem the following phase transition can occur: For a certain family of training sets Ω , one can prove the existence of optimal NNs for solving the inverse problem for each 𝒯∈Ω , however, these optimal neural networks can only be computed to a certain accuracy ϵ _1 > 0 . Below the approximation threshold ϵ _1 , not only does it become intractable to compute the NNs, it becomes impossible regardless of computing power, and no randomised algorithm can solve the problem with probability better than 1/2. Moreover, despite the existence of a stable optimal NN, any attempts of computing it below two times the approximation threshold 2ϵ _1 will yield an unstable NN. Our results use and extend the current mathematical framework of the Solvability Complexity Index (SCI) hierarchy and initiate a program for analysing the GHA phenomenon throughout computational mathematics and AI. GHA generalises the phenomenon of hardness of approximation in discrete computations to arbitrary computational problems.
