cs.NA — Pith

Top Pith

4

math.NA 2026-05-20 2 theorems

Discretization preserves area and volume in curvature flows at any order

by Ganghui Zhang, Boris D. Andrews +1 more

Arbitrary-order structure-preserving discretizations for geometric curvature flows

Auxiliary variables let the continuous area dissipation law hold exactly after arbitrary-order discretization.

abstract click to expand

Geometric flows, where an immersed manifold evolves in time according to its own geometry, exhibit important structural properties. For example, surface diffusion dissipates surface area while conserving volume; it is desirable to preserve these properties on discretization. This has motivated a substantial body of research on structure-preserving discretizations for these flows, albeit at low order in time. In this work, we present the first discretization of geometric curvature flows (curve shortening/mean curvature flow and curve/surface diffusion) that preserves the evolution of area and volume at arbitrary order in space and time. The key idea is to introduce auxiliary variables in a particular way so that the derivation of the area dissipation law can be replicated after discretization with continuous Petrov--Galerkin in time. These auxiliary variables are indicated by a general strategy for structure-preservation in time that applies to many other problems. The proposed scheme also preserves mesh quality in the same manner as the minimal deformation rate strategy. We demonstrate its structure-preserving properties and high-order convergence on several benchmark examples.

0

Top Pith

2

cs.LG 2026-05-15

Operator networks approximate nonlinear maps and their derivatives

by Filippo de Feo

Universal Approximation of Nonlinear Operators and Their Derivatives

Universal approximation theorems cover k-times differentiable operators in weighted Bastiani-Sobolev spaces on general Banach spaces.

abstract click to expand

