Papers

We develop a framework for offline full-window branch discovery, optionally followed by online continuation with an ensemble Kalman filter (EnKF). Three mechanisms drive the branch search: adjoint path-kernel (APK) differentiation balances kernel differentiation and correction-stabilized path perturbation, shifting the optimization from exploration to exploitation; an optimized Gaussian initial law broadens the search over initial-state basins; and loss-weighted mixing across independent runs recombines successful path components. We may then select an interior state using a local loss and continue online with an EnKF. In 40-dimensional Lorenz-96 experiments, the mean offline path RMSE of APK is 4.3 times smaller than that of population weak-$\mathrm{4D\text{-}Var}_x$. The resulting APK-EnKF method has a mean online RMSE 64 times smaller than that of ordinary EnKF.

In the Navier–Stokes equations, incompressibility allows rewriting the viscous term in various forms leading to distinct numerical properties and flow descriptions. Furthermore, models accounting for non-Newtonian, thermal or turbulent effects often break the constant-viscosity assumption, thereby producing additional consistency terms. In this context, the present work compares the classical symmetric-gradient diffusion term with more recent variable-viscosity generalizations of the Laplacian and rotational forms. We discuss, analyze and test their differences with respect to implementation, efficiency, numerical stability and outflow boundary conditions. With a focus on time-dependent flows, we consider second-order implicit-explicit (IMEX) temporal discretizations aimed at improving efficiency and numerical stability. Through a rigorous stability analysis, we show how selected explicit treatments can bypass algorithmic nonlinearities without inducing CFL conditions. Our numerical results highlight important differences between the three viscous formulations—especially in the presence of outflow boundaries, for which the generalized Laplacian form proves more suitable in diffusion-dominated regimes. %(as widely known for constant viscosity).

Considering the computational advantages of virtual element methods (VEM), this work employs a conforming VEM for the numerical approximation of the Kelvin–Voigt model of viscoelasticity. For clarity and simplicity, we focus on the primal formulation. The spatial discretization is carried out using the virtual element method, while the temporal discretization is handled via the {second-order Crank–Nicolson} scheme. We establish the well-posedness of both the semi-discrete and fully discrete problems and derive a priori error estimates. Several representative numerical examples are presented to validate the theoretical results and to demonstrate the effectiveness of the proposed formulation.

Hyperbolic conservation laws pose a challenging setting for variational data assimilation, since the Gaussian assumption imposes a smoothness that smears jump discontinuities. A sparsity promoting regularization can help to mitigate this problem by incorporating a structural prior into the variational objective. Care must be taken when choosing this prior because it serves to counterbalance both a structurally misspecified background-error covariance, which introduces spurious oscillations over smooth regions, as well as spatially sparse observations, which leave too little data to constrain the state. Moreover, state variables of hyperbolic conservation laws typically do not have piecewise constant structure, which is an assumption made when using standard sparsity-promoting operators, such as total variation. Indeed the underlying variability of hyperbolic conservation law solutions are neither known in advance nor fixed. Hence using higher order total variation is also not suitable. Here we introduce a new regularization term for data assimilation, built on the sparsity promoting residual transform that neither commits to a fixed smoothness order nor requires the underlying variability to be known a priori. We solve the resulting nonconvex objective within a three-dimensional variational framework using generalized sparse Bayesian learning. This approach replaces a global \(\ell _{1}\) penalty with a hierarchical conditional Gaussian prior on the residual transform, where a learned hyper-prior yields location-specific weights attached to each spatial point. Numerical experiments across scalar, shallow-water, and Euler test problems show that this regularization term resolves discontinuities more accurately than total variation, and that its advantage grows as the background-error covariance becomes misspecified and the observations become sparser.

In this note, we provide a short and direct proof that approximately solving $Ax=b$ to relative error $\varepsilon$, where $A$ has condition number $\kappa$ and unrestricted dimension, requires $\Omega(\kappa\log(1/\varepsilon))$ matrix-vector multiplications in the worst case, even for randomized algorithms. This essentially recovers the lower bound of Dereziński, Epperly and Meyer [2026] for this setting, whose elegant and more general approach inspired us to seek a short direct proof. A straightforward reduction implies the classical $\Omega(\sqrt{\kappa}\log(1/\varepsilon))$ lower bound for optimizing strongly convex quadratic functions, applicable to randomized algorithms.

