cs.NA

arXiv:2503.22646v1 Announce Type: new Abstract: Simulation-based approaches are among the most practical means to search for safety violations, bugs, and other unexpected events in cyber-physical systems (CPS). Where existing approaches search for simulations violating a formal specification or maximizing a notion of coverage, in this work we propose a new goal for testing: to discover unknown rare behaviors by examining discrete mode sequences. We assume a CPS simulator outputs mode information, and strive to explore the sequences of modes produced by varying the initial state or time-varying uncertainties. We hypothesize that rare mode sequences are often the most interesting to a designer, and we develop two accelerated sampling algorithms that speed up the process of finding such sequences. We evaluate our approach on several benchmarks, ranging from synthetic examples to Simulink diagrams of a CPS, demonstrating in some cases a speedup of over 100x compared with a random sampling strategy.

arXiv:2503.20742v2 Announce Type: replace Abstract: Despite the many challenges in exploratory data analysis, artificial neural networks have motivated strong interests in scientists and researchers both in theoretical as well as practical applications. Among sources of such popularity of artificial neural networks the ability of modeling non-linear dynamical systems, generalization, and adaptation possibilities should be mentioned. Despite this, there is still significant debate about the role of various underlying stochastic processes in stabilizing a unique structure for data learning and prediction. One of such obstacles to the theoretical and numerical study of machine intelligent systems is the curse of dimensionality and the sampling from high-dimensional probability distributions. In general, this curse prevents efficient description of states, providing a significant complexity barrier for the system to be efficiently described and studied. In this strand of research, direct treatment and description of such abstract notions of learning theory in terms of quantum information be one of the most favorable candidates. Hence, the subject matter of these articles is devoted to problems of design, adaptation and the formulations of computationally hard problems in terms of quantum mechanical systems. In order to characterize the microscopic description of such dynamics in the language of inferential statistics, covariance matrix estimation of d-dimensional Gaussian densities and Bayesian interpretation of eigenvalue problem for dynamical systems is assessed.

arXiv:2503.22360v1 Announce Type: new Abstract: F. Stenger proposed efficient approximation formulas for derivatives over infinite intervals. Those formulas were derived by the combination of the Sinc approximation and appropriate conformal maps. It has been shown that those formulas can attain root-exponential convergence. In this study, we enhance the convergence rate by improving the conformal maps employed in those formulas. We provide theoretical error analysis and numerical experiments that confirm the effectiveness of our new formulas.

arXiv:2503.22631v1 Announce Type: new Abstract: Many scientific applications require the evaluation of the action of the matrix function over a vector and the most common methods for this task are those based on the Krylov subspace. Since the orthogonalization cost and memory requirement can quickly become overwhelming as the basis grows, the Krylov method is often restarted after a few iterations. This paper proposes a new acceleration technique for restarted Krylov methods based on randomization. The numerical experiments show that the randomized method greatly outperforms the classical approach with the same level of accuracy. In fact, randomization can actually improve the convergence rate of restarted methods in some cases. The paper also compares the performance and stability of the randomized methods proposed so far for solving very large finite element problems, complementing the numerical analyses from previous studies.