Establishing Universal Approximation Theorems (UATs) for nonlinear operators and their derivatives is a foundational open problem in Operator Learning (OL) and raises delicate questions in Nonlinear Functional Analysis. We prove the first UATs for $k$-times differentiable nonlinear operators and their derivatives via OL architectures, uniformly on compact sets and in weighted Bastiani--Sobolev spaces for general finite input measures. In full Banach-space generality, these are the first complete generalizations of the corresponding influential classical UATs in [Hornik, 1991] to infinite-dimensional spaces and OL and they launch Derivative-Informed Operator Learning (DIOL)-learning nonlinear operators and their derivatives-on general Banach spaces. Based on our UATs, we formulate Bastiani--Sobolev training in DIOL. We present open frontiers where DIOL and our UATs find applications: high-order accuracy in OL; fast constrained optimization in Banach spaces (e.g. optimal control of PDEs, inverse problems) via Learn-Then-Optimize; numerical methods for infinite-dimensional PDEs (e.g. HJB PDEs on Banach spaces from infinite-dimensional optimal control via Optimize-Then-Learn, such as optimal control of PDEs, SPDEs, path-dependent systems, partially observed systems, mean-field control). We parameterize nonlinear operators via Encoder-Decoder Architectures, classical OL architectures. These include DeepONets, Deep-H-ONets, and PCA-Nets, which our UATs cover. Our UATs are based on (i) Approximation Properties of Banach spaces; (ii) continuous Bastiani differentiability (weaker than continuous Fr\'echet differentiability); (iii) $C^k_B$ (Bastiani) compact-open topologies; indeed, UA in $C^k$ (Fr\'echet) compact-open topologies (induced by operator norms) fails; (iv) construction of weighted Bastiani--Sobolev spaces, generalizing classical Gaussian Sobolev spaces on Banach spaces.

0

Top Pith

1

math.NA 2026-05-06 3 theorems

Reduced model speeds up crowd-contact simulations 30 to 60 times

by Giulia Sambataro, Virginie Ehrlacher

Model order reduction for parametrized variational inequalities: application to crowd motion

A greedy index-picker for contact forces plus a learned correction keeps agents from overlapping in dense crowds

abstract click to expand

This work investigates model order reduction for time-dependent parametrized variational inequalities, with a focus on discrete contact problems. As a prototypical example, we consider an agent-based crowd model [Maury et al., 2011] in which agent velocities are obtained at each time step from a constrained least-squares problem. Geometric parameter variations induce significant variability in both agent positions and contact forces, leading to a slowly decaying Kolmogorov $n$-width of the solution manifold. We propose a nonlinear approach that combines a linear reduced-order model with a deep-learning-based correction. The method utilizes a greedy index selection (gIS) algorithm for compressing Lagrange multipliers and Proper Orthogonal Decomposition (POD) applied to velocity snapshots. Additionally, we explore hyper-reduction techniques, comparing the Empirical Interpolation Method (EIM) and the Empirical Quadrature (EQ) procedure from both computational complexity and accuracy perspectives. Finally, we demonstrate the applicability of the methodology in a complex scenario involving many agents in a highly congested geometric configuration. This work represents the first attempt to apply model order reduction to a discrete contact problem of the type introduced in [Maury et al., 2011] and paves the way for future advancements in nonlinear MOR specifically for this class of problems.

0

math.NA 2026-07-03

L²-projection stays W¹,²-stable on hybrid meshes for K≥2

by Lars Diening, Viktoria Lingert +1 more

Sobolev stability of the L²-projection on hybrid meshes

The bound holds for all polynomial degrees starting at two on meshes mixing triangles with convex quadrilaterals from adaptive refinement.

abstract click to expand

We establish $L^p$- and $W^{1,p}$-stability of the $L^2$-projection onto mapped Lagrange finite elements on hybrid meshes consisting of triangles and convex quadrilaterals arising from adaptive mesh refinement. If $K$ is the (tensor product) degree of polynomials of the discretisation, then we show, in particular, $W^{1,2}$-stability for all $K\geq 2$ for the Q-RG and Q-RB refinements. This extends results by Ali, Funken, and Schmidt (2022) which hold for the range $2 \leq K \leq 9$ for initial meshes consisting of parallelograms. Our proof relies on an extension of the technique by Diening, Storn and Tscherpel (2021) to general convex quadrilaterals.

0

math.NA 2026-07-03

Nonlinear drift yields reverse generative PDE beyond score matching

by Horacio Tettamanti, Michael Herty

A PDE-Based Framework for Generative Modeling Beyond Classical Score-Based Diffusion

A superlinear term in the Fokker-Planck forward process creates condensation from which a stabilized reverse equation recovers the original

abstract click to expand

We introduce an alternative generative framework based on a nonlinear modification of the classical Ornstein--Uhlenbeck dynamics. The proposed dynamics admits both a microscopic description through an interacting particle system and, in the mean-field limit, a macroscopic formulation given by a nonlinear Fokker--Planck equation with a superlinear drift term. We show that, for suitable choices of the model parameters and sufficiently large initial mass, the forward dynamics exhibits condensation phenomena by proving the loss of $L^2$ regularity of the solution in finite time. Building upon this formulation, we derive a stabilized reverse-time partial differential equation that reconstructs the initial distribution from the asymptotic state of the forward dynamics, thereby extending the generative paradigm beyond the classical score-based framework. Furthermore, we introduce numerical discretizations of both the forward and reverse processes that accurately capture the asymptotic behavior of the continuous model while successfully reconstructing the initial distribution. Numerical experiments in one and two spatial dimensions validate the proposed methodology and illustrate its application to density filtering through successive iterations of the generative process.

0

math.NA 2026-07-03

Block preconditioner yields mesh-independent GMRES for SBM Stokes

by Micha{l} Wichrowski, Ajay Ajith

Block Preconditioning for Shifted Boundary Method Discretisations of the Stokes Problem

Field-of-values analysis shows non-symmetric SBM terms become small perturbations of a standard saddle-point operator on fine meshes.

abstract click to expand

The Shifted Boundary Method (SBM) sidesteps body-fitted meshing by shifting boundary conditions onto a surrogate boundary and correcting for the displacement through Taylor expansions. Despite its broad analysis and application, scalable iterative solvers for the incompressible Stokes equations remain underdeveloped. We present a block preconditioner for SBM--Stokes discretisations that uses the velocity block together with a pressure mass matrix as a Schur complement approximation. Because the SBM system is non-symmetric, classical operator preconditioning does not apply directly; a field-of-values analysis instead shows that the non-symmetric SBM contributions act as asymptotically small perturbations of a standard saddle-point operator, yielding mesh-independent GMRES convergence on sufficiently fine meshes. Numerical experiments demonstrate iteration counts under refinement across geometries of increasing complexity. We expose a coarse-mesh regime in which an under-resolved grid produces elevated iteration counts, an artefact of insufficient resolution that vanishes once the mesh captures the geometry.

0

math.NA 2026-07-03

Flexible CG runs exactly on cell-local vectors with continuous preconditioners

by Micha{l} Wichrowski

Coalesced Matrix-Free Finite Elements in Cell-Wise Storage

Primal-dual pairing keeps every Krylov scalar identical to the assembled solve and confines communication to the preconditioner step alone.

abstract click to expand

We present a GPU-oriented formulation of continuous high-order finite elements in which the redundant, cell-wise (element-local) vector is the persistent primary representation of all field data, rather than a transient stage of matrix-free operator evaluation. We prove that, given a preconditioner whose image is continuous, the entire flexible conjugate gradient iteration can be carried out exactly on this unassembled representation: a simple primal-dual pairing identity shows that all Krylov scalars computed from local data coincide with those of the assembled solve, so inter-element communication is confined entirely to the preconditioner. The required direct stiffness summation (DSS) is then realized without indirect gather-scatter, atomics, or coloring, by a dimensionally-split cascade of one-to-one face exchanges that provably accumulates edge and vertex contributions as a byproduct of sequential axis passes; unstructured macro-block interfaces and $h$-adaptive hanging nodes are handled by disjoint topological kernels and a shadow-cell wrapper that leaves the high-throughput sweeps untouched. The cell-wise storage decouples the memory layout from the mesh topology, and we exploit this freedom to benchmark blocked layouts that trade memory coalescing against element contiguity. Numerical experiments on modern GPUs demonstrate that the resulting operator evaluation and solver outperform state-of-the-art matrix-free implementations, signifficantly exceeding throughput of existing implementations.

0

math.NA 2026-07-03

Energy stationarity defines CutFEM for finite-strain elasticity

by Micha{l} Tomasz Wichrowski, Ella Godiva Noomen

A Unified CutFEM Formulation for Finite-Strain Elasticity: Energy Minimisation and Corner Singularities

Automatic differentiation supplies model independence while analysis proves cut-independent stability and shows singularity limits match tho

abstract click to expand

We present a fully variational, model-independent formulation of the Cut Finite Element Method (CutFEM) for finite-strain elasticity. The discrete problem is the stationarity condition of a augmented energy functional consisting of the bulk hyperelastic energy, the Nitsche terms that impose the boundary conditions weakly, and the ghost-penalty stabilisation. The residual and the (symmetrised) tangent follow from this functional by successive variations. Automatic differentiation (AD) generates the first Piola--Kirchhoff stress tensor and the elasticity tensor directly from the scalar energy density, avoiding manual re-derivation when exchanging hyperelastic models. To our knowledge, this is the first unfitted finite-strain scheme combining an energy-only, model-independent construction with AD and an accuracy analysis at unfitted boundaries. Analysis of the linearised problem solved at each Newton step establishes cut-independent coercivity, continuity, and an $O(h^{-2})$ condition number bound, yielding a quasi-optimal convergence theorem for regular solutions through the Brezzi--Rappaz--Raviart framework. Numerically, the method attains optimal $h$-convergence for linear, quadratic, and cubic elements on a smooth test case. Furthermore, we quantify the method's accuracy limit at mixed Dirichlet--Neumann junctions using the Kolosov--Muskhelishvili characteristic equation. The exact solution's corner singularity caps the convergence rate identically for fitted and unfitted methods. We demonstrate that local mesh refinement removes this bound, with the unfitted discretisation inheriting the recovered optimal rates and cut-independent constants.

0

math.NA 2026-07-03

Space-time method keeps acoustics stable for long industrial runs

by Simon Schneider, Ceyhun Özdemir +3 more

A Stable Boundary Element Method for Reliable Long-Time Industrial Sound Emission

The boundary element scheme stays accurate over extended times where standard methods diverge, matching real measurements.

abstract click to expand

In this paper we investigate a stable space-time formulation for long-time industrial sound emission problems. To this end, we use a well-posed Galerkin formulation in space and time of the acoustic wave equation in $\mathbb{R}^3$, involving a hypersingular boundary integral operator. Our numerical experiments confirm that the resulting time stepping scheme is stable and accurate for complex acoustic problems in industrial geometries, in contrast to alternative well-known schemes. The proposed method is shown to be efficient for real-world problems, and we obtain very good agreement with physical acoustic measurements.

0

math.NA 2026-07-03

Homothetic surfaces enable time-domain reconstruction of 3D obstacles

by Lu Zhao, Heping Dong +1 more

A novel time-domain iterative method for a three-dimensional inverse acoustic obstacle scattering problem

Retarded integrals on scaled surfaces converge to exact fields and support stable iterative recovery of shape and location.

abstract click to expand

This paper concerns the three-dimensional forward and inverse acoustic obstacle scattering problem in the time domain. For the forward problem, a retarded potential formulation discretized by convolution quadrature and Galerkin methods is introduced. By introducing the retarded boundary integral defined on a homothetic surface, we propose a novel time-domain convolution quadrature based iterative method to reconstruct both the shape and location of a rigid obstacle. The retarded integral in the time domain is reformulated into a system of integrals in the s-domain. The resulting s-domain integrals are very fast to compute, as they only involve non-singular integrals over the homothetic surfaces. Moreover, the Fr\'echet derivative with respect to the boundary can be derived straightforwardly. We also prove that the scattered field generated by the homothetic surface converges to the exact field in the time domain. To improve the stability of the inversion algorithm, an incremental truncation technique is proposed, and numerical experiments confirm the effectiveness and robustness of our method.

0

cs.LG 2026-07-03

DSGNAR hits 3e-16 error on PINN PDE problems

by Joseph Webb, Sadok Jerad +1 more

An Optimisation Framework for the Well-Conditioned Training of Physics-Informed Neural Networks

Doubly-sketched Gauss-Newton with adaptive ratio improves accuracy five to eight orders of magnitude over prior methods while staying faster

abstract click to expand

Physics-informed neural networks (PINNs) have emerged as a promising route to solve partial differential equations, yet they have struggled to reach the precision of classical solvers. The obstacle is increasingly understood to be one of optimisation, owing to the severely ill-conditioned loss landscape. We present $\textbf{DSGNAR}$: Doubly-Sketched Gauss-Newton with Adaptive Ratio, a scalable second-order optimisation framework that confronts this ill-conditioning and, in doing so, obtains unprecedented accuracy and speed. $\textbf{DSGNAR}$ couples a doubly-sketched Gauss-Newton model with a novel strategy that carefully controls both regularisation and step length. Across a suite of problems spanning nonlinear, chaotic, multi-scale, high-dimensional, and Navier-Stokes, the framework greatly improves on the state of the art: able to attain relative $\ell_2$ errors as low as $3\times10^{-16}$ in double precision, improve contemporary results by five orders of magnitude on the canonical Burgers' equation, and as much as eight orders on a high-dimensional Poisson problem, while remaining markedly faster. We further show that, in single precision, solutions at the limit of round-off error can be obtained very quickly: Burgers' equation to $\ell_2^{\text{rel}} = 4.75 \times 10^{-7}$ in under ten seconds. The framework is also robust to the choice of architecture, arithmetic precision, and initial hyperparameters. The code is available at https://www.github.com/wephy/physics-informed-neural-networks

0

math.NA 2026-07-03

Iteration reconstructs 3D obstacles from far-field data without singularities

by Zhiyong Cheng, Heping Dong +1 more

A highly efficient iterative approach for inverse acoustic obstacle scattering problems in three dimensions

Boundary integrals on a homothetic surface let the method avoid singular kernels while recovering the shape from phased or phaseless measure

abstract click to expand

This paper concerns a three-dimensional inverse acoustic obstacle scattering problem from scattered field or phased/phaseless far-field data. Based on the boundary integral defined on a homothetic surface, we propose a highly efficient iterative approach for obstacle reconstruction that completely avoids dealing with any singularity. Here, the injectivity and dense-range property of the Fr\'echet derivative have been proved to ensure the solvability of the linearized equivalent data equation. We also prove that the scattered field generated by the homothetic surface can arbitrarily approximate the exact one. Numerical experiments are presented to verify the superiority and robustness of the proposed approach.

0

math.NA 2026-07-03

Projection method gives unconditional optimal errors for harmonic maps

by Yongyong Cai, Xingwei Yang

Unconditional Optimal Error Estimates and Energy Stability for a Linearly Implicit Mass-Lumped Projection Finite Element Method for the Harmonic Map Flow

Converges at O(Δt + h²) in discrete L² and O(Δt + h) in H¹ on rectangular meshes with no time-step restriction

abstract click to expand

We propose and analyze a linearly implicit mass-lumped finite element method for the heat flow of harmonic maps into the unit sphere. The method consists of a linear predictor followed by a nodal projection and therefore preserves the unit-length constraint exactly at all finite element nodes. The predictor is derived from a cross-product reformulation of the equation and is shown to be equivalent to a mass-lumped discretization of the original formulation with a correction term enforcing nodal orthogonality, as well as to a tangent plane scheme. A key ingredient is the consistent use of the discrete inner product in both the mass and stiffness terms. This yields a nodal orthogonality relation implying that the auxiliary solution lies on or outside the unit sphere at every node. Consequently, the projection is well defined and the projected error satisfies a contraction property in the discrete $L^2$-norm. On Cartesian rectangular and cuboidal tensor-product meshes, the nodal projection is also nonexpansive in a discrete Dirichlet energy, which gives an unconditional discrete energy dissipation law. For sufficiently smooth solutions, we prove optimal error estimates without any coupling condition between the time step and the mesh size: the method converges with order $O(\Delta t+h^2)$ in $\ell^\infty(0,T;L^2)$ and order $O(\Delta t+h)$ in $\ell^2(0,T;H^1)$. The proof combines the projected-error contraction, quadrature consistency estimates, edge-based cancellation identities, and a bootstrap argument for controlling nonlinear terms. Numerical experiments confirm the predicted convergence rates and the discrete energy decay.

0

math.NA 2026-07-03

Derivative moments on intervals fix polynomials up to constant

by Allal Guessab, Federico Nudo

Weighted Derivative Histopolation on Arbitrary Grids: Admissibility and Exact Factorizations

Separation of lower-degree polynomials via weighted moments plus nonzero normalization on constants yields a well-posed scheme on degree-N p

abstract click to expand

In this paper, we introduce a weighted derivative histopolation framework on families of intervals. The degrees of freedom consist of one scalar normalization and weighted integral moments of the derivative over a prescribed family of subintervals. We prove that the resulting scheme is unisolvent on $\Pi_N$ when the interval family separates polynomials of degree at most $N-1$ through weighted moments and the normalization is nonzero on constants. Thus, the derivative moments determine the polynomial up to an additive constant, and the scalar normalization fixes this remaining degree of freedom. This gives a sharp criterion for the well-posedness of the interpolation problem and a complete characterization of the admissible scalar normalizations. We then show how admissible families of intervals can be constructed from a fixed grid. When the endpoints of the intervals belong to the grid, admissibility is reduced to the nonsingularity of an interval matrix associated with the family, which depends only on the representation of the intervals in terms of consecutive cells. For Jacobi weights, the associated data matrices have a natural block structure in Jacobi polynomial bases, and the reduced derivative matrix can be expressed in terms of shifted Jacobi moment matrices. We next study Chebyshev configurations in which this structure becomes explicit. For the four classical Chebyshev families, suitable polynomial bases lead to diagonal Gram matrices for the reduced derivative matrices. We show that this diagonal structure depends on the simultaneous choice of the weight, the basis, and the grid. Numerical experiments on equispaced and Chebyshev--Lobatto nodes show the behaviour of the method for different interval families and for different Jacobi parameters.

0

math.NA 2026-07-03

Nonlinear Laplacian solves network flows after 2-4 linear solves

by Oren E. Livne

NLF: A Resistor-Network Framework and Linear-Time Solver for Convex Network-Flow Equilibria

Damped chord-Newton iteration on Bρ(B^T φ)=α d finishes every SuiteSparse graph up to 18 million edges and beats interior-point and L-BFGS b

abstract click to expand

We present NLF (Nonlinear Laplacian Flow), a unified framework and linear-time solver for convex network-flow equilibria. Congestion routing, minimum-delay routing, and maximum flow share one form: the nonlinear graph Laplacian $B\rho(B^T\phi)=\alpha d$, where a monotone edge law $\rho_e$ encodes the physics (undirected graphs; directed variants are future work). NLF solves it by a damped chord-Newton iteration whose frozen linearization -- a weighted graph Laplacian -- is inverted by a near-linear Laplacian solver (default: approximate Cholesky, LAMG+ interchangeable). The nonlinear solve costs $2$--$4$ linear Laplacian solves, making the wall-clock empirically $O(m)$ in the edge count $m$ (not a proved bound). On single-commodity congestion (BPR cost), NLF converges on all 2,003 SuiteSparse corpus graphs up to $1.8\times10^7$ edges. Against a state-of-the-art interior-point method, NLF is a median $2.6\times$ faster where both converge and $>45\times$ on poorly-separable graphs where the IPM's direct core is superlinear; against L-BFGS, a median $4.2\times$ faster and the only solver to finish on the 90 hardest instances. A multicommodity extension routes $K$ commodities through one shared hierarchy at $O(Km)$ per step. The same machinery recovers the exact max-flow as a short sequence of Laplacian solves, with the cut potential as a by-product. Code: https://github.com/orenlivne/nlf

0

math.NA 2026-07-03

Biot-contact scheme shown to have unique solutions and optimal errors

by Marius Nevland, Kundan Kumar +3 more

Numerical analysis of the Biot equations coupled to frictional contact mechanics

Conforming elements plus backward Euler yield stability and convergence rates for the coupled poro-visco-elastic frictional model.

abstract click to expand

We consider a mathematical model of a poro-visco-elastic medium subject to frictional contact with a rigid obstacle, and study its numerical approximation. This model couples the Biot equations and contact conditions in the form of normal compliance and Coulomb friction. The resulting variational problem consists of a linear partial differential equation coupled to a nonlinear variational inequality. We propose and analyze a fully discrete numerical scheme for this problem, using conformal finite elements in space and the implicit Euler method in time. Existence and uniqueness of the discrete solution is established, and stability and a priori error estimates are derived. A numerical experiment is performed in which numerical error estimates are computed and compared to the theoretical results.

0

math.NA 2026-07-03

Alfeld splits yield low-regularity elasticity finite element complexes

by Johnny Guzmán, Xuehai Huang

Low-regularity finite element elasticity complexes with hybridizable stresses on tetrahedral Alfeld splits

The complexes reduce Sobolev regularity and polynomial degree while supplying hybridizable H(div) symmetric stresses without vertex degrees

abstract click to expand

Finite element elasticity complexes of low regularity are constructed on tetrahedral Alfeld splits. In comparison with existing three-dimensional elasticity complexes on such splits, the complexes constructed here lower both the Sobolev regularity and the polynomial degrees, while ending in a hybridizable $H({\rm div};\mathbb S)$-conforming symmetric stress space with no vertex degrees of freedom. The construction is obtained from local Bernstein-Gelfand-Gelfand arguments applied to polynomial de Rham complexes on the Alfeld split. Two local polynomial elasticity complexes are proved: an $H^2$-$H^1({\rm inc})$ complex and a lower-regularity $H^1({\rm curl})$-$H({\rm inc}^+)$ complex. Their bubble subcomplexes and dimension formulas are derived. These local exact sequences lead to unisolvent finite elements for the displacement and incompatibility spaces and to global finite element subcomplexes of the corresponding elasticity sequences. In the lowest-order $H^1({\rm curl})$-$H({\rm inc}^+)$ finite element complex, the $H({\rm inc}^+;\mathbb S)$-conforming tensor space is piecewise cubic. At the same order, the terminal stress-displacement pair recovers the Johnson-Mercier-K\v{r}\'{i}\v{z}ek element, while the construction covers higher-order hybridizable symmetric stresses for all $k\ge1$. A second family gives a low-regularity $H^1$-$H({\rm inc})$ finite element complex for the standard elasticity sequence for all $k\ge2$. Commuting interpolation diagrams are established for both global complexes.

0

math.NA 2026-07-03

Nonlocal model captures dynamic crack branching below 0.6 wave speed

by Qing Cheng, Yuqi Sun +1 more

A Variational Nonlocal Phase-Field Model for Dynamic Fracture in Elastic Solids

Implicit crack length from interaction domain reproduces 48-degree shear paths and branching without explicit tracking.

abstract click to expand

We develop a variational nonlocal phase-field model for dynamic fracture in elastic solids. The proposed formulation is distinguished by three main features. First, the model is formulated through nonlocal kinematics and kernel-dependent function spaces, allowing weaker regularity requirements while recovering the classical local theory as the nonlocal interaction domain vanishes. Second, a nonlocal crack-surface functional is introduced as an integral counterpart of the Ambrosio--Tortorelli regularization, so that the characteristic length of the diffusive crack is implicitly determined by the nonlocal interaction domain rather than by a prescribed length scale. Third, the degraded nonlocal elastic energy and the nonlocal crack-surface functional are combined into a variationally consistent dynamic fracture system, consisting of a nonlocal momentum balance and an irreversible nonlocal gradient-flow evolution law for the phase field. The coupled system is solved using two temporal discretization strategies: a structure-preserving scalar auxiliary-variable scheme and a staggered alternating scheme, both combined with finite element discretization in space. Numerical examples involving Mode-I fracture, dynamic crack branching, Kalthoff--Winkler-type shear fracture, and fragmentation show that the proposed model captures complex crack initiation, propagation, branching, and interaction without explicit crack tracking. Quantitatively, the predicted crack-tip velocities remain below $0.6c_R$ in the dynamic branching and shear-loading tests, and the shear-loading benchmark gives an inclined crack path of approximately $48^\circ$, consistent with the characteristic Kalthoff--Winkler fracture pattern.

0

cs.DC 2026-07-03

Task-based Gauss-Seidel matches divide-and-conquer speed on 2D grids

by Shiting Long, Gustavo Ramirez-Hidalgo +2 more

Exploiting Task-Based Parallelism for the Red-Black Gauss-Seidel Method on 2D Grids

Task versions run at comparable speed yet recover better when processor timing varies and extend to irregular grids.

abstract click to expand

Gauss-Seidel is a well-established iterative method for the solution of linear systems, and multicoloring has been widely used to increase parallelism in iterative solution techniques. Implementing multi-color Gauss-Seidel with conventional divide-and-conquer parallelization strategies, however, may be inefficient due to global synchronization requirements and load imbalances. Task-based programming models can mitigate these issues by enabling fine-grained parallelism, removing global barriers and allowing updates of different colors to partially overlap in time. In this work, we implement the red-black Gauss-Seidel method using two task-based programming models and compare them with a classical divide-and-conquer parallel implementation to evaluate the impact of fine-grained parallelism on execution efficiency. The red-black scheme serves as a representative example, as task-based approaches naturally extend to more general multi-color schemes arising from unstructured grids and wider stencils. Using the solve of the 2D Poisson equation as benchmark, our results show that task-based implementations can achieve performance comparable to conventional divide-and-conquer parallelization while providing greater resilience to hardware-level asynchronicity.

0

math.NA 2026-07-03

Transport signatures compress density PDE solutions to lower rank

by Jiajia Yu, Jingwei Hu +4 more

Structure-Preserving Reduced-Order Modeling via Low-Rank Transport Signatures

Kantorovich potentials from one reference density map to signatures whose low-rank form preserves mass and bounds Wasserstein error.

abstract click to expand

Parametrized PDEs with density-valued solutions are often difficult to approximate with classical linear reduced-order models, especially in transport-dominated regimes. We introduce an optimal-transport-based reduced-order modeling that represents each density by the Kantorovich potential transporting a fixed reference density to the target density, and then maps these potentials to transport signatures using a weighted Laplacian associated with the reference measure. This embeds the density-valued solution map in a Hilbert space while preserving control of the induced transport maps and Wasserstein error. We treat the signature map as a continuous matrix indexed by parameters and space, construct a low-rank skeleton decomposition using a maximal-volume criterion, and learn the parameter-to-coefficient map with a neural network for efficient non-intrusive online evaluation. The reconstructed solution is obtained by pushing forward the reference density, so mass preservation is built into the method. We prove a mean-squared Wasserstein error bound separating low-rank approximation, discretization, sampling, and learning errors, and demonstrate the method on a two-dimensional continuity equation, where transport signatures yield substantially lower-rank structure than the original density snapshots.

0

physics.flu-dyn 2026-07-03

RIS matches IB aortic valve pressure gradients within 15%

by Han Zhao, Alexander D. Kaiser +6 more

Comparative analysis of resistive immersed surface and immersed boundary methods for aortic valve simulation

Prescribed-motion method recovers bulk flow features of fully coupled simulation at 60% lower computational cost.

abstract click to expand

Numerical modeling of aortic valve dynamics is essential for understanding the complex fluid-structure interaction (FSI) governing valve biomechanics in health and disease. Immersed methods provide a flexible computational framework for simulating the large deformations of valve leaflets and associated blood flow without requiring body-fitted meshes. Among these approaches, the Resistive Immersed Surface (RIS) and Immersed Boundary (IB) methods are widely used. However, systematic comparative analysis of these methods for realistic aortic valve simulations has not been performed. In this work, we compare a prescribed-kinematics RIS workflow implemented in SimVascular's svMultiPhysics solver with a fully coupled IB workflow using IBAMR for trileaflet and bicuspid aortic valve configurations. The RIS method represents the valve as a surface with prescribed kinematics embedded in the fluid domain and introduces a penalty force that drives the surrounding fluid velocity toward the prescribed leaflet velocity. This formulation reduces modeling complexity and provides useful hemodynamic predictions when representative leaflet kinematics are available. In contrast, the IB method models the leaflets as elastic structures fully immersed in the fluid domain and resolves leaflet deformation through fully coupled two-way FSI. The study focuses on the extent to which RIS reproduces bulk hemodynamic features and transvalvular pressure gradients. Results show that the RIS method captures the large-scale flow structures and predicts the mean transvalvular pressure gradient with a relative error within 15% of the fully coupled IB simulation, improving to within 5% when inlet boundary conditions are matched, while reducing computational cost by approximately 60%.

0

cs.CV 2026-07-03

Plug-and-play recovers cellular structures from compressed light-sheet data

by Jianqing Jia, Yi Gong +4 more

Plug-and-Play Volumetric Reconstruction for Compressive Sensing Light-Sheet Microscopy

Axial coupling and Woodbury updates let any denoiser restore 3D continuity from multiplexed single-exposure measurements.

abstract click to expand

We investigate volumetric reconstruction for compressive sensing light-sheet microscopy (CS-LSM), where fast volumetric imaging is achieved by encoding multiple axial planes into each camera exposure. To recover the underlying volume from highly multiplexed measurements, we propose a plug-and-play (PnP) framework that flexibly incorporates any user-specified denoiser into the reconstruction process. Building on a slice-based formulation, we further introduce an axial-coupled model that exploits correlations between adjacent slices to improve volumetric continuity. For efficient computation, we derive a Woodbury-based update for the data-consistency step in both the slice-based and axial-coupled formulations, and employ a Gauss-Seidel sweep for the denoising step in the axial-coupled model. Under a weakly convex regularization assumption, we establish subsequential convergence of the proposed algorithm. Experiments on synthetic and real zebrafish-heart data demonstrate that the proposed framework successfully recovers cellular structures from compressed measurements, and provide practical insights into the comparative performance of commonly used denoisers within the PnP framework under the CS-LSM setup.

0

math.NA 2026-07-03

Fourier symbols make natural gradients identical to mirror descent

by Jeongbin Jo

Fourier-Diagonalized Natural Gradients and Sobolev Mirror Descent

When the spectral symbols of the Fisher metric and Sobolev norm coincide, the optimization updates match exactly; otherwise one precondition

abstract click to expand

We study natural-gradient updates whose metric operators are diagonalized by the Fourier transform and relate them to Sobolev mirror descent. Translation-invariant Fisher geometries and Sobolev mirror geometries share a common inverse-map structure in the spectral domain. The Fisher metric is represented by a positive Fourier symbol, while Sobolev mirror geometry corresponds to the specific Bessel-potential symbol associated with the Sobolev norm. When these symbols coincide, the natural-gradient and mirror-descent updates are identical; otherwise, Sobolev mirror descent provides a canonical spectral preconditioner for the Fisher inverse geometry. This gives a mathematical lens through which spectral filtering and truncation techniques in PDE and operator learning can be viewed as natural actions of inverse metric geometry. We introduce Spectral Natural Gradient, an FFT-based implementation of these geometric updates.

0

math.NA 2026-07-03

Stream function VEM solves Stokes-Darcy on polygonal meshes

by Franco Dassi, Rekha Khot +3 more

Stream function -- pressure virtual element methods for the Stokes--Darcy interface problem

Enforces mass conservation, normal stress balance and slip condition at the interface without remeshing.

abstract click to expand

This paper introduces a novel Virtual Element Method (VEM) for the coupled Stokes--Darcy system in primal-primal form. In the free-flow Stokes domain, we implement a stream function formulation that inherently satisfies the incompressibility constraint and reduces computational cost. Across the interface, mass conservation, normal stress balance, and the Beavers--Joseph--Saffman slip condition are enforced to couple the biharmonic stream function equation with the Darcy's pressure equation. Leveraging VEM's ability to handle general polygonal meshes, the proposed method naturally accommodates irregular interface geometries without requiring remeshing or adaptive refinement. The accuracy of the method is validated through several numerical simulations that include applications to dead-end filtration, and network flow in bioartificial organs.

0

math.NA 2026-07-03

Gegenbauer post-processing removes ROM oscillations for transport problems

by Lei Yan, Yan Jiang +1 more

Post-Processing Reduced-Order Models for Transport-Dominated Problems by Gegenbauer Reconstruction

Re-projection onto local Gegenbauer bases cuts errors by up to two orders and sharpens discontinuities better than total variation regulariz

abstract click to expand

In this paper, we develop a physics-based post-processing technique for data-driven reduced-order models (ROMs) of transport-dominated problems. Besides the slow decay of the Kolmogorov n-width, ROMs based on globally supported bases often produce unphysical oscillations when approximating solutions with shocks or sharp gradients, a phenomenon analogous to Gibbs oscillations in spectral approximations. To address this issue, we introduce a post-processing framework based on Gegenbauer polynomial reconstruction. The key idea is to re-project the ROM solution onto a Gegenbauer polynomial basis over each interval of analyticity. Originally developed for spectral approximations, Gegenbauer reconstruction achieves spectral accuracy while effectively suppressing Gibbs oscillations. We extend this technique to data-driven ROMs and consider three representative approaches: Proper Orthogonal Decomposition (POD)-Galerkin ROM, Operator Inference (OpInf), and nonlinear manifold ROMs based on convolutional autoencoders (CAE). Numerical results show that the proposed post-processing consistently removes spurious oscillations and substantially improves solution quality for all three ROMs. For one-dimensional problems, the method is straightforward to implement once discontinuities are detected. We further develop a practical extension to two-dimensional problems using line-by-line reconstruction in each coordinate direction. Extensive numerical experiments demonstrate that the proposed method reduces errors by up to one or two orders of magnitude for inviscid transport problems and significantly outperforms total variation regularization in both numerical accuracy and the sharp resolution of discontinuities.

0

math.NA 2026-07-02

SS-POD reaches 3.9e-8 error with 91 functions from 5 snapshots

by Tianhao Hu, Zecheng Gan

A spectral-subspace-augmented POD-Galerkin method for parametrized PDEs with limited snapshot data

Energy-balanced local POD on spectral subspaces beats standard POD saturation and full spectral expense.

abstract click to expand

Parametrized partial differential equations (PDEs) arise in many-query simulation, optimization, control, and uncertainty quantification, where repeated full-order solves restrict the number of high-fidelity snapshots that can be generated. This limitation is particularly pronounced in computational energy science, where multiscale models of porous-media flow, transport, and energy materials often make large snapshot datasets impractical. Proper orthogonal decomposition (POD) constructs compact reduced bases from solution snapshots, but it may exhibit limited out-of-sample predictive capability when the snapshots insufficiently sample the solution manifold. To address this limitation, we propose a spectral-subspace-augmented POD-Galerkin method (SS-POD) tailored to limited-data regimes. SS-POD starts from a problem-adapted spectral approximation space, partitions it into orthogonal subspaces, and performs POD locally on the projected snapshot matrices. An energy-balancing rule determines the spectral partition so that the resulting local POD problems are assigned comparable amounts of snapshot energy. For nonlinear parametrized PDEs, SS-POD is coupled with the discrete empirical interpolation method (DEIM). Numerical experiments show that SS-POD improves out-of-sample accuracy over standard POD-Galerkin while retaining compact reduced bases in limited-snapshot regimes. In particular, for a Laplace-Beltrami problem on the unit sphere with only 5 snapshots, SS-POD achieves a relative error of $3.9*10^{-8}$ using 91 basis functions, whereas the standard POD error saturates at $7.8*10^{-4}$ and the spectral-Galerkin method requires 256 basis functions for comparable accuracy. These results indicate that SS-POD provides an effective strategy for high-fidelity reduced-order modeling from limited snapshot data.

0

cs.LG 2026-07-02

KANs with analytic R-functions cut test RMSE by up to 67%

by Sergei Kucherenko, Nilay Shah

Geometry-Aware R-Structured Kolmogorov-Arnold Networks

The models encode known geometric boundaries explicitly while KAN branches learn the smooth nonlinear parts, sharpening both accuracy and bo

abstract click to expand

We propose a novel hybrid neural architecture, the Geometry-aware R-Structured Kolmogorov-Arnold Network (GRS-KAN), which integrates V.L.Rvachev's R-functions into the Kolmogorov-Arnold Network (KAN) framework. The proposed approach combines two complementary modeling mechanisms: smooth nonlinear structure is learned by KAN branches, while known geometric or logical constraints are encoded analytically using differentiable R-functions. This enables explicit representation of discontinuities, feasible regions, and implicit geometric boundaries within a trainable neural architecture. The framework implements differentiable logical operations through R-conjunctions and R-disjunctions, allowing complex geometric supports to be represented analytically and incorporated directly into regression models. Several GRS-KAN variants are introduced, including additive, multiplicative, and agnostic branch-weighted architectures. The method is demonstrated on regression problems involving discontinuities with circular and rectangular supports. Numerical experiments show that explicit geometric encoding substantially improves predictive accuracy and boundary localization compared with standard KANs. In the considered benchmarks, geometry-aware GRS-KAN models reduce test RMSE by up to 67% while simultaneously improving interpretability through explicit analytical representation of the learned geometric structure. The agnostic variant further demonstrates the ability to automatically determine whether geometric priors are beneficial for a given learning task.

0

math.NA 2026-07-02

Three equivalent ways fix the unsteady Kutta wake amplitude

by Jiguang Yu, Louis Shuo Wang

Fredholm--residue selection of the unsteady Kutta amplitude

Singularity cancellation, Fredholm solvability and pole residue agree under structural hypotheses on the edge operator.

abstract click to expand

We give an operator-theoretic interpretation of unsteady Kutta selection in trailing-edge acoustic receptivity. The inviscid acoustic--wake problem leaves one outgoing wake amplitude undetermined. We show that, under explicit structural hypotheses, this amplitude is the same scalar obtained from three representations: cancellation of the inverse-square-root edge singularity, Fredholm compatibility of the viscous lower-deck problem, and the residue of the Kutta-normalized transform solution at the downstream wake pole: $\displaystyle A = -\frac{C_-^{(0)}}{C_-^{(KH)}} = -\frac{\langle \mathbf F_{\rm inc},\Psi^\ast\rangle}{\langle \mathbf F_{KH},\Psi^\ast\rangle} = i\operatorname*{Res}_{\alpha=\alpha_{KH}}\mathcal M(\alpha)$. The inner Fredholm--edge mechanism is verified exactly in a linear-shear lower-deck model, where the primal shear and adjoint velocity are Airy fields and the edge concomitant is nonzero outside a discrete resonance set.

0

math.NA 2026-07-02

Projection forces low-rank VM solution to conserve mass

by Shadi Heenatigala, Wei Guo

A Local Macroscopic Conservative (LoMaC) low rank tensor method for the Vlasov-Maxwell system

Hierarchical Tucker approximation is adjusted after each step to match moments from separately evolved conservation laws

abstract click to expand

The main computational challenges of solving the Vlasov-Maxwell (VM) system include the high dimensionality of the phase space, nonlinearity, inherent conservation properties, among others. In this paper, we develop a novel Local Macroscopic Conservative (LoMaC) low rank tensor method for the VM system, as a continuation of our previous work (arXiv:2207.00518). The method takes advantage of the tensor friendly structure of the Vlasov equation and employs the low rank hierarchical Tucker decomposition to approximate the Vlasov solution in high dimensions. Hence, the curse of dimensionality can be mitigated. Furthermore, to realize the LoMaC property, the algorithm simultaneously evolves the conservation laws of mass, momentum and energy alongside the Vlasov equation using a high order conservative method with the kinetic flux vector splitting. By a conservative orthogonal projection, the low rank solution is guaranteed to have the same macroscopic observables updated from the conservation laws. A collection of numerical tests on the VM system are presented to demonstrate the efficiency and efficacy of the proposed algorithm.

0

math.NA 2026-07-02

Macro-micro split reduces Monte Carlo variance in transport

by Caleb A. Shaw, Dmitriy Y. Anistratov

Hybrid Two-Level Transport Method with Solution Decomposition in Macro and Micro Components

P1 macro moments with exact closures plus MC micro component solved by fixed-point iteration for the steady-state Boltzmann equation.

abstract click to expand

This paper presents a new hybrid MC/deterministic method for solving the one-group steady-state Boltzmann transport equation based on decomposition of solution in macro and micro components. The macro component captures the large-scale structure of the solution. It is represented by angular moments of the high-order transport solution. The $P_1$ approximation is applied to define the macro component. The first two angular moments are obtained as a solution of hybrid low-order moment equations with exact closures. The equation for the micro component is solved using a MC simulation. The hybrid two-level system of equations for macro and micro components is solved by fixed-point iteration scheme. Numerical results are presented to demonstrate variance reduction of stochastic numerical solution and improvement in computational efficiency.

0

cs.LG 2026-07-02

One model solves forward and inverse PDEs on any geometry

by Meenakshi Krishnan, Pranav Pulijala +3 more

GAIA: Geometry-Adaptive Operator Learning for Forward and Inverse Problems

Geometry tokens condition cross-attention layers so the same architecture handles BVPs and inverses without retraining or graphs.

abstract click to expand

Operator learning for partial differential equations (PDEs) on arbitrary geometries builds fast neural surrogates for large-scale simulation. Although recent geometry-adaptive neural operators have made substantial progress, they are mainly designed for forward problems in which inputs and outputs share the same spatial domain. This limits their applicability for boundary value problems (BVPs) and inverse problems, where inputs and outputs may live on different domains. We introduce the Geometry-Adaptive Integral Autoencoder (GAIA), an operator learning model that encodes the domain boundary and the interior field distribution into geometry tokens, and conditions integral transform layers on these tokens via cross-attention, allowing the kernel to adapt locally to geometric features. This yields a single architecture for forward (including BVPs) and inverse problems on arbitrary domains in one pass, without retraining, iterative optimization, or graph construction. We evaluate GAIA on seven 2D and 3D benchmarks, four of which are new or substantially extended benchmarks for inverse problems and BVP: electrical impedance tomography, optical tomography, 3D Darcy flow on varying geometries, and a modified setting of Poisson BVP on mechanical components benchmark (MCB). GAIA sets new state-of-the-art results on every inverse and BVP task, reducing median relative $L^2$ error by 64% on airfoil flow reconstruction and 27% on EIT relative to the next best amortized method, and outperforming all baselines on every shape category of MCB. On other forward problems, GAIA is competitive with specialized solvers while maintaining stable accuracy across point resolutions on which transformer-based baselines degrade.

0

physics.comp-ph 2026-07-02

Lanczos method reduces cost for nuclear QRPA strength functions

by Dong Min Roh, Chao Yang +2 more

Lanczos Method for QRPA Strength Functions in Atomic Nuclei

Single Krylov run matches GMRES accuracy for two nuclei over broad energy range

abstract click to expand

We present a symmetric Lanczos method for computing charge-changing QRPA strength functions in atomic nuclei. Starting from the finite-amplitude-method formulation of the QRPA linear-response problem, we derive equivalent spectral representations and, in the real case, a reduced eigenvalue problem involving the matrix products $MK$ and $KM$, where $M\equiv A+B$ and $K\equiv A-B$ are formed from the usual QRPA matrices $A$ and $B$. The resulting formulation enables a matrix-free Lanczos approximation of the Lorentzian-smeared strength function over a broad energy interval from a single Krylov run, in contrast to conventional frequency-by-frequency response calculations. Numerical tests for $^{112}$Sn and $^{150}$Nd first show that GMRES reproduces the converged iterative FAM strength profiles while requiring fewer iterations. Using GMRES as the frequency-by-frequency reference, we then show that the Lanczos approximation reproduces the same strength profiles with reduced overall cost. These results indicate that symmetric Lanczos projection provides an efficient and accurate approach for QRPA strength-function calculations when spectral information is required over an extended frequency range.

0

math.NA 2026-07-02

Monotonic ICNNs approximate polyconvex envelopes with lower cost

by Timo Neumeier, Julian Salmon

Compression of Polyconvex Envelopes of Isotropic Functions via Monotonic Input Convex Neural Networks

A sufficient criterion and positive-octant reduction cut computation while matching the Saint Venant-Kirchhoff envelope closely.

abstract click to expand

This work presents a novel neural-network compression approach for polyconvex envelopes of isotropic functions. The approach relies on a classical sufficient criterion for polyconvexity and is particularly suited for the representation of determinant-constrained energy densities arising in non-linear elasticity. Compared with existing compression methods based on the necessary and sufficient characterisation of polyconvex isotropic functions, the proposed framework reduces computational costs, due to the domain reduction through the restriction to the positive octant in the singed singular value space. The underlying neural-network architecture employs input-convex neural networks (ICNNs) with non-negative weight constraints to enforce the required convexity and monotonicity properties. The additional symmetry and inequality conditions characterising the polyconvex envelope are incorporated weakly through the loss function during training. Although the employed criterion is only sufficient and thus generally yields only a lower bound on the polyconvex envelope, numerical experiments based on the classical Saint Venant--Kirchhoff energy demonstrate that the proposed approach produces accurate approximations in practice while offering a computationally more efficient alternative to existing methods.

0

physics.flu-dyn 2026-07-02

PICNN replaces traditional solver in two-phase porous flow

by Yuanshuo Kong, Xue Wang +1 more

The PICNN-Assisted Physics-Preserving Scheme for Thermodynamically Consistent Two-Phase Flow in Porous Media

Network trained on finite-volume residuals with TPFA yields energy-stable, mass-conservative solutions after correction.

abstract click to expand

In this paper, we develop a physics-informed convolutional neural network (PICNN) assisted physics-preserving method for a thermodynamically consistent model of incompressible and immiscible two-phase flow in porous media. Following the physics-preserving prediction-correction scheme of Li et al. \cite{li2025class}, the prediction step is performed by a PICNN trained with finite-volume residuals, where the interfacial fluxes are evaluated by the two-point flux approximation (TPFA) using two-point difference quotients of neighboring cell-centered unknowns to approximate interfacial normal gradients. The PICNN output is further corrected by a post-processing procedure to obtain energy-stable, mass-conservative, and bounds-preserving solutions. Numerical results show that the finite-volume residuals trained PICNN can replace the traditional prediction solver within the physics-preserving framework. Compared with conventional physics-informed neural networks (PINNs), the PICNN better captures local spatial interactions between each control volume and its neighboring cells, while the finite-volume residuals accommodate discontinuous permeability fields and interfacial flux continuity.

0

math.NA 2026-07-02

HSS approximation enables direct 2D inverse NUDFT solve

by Yingzhou Li, Jingyu Liu

A Superfast Direct Solver for 2D Type-II Inverse Nonuniform Discrete Fourier Transform Based on Hierarchically Semiseparable Matrix

Factoring A as G F and approximating the kernel G yields O(M + N^{3/2} log^3 N) offline and O(M + N log^3 N) online costs plus preconditioni

abstract click to expand

This paper proposes a direct inversion method for the 2D type-II nonuniform discrete Fourier transform~(NUDFT). The NUDFT matrix $A$ is factored as $A = G F$, where $G$ can be expressed as a kernel matrix and $F$ is the 2D DFT matrix. We show that $G$ can be approximated by a hierarchically semiseparable~(HSS) matrix and give an estimate of the HSS rank. Then, using the least-squares solver for HSS matrix and the two-dimensional inverse fast Fourier transform, the inverse NUDFT problem can be solved efficiently. Our algorithm has an offline complexity of $O\bigl(M+ N^{3 / 2} \log^{3} N\bigr)$ where $M$ and $N$ are the size of rows and columns of the NUDFT matrix, respectively. Once the direct solver is built, it can be applied to a vector with an online complexity of $O\bigl(M+ N \log^{3} N\bigr)$. The proposed method can be used as a preconditioner for iterative methods, especially when the sample points are distributed on a grid such that $A$ is ill-conditioned. Numerical results are provided to show the scaling performance of the inversion method and demonstrate the efficiency and robustness of it as a preconditioner.

0

math.NA 2026-07-02

CRE method gives guaranteed bounds for PGD plate models

by Jean Ruel, Ludovic Chamoin +2 more

Certification of PGD reduced-order models with separated spatial variables

Dedicated flux construction supports adaptive control of discretization error and PGD mode count.

abstract click to expand

Model order reduction techniques have become an attractive approach for obtaining fast approximations of multidimensional problems. Besides computational efficiency, ensuring the reliability of the resulting approximations is of primary importance. This work focuses on the certification of PGD-based reduced-order models based on the separation of spatial variables, which are particularly well suited to plate and shell geometries. Considering diffusion problems defined in plate-like domains, we introduce a guaranteed global error estimate associated with the PGD approximation. To this end, the error bounds are derived from the Constitutive Relation Error (CRE) method. The main difficulty of this approach lies in the construction of equilibrated fluxes, for which a dedicated procedure is proposed. Based on the resulting estimator, an adaptive strategy is developed to control both the discretization error and the number of PGD modes. This certification procedure is further extended to the error control in quantities of interest. We provide several numerical examples illustrating the reliability and efficiency of our procedure.

0

math.NA 2026-07-02

Small overlaps cut iterations in meshless domain decomposition

by Alexander Westermann, Oleg Davydov +1 more

Overlapping Domain Decomposition for Meshless Finite Difference Methods

Numerical tests on Poisson and Stokes show RBF-FD Schwarz methods need no disjoint subdomains for good performance.

abstract click to expand

Schwarz type domain decomposition methods generally require a partition of unity to combine solutions on subdomains. However, in mesh-based methods it is common to organize subdomains with minimal overlap, if any, which is facilitated by the availability of a mesh. This study analyzes how the continuity of the partition of unity affects the algebraic Schwarz method for Poisson and Stokes equations from a meshless point of view, whereby the underlying differential operators are discretized using the radial basis function finite difference (RBF-FD) method. We demonstrate numerically that, in this setting, small overlaps improve the performance of the domain decomposition, leading to smaller iteration counts, and therefore no disjoint partitioning technique is required.

0

math.NA 2026-07-02

Staggered DG methods solve any-order polyharmonic equations on polytopes

by Long Chen, Xuehai Huang +2 more

Hybridizable Staggered Discontinuous Galerkin Methods for Polyharmonic Equations on Polytopes

Hybridizable formulation on shape-regular meshes delivers optimal energy errors without dimension-order restrictions.

abstract click to expand

Hybridizable staggered discontinuous Galerkin methods are developed for arbitrary-order polyharmonic equations $(-\Delta)^m u=f$ on shape-regular polytopal meshes in $\mathbb R^d$, for any $m\ge1$, $d\ge2$, and polynomial degree $k\ge0$. The method uses the mixed variable $\sigma=\nabla^m u$ and a staggered primal--dual mesh to impose complementary continuity on scalar and tensor unknowns, without restrictions such as $d\ge m$. Local trace and bubble enrichments stabilize low-order tensor spaces without adding global unknowns. Hybridization localizes the tensor variable and yields an equivalent stabilization-free weak Galerkin formulation. Well-posedness and optimal energy error estimates are proved, and numerical experiments on polygonal and tetrahedral meshes confirm the predicted rates.

0

math.NA 2026-07-02

One-point correction lifts HSS solver to cubic convergence for scattering

by Yasuhiro Matsumoto, Rio Yokota

Proxy-surface-based fast direct solver for TE-mode scattering problems on distributed memory systems

Weakly singular Burton-Miller BIE plus Nyström discretization yields O(h^3) accuracy while preserving low-rank blocks for parallel direct so

abstract click to expand

This paper describes an MPI/OpenMP hybrid parallelized fast direct solver for the scattering problem of transverse electric (TE)-mode electromagnetic waves. Because TE-mode scattering can be reduced to the two-dimensional Helmholtz equation, solvers based on the hierarchically semiseparable (HSS) representation are highly attractive due to their high parallel efficiency. However, as the HSS representation applies low-rank approximations to all off-diagonal blocks, it exhibits poor compatibility with high-order discretization methods. We developed a fast direct solver with $O(h^3)$ convergence for Helmholtz transmission problems, whereas conventional HSS solvers typically yield only $O(h)$ convergence (where $h$ represents intervals between the quadrature nodes). It is based on the weakly singular Burton-Miller boundary integral equation and the Nystr\"om method with a one-point correction. Furthermore, recognizing that matrix component calculation, rather than matrix factorization, dominates the total computational time of HSS-type boundary integral solvers, we introduced a load-balancing method to maximize parallel efficiency. Numerical results demonstrate that the direct solver achieves high-accuracy convergence and nearly ideal strong and weak scalabilities.

0

cs.LG 2026-07-02

Parallel sampling cuts discrete diffusion complexity to log scale

by Yu Yao, Huanjian Zhou +3 more

Accelerating Discrete Diffusion Models with Parallel-In-Time Sampling

Reformulates τ-leaping with Picard iteration to reach O(log(d log S) log d) NFE while matching original sample quality.

abstract click to expand

Discrete diffusion models are widely used for learning and generating discrete distributions. As the generation process is inherently sequential, the acceleration of sampling is of significant importance. In this work, we parallelize the mainstream $\tau$-leaping algorithm for absorbing discrete diffusion in a Continuous-Time Markov Chain (CTMC) framework. By leveraging the continuous-time stochastic integral form of the $\tau$-leaping algorithm and the Picard iteration method, we achieve parallel-in-time sampling acceleration and provide a proof of exponential-factorial convergence for our algorithm. We improve the overall time complexity of $\tau$-leaping under absorbing settings from ${\mathcal{O}}(d \log S)$ to ${\mathcal{O}}(\log (d\log S)\cdot \log d)$ with respect to NFE. Empirically, our method shows consistent acceleration across synthetic and real-data settings. The new sampler achieves at most $7$--$9\times$ runtime speedup for synthetic distribution, and maintains the same quality with $50\%$ fewer NFE and $1.45$--$1.86\times$ runtime speedups in image/text tasks on a single GPU. Our research expands the potential of discrete diffusion models for efficient parallel inference, with broader implications for applications such as molecular structure and language generation.

0

math.NA 2026-07-02

Scheme for multi-species chemotaxis preserves positivity and mass

by Ao Zhang, Bingyin Zhang +1 more

A linear, decoupled and positivity-preserving time-staggered block-centered finite difference method for the multi-species Keller-Segel chemotaxis system

The decoupled method achieves second-order accuracy on non-uniform grids and supports adaptive resolution near blow-up.

abstract click to expand

In this paper, we present a linearly implicit, second-order block-centered finite difference (BCFD) prediction-then-projection scheme for the multi-species Keller-Segel chemotaxis system on non-uniform spatio-temporal grids. The proposed scheme integrates a standard Crank-Nicolson time-marching algorithm with an $L^2$ projection step to enforce positivity and mass conservation. The use of variable time stepsize and time-staggered discretization fully decouples the solutions of the multi-species cell density variables and the chemoattractant concentration variable while facilitating linearization, thereby greatly enhancing computational efficiency. Notably, the variable time-stepping algorithm and non-uniform grid BCFD discretization jointly enable adaptive resolution and local refinement near blow-up, thereby improving efficiency and accuracy without compromising the desired physical property-preserving in the simulation. Furthermore, using the mathematical induction method and the energy analysis approach, the unique solvability of the proposed scheme is rigorously proved, and we show that cell densities achieve second-order convergence in both time and space in the discrete $L^2$ norm, while the chemoattractant concentration achieves second-order convergence in the discrete $H^1$ norm. Representative numerical experiments are presented to validate the theoretical findings and demonstrate the reliability of the proposed scheme in simulating the blow-up phenomenon.

0

math.NA 2026-07-02

DWR estimator drives goal-oriented space-time meshes for Navier-Stokes

by Marius Paul Bruchhäuser, Nils Margenberg +1 more

Goal-oriented space-time adaptivity for the Navier--Stokes equations based on the dual weighted residual method

The approach controls errors in chosen flow quantities while keeping meshes computationally feasible and algebraic solves robust.

abstract click to expand

This work presents a goal-oriented a posteriori error estimator based on the Dual Weighted Residual (DWR) method together with space-time mesh adaptivity for the Navier--Stokes equations. The resulting nonlinear algebraic systems on the space-time slabs are solved by Newton's method with GMRES, preconditioned by a slab-wise geometric multigrid method. This combination yields reliable control of target quantities on computationally feasible space-time meshes together with a robust and efficient solution of the algebraic systems. The implementation is based on a MPI-parallel programming model in the deal.II library. Further ingredients are a discontinuous Galerkin discretization in time and inf-sup stable finite element pairs with discontinuous pressure on tensor-product meshes. The performance of the approach is investigated in benchmark computations with regard to accuracy, efficiency, and stability.

0

math.NA 2026-07-02

De Rham complex removes penalty tuning for high-dimensional Dirichlet problems

by Jiarong Chen, Xia Ji +2 more

Penalty-Free Natural Deep Ritz Method Based on de Rham Complex for High-Dimensional Dirichlet Boundary Value Problems

Boundary constraints become three coupled natural subproblems whose errors decay together up to six dimensions without any boundary weight.

abstract click to expand

Deep neural networks show great promise for high-dimensional PDEs, yet enforcing essential boundary conditions remains challenging, especially as penalty parameters require problem-specific retuning with increasing dimensionality. In this work, we extend the Natural Deep Ritz Method (NatDRM) [H. Yu and S. Zhang, J. Comput. Phys., 537 (2025)] to a unified framework for all dimensions $d \geq 2$ based on the de Rham complex and its penalty-free boundary decomposition: curl-type operators act on scalar potentials in 2D, vector potentials in 3D, and antisymmetric second-order tensor potentials in $d \geq 4$, respectively. This method converts Dirichlet constraints into three coupled natural (Neumann-type) subproblems with corresponding Ritz-type losses, eliminating the need for a boundary penalty parameter $\beta$. We derive dimension-unified discrete losses, lightweight boundary-based gauge-fixing regularizations to resolve curl-kernel non-uniqueness, and a joint training procedure; extensions to variable-coefficient elliptic and semilinear Poisson problems are formulated at the first subproblem level. Numerical experiments on smooth benchmarks up to 6D show that NatDRM, without any penalty tuning, matches or exceeds the accuracy of optimally tuned DRM and PINN in most cases. It converges stably in 6D where penalized DRM fails for most penalty values, and exhibits synchronous decay of interior and boundary errors, resolving the inherent imbalance of penalty-based methods.

0

math.NA 2026-07-02

Symmetric CAEs yield more accurate latent trajectories in PDE models

by G. Li Causi, N.Tonicello +2 more

Convolutional Symmetric AutoEncoders: enhancing latent stability via differential geometry

Extending representation consistency to convolutional layers cuts reconstruction errors and boosts robustness on advection, Burgers and Kura

abstract click to expand

Autoencoders (AEs) have emerged as powerful tools for non-linear dimensionality reduction, often surpassing traditional linear methods such as Proper Orthogonal Decomposition (POD) in scenarios characterized by slowly decaying Kolmogorov $n$-widths. In the realm of Reduced-Order Modelling (ROM), these models are increasingly utilized to learn low-dimensional representations of solution manifolds associated with parametric Partial Differential Equations (PDEs). However, the high expressivity of AEs presents a challenge: although trained networks typically minimize reconstruction error, they often struggle to capture the essential properties necessary for building accurate and robust ROMs. Recent works by arXiv:2307.15288v2 and arXiv:2506.11641v1 have tackled this challenge in fully connected AEs by proposing representation-consistent architectures, which preserve some of the properties belonging to POD. This study builds upon that concept by extending representation consistency for convolutional layers. We introduce a novel class of symmetric Convolutional AutoEncoders (CAEs) designed to embody the primary properties of manifold parametrization mappings. When integrated into a ROM framework, this architecture demonstrates significantly improved predictive capabilities. Specifically, we compared the performance of the ROMs based on classical and symmetric CAEs on three one dimensional academic test cases, namely the Linear Advection, the Viscous Burger and the Kuramoto Sivashinsky equation. Numerical results demonstrate that our proposed symmetric approach consistently yields more accurate latent trajectories, lower reconstruction errors, and enhanced model robustness.

0

math.NA 2026-07-02

Multilevel neural training solves BIEs with work linear in finest grid

by Bing-Ze Lu, Richard Tsai

A multilevel stochastic-gradient neural solver for boundary integral equations

Refining quadrature grids turn spectral bias into a smoother and convert the residual into an a posteriori error bound.

abstract click to expand

We develop a multilevel stochastic-gradient neural solver for boundary integral equations of the second kind. The unknown density is represented by a multilayer perceptron, trained by minimizing the Nystr\"om-discretized residual on a ladder of refining quadrature grids, each level warm-started from the parameters of the previous one. Each step requires only dense matrix-vector products on mini-batches of collocation rows and network passes, operations that map directly onto GPU hardware. The residual contraction is governed by the empirical neural tangent kernel (NTK), the discrete sample of a single continuum kernel. On a fixed grid, training stalls once the residual concentrates in modes the network contracts slowly, the plateau described by the frequency principle; a spectral analysis explains, and experiments confirm, how refining the quadrature resolves more of the continuum kernel's spectrum and returns these modes to the optimizer's reach. Spectral bias, elsewhere an obstruction to neural network solvers, thus serves as the smoother of a multigrid-type iteration, with quadrature refinement in place of coarse-grid correction. Under a uniform regularity bound on the network, the total work is a constant multiple of the work on the finest grid, and the uniform conditioning of the discrete second-kind operator leaves the NTK as the sole rate-determining spectrum while converting the training residual into an a posteriori error bound. Experiments on interior Dirichlet Laplace/Poisson problems and exterior Neumann Helmholtz problems, using both parametric and signed-distance surface representations, demonstrate the effectiveness and efficiency of the proposed method compared with GMRES at comparable tolerances.

0

math.NA 2026-07-02

NSFD scheme keeps p-Laplacian solutions positive at large time steps

by Achraf Zinihi, Matthias Ehrhardt +1 more

A Nonstandard Finite Difference Scheme for a Nonlinear Parabolic Equation with p-Laplacian-Type Diffusion

Nonlinear denominator and nonlocal diffusion approximation prevent negative values and oscillations that appear in standard methods.

abstract click to expand

We propose and analyze a nonstandard finite difference (NSFD) scheme for nonlinear parabolic equations involving a p-Laplacian-type diffusion operator in one- and two-dimensional spatial domains. Following Mickens' design principles, the proposed discretization employs a nonlinear denominator function phi(.) together with a nonlocal approximation of the nonlinear diffusion term Delta_p, yielding a structure-preserving discrete model. The scheme is designed to retain key qualitative properties of the continuous problem, including positivity, boundedness, and stability, which may be lost by standard finite difference methods (FDMs). We establish the well-posedness of the continuous model, derive the NSFD scheme, and investigate its consistency, convergence, and local truncation error. Numerical experiments confirm the theoretical results and demonstrate that, unlike the standard explicit FDM, the proposed NSFD scheme avoids spurious oscillations and nonphysical negative solutions even for relatively large time-step sizes.

0

math.NA 2026-07-02

Inner-outer iteration with optimal tuning solves stochastic Lyapunov equations

by Donghuan He, Xiaowen Su +2 more

An Inner-Outer Iteration Algorithm with Optimal Parameters for Stochastic Lyapunov Matrix Equation

Convergence from any start follows from the spectral radius of the iteration matrix when the system is mean-square stable.

abstract click to expand

This paper proposes an inner--outer (IO) iterative algorithm with optimal parameters for solving stochastic Lyapunov matrix equation associated with discrete-time stochastic linear system. First, under the assumption that the underlying stochastic linear system is asymptotically mean-square stable, the monotonicity and boundedness of the iterative sequence generated by the proposed algorithm are analyzed. On this basis, a sufficient convergence result is established for the zero initial condition. Second, by deriving the spectral radius of the corresponding iteration matrix, several necessary and sufficient convergence conditions are obtained for arbitrary initial conditions. In addition, the optimal parameter-selection strategies are developed to improve the convergence performance of the algorithm. Finally, numerical examples are presented to verify the theoretical results and demonstrate the advantages of the proposed algorithm over several existing iterative methods.

0

math.NA 2026-07-02

Operator matching yields reliable long-time dynamics

by Yinong Huang, Jonah Botvinick-Greenhouse +1 more

Data-Adaptive Learning of Dynamical Systems by Matching Transfer Operators and Invariant Measures

Transition-statistics learning on adaptive meshes outperforms pointwise fitting under noise and sparse sampling in chaotic system identifica

abstract click to expand

Trajectory-based learning of dynamical systems is often fragile in the presence of noise, chaos, or sparse observations, as small pointwise errors can rapidly amplify. We introduce a transition-statistics approach to system identification that learns dynamics by matching the induced motion of probability mass across a data-adaptive mesh. Given trajectory data, we build an unstructured partition of state space and approximate the Perron--Frobenius operator with a regularized Ulam transition matrix. We replace hard cell indicators with continuous, piecewise-smooth partition-of-unity weights, yielding a Markov matrix supporting gradient-based optimization with respect to the parameters of a learned vector field. This enables two related training objectives: matching invariant measures through the stationary eigenvectors of the transition matrices, and matching the full transition matrices to capture transport between regions of state space. Numerical experiments on Lorenz-63, Lorenz-96, and a reduced-order NOAA sea surface temperature forecasting problem show that transition-statistics matching gives more reliable long-time dynamics than pointwise trajectory matching, particularly under measurement noise and sparse sampling. The approach provides a robust operator-theoretic alternative to trajectory-level losses for learning chaotic and partially observed dynamical systems.

0

math.NA 2026-07-02

Low-regularity integrators give linear parareal convergence for NLS

by Qingle Lin, Zhi Zhou

Exponential Low-Regularity Parareal Algorithms for Nonlinear Schr\"odinger Equations

Contraction factor scales with coarse step size for solutions of limited regularity.

abstract click to expand

The parareal algorithm is one of the most widely studied parallel-in-time methods for the numerical approximation of time-dependent problems. For non-diffusive equations, however, standard parareal methods may converge slowly or even become unstable due to the absence of damping, while nonlinear interactions can transfer and amplify phase errors across Fourier modes. In this work, we consider the nonlinear Schr\"odinger equation (NLS) as a representative non-diffusive model and analyze parareal algorithms with an exact fine propagator, with particular emphasis on the design of suitable coarse propagators. We establish a general convergence framework, valid for solutions with limited regularity, under stability and local truncation error assumptions on the coarse propagator. These assumptions are verified for selected exponential low-regularity integrators designed for one-dimensional quadratic and cubic NLS equations, which achieve optimal approximation orders without derivative loss. To the best of our knowledge, this is the first construction of parareal algorithms for NLS equations that are provably linearly convergent, with a contraction factor proportional to the coarse time-step size even for solutions of limited regularity. Numerical experiments on quadratic, cubic, and quintic NLS equations demonstrate rapid convergence and improved performance over parareal variants using classical coarse propagators, including Lie and Strang splitting methods and first- and third-order exponential Runge--Kutta integrators.

0

math.NA 2026-07-02

PINNs solve delay differential equations with history switch and segments

by Stone Yao, Vipin Kumar +1 more

Physics Informed Neural Networks for Nonlinear Delay Differential Equations

Differentiable history tracking plus segmented collocation keeps gradients stable and solutions continuous over long intervals.

abstract click to expand

In this paper we propose a novel physics-informed neural network framework for solving general first-order delay differential equations. Our approach combines a differentiable history switch, a trial-solution formulation that explicitly enforces history constraints, and a segmented collocation strategy to stabilize gradient propagation across large temporal domains. The method enables a scalable and physics-consistent approximation of delay differential equation solutions while maintaining continuity across subintervals. Numerical experiments demonstrate the effectiveness of the proposed method.

0

math.NA 2026-07-02

Two-stage PINN matches FEM warpage at 1000x speed

by Xinyu Li, Min Tang +4 more

WarpagePINN: Thermal Warpage Prediction in Advanced Packaging via a Two-Stage Physics-Informed Neural Networks

WarpagePINN predicts 0.2-micron chiplet deformation using only physical equations, without simulation training data.

abstract click to expand

Thermal warpage has become a critical issue in advanced packaging, primarily caused by the mismatch in coefficients of thermal expansion (CTE) among heterogeneously integrated materials. However, only a limited number of studies have focused on developing computational methods for coupled thermal-warpage prediction in the chiplet. This paper proposes a two-stage physics-informed neural network (WarpagePINN) framework to compute both temperature profile and warpage deformation of chiplets. The neural networks are trained without relying on labeled datasets generated by conventional simulators. In the first stage, the temperature field is modeled using a Fourier series representation that inherently satisfies boundary conditions, and the network is trained solely through a loss function derived from the governing equation. In the second stage, a multilayer perceptron (MLP) is employed for warpage prediction, utilizing a novel hybrid supervisory strategy to optimize the energy-based loss function instead of residual loss. A parametric WarpagePINN is also developed to quantify uncertainties associated with the CTE. Numerical results show that the proposed WarpagePINN framework achieves excellent agreement with conventional finite element methods, with a mean absolute error (MAE) of 0.2 {\mu}m, while achieving a speedup of approximately 1000 {\times} in CTE parameterization studies.

0

math.NA 2026-07-02

RLM schemes dissipate relaxed original energy in phase fields

by Xiaobo Jing, Jia Zhao

Relaxed Lagrange Multiplier (RLM) Schemes for Phase Field Models Preserving the Relaxed Original Energy Dissipation Law

They track the true dissipation rate and stay uniquely solvable over wide time steps at SAV-like cost.

abstract click to expand

Phase-field models are typically derived from variational principles for a free-energy functional and are widely used to simulate complex multiphase phenomena in science and engineering. A central goal in designing numerical schemes for these models is to preserve the underlying energy-dissipation law. In this paper, we propose a class of relaxed Lagrange multiplier (RLM) schemes for phase field models. In contrast to popular scalar auxiliary variable (SAV) and invariant energy quadratization (IEQ) methods, which dissipate a modified energy involving auxiliary variables, the RLM schemes dissipate a relaxed version of the original energy and closely track the original energy dissipation rate. Compared with the classical Lagrange multiplier (LM) approach, the RLM schemes ensure that the resulting discrete system is uniquely solvable over a broad range of time steps. The key idea is to augment the LM formulation with a relaxation term, yielding a scalar quadratic equation for the multiplier with an explicit closed-form solution. The resulting schemes are linear and efficient because each time step requires solving only two linear systems with constant coefficients, at a cost comparable to that of SAV schemes. We construct both first-order and second-order variants and prove their energy stability. Numerical experiments verify the expected convergence rates and demonstrate that the RLM schemes accurately capture interface dynamics.

0

math.NA 2026-07-02

Dirac equation splitting errors scale as τ/ε² and τ²/ε³

by He Wang, Jia Yin

Rigorous analysis of the time-splitting methods for the semiclassical Dirac equation

Rigorous analysis gives exact dependence of approximation errors on time step, mesh size and small parameter ε.

abstract click to expand

We provide rigorous error analysis of the mass-preserving time-splitting methods for solving the semiclassical Dirac equation. The scaled Planck constant $\epsilon$ in the equation gives rise to rapid oscillations in both space and time when $0<\epsilon\ll 1$ with wavelengths of order $O(\epsilon)$. %We prove that the first-order splitting $S_1$ and the second-order splitting $S_2$ schemes preserve the total discretized mass. Rigorous error estimates reveal the precise dependence of the approximation errors on the time step $\tau$, the spatial mesh size $h$, and the parameter $\epsilon$. Specifically, the temporal error scales as $O\left(\tau/\epsilon^2\right)$ for the first-order splitting $S_1$ and as $O\left(\tau^2/\epsilon^3\right)$ for the second-order splitting $S_2$, while the spatial error scales as $O(h^m/\epsilon^m)$ for both methods, where $m$ is related to the regularity of the solution. In addition, we obtain error bounds for key physical observables, including the total probability density $\rho$ and the current density $\mathbf{J}$. Compared with finite difference time domain (FDTD) methods, time-splitting approaches exhibit spectral accuracy in space and retain a relatively low computational cost. Furthermore, we demonstrate that higher accuracy can be achieved by employing the fourth-order compact time-splitting ($S_\text{4c}$) method. Numerical experiments are conducted to verify the reliability of the error estimates.

0

stat.ML 2026-07-02

FNOs achieve polynomial sample complexity on dissipative PDE operators

by Nisha Chandramoorthy, Daniel Sanz-Alonso +1 more

From Spectral Methods to Sample Complexity Bounds for Fourier Neural Operators

Bounds hold uniformly over families of equations when operators admit spectral discretizations, with rates set by smoothness and dimension.

abstract click to expand

We establish approximation and learning guarantees for Fourier neural operators (FNOs) applied to time-$T$ solution operators of dissipative evolution equations. The analysis builds on the premise that FNOs can efficiently approximate and learn solution operators whenever these operators admit stable and accurate spectral discretizations. To formalize this idea, we introduce classes of evolution operators defined through spectral methods and derive FNO approximation bounds and polynomial sample complexity guarantees for these classes. For equations with polynomial nonlinearities, the learning rates depend primarily on the smoothness of the input space and the dimension of the physical domain. Our results hold uniformly over broad families of dissipative equations, rather than for a single fixed PDE, and apply in particular to the Navier--Stokes, Allen--Cahn, and Cahn--Hilliard equations. For equations with non-polynomial smooth nonlinearities, we prove that polynomial sample complexity still holds with rates that now additionally depend on the smoothness of the nonlinear terms and the dissipation strength. Overall, we connect classical spectral approximation theory with modern operator learning and explain when FNOs can learn nonlinear evolution operators efficiently.

0

math.NA 2026-07-02

Quantum algorithm speeds Stokes flow simulation exponentially in dimensions

by Shi Jin, Jiaqi Tang +2 more

Quantum Simulation of Stokes Flow via Schr\"odingerisation and Artificial Compressibility

Artificial compressibility and Schrödingerisation turn the saddle-point problem into an explicit quantum circuit with exponential dimensiona

abstract click to expand

Simulating incompressible Stokes flow is essential for studies in microfluidics and low-Reynolds number hydrodynamics. However, the computational cost of resolving the associated saddle-point problem grows prohibitively with the dimensionality of the problem. In this work, we present a quantum algorithm based on the Schr\"odingerisation technique for the Stokes equations, incorporating an artificial compressibility regularization. The core of our approach is the design of an explicit quantum circuit that encodes the resulting regularized system. The artificial compressibility formulation provides a unified framework for the system, which is then efficiently mapped to a quantum circuit via the Schr\"odingerisation procedure. A rigorous complexity analysis demonstrates the quantum computational advantage of our algorithms in high-dimensional settings, notably an exponential speedup in problem dimensionality. The validity and scalability of the proposed method are corroborated by numerical simulations performed on Qiskit.

0

math.NA 2026-07-01

Floquet-Bloch transform locates defects in infinite periodic arrays

by Dinh-Liem Nguyen, Nhung H. Nguyen +1 more

Identifying defective units in infinite periodic arrays of point sources

Boundary data and a sampling-algebraic step recover the number, positions and intensities of defective point sources after the periodicity-b

abstract click to expand

This paper focuses on identifying defective units in unbounded periodic arrays of point sources using boundary data. The study is motivated by the noninvasive evaluation of large-scale periodic source systems. Unlike classical inverse source problems in free space, the key challenge here lies in the disruption of periodicity caused by defective sources in the infinite array. To address this, we employ the Floquet - Bloch transform to reformulate the original inverse source problem as a quasi-periodic inverse source problem. We first establish uniqueness theorems for both the original and the quasi-periodic formulations. Then, we develop a new numerical method for identifying defective sources. This method combines a sampling indicator function with an algebraic technique to determine not only the number of defective sources, but also their locations and intensities. Numerical experiments are presented to validate the effectiveness of the proposed method.

0

math.DS 2026-07-01

Filters split antagonistic dynamics into three delay-stability regimes

by Alexander Omelchenko

Implementation Filters and Delay-Budget Instability in Coupled Replicator--Mutator Dynamics

Hard lags sum while implementation rates set separate margins, so only certain reductions restore stability in the intermediate window.

abstract click to expand

We model an adaptive contest in which two antagonistically coupled populations continually reallocate effort among competing methods, but decisions are not fielded instantly. Each side has an intended portfolio and a deployed portfolio: intended reallocations follow delayed observations of the opponent, while deployment follows intent through a first-order implementation filter. Under barycentric balance and uniform exploration, the linearized scalar branches have a characteristic factor in which hard observation and deployment lags enter only through their total sum, whereas implementation rates enter through real filter factors that cannot be absorbed into selection or exploration. In the strictly antagonistic class, negative spectral branches split into three regimes: weak branches have no positive-frequency crossing, intermediate branches lose stability through a delay-induced Hopf bifurcation, and strong branches are at or beyond the implementation-filter instability margin already at zero hard delay. This gives an operational delay-budget rule: in the delay-induced window, reducing any hard lag has the same first-order stabilizing leverage at onset; in the filter-induced regime, hard-lag reduction alone cannot restore stability. Balanced scalar performance observables generically show a mean shift and a second harmonic at twice the compositional frequency, and under strict antagonism the two performance signals are locked in antiphase with fixed amplitude ratio. For a baseline branch, a finite-dimensional Hopf normal-form calculation gives a negative cubic coefficient, and direct simulations reproduce the predicted threshold, amplitude scaling, and observable signatures. Motivating applications include cybersecurity and rapid technological countermeasure adaptation.

0

math.NA 2026-07-01

GQL reformulation yields non-iterative PCP limiters for relativistic hydro

by Linfeng Xu, Shengrong Ding +1 more

GQL-Based Physical-Constraint-Preserving High-Order Finite Difference Schemes for Special Relativistic Hydrodynamics in Arbitrary Dimensions

Linear inequalities and small eigenvalue solves enforce positive density and subluminal velocity in high-order WENO schemes up to 3D

abstract click to expand

High-order accurate simulations of special relativistic hydrodynamics (RHD) are prone to numerical breakdown if intrinsic physical constraints (positive rest-mass density/pressure and subluminal velocity) are violated near strong discontinuities. In this work, we develop a robust and efficient physical-constraint-preserving (PCP) flux-limiting framework for high-order schemes, using finite-difference WENO as a representative example. By leveraging the geometric quasilinearization (GQL) representation, which equivalently reformulates the nonlinear RHD constraints into a family of linear inequalities, we integrate a Zalesak-type Flux-Corrected Transport (FCT) update into a scalar-style limiter that acts directly on conservative variables. A critical innovation is the explicit, non-iterative determination of limiting parameters via a rational stereographic parameterization of the GQL normal vector. This technique transforms the required worst-case minimization over auxiliary variables into a generalized Rayleigh-quotient formulation, allowing the optimal parameters to be obtained by solving small symmetric eigenvalue problems ($2\times2$ in 1D; $(d+1)\times(d+1)$ in $d$ dimensions). Relaxed variants are further introduced to reduce computational costs in multidimensions while retaining the PCP guarantee. Extensive numerical benchmarks ranging from 1D to 3D, including ultra-relativistic Riemann problems and astrophysical jets, demonstrate that the proposed method robustly enforces physical admissibility, sharply resolves discontinuities, and maintains design-order accuracy for smooth solutions.

0

math.NA 2026-07-01

CDG on overlapping meshes reaches proven O(h^{k+2}) superconvergence

by Manting Peng, Kailiang Wu

Hidden Accuracy and Superconvergence Analysis of Central Discontinuous Galerkin Methods on Overlapping Meshes

Projection-correction analysis recovers hidden residual cancellations to prove pointwise accuracy and stronger cell averages under corrected

abstract click to expand

This paper establishes the first rigorous superconvergence theory for semidiscrete and fully discrete central discontinuous Galerkin (CDG) methods for linear hyperbolic equations on overlapping meshes. While the optimal $L^2$ convergence of $\mathbb{Q}^k$ CDG schemes was established on uniform Cartesian meshes by Liu, Shu, and Zhang [ SIAM J. Numer. Anal.}, 56 (2018), pp. 520--541], their observed $\mathcal{O}(h^{k+2})$ pointwise superconvergence has remained unproven, due to the loss of standard single-mesh Galerkin orthogonality inherent in the CDG overlapping structure. To overcome this fundamental barrier, we introduce a projection-correction framework that identifies a hidden superconvergent mechanism: an asymptotic weak residual cancellation in one dimension, and a high-order cancellation-by-aggregation (HOCA) mechanism in multiple dimensions. This HOCA approach overcomes the analytical challenge posed by coupled primal-dual directional residuals, recovering critical error cancellation properties absent from the standard variational formulation. Consequently, we provide the rigorous proof of the conjectured $\mathcal{O}(h^{k+2})$ pointwise superconvergence in the discrete $\ell^{\infty}$ norm across all superconvergent points. Furthermore, we reveal that under a systematically corrected initialization, this framework yields a previously undiscovered, stronger cell-average superconvergence estimate of order $\mathcal{O}(h^{\min\{2k+1,k+3\}})$. The theory is extended to fully discrete explicit Runge--Kutta CDG schemes, where stagewise corrected errors are constructed to preserve spatial superconvergence up to temporal truncation errors, yielding a stable reconstruction-based postprocessing estimate. Numerical experiments in one and two spatial dimensions confirm the sharpness of the theoretical rates.