In a recent numerical study, we found a high-precision singular profile for the Euler equations on the unbounded domain $\mathbb{R}^3$. The present manuscript complements that study by establishing in detail a preliminary framework for proving (nonlinear) stability of the approximate self-similar profile, reducing the analysis to a large but finite collection of explicit estimates and computable constants. Conditional on rigorous certification of the estimates and constants appearing in the argument, and on the candidate profile satisfying the required nonlinear stability conditions, the framework closes the stability proof and, crucially, allows the resulting stable rescaled profile to be reconstructed as an admissible solution in the original variables that becomes singular in finite physical time. With the overall stability and reconstruction mechanisms formulated, the remaining work within this approach is largely quantitative: determining whether the explicit constants and margins can be rigorously certified with sufficient positive margin and, where necessary, sharpening selected analytic estimates.

Evaluating layer potentials with rapidly decaying kernels becomes computationally challenging in the presence of boundary layers, as occurs for Helmholtz layer potentials with large imaginary wavenumbers. We introduce Quadrature by Matched Asymptotic Expansion (QBMAX), a high-order method that extends the Quadrature by Expansion (QBX) framework by incorporating the kernel's asymptotic decay behavior directly into the Taylor expansion of the potential. Numerical results in both 2D and 3D show improvements of up to four digits over QBX for the tested large-parameter cases, while the methods can perform similarly for smaller parameters or more strongly curved geometries. A symbolic weighted operation-count model indicates comparable costs for forming the two expansion kernels. The Taylor expansions have the same formal order, while QBMAX reduces the truncation-error constants in the flat-boundary model and generally exhibits less convergence degradation in the reported experiments.

We develop a high-order, fully explicit operator-splitting method for 2.5D ideal magnetohydrodynamics on Cartesian meshes. The hydrodynamic subflow is advanced by an entropy-stable discontinuous Galerkin spectral element method equipped with stagewise oscillation elimination and a positivity-preserving limiter. The magnetic–velocity subflow uses compatible finite elements and mass-lumped reconstructions of electric field and current density. Its discrete-curl update exactly preserves the global $H(\mathrm{div})$ divergence-free subspace, while the ideal semidiscretization balances kinetic, magnetic, and internal energy. Curl-form artificial resistivity and a direction-resolved velocity filter return removed magnetic and kinetic energy to internal energy; consequently, the stabilized magnetic stage preserves positive internal energy and satisfies a discrete entropy inequality. The two solvers are composed by a second-order hydrodynamic–magnetic–hydrodynamic Strang splitting, yielding a matrix-free scheme that retains global mass, nodal admissibility, and magnetic divergence on accepted steps. Smooth Alfvén-wave and advected-vortex tests reveal an even–odd convergence pattern in the magnetic polynomial degree, verify second-order temporal accuracy, and show that stabilization preserves high-order accuracy. Field-loop, Orszag–Tang, rotor, MHD blast-wave, and Kelvin–Helmholtz calculations demonstrate robust performance for nonsmooth multidimensional flows and show that artificial resistivity suppresses grid-scale magnetic oscillations while retaining the resolved structures.

A robust finite-volume framework is presented for the simulation of compressible multiphase flows with surface tension across a wide range of Mach numbers. The method is based on a two-pressure, six-equation diffuse interface model incorporating viscous, gravitational, and capillary effects through the continuum surface force formulation. To consistently account for capillary-induced pressure jumps, an HLLC Riemann solver is developed using generalized Riemann invariant analysis, and the instantaneous pressure relaxation procedure is modified to preserve the Laplace pressure jump during phase equilibration. To overcome the excessive numerical diffusion of conventional approximate Riemann solvers in the low-Mach regime, a robust low-Mach correction is proposed by extending our previous formulation with a modified scaling strategy that remains stable in regions of strong pressure variation. The resulting method retains accuracy from nearly incompressible flows to compressible regimes while preserving the robustness of the six-equation formulation. The numerical framework is validated using a series of benchmark problems involving surface tension, viscosity, gravity, and compressibility. The results demonstrate accurate prediction of interface dynamics, capillary pressure, and low-Mach flow features, while significantly reducing numerical dissipation without compromising stability. The proposed methodology provides an efficient and reliable approach for the simulation of complex multiphase flows spanning a broad range of flow regimes.