arXiv:2407.02009v3 Announce Type: replace Abstract: We theoretically explore boundary conditions for lattice Boltzmann methods, focusing on a toy two-velocities scheme to tackle a linear one-dimensional advection equation. By mapping lattice Boltzmann schemes to Finite Difference schemes, we facilitate rigorous consistency and stability analyses. We develop kinetic boundary conditions for inflows and outflows, highlighting the trade-off between accuracy and stability, which we successfully overcome. Consistency analysis relies on modified equations, whereas stability is assessed using GKS (Gustafsson, Kreiss, and Sundstr{\"o}m) theory and -- when this approach fails on coarse meshes -- spectral and pseudo-spectral analyses of the scheme's matrix that explain effects germane to low resolutions.

arXiv:2503.17807v2 Announce Type: replace Abstract: In this paper we consider a new probability sampling methods based on Langevin diffusion dynamics to resolve the problem of existing Monte Carlo algorithms when draw samples from high dimensional target densities. We extent Metropolis-Adjusted Langevin Diffusion algorithm by modelling the stochasticity of precondition matrix as a random matrix. An advantage compared to other proposal method is that it only requires the gradient of log-posterior. The proposed method provides fully adaptation mechanisms to tune proposal densities to exploits and adapts the geometry of local structures of statistical models. We clarify the benefits of the new proposal by modelling a Quantum Probability Density Functions of a free particle in a plane (energy Eigen-functions). The proposed model represents a remarkable improvement in terms of performance accuracy and computational time over standard MCMC method.

arXiv:2503.22286v1 Announce Type: new Abstract: This paper presents some theoretical results relating the Bregman log determinant matrix divergence to Kaporin's condition number. These can be viewed as nearness measures between a preconditioner and a given matrix, and we show under which conditions these two functions coincide. We also give examples of constraint sets over which it is equivalent to minimise these two objectives. We focus on preconditioners that are the sum of a positive definite and low-rank matrix, which were developed in a previous work. These were constructed as minimisers of the aforementioned divergence, and we show that they are only a constant scaling from also minimising Kaporin's condition number. We highlight connections to information geometry and comment on future directions.

arXiv:2503.22335v1 Announce Type: new Abstract: The Hilfer fractional derivative interpolates the commonly used Riemann-Liouville and Caputo fractional derivative. In general, solutions to Hilfer fractional differential equations are singular for $t \downarrow 0$ and are difficult to approximate with high numerical accuracy. We propose a numerical Bernstein splines technique to approximate solutions to generalized nonlinear initial values problems with Hilfer fractional derivatives. Convergent approximations are obtained using an efficient vectorized solution setup with few convergence requirements for a wide range of nonlinear fractional differential equations. We demonstrate efficiency of the developed method by applying it to the fractional Van der Pol oscillator, a system with applications in control systems and electronic circuits.

arXiv:2503.22455v1 Announce Type: new Abstract: This work presents a multigrid preconditioned high order immersed finite difference solver to accurately and efficiently solve the Poisson equation on complex 2D and 3D domains. The solver employs a low order Shortley-Weller multigrid method to precondition a high order matrix-free Krylov subspace solver. The matrix-free approach enables full compatibility with high order IIM discretizations of boundary and interface conditions, as well as high order wavelet-adapted multiresolution grids. Through verification and analysis on 2D domains, we demonstrate the ability of the algorithm to provide high order accurate results to Laplace and Poisson problems with Dirichlet, Neumann, and/or interface jump boundary conditions, all effectively preconditioned using the multigrid method. We further show that the proposed method is able to efficiently solve high order discretizations of Laplace and Poisson problems on complex 3D domains using thousands of compute cores and on multiresolution grids. To our knowledge, this work presents the largest problem sizes tackled with high order immersed methods applied to elliptic partial differential equations, and the first high order results on 3D multiresolution adaptive grids. Together, this work paves the way for employing high order immersed methods to a variety of 3D partial differential equations with boundary or inter-face conditions, including linear and non-linear elasticity problems, the incompressible Navier-Stokes equations, and fluid-structure interactions.

arXiv:2503.22621v1 Announce Type: new Abstract: We prove optimal convergence rates for certain low-regularity integrators applied to the one-dimensional periodic nonlinear Schr\"odinger and wave equations under the assumption of $H^1$ solutions. For the Schr\"odinger equation we analyze the exponential-type scheme proposed by Ostermann and Schratz in 2018, whereas in the wave case we treat the corrected Lie splitting proposed by Li, Schratz, and Zivcovich in 2023. We show that the integrators converge with their full order of one and two, respectively. In this situation only fractional convergence rates were previously known. The crucial ingredients in the proofs are known space-time bounds for the solutions to the corresponding linear problems. More precisely, in the Schr\"odinger case we use the $L^4$ Strichartz inequality, and for the wave equation a null form estimate. To our knowledge, this is the first time that a null form estimate is exploited in numerical analysis. We apply the estimates for continuous time, thus avoiding potential losses resulting from discrete-time estimates.

arXiv:2503.22652v1 Announce Type: cross Abstract: Chebyshev Filtered Subspace Iteration (ChFSI) has been widely adopted for computing a small subset of extreme eigenvalues in large sparse matrices. This work introduces a residual-based reformulation of ChFSI, referred to as R-ChFSI, designed to accommodate inexact matrix-vector products while maintaining robust convergence properties. By reformulating the traditional Chebyshev recurrence to operate on residuals rather than eigenvector estimates, the R-ChFSI approach effectively suppresses the errors made in matrix-vector products, improving the convergence behaviour for both standard and generalized eigenproblems. This ability of R-ChFSI to be tolerant to inexact matrix-vector products allows one to incorporate approximate inverses for large-scale generalized eigenproblems, making the method particularly attractive where exact matrix factorizations or iterative methods become computationally expensive for evaluating inverses. It also allows us to compute the matrix-vector products in lower-precision arithmetic allowing us to leverage modern hardware accelerators. Through extensive benchmarking, we demonstrate that R-ChFSI achieves desired residual tolerances while leveraging low-precision arithmetic. For problems with millions of degrees of freedom and thousands of eigenvalues, R-ChFSI attains final residual norms in the range of 10$^{-12}$ to 10$^{-14}$, even with FP32 and TF32 arithmetic, significantly outperforming standard ChFSI in similar settings. In generalized eigenproblems, where approximate inverses are used, R-ChFSI achieves residual tolerances up to ten orders of magnitude lower, demonstrating its robustness to approximation errors. Finally, R-ChFSI provides a scalable and computationally efficient alternative for solving large-scale eigenproblems in high-performance computing environments.

arXiv:2310.15457v3 Announce Type: replace Abstract: In this work, we present an iteratively decoupled algorithm for solving the quasi-static multiple-network poroelastic model. Our approach employs a total-pressure-based formulation with solid displacement, total pressure, and network pressures as primary unknowns. This reformulation decomposes the original problem into a generalized Stokes problem and a parabolic problem, offering key advantages such as reduced elastic locking effects and simplified discretization. The algorithm guarantees unconditional convergence to the solution of the fully coupled system. Numerical experiments demonstrate the accuracy, efficiency, and robustness of the method with respect to physical parameters and discretization. We further apply the algorithm to simulate brain flow dynamics, showcasing its practical utility in biomechanical modeling.

arXiv:2409.08734v3 Announce Type: replace Abstract: The blind image deconvolution is a challenging, highly ill-posed nonlinear inverse problem. We introduce a Multiscale Hierarchical Decomposition Method (MHDM) that is iteratively solving variational problems with adaptive data and regularization parameters, towards obtaining finer and finer details of the unknown kernel and image. We establish convergence of the residual in the noise-free data case, and then in the noisy data case when the algorithm is stopped early by means of a discrepancy principle. Fractional Sobolev norms are employed as regularizers for both kernel and image, with the advantage of computing the minimizers explicitly in a pointwise manner. In order to break the notorious symmetry occurring during each minimization step, we enforce a positivity constraint on the Fourier transform of the kernels. Numerical comparisons with a single-step variational method and a non-blind MHDM show that our approach produces comparable results, while less laborious parameter tuning is necessary at the price of more computations. Additionally, the scale decomposition of both reconstructed kernel and image provides a meaningful interpretation of the involved iteration steps.

arXiv:2409.16796v2 Announce Type: replace Abstract: As computational machines become larger and more complex, the probability of hardware failure rises. ``Silent errors'', or bit flips, may not be immediately apparent but can cause detrimental effects to algorithm behavior. In this work, we examine an algorithm-based approach to silent error detection in the context of pipelined Krylov subspace methods, in particular, Pipe-PR-CG, for the solution of linear systems. Our approach is based on using finite precision error analysis to bound the differences between quantities which should be equal in exact arithmetic. By monitoring select quantities during the iteration, we can detect when these bounds are violated, which indicates that a silent error has occurred. We use this approach to develop a fault-tolerant variant and also suggest a strategy for dynamically adapting the detection criteria. Our numerical experiments demonstrate the effectiveness of our approach.

arXiv:2503.21855v1 Announce Type: new Abstract: We present a new class of numerical methods for solving stochastic differential equations with additive noise on general Riemannian manifolds with high weak order of accuracy. In opposition to the popular approach with projection methods, the proposed methods are intrinsic: they only rely on geometric operations and avoid coordinates and embeddings. We provide a robust and general convergence analysis and an algebraic formalism of exotic planar Butcher series for the computation of order conditions at any high order. To illustrate the methodology, an explicit method of second weak order is introduced, and several numerical experiments confirm the theoretical findings and extend the approach for the sampling of the invariant measure of Riemannian Langevin dynamics.

arXiv:2503.22273v1 Announce Type: new Abstract: The Gauss-Seidel method has been used for more than 100 years as the standard method for the solution of linear systems of equations under certain restrictions. This method, as well as Cramer and Jacobi, is widely used in education and engineering, but there is a theoretical gap when we want to solve less restricted systems, or even non-square or non-exact systems of equation. Here, the solution goes through the use of numerical systems, such as the minimization theories or the Moore-Penrose pseudoinverse. In this paper we fill this gap with a global analytical iterative formulation that is capable to reach the solutions obtained with the Moore-Penrose pseudoinverse and the minimization methodologies, but that analytically lies to the solutions of Gauss-Seidel, Jacobi, or Cramer when the system is simplified.

arXiv:2503.22297v1 Announce Type: new Abstract: We consider second-order PDE problems set in unbounded domains and discretized by Lagrange finite elements on a finite mesh, thus introducing an artificial boundary in the discretization. Specifically, we consider the reaction diffusion equation as well as Helmholtz problems in waveguides with perfectly matched layers. The usual procedure to deal with such problems is to first consider a modeling error due to the introduction of the artificial boundary, and estimate the remaining discretization error with a standard a posteriori technique. A shortcoming of this method, however, is that it is typically hard to obtain sharp bounds on the modeling error. In this work, we propose a new technique that allows to control the whole error by an a posteriori error estimator. Specifically, we propose a flux-equilibrated estimator that is slightly modified to handle the truncation boundary. For the reaction diffusion equation, we obtain fully-computable guaranteed error bounds, and the estimator is locally efficient and polynomial-degree-robust provided that the elements touching the truncation boundary are not too refined. This last condition may be seen as an extension of the notion of shape-regularity of the mesh, and does not prevent the design of efficient adaptive algorithms. For the Helmholtz problem, as usual, these statements remain valid if the mesh is sufficiently refined. Our theoretical findings are completed with numerical examples which indicate that the estimator is suited to drive optimal adaptive mesh refinements.

arXiv:2503.22301v1 Announce Type: new Abstract: In the present paper, we introduce three neural network operators of convolution type activated by symmetrized, deformed and parametrized B-generalized logistic function. We deal with the approximation properties of these operators to the identity by using modulus of continuity. Furthermore, we show that our operators preserve global smoothness and consider the iterated versions of them. Here, we find it is worthy to mention that these operators play important roles in neural network approximation since most of the basic network models are activated by logistic functions.

arXiv:2503.22386v1 Announce Type: new Abstract: The study of parametric differential equations plays a crucial role in weather forecasting and epidemiological modeling. These phenomena are better represented using fractional derivatives due to their inherent memory or hereditary effects. This paper introduces a novel scientific machine learning approach for solving parametric time-fractional differential equations by combining traditional spectral methods with neural networks. Instead of relying on automatic differentiation techniques, commonly used in traditional Physics-Informed Neural Networks (PINNs), we propose a more efficient global discretization method based on Legendre polynomials. This approach eliminates the need to simulate the parametric fractional differential equations across multiple parameter values. By applying the Legendre-Galerkin weak formulation to the differential equation, we construct a loss function for training the neural network. The trial solutions are represented as linear combinations of Legendre polynomials, with the coefficients learned by the neural network. The convergence of this method is theoretically established, and the theoretical results are validated through numerical experiments on several well-known differential equations.

arXiv:2301.06765v4 Announce Type: replace-cross Abstract: The Chernoff approximation method is a powerful and flexible tool of functional analysis, which allows in many cases to express exp(tL) in terms of variable coefficients of a linear differential operator L. In this paper, we prove a theorem that allows us to apply this method to find the resolvent of L. Our theorem states that the Laplace transforms of Chernoff approximations of a $C_0$-semigroup converge to the resolvent of the generator of this semigroup. We demonstrate the proposed method on a second-order differential operator with variable coefficients. As a consequence, we obtain a new representation of the solution of a nonhomogeneous linear ordinary differential equation of the second order in terms of functions that are coefficients of this equation, playing the role of parameters of the problem. For the Chernoff function, based on the shift operator, we give an estimate for the rate of convergence of approximations to the solution.