0

cs.LG 2026-07-01

Neural preconditioner cuts Newton failures in CZM simulations

by Zhangyong Liang, Huanhuan Gao

Interface-Aware Neural Newton Preconditioning for Robust Cohesive Zone Model Simulations

Learned interface corrections recover branches and reduce step cuts while leaving the force-displacement response unchanged.

abstract click to expand

Cohesive Zone Models (CZMs) are widely used to simulate interface fracture, delamination, adhesive failure, and fiber--matrix debonding in aerospace composite structures. In implicit quasi-static finite element analyses, cohesive softening may introduce negative interface tangents, solution jumps, and Newton-basin mismatch, so the previous converged state can become a poor initial guess for the next increment. This may lead to stagnation, wrong-branch convergence, or repeated step cuts. Existing remedies, including viscous regularization, path following, dynamic relaxation, and manual Newton--Raphson (NR) modification, either alter the effective response, increase cost, or rely on hand-crafted interface rules. This work proposes an Interface-Aware Neural Newton Preconditioner (IA-NNP) for difficult CZM increments. IA-NNP recasts manual NR modification as rule-based interface lifting and generalizes it into a learned, state-dependent interface correction. The method acts only on active interface variables and preserves the original traction--separation law, residual assembly, tangent evaluation, history update, and dissipation checks. Two realizations are developed: IA-NNP-Init for learned initial-guess lifting and IA-NNP-NL for iteration-level nonlinear right preconditioning. Interface graph features encode opening, traction, tangent, damage/history variables, mode mixity, residuals, and neighboring states. The correction is bounded, confidence-gated, and accepted only through the original CZM Newton solve. A root-equivalence property shows that IA-NNP changes the path to convergence but not the discrete CZM solution set. Tests on horizontal, circular, two-interface, and active-front benchmarks show improved difficult-increment convergence, better branch recovery, and fewer failures than standard NR and manual NR modification, while preserving the force--displacement response.