We develop a vectorial lattice Boltzmann model with a space-time adaptive relaxation coefficient, combined with adaptive time-stepping for compressible Euler dynamics with generic equation of state. The model, due to special form of the relaxation coefficient devised here is shown to converge to the Euler limit with second-order accuracy in the absence of shocks under acoustic scaling. The solver is robust and able to properly capture compressible gas dynamics for both ideal and non-ideal equations of state, including the Bethe–Zel'dovich–Thompson regime. This is verified and demonstrated through a variety of configurations of increasing complexity.

This paper presents a unified invariant-domain-preserving (IDP) framework for hybrid discretizations of hyperbolic conservation laws, including active flux and PAMPA methods, in which cell averages are updated conservatively while cell-boundary point values evolve under a possibly non-conservative operator. The main challenge is to preserve admissibility for these two coupled sets of states under a single local stability condition, without adding evolved degrees of freedom or relying on post-update repairs. For the point-value update, we introduce an admissibility transform based on a barrier–Legendre map for convex admissible interiors described by concave constraints, and prove that the inverse map is globally defined and Lipschitz continuous on the relevant sets. For the conservative cell-average update, we establish a structural obstruction theorem: the single-state continuous physical flux built from admissible boundary traces alone cannot provide a conservative IDP guarantee, for any prescribed CFL number, when internal reconstruction values are uncontrolled. This identifies the missing local control that must be supplied by an additional IDP flux mechanism. To realize this mechanism explicitly, we combine cell average decompositions (CAD), geometric quasilinearization, and local a priori scaling to construct admissible, generally discontinuous trace states before flux evaluation. Under an explicit trace-based CFL condition, with constants determined by local CAD weights and trace-state wave-speed bounds, the coupled hybrid update preserves the prescribed invariant domain. Concrete third-order schemes are developed on triangular, Cartesian, convex quadrilateral, and general convex polygonal meshes. Numerical results demonstrate the designed order of accuracy for smooth solutions and the strict preservation of physical admissibility.

A combination of Taylor polynomial approximation and averaging results in a powerful quasi-interpolation technique that can be used to analyze the approximation properties of finite element spaces. With this particular approach, the notion of affine equivalent families is completely bypassed in the sense that the degrees of freedom are defined locally rather than being exported from a reference element. In particular, domains with curved boundaries can be handled in a natural manner. Finally, the ease of its application to a particular finite element makes it a valuable pedagogical tool.

Mean field planning on a surface prescribes initial and terminal densities and minimizes a transport energy subject to the continuity equation. Proximal algorithms for this problem repeatedly solve a time–space Poisson equation, whose temporal derivative and surface gradient determine the density and momentum corrections. We study a gradient-enhanced approximation of this constraint projection using finite differences in time, surface finite elements in space, temporal polynomial preserving recovery, and spatial parametric polynomial preserving recovery. The same construction is incorporated into ISTA, FISTA, and Douglas–Rachford splitting. We distinguish the recovered update from an exact discrete projection and derive residual identities and conditional finite-iteration perturbation bounds that retain data, boundary, and linear-solver errors. Existing derivative-recovery estimates identify a higher-order contribution under suitable regularity and mesh assumptions; they do not by themselves establish convergence of the outer optimization iteration. Available numerical illustrations on the sphere and a more complicated algebraic surface are discussed together with the limits of the recorded refinement data.

We develop a family of weakly symmetric and pointwise traceless tangential-normal tensor finite elements in arbitrary space dimension and for all polynomial orders. Symmetry is imposed through local cell moments, while the only globally coupled stress degrees of freedom are tangential-normal facet moments; no vertex degrees of freedom are required. In dimensions three and higher, a lowest-order linear enrichment restores the rigid-motion facet control required for discrete Korn stability. As a principal application, we construct a distributional mixed method for the incompressible Brinkman equations using the physical viscous stress. Coupled with divergence-conforming BDM velocities and discontinuous pressures, the method is stabilization-free, uniformly stable with respect to the viscosity parameter, exactly divergence-free, and pressure-robust. We establish optimal-order error estimates in the natural norms. Under suitable parameter-explicit regularity assumptions, we also obtain a parameter-uniform boundary-layer estimate with optimal Darcy approximation order. Relaxing tangential-normal continuity yields an algebraically equivalent stress-hybridized formulation and a stabilization-free virtual element realization.

We present and analyze an iterative refinement (IR) framework for improving the precision dependence of algorithms that combine a quantum linear system algorithm (QLSA) with quantum state tomography. Existing QLSAs achieve polylogarithmic dependence on the inverse error tolerance to prepare a quantum state encoding the solution, but extracting a classical description of such a state via tomography typically introduces a polynomial dependence on the target precision. To retain polylogarithmic dependence in the inverse error tolerance throughout the entire process, we develop an IR scheme that solves a sequence of related linear systems to progressively refine the solution while requiring only fixed-precision quantum subroutines, finally obtaining a classical description of a high-precision solution. We analyze the proposed framework under three input models: quantum-read/classical-write RAM (QRAM), linear combinations of unitaries (LCU), and sparse-access oracles. We evaluate the proposed scheme through numerical experiments on both quantum simulators and real quantum hardware. The results demonstrate that iterative refinement efficiently improves the precision of the solution and exhibits robustness to hardware noise.

Directly self-starting implicit methods are attractive for transient analysis because they avoid auxiliary starting procedures while retaining the original first- or second-order governing equations. However, most existing formulations usually fix the last sub-step at the end of each time interval, which restricts the attainable order. This study develops a generalized $s$-sub-step implicit framework by releasing this constraint and treating all sub-step locations as design variables. The resulting methods admit a unified Runge–Kutta representation for both first- and second-order transient systems and preserve identical effective matrices over all sub-steps. Accuracy conditions are derived by simultaneously matching the numerical amplification factor and load operator, thereby accounting for both homogeneous and forced responses. For $s=1,~\cdots,~6$, two complementary families are obtained: $s$th-order members with user-controllable high-frequency numerical dissipation and adjustable sub-step locations, and $(s+1)$th-order members obtained by selecting the sub-step locations, with fixed dissipation. The latter reach up to seventh-order accuracy without increasing the number of sub-steps, although some high-order members are $A(\alpha)$-stable with stability angles extremely close to $90^\circ$. Analytical amplitude and phase errors further reveal parity-dependent superconvergence in undamped systems, and appropriate parameter selections can substantially increase either phase or amplitude accuracy beyond the formal order. Numerical benchmarks confirm the predicted convergence orders and the controllable suppression of spurious high-frequency responses.

Computing the closest point on a clothoid is a recurring task in geometric design, road and railway alignment, and path planning. The efficient algorithm of Frego and Bertolazzi addresses this problem, but its candidate-selection analysis assumes at most one local minimum per search interval. We exhibit admissible configurations with two local minima, raising the question of whether the existing strategy accounts for every possible minimum. Using the geometry of the clothoid evolute, we prove that, for any query point and any proper no-inflection planar clothoid segment with tangent-angle variation at most $2\pi$, the squared-distance function has at most three stationary points; if all three are local extrema, their order is min-max-min. This establishes the completeness of the original candidate-selection logic beyond the one-minimum premise. It also shows that no interior search is needed when neither endpoint derivative test is active, allowing unnecessary midpoint searches to be omitted while retaining numerical fallback. Numerical experiments demonstrate reductions in iteration count and evaluation time.