0

math.NA 2026-07-01

Closed-form eigenvalues make spline bases outlier-free for Laplace spectra

by Damiano Ricci, Carla Manni +1 more

Isogeometric discretizations for the spectrum of the Laplace operator: outlier-free spline bases

Optimal spline subspaces get bases with B-spline support whose matrices diagonalize simultaneously, yielding explicit eigenvalues.

abstract click to expand

Optimal spline subspaces are an elegant and efficient tool to remove spurious outliers in isogeometric Galerkin discretizations for the approximation of the spectrum of the Laplace operator. For practical purposes, it is valuable to have a basis construction for such spaces with good computational and spectral properties. We provide a characterization of the bases that enjoy a B-spline-like support structure and whose mass and stiffness matrices are simultaneously diagonalizable. It turns out that these mass and stiffness matrices admit explicitly known closed-form expressions for their eigenvalues, implying that the considered bases are outlier-free. A numerical procedure for constructing such bases is also presented.

0

math.NA 2026-07-01

Deflating eigenvectors keeps WaveHoltz rate constant at high frequency

by Daniel Appelo, William D. Henshaw +1 more

Numerical Study of Eigenvector Deflation to Accelerate the WaveHoltz Method

In 2D the iteration count stops growing with omega once the deflation set size scales quadratically with frequency.

abstract click to expand

We present a numerical study of eigenvector deflation as a means of accelerating the WaveHoltz method for solving the Helmholtz equation. For energy-conserving (Dirichlet or Neumann) boundary conditions the WaveHoltz fixed-point iteration converges slowly at high frequency, requiring approximately $\mathcal{O}(\omega^{2d})$ iterations in $d$ dimensions. We show that deflating the eigenvectors whose eigenvalues lie nearest the driving frequency substantially reduces iteration counts, and we examine two ways of incorporating the eigenvectors: direct eigenvector deflation (DEVD), in which the forcing and iterate are projected against the deflation set, and augmented-Krylov eigenvector deflation (AUKED) using deflated conjugate gradient (DCG), augmented GMRES (AGMRES), and augmented (recycled) BICGSTAB (ABICGSTAB). The required eigenpairs can be computed efficiently with the EigenWave approach, and we demonstrate, in two dimensions, that when the number of deflation vectors grows quadratically with $\omega$ the asymptotic convergence rate remains essentially constant. Because the eigenvectors on structured grids are naturally represented as matrices, we further apply SVD-based compression to reduce their storage. Numerical experiments on single curvilinear grids discretized with summation-by-parts operators, and on overset grids illustrate the robustness and efficiency of the approach, with the deflated solver breaking even against the undeflated solver after as few as two right-hand sides, when accounting for the cost of precomputing the eigenvectors.