We develop a thermodynamically consistent diffuse-interface model for multicomponent membranes. The proposed model is derived from a coupled free energy functional that incorporates protein-dependent bending elasticity, diffuse surface tension, a volume penalty, and a membrane-associated Ohta–Kawasaki energy. Applying the Onsager variational principle, we derive a coupled $L^2$ gradient flow system and its energy dissipation law. We then construct a stabilized alternating ETD1 scheme and a stabilized alternating ETDRK2 scheme with Strang-type composition (alternating Strang-ETDRK2). To the best of our knowledge, the proposed alternating Strang-ETDRK2 scheme has not previously been developed and analyzed for coupled phase field systems. We further prove the discrete energy dissipation for both schemes under some regularity assumptions on the numerical solutions. Numerical experiments in two and three dimensions validate the discrete energy dissipation law, and present the protein segregation and membrane deformation produced by the proposed model.

In explicit dynamics, the critical time step is governed by the maximum eigenfrequency of the semi-discrete system. In Isogeometric Analysis (IGA), open knot vectors introduce a characteristic spectral boundary effect that can lead to high-frequency outliers and restrict the admissible time step. In the row-sum-lumped setting considered here, extending the computational patch and trimming away the exterior boundary functions mitigates this time-step penalty, but does not necessarily prevent boundary-adjacent basis functions from governing the maximum eigenfrequency. We show that reduced-support basis functions adjacent to the trimming boundary can remain critical spectral contributors even for knot-exact trimming, i.e., in the absence of small cut cells. Hence, small cut cells are not required for trimming-induced time-step penalties. Knot-exact and arbitrary trimming represent different degrees of support reduction within the same underlying mechanism. To control this effect, we propose the Boundary-Level-Constrained Refinement (BLCR) strategy, a local refinement constraint for LR- and THB-splines. BLCR constrains the refinement level of basis functions whose support intersects the trimming boundary relative to that of the refined interior. In all configurations investigated in this work, this constraint ensures that the maximum eigenfrequency is governed by untrimmed interior basis functions rather than by refined trimmed functions. Consequently, the trimming-induced high-frequency outliers governing the critical time step are suppressed, yielding a larger admissible time step than globally refined trimmed B-splines with the same interior resolution.[...]

We develop and analyse a new arbitrary-order method for the Kirchhoff–Love plate problem on general polygonal meshes. The construction relies on a discrete Hessian complex obtained through the Bernstein–Gelfand–Gelfand construction, which stacks a Stokes complex on top of a tensorised de Rham complex. A key ingredient of the analysis is a detailed study of the former. We specifically establish primal and adjoint consistency estimates, as well as uniform Poincaré inequalities, and use these results to derive the analytical properties of the resulting discrete Hessian complex. The proposed method is then shown to be coercive and to converge with order $k+1$ with respect to the mesh size, where $k\ge 0$ is the polynomial degree of the complex. Numerical experiments on polygonal meshes confirm the theoretical convergence rates.

A shell finite-element model of a stiffened aluminium wingbox panel, a laboratory structure shared by several experimental studies, is built with its plate-to-stiffener joints at the $78$ physical fastener positions and calibrated against measured modal data, a six-parameter sensitivity update of the substructure moduli within $\pm20\%$ reproducing the first five measured modes to $1.1$–$4.3\%$ with modal-assurance values of $0.88$–$1.00$. Making the fasteners nonlinear, at equal mesh, mass and damping, shows that joint nonlinearity enters the response strongly asymmetrically. Hardening is almost invisible, at most $+1\%$ in frequency, whereas softening or slipping the joints moves the first three resonances by $-2.1$, $-4.9$ and $-9.2\%$, and friction removes up to $80\%$ of the resonant peak before it recovers. Newmark integration of all $18{,}804$ degrees of freedom and a harmonic balance condensed exactly onto the joints and continued in arclength agree to $5.5\times10^{-4}$. A projection-based nonlinear model order reduction then locates the criterion that a jointed structure imposes on a reduced basis. The binding quantity is not the linear response, which eleven vectors reproduce to better than $0.01\%$, but the receptance of the structure between the joints, of which $126$ global eigenvectors carry about $2\%$, and without which the reduced model overpredicts the hardening shift of the fundamental by a factor of thirty and returns a value beyond the rigid-joint limit of the panel.