0

math.NA 2026-07-01

Hermite splines plus Legendre derivatives form C^m bases for any m

by Peter K Moore

Complementary families of approximating polynomials with applications to finite element methods applied to differential equations of arbitrary even spatial order

Explicit nodal and bubble polynomials yield p >= 2m+2 elements with proven interpolation error and H^{m+1} superconvergence for linear probl

abstract click to expand

Complementary families of polynomials are introduced to generate $C^m$ finite element basis functions of order $p \geq 2m+2$ for arbitrary $m \ge 0$. One family consists of the Hermite splines that serve as the nodal basis functions by ensuring $C^m$ continuity across element boundaries. Explicit formulas for these splines for any $m \ge 0$ are presented on the canonical interval $[0,1]$. The second family is derived on the interval $[-1,1]$ from derivatives of order $m+1$ of the Legendre polynomials of degree $p-m-1$ multiplied by binomial powers of degree $m+1$ at -1 and 1, respectively. These polynomials, related to the ultraspherical polynomials, serve as the interior or bubble basis functions. A relationship between the two families of polynomials is demonstrated. For a particular $m$ and $p$, an interpolant is constructed using these basis pairs together with the roots of the related ultraspherical polynomial and the interval endpoints. A formula for the interpolation error that extends the results for $m=0$ and $m=1$ is given. To prove the formula extensions of the Lagrange interpolants are introduced. A superconvergence result along with the related asymptotic equivalence of the interpolant and finite element solution is proved in the linear case in $H^{m+1}$. Computational results demonstrate the theory for a model problem.

0

math.NA 2026-07-01

Adaptive differentiation recovers nonlinear PDEs from noisy data in strong form

by Roy Y. He, Sung Ha Kang

PDE Identification Using Noise Adaptive Differentiation in Strong Form (S-IDENT)

S-IDENT matches weak-form accuracy at noise levels where earlier strong-form methods fail by choosing SG windows via SURE.

abstract click to expand

We explore identifying partial differential equations (PDEs) from noisy observations of single time-space trajectories. Recent developments show the benefits of identifying PDEs in their weak forms. We investigate the use of differential Strong-form dictionaries for PDE IDENTification (S-IDENT), which enables finding more general linear and nonlinear PDEs. Building on an extensive exploration of integral-type denoised differentiation approaches, we propose to use Savitzky--Golay (SG) differentiation with an adaptive window length chosen based on Stein's Unbiased Risk Estimate (SURE). This offers a guaranteed order of accuracy while producing estimators with minimal variance. The identification process is further refined and stabilized through trimming and reduction-in-residual model selection. Numerical evidence shows that S-IDENT can successfully identify nonlinear PDEs at higher levels of noise than existing strong-form methods, while also yielding results comparable to weak-form approaches. We further verify the effectiveness of S-IDENT through comparisons with various strategies to approximate differential features. We provide numerical evidence that general differential-form dictionaries are larger and more ill-conditioned than those used for weak-form identification, yet S-IDENT does not significantly suffer from this combinatorial increase in dictionary size.

0

math.NA 2026-07-01

Clustering method dictates FMM stability for Helmholtz BEM

by Wolfgang Kreuzer, Christian Kasess

Effect of different clustering approaches on the multilevel fast multipole method for the Helmholtz equation

Tests on uniform and irregular meshes show cluster formation affects speed and whether the computation runs without numerical breakdown.

abstract click to expand

The fast multipole method (FMM) is an important component for the boundary element method (BEM), because with the FMM the efficiency and feasibility of the BEM can be enhanced to a large degree. Part of the FMM is grouping the elements of the boundary element mesh into different clusters. The size of these clusters in terms of number of elements and spatial expansion has a huge impact on the efficiency and stability of the method. However, while the theory behind the multipole expansion has been broadly researched, the clustering process itself and its effect on the FMM has been neglected in comparison. Most of the time, for example, it is implicitly assumed that the elements of the mesh have about the same size, which is often not the case in practical applications, e.g., when calculating the sound field around the human head. In this study we compare different types of clustering approaches with respect to stability and efficiency of the underlying FMM applied to meshes that have uniform as well as non-uniform element sizes. Also, some examples are provided for cases where a wrong clustering can lead to numerical problems and instabilities of the FMM-BEM.

0

math.OC 2026-07-01

Regularized Newton converges linearly to manifold-affine intersections

by Dengyu Zheng, Shixiang Chen

A Geometry-Adaptive Regularized Newton-Type Method for Manifold-Affine Intersection Problems

A residual-dependent regularization parameter restores quadratic convergence when full transversality holds.

abstract click to expand

We propose Regularized Newton-SLRA (RN-SLRA), a regularized Newton-type method for local manifold--affine intersection problems motivated by structured low-rank approximation. Classical Newton-SLRA achieves fast local convergence under transversality, but its tangent-space intersection step may become ill-defined, singular, or severely ill-conditioned when transversality fails. RN-SLRA overcomes this difficulty by replacing the exact tangent-space intersection step with a regularized quadratic subproblem over the affine space. Under intrinsic transversality, we prove local linear convergence to the intersection. Under transversality, we show that a residual-dependent choice of the regularization parameter yields higher-order local convergence; in particular, the method converges quadratically for the linear residual rule. We also analyze an inexact variant based on quasioptimal manifold projections. When the quasioptimality constant is sufficiently accurate, the inexact method retains local residual convergence. Numerical experiments on constructed degenerate SLRA instances and Hankel-structured examples illustrate the robustness of RN-SLRA in settings where Newton-SLRA may fail, and show that the inexact variant can reduce the projection cost in large-scale problems.

0

quant-ph 2026-07-01

Joint oracle gives quadratic speedup for quantum boundary-value search

by Daniel Jaroszewski, Bastian Harrach

A Quantum Collocation Approach to One-Dimensional Boundary Value Problems with Coherent Amplitude Amplification

Residual-threshold construction on spatial and parameter registers achieves poly-log gate cost while keeping amplitude-amplification gains.

abstract click to expand

We propose a quantum collocation framework for approximating solutions of one-dimensional linear and nonlinear boundary value problems. The method formulates the search for admissible solutions as a residual-based quantum search over a discretized ansatz space, where candidate solutions are evaluated through residual conditions imposed at collocation points. A residual-threshold oracle is constructed that acts jointly on spatial and parameter registers. This joint oracle structure leads to amplification dynamics that decompose into a coherent superposition of spatially conditioned amplitude-amplification processes rather than a single global amplification mechanism. We derive the corresponding amplification geometry and show that the success probability is governed by a weighted combination of spatially dependent amplification angles. Furthermore, we prove that the reversible residual oracle can be implemented with gate complexity polynomial in the logarithm of the number of collocation points, while retaining the quadratic search acceleration associated with amplitude amplification in the parameter space. We analyze how the spatially dependent oracle structure influences the amplification dynamics and corresponding success probabilities. Furthermore, we investigate how discretization, ansatz expressivity, oracle tolerance, and finite-precision effects influence both approximation quality and amplification behavior. Numerical experiments validate the theoretical predictions and illustrate the resulting search dynamics across different discretization and precision regimes.

0

math.NA 2026-07-01

Quadratic correction restores positivity to low-rank Vlasov solutions

by Katharina Kormann, Murtazo Nazarov +1 more

Positivity-preserving dynamical low-rank methods for the Vlasov equation

Minimal adjustment term enforces physical constraints while preserving the efficiency of the low-rank form in plasma kinetic models

abstract click to expand

In this manuscript, we introduce positivity-preserving correction methods for low-rank approximations of the Vlasov equation. The key idea is to formulate structural properties, including positivity-preservation, as constraints and to seek a minimal correction term that is added to the low-rank solution, by solving a quadratic programming problem. As a result, the corrected solution satisfies the constraints and preserve these properties, while remaining close to the original low-rank solution. Two positivity-preserving schemes are proposed in this work, and one of them also preserves the total mass and momentum of the system. We apply the proposed methods to a Vlasov--Poisson and Vlasov--Poisson-BGK employing a spectral discretization in space and an explicit Runge--Kutta scheme in time. Numerical experiments demonstrate the effectiveness of the proposed methods.

0

math.NA 2026-07-01

Kernel of V-line tensor transforms characterized inside disk

by Rahul Bhardwaj, Madhu Gupta

V-Line Tensor Tomography in a Disk: Theoretical and Numerical Reconstruction

Decomposition yields explicit inversion formula; numerical tests recover vectors and 2-tensors with noise.

abstract click to expand

In this article, we investigate V-line transforms for symmetric $m$-tensor fields whose support lies inside a disk of radius $R$ and centered at the origin. We provide an explicit characterization of the kernel of the V-line transforms acting on a symmetric $m$-tensor field and derive a new inversion formula using a decomposition result. In addition, we present a comprehensive numerical verification and validation of the inversion algorithms for these V-line transforms for vector fields and symmetric $2$-tensor fields, which were recently developed in \cite{bhardwaj_2024,bhardwaj2025tensor}. The reconstruction results obtained for various phantoms demonstrate the effectiveness and robustness of the proposed numerical methods, including in the presence of noise.

0

cs.AR 2026-07-01

PEERS computes exact resistances on 1M-node graphs in 18.8 seconds

by Baiyu Chen, Lin Gan +2 more

PEERS: A Parallel and Exact Effective Resistance Solver via Implicit Inversion and Augmented Symbolic Analysis

Implicit Cholesky inversion lets the solver answer all-edge queries exactly in one parallel sweep and scale to 17M nodes in under an hour.

abstract click to expand