Applying a matrix function to a vector is a common operation in large-scale scientific computing, such as Bayesian inference. When the matrix is symmetric positive definite and the function is its inverse, the Lanczos algorithm gives a matrix-free approximation from matrix-vector products alone. Stopping after a fixed number of steps keeps the approximation on the explored Krylov subspace but assigns nothing to its orthogonal complement, which holds the smallest eigenvalues and therefore the largest contributions to the inverse. This truncation is usually treated as a silent numerical error rather than as measurable uncertainty. We show that a single Lanczos run yields both a matrix-free approximation of the inverse and a Gaussian distribution over the exact result, whose covariance measures what the truncation discarded and concentrates as the number of steps grows. The covariance restores the complement in two parts: a rank-one coupling at the truncation boundary, given exactly by the Lanczos residual, and an isotropic term on the remaining bulk, whose size we estimate by Projected Stochastic Lanczos Quadrature, reusing one Krylov basis across all probes. We prove consistency, with bias decaying exponentially in the quadrature depth and variance at the Monte Carlo rate

Numerically solving McKean-Vlasov stochastic differential equations is computationally challenging due to the compounding costs of discretizing in time and in the distribution of the solution. Multilevel ideas have been proposed to provide speed-ups. In this work, we study the multilevel Monte Carlo method proposed by Ricketson (2015) for equations whose drift and diffusion terms depend on the law of the solution $X_t$ through the expectation $E[R(X_t)]$. The scheme follows the single-ensemble paradigm, where particles interact across levels at each timestep. While cross-level feedback makes this scheme attractive in practice, the correlations it introduces have so far confined its cost-error analysis to a model problem with linear drift, deterministic diffusion, and $R$ the identity. We use additional coarse particles to enforce geometrically decaying coupling errors towards the coarser levels. This allows us to prove our main contribution, an $L^p$-error of $O(\epsilon)$ at cost $O(\epsilon^{-2-\delta})$ for any $p\ge2$ and $\delta>0$ (with a constant that grows as $\delta\to0$), assuming only global Lipschitz bounds on the drift, the diffusion, and $R$. An exploratory experiment is consistent with the derived rates and suggests that in practice the constant does not grow significantly for small $\delta$. Our methodology and proof strategy may also be useful for other single-ensemble multilevel schemes, such as multilevel ensemble Kalman filters.

We consider the problem of estimating a high-dimensional covariance matrix from a very limited number of samples. This problem is ubiquitous in computational fluid dynamics, where a small number of fluid snapshots must be used to construct a Gramian matrix determining a reduced-order model, as well as in computational geoscience, where a small ensemble of Earth system forecasts must be used to estimate the covariance matrix associated with the forecast uncertainty. It is common practice to regularize the small-sample covariance by imposing a "localization" structure that enforces a physically realistic correlation length scale, imposing a sparsity constraint, "shrinking" towards a prescribed target, or attenuating small correlations. We propose an alternate technique that regularizes the small-sample covariance by imposing hierarchical rank structure. Compared to regularization methods that assume sparsity such as spatial localization, hierarchical rank structure accommodates a wider range of covariance matrices, roughly corresponding to situations where long-range correlations vary more smoothly than short-range ones. It also results in a data-sparse matrix format that permits highly efficient matrix-vector products. We present theory and algorithms which show how to efficiently estimate a high-dimensional, hierarchically rank structured covariance matrix from limited samples. Through an error analysis and numerical experiments with a variety of model problems, we demonstrate that these techniques are effective at reducing sampling errors, and that in many cases they achieve smaller estimation error than conventional techniques.

This paper extends the convergence analysis of explicit exponential Runge–Kutta methods for linear parabolic problems $u'(t) + Au(t) = Bu(t)$, where $A$ generates an analytic semigroup and $B$ is relatively bounded with respect to $A$, from the third-order case to fourth-order schemes. By establishing the global error recursion relation and extending the defect-based analytical framework, we identify the terms responsible for stiff order reduction when $A$ and $B$ do not commute. Numerical experiments are performed on a non-commuting advection-diffusion problem to validate the theoretical results. Numerical tests using the classical four-stage, fourth-order schemes of Krogstad and Strehmel & Weiner exhibit an observed convergence order of approximately 2.75, which matches the theoretical prediction from the convergence analysis.