High-precision effective resistance computation is a cornerstone of Electronic Design Automation (EDA) sign-off, yet it remains a fundamental bottleneck in large-scale power grid analysis, spectral sparsification, and circuit reliability. Existing approaches face a prohibitive "precision-memory impasse": approximate methods lack the stringent accuracy required for high-stakes industrial sign-off, while exact methods either suffer from redundant query overheads or trigger $O(n^2)$ memory explosions. To resolve this, we propose PEERS, a Parallel and Exact Effective Resistance Solver powered by an implicit inverse computing model of the Cholesky factor. By integrating a state-inherited augmented depth-first search (DFS) with a dynamic query update mechanism, PEERS eliminates numerical redundancy and evaluates all-edge resistance queries in a single parallel sweep. We provide a rigorous Work-Span analysis, proving that for graphs satisfying an $O(n^\alpha)$ separator theorem, PEERS achieves a theoretically optimal parallel span of $O(n^\alpha)$ while strictly maintaining $O(nnz(L))$ space complexity. Numerical evaluations on industrial benchmarks demonstrate that PEERS achieves an average speedup of 83.3x over state-of-the-art parallel solvers under identical memory constraints. Notably, PEERS processes a 1-million-node industrial graph in just 18.8 seconds and scales to 17 million nodes in under an hour, providing the first computationally feasible path for exact all-edge resistance analysis in multi-million-gate designs.

0

math.NA 2026-07-01

DKT error for biharmonic equation bounded by Hessian approximation

by Dietmar Gallistl, Ngoc Tien Tran

An error analysis of discrete Kirchhoff elements

Bound uses only piecewise constant approximation and right-hand side oscillation, implying first-order convergence without extra regularity

abstract click to expand

The Discrete Kirchhoff Triangle (DKT) method for the biharmonic equation is analyzed in the discrete energy norm. The error is bounded by the best approximation of the Hessian by piecewise constants and the oscillation of the right-hand side, without additional regularity assumptions on the exact solution. This result implies first-order convergence of the classical DKT element and the analysis yields a canonical extension to three space dimensions with the same approximation properties. Residual-based a posteriori error estimates are derived. The analysis is formulated within a general framework for low-order nonconforming methods, which also applies to various classical elements and yields best-approximation results by constants. It is furthermore shown how known stable pairs for the planar Stokes system have discrete stream functions in discrete Kirchhoff spaces. This yields variants of the known schemes with positive definite formulations and pressure-robust error bounds.

0

math.NA 2026-07-01

Fractal-fractional operator yields stable HIV trajectories in simulations

by Niaz Ali Shah, Samad Noeiaghdam

Fractal-Fractional HIV Dynamics with Mittag-Leffler Kernel: Analysis, Stability, and Numerical Simulations

Existence, Hyers-Ulam stability, and parameter-sensitivity plots confirm the Mittag-Leffler kernel captures both transient outbreaks and lon

abstract click to expand

In this paper, a fractal--fractional HIV model with the Mittag--Leffler kernel is proposed using the Atangana--Baleanu--Caputo operator to capture the memory and hereditary properties of the disease dynamics. The existence and uniqueness of the solutions are investigated using suitable analytical techniques, and the Hyers--Ulam stability analysis is carried out to verify the stability behavior of the proposed system. For the numerical simulations, the Newton polynomial approximation method together with the Atangana--Toufik numerical scheme is employed to obtain approximate solutions for different parameter settings. Furthermore, several visualization techniques, including sensitivity heatmap representation and tornado diagram analysis, are utilized to study the influence of model parameters on the HIV dynamics. The obtained numerical results demonstrate that the proposed fractal--fractional framework provides an effective and reliable approach for analyzing the transient and long-term behavior of HIV transmission dynamics.

0

math.NA 2026-07-01

FFT Nyström solver for 1000 acoustic scatterers

by Jun Lai, Yuxin Li

A Spectral Solver for Acoustic Scattering by Multiple Quasi-Axisymmetric Structures

Reformulating singular integrals as modal Green's functions computed by FFT enables fast, accurate solutions for large collections of quasi-

abstract click to expand

Acoustic scattering arises in a wide range of applications, including medical imaging, geophysical exploration, acoustic metamaterials, etc. In this paper, we develop a fast and highly accurate algorithm for acoustic scattering by multiple quasi-axisymmetric objects, whose axis of rotation is an arbitrary curve. The method is based on a Nystr\"om discretization that combines Gauss-Legendre quadrature with the trapezoidal rule. To treat the singular integrals that occur when target points are close to or coincide with source points, we reformulate them as evaluations of the modal Green's function and its derivatives, which are computed efficiently using the fast Fourier transform and convolution. The multiple scattering solver is then constructed by coupling the single scatterer discretizations through inter-body boundary integral interactions. We also present a convergence analysis for scattering problems with smooth geometries. Numerical examples demonstrate the efficiency and accuracy of the proposed method for solving multiple scattering problems involving up to 1000 quasi-axisymmetric structures.

0

math.NA 2026-07-01

Domain split lets randomized nets solve unbounded PDEs

by Haixin Wang, Haoning Dang +2 more

Domain-Decomposed Randomized Neural Networks for Partial Differential Equations in Unbounded Domains

Near-field and far-field subnetworks coupled by least-squares yield conditional error bounds without artificial boundaries.

abstract click to expand

Partial differential equations on unbounded domains are challenging because the exterior region must be represented without excessive truncation error. Truncation-based methods often require problem-dependent artificial boundary conditions, while global spectral bases may be inefficient for localized structures, irregular geometries, or solutions with different near-field and far-field behaviors. We propose a domain-decomposed randomized neural network framework for such problems. Different randomized subnetworks are assigned to different spatial regimes: a near-field subnetwork captures local and geometric features, whereas a far-field subnetwork represents exterior decay. The subnetworks are coupled by boundary and interface conditions, and only the output-layer coefficients are solved from linear least-squares systems arising from Petrov--Galerkin or collocation formulations. We develop a Petrov--Galerkin method for semi-unbounded elliptic problems and a collocation method for fully unbounded, perforated, and time-dependent problems. A conditional bounded-parameter approximation result is proved in a broken Sobolev norm, together with an error decomposition covering approximation, empirical-consistency/quadrature, and least-squares optimization errors. Numerical experiments for Poisson and time-dependent Schr\"odinger equations demonstrate the accuracy and flexibility of the proposed method.

0

math.NA 2026-07-01

EERK time stepping yields rate min(1+γ/2+ρ1(γ)/2,p) on rough parabolic data

by Shuo Yang, Runjie Zhang +2 more

Higher-order exponential Runge-Kutta Galerkin finite element method for semilinear parabolic problems with nonsmooth data

Semigroup and fractional-power bounds replace missing a priori estimates, with experiments confirming the rate is sharp.

abstract click to expand

We develop a rigorous numerical analysis framework for a class of semilinear parabolic problems with nonsmooth initial data. We employ a linear Galerkin finite element method for spatial discretization coupled with a high-order explicit exponential Runge-Kutta (EERK) temporal integration scheme. In contrast to conventional smooth error analysis, the nonsmooth case lacks a priori estimates for the higher-order derivatives of both the nonlinear term and the exact solution. By combining analytic semigroup techniques with fractional power space theory, we establish rigorous bounds for these derivatives. Finally, our analysis proves that the $p$th-order EERK method achieves a convergence rate of $\min(1 + \gamma/2 + \rho_1(\gamma)/2,\:p)$, where $\gamma$ characterizes the initial data regularity and $\rho_1(\gamma)$ quantifies the boundedness of the nonlinearity's first Fr\'echet derivative. Numerical experiments confirm the sharpness of these estimates.

0

math.NA 2026-07-01

AB-BA GKB variants stabilize CT reconstructions with unmatched projectors

by Abdulmajeed Alsubhi

GKB Methods for X-Ray Computed Tomography with an Unmatched Back Projector

Preconditioned bidiagonalization resolves inconsistent equations from mismatched forward and back projectors while resisting semiconvergence

abstract click to expand

In large scale X ray Computed Tomography (CT) inverse problems, the forward and back projectors are often generated using different discretizations. This discrepancy leads to unmatched pairs of projections, resulting in inconsistent normal equations. Consequently, employing the Conjugate Gradient method does not produce a useful solution. For matched operator pairs, the Golub Kahan bidiagonalization (GKB) method provides an efficient solution strategy. It works by projecting the original large-scale problem onto a lower-dimensional subspace, enabling the solution to be computed via a singular value decomposition of a sparse lower bidiagonal matrix. To address unmatched-pair problems in CT, we propose the AB and BA GKB algorithms as preconditioned forms of the GKB. These methods are straightforward to implement and allow for parameter tuning. We provide a discussion on the theoretical computational costs of our proposed algorithms in terms of floating point operations and compare with existing methods. While many Krylov methods tend to amplify noise in solutions, leading to semiconvergence, our proposed algorithms demonstrate greater resilience against this effect. We validate the effectiveness of our approach through numerical examples across various CT problems, showcasing its ability to deliver more stable solutions.

0

math.NA 2026-07-01

GSAV scheme unconditionally stable for perturbed Boussinesq

by M Nader Alhomsi, Jiahong Wu +1 more

A consistent-splitting generalized scalar auxiliary variable scheme for the perturbed Boussinesq system

Consistent splitting reduces each step to decoupled linear solves while delivering second-order accuracy for velocity, pressure and temperat

abstract click to expand

We propose and analyze a second-order consistent-splitting scheme, based on the generalized scalar auxiliary variable (GSAV) approach, for the two-dimensional perturbed Boussinesq system. The system is obtained by subtracting a stable, linearly stratified hydrostatic equilibrium from the standard Boussinesq system. The time discretization extends the consistent-splitting generalized BDF2 framework of Huang and Shen [17] for the Navier-Stokes equations, treating the nonlinear convection and advection together with the linear buoyancy and stratification couplings explicitly, so that each time step reduces to a small number of decoupled linear systems. We prove an unconditional weak stability theorem for the GSAV scheme and derive optimal second-order error estimates for the velocity, pressure, and temperature. A careful tracing reveals that the error constant depends on the inverse viscosity and inverse thermal diffusivity through a quadruply-nested exponential, so the scheme is not robust as either tends to zero. Numerical experiments confirm the second-order convergence and reproduce the expected internal-wave dynamics and exponential relaxation toward hydrostatic balance in a long-time stratified-flow simulation.

0

quant-ph 2026-07-01

Quantum methods cut SPDE pricing samples from O(1/ε²) to O(1/ε)

by Xinmiao Li, Yanqiao Wang +3 more

Quantum Derivative Pricing for SPDEs via BDSDE Representation

QA-MLMC estimators with order-one Forward-Backward Taylor schemes deliver the improved complexity for prices and first- and second-order Gre

abstract click to expand

We study quantum speedups of derivative pricing for stochastic partial differential equation (SPDE) models through their backward doubly stochastic differential equation (BDSDE) representations. We develop conditional and nested quantum-accelerated multilevel Monte Carlo (QA-MLMC) methods for estimating the resulting conditional and nested expectations, improving the sampling complexity of classical Monte Carlo methods from $\widetilde{O}(\epsilon^{-2})$ to $\widetilde{O}(\epsilon^{-1})$ within additive error $\epsilon$. We apply the framework to derivative pricing and sensitivity analysis, providing quantum-accelerated estimators for prices as well as first-order and second-order Greeks, likelihood-ratio and Malliavin-weight representations for Greeks, and Heston-type stochastic-volatility models. To enable efficient multilevel coupling, we construct a family of Forward--Backward Taylor discretization schemes for the stochastic integrals arising in the BDSDE representations and establish global strong-error order one convergence for pricing and Greek estimators. Numerical experiments showcase our schemes for first-order and second-order Greeks can reach the required orders for the full quadratic quantum speedups.

0

math.NA 2026-07-01

Online TT-ALS cuts streaming tensor updates to linear rank cost

by Hiroki Takeda, Yuto Miyatake +1 more

Online TT-ALS for Streaming Tensor Decomposition with Incremental Orthogonalization

Sequential orthogonal constraints keep the objective falling while delivering monotonic decrease and O(I^{n-1} r) complexity.

abstract click to expand

Tensor Train (TT) decomposition is a powerful technique for analyzing high-dimensional data. Existing algorithms for computing TT decompositions can be categorized into two main types: conventional batch-based approaches and recursive online methods. In the context of streaming data, batch methods typically achieve higher reconstruction accuracy but often suffer from memory exhaustion, while online methods provide greater computational efficiency. In this work, we introduce Online TT-ALS (Alternating Least Squares), an algorithm that sequentially enforces orthogonality constraints. This approach allows for efficient and exact updates of the core tensor while maintaining high reconstruction accuracy. Theoretically, we prove that enforcing these orthogonal gauge constraints guarantees monotonic decrease of the local objective function and temporal smoothness. Computationally, our deterministic single-sweep update reduces the rank dependence from quadratic to linear, achieving an overall complexity of $\mathcal{O}(I^{n-1} r)$. Experimental results demonstrate that the proposed method outperforms existing online techniques not only in terms of mathematical approximation accuracy but also in human perception-based video quality metrics. Furthermore, compared to recent deep learning-based paradigms, our algebraic approach achieves speedups of several orders of magnitude. Consequently, our method exhibits high computational efficiency and is suitable for low-latency real-time processing applications.

0

cond-mat.soft 2026-07-01

Three-species model links micellar relaxation to shear bands

by Rongxin Lu, Jiwei Jia +1 more

Beyond binary scission: a generalized three-species cascade breakage model for wormlike micellar solutions

Cascade breakage connects short-fragment relaxation to gel-network recovery, producing three-mode response and stress plateau.

abstract click to expand

Wormlike micellar fluids exhibit complex rheological behavior driven by the continuous breakage and recombination of self-assembled micellar networks. Existing two-species models provide a coarse binary representation of the micellar population, limiting their ability to resolve intermediate structural states and broad relaxation spectra. To address this limitation, we develop a three-species cascade breakage model consisting of gel-network, long chains, and short chains. By introducing an intermediate micellar state, the model links the rapid relaxation of short fragments to the slow recovery of the gel-network within a unified kinetic framework. This additional structural pathway gives rise to a three-mode viscoelastic response, improves the high-frequency description of the dynamic moduli, and produces a non-monotone constitutive curve that evolves into a stress plateau with coexisting shear bands in Couette flow. This cascade mechanism also governs the transient response, including stress overshoot, hysteresis, and multistep relaxation after shear cessation. Overall, the proposed three-species model provides a physically interpretable framework for worm-like micellar shear banding, capturing the connection between cascade microstructural evolution, broad relaxation dynamics, and macroscopic flow localization.

0