A nonlocal diffusion system describing the dynamics of multi-species populations and approximating the Shigesada-Kawasaki-Teramoto (SKT) model is analyzed in a bounded domain. The model consists of a system of integro-differential equations of Andreu-Mazón-Rossi-Toledo type, in which a nonlocal cross-diffusion operator acts on nonlinear diffusion potentials. We establish the global existence and uniqueness of strong solutions. The existence proof relies on suitable a priori estimates derived from an entropy inequality, whose proof requires the development of novel discrete chain-rule inequalities. As a by-product, these inequalities provide the basis for the construction of structure-preserving finite-volume schemes for the corresponding local SKT system. Uniqueness of solutions is established by means of a duality argument. Finally, one-dimensional numerical simulations illustrate the nonlocal-to-local limit and the qualitative behavior of solutions for superlinear and sublinear diffusion potentials.

This paper establishes a convergence theory for a continuous domain decomposition method for elliptic equations on manifolds. This method originated in \cite{lions2} in the setting of Euclidean domains and was later adapted and generalized to manifolds by \cite{qin_wang_wang}. Although its convergence was known, whether the convergence is geometric has remained open. We prove its geometric convergence and provide various estimates on the convergence rate.

Lewis and Overton conjectured in 2008 that, for BFGS applied to the one-dimensional tilted absolute value with their Armijo–Wolfe line search, every nonterminating execution converges at a trial-normalized rate determined only by the tilt parameter and independent of the initial data. We prove that this initialization-independent rate assertion is false. In fact, for every tilt parameter in an explicit open interval, we construct two nonterminating executions with different sharp trial-normalized convergence rates. The proof reduces the line search to a scale-free state consisting of the iterate sign and the zero-crossing stepsize, and verifies two periodic state cycles by elementary rational inequalities. This yields an explicit family of counterexamples and shows that the asymptotic trial-normalized behavior of nonsmooth BFGS can depend essentially on the initialization even for this canonical one-dimensional model.

The DIRECT algorithm is a deterministic global optimization method known for its versatility and balanced exploration-exploitation strategy. However, DIRECT-type algorithms are primarily effective for low-dimensional problems and often exhibit slow convergence as dimensionality increases, limiting their applicability to more complex optimization tasks. To address this limitation, this paper introduces X-DTC-GL, a novel DIRECT-type algorithm that incorporates dynamic partitioning and hybridization techniques. The dynamic partitioning approach adaptively refines the search space based on local one-dimensional surrogate models, enabling rapid subdivision of promising hyper-rectangles. The hybridization strategy selectively employs a hill-climbing method to exploit promising regions identified by the surrogate models. Extensive experiments on four diverse benchmark suites demonstrate that X-DTC-GL significantly outperforms existing DIRECT-type baselines, achieving improvements of ~12% in solvability and ~27% in solution quality. Performance-profile analyses indicate the fastest convergence on up to ~40% of instances, the best runtime performance on ~17% of problems, and competitive overall execution times. By improving performance within the partition-based framework, these advances strengthen the algorithm's competitiveness in state-of-the-art black-box optimization.

We prove uniform approximation theorems for Chebyshev expansions of continuous functions of two variables. Under a variation-type condition on the square $[-1,1]^2$, a continuous function admits a uniformly convergent Chebyshev expansion in the first variable whose coefficients are continuous functions of the second variable. These results are applied to the spectral graph wavelet transform of time-varying signals on finite weighted graphs: the scaling and wavelet kernels are approximated by expansions with time-varying coefficients whose degrees do not depend on time, the composition of the approximate transform with its adjoint admits the same explicit coefficient formulas as in the time-independent case, and reconstruction by the pseudoinverse is stable, with explicit bounds in terms of the uniform kernel errors. Numerical experiments on a sensor network confirm the convergence and stability estimates. In a denoising problem, soft thresholding of the graph wavelet coefficients with a time-varying transform parameter improves over its time-independent, fixed-parameter counterpart.