Search

Surface Stokes and Navier--Stokes equations are used to model fluid flow on surfaces. They have attracted significant recent attention in the numerical analysis literature because approximation of their solutions poses significant challenges not encountered in the Euclidean context. One challenge comes from the need to simultaneously enforce tangentiality and H1 conformity (continuity) of discrete vector fields used to approximate solutions in the velocity-pressure formulation. Existing methods in the literature all enforce one of these two constraints weakly either by penalization or by use of Lagrange multipliers. Missing so far is a robust and systematic construction of surface Stokes finite element spaces which employ nodal degrees of freedom, including MINI, Taylor--Hood, Scott--Vogelius, and other composite elements which can lead to divergence-conforming or pressure-robust discretizations. In this paper we construct surface MINI spaces whose velocity fields are tangential. They are not H1-conforming, but do lie in H(div) and do not require penalization to achieve optimal convergence rates. We prove stability and optimal-order energy-norm convergence of the method and demonstrate optimal-order convergence of the velocity field in L2 via numerical experiments. The core advance in the paper is the construction of nodal degrees of freedom for the velocity field. This technique also may be used to construct surface counterparts to many other standard Euclidean Stokes spaces, and we accordingly present numerical experiments indicating optimal-order convergence of nonconforming tangential surface Taylor--Hood P2 P1 elements. [ABSTRACT FROM AUTHOR]

In this paper, we consider an approach to improve the performance of exponential Runge–Kutta integrators and Lawson schemes in cases where the solution of a related, but usually much simpler, problem can be computed efficiently. While for implicit methods such an approach is common (e.g., by using preconditioners), for exponential integrators this has proven more challenging. Here we propose to extract a constant coefficient differential operator from the semilinear advection-diffusion-reaction equation for which, in many situations, efficient methods are known to compute the required matrix functions. Both a linear stability analysis and extensive numerical experiments show that the resulting schemes can be unconditionally stable. In fact, we find that exponential integrators of Runge–Kutta type and Lawson schemes can have better stability properties than similarly constructed implicit-explicit schemes. We also derive two new Lawson-type integrators that further improve on these stability properties. The overall effectiveness of the approach is highlighted by a number of performance comparisons on examples in two and three space dimensions. [ABSTRACT FROM AUTHOR]

Partial differential equations (PDEs) with uncertain or random inputs have been considered in many studies of uncertainty quantification. In forward uncertainty quantification, one is interested in analyzing the stochastic response of the PDE subject to input uncertainty, which usually involves solving high-dimensional integrals of the PDE output over a sequence of stochastic variables. In practical computations, one typically needs to discretize the problem in several ways: approximating an infinite-dimensional input random field with a finite-dimensional random field, spatial discretization of the PDE using, e.g., finite elements, and approximating high-dimensional integrals using cubatures such as quasi--Monte Carlo methods. In this paper, we focus on the error resulting from dimension truncation of an input random field. We show how Taylor series can be used to derive theoretical dimension truncation rates for a wide class of problems and we provide a simple checklist of conditions that a parametric mathematical model needs to satisfy in order for our dimension truncation error bound to hold. Some of the novel features of our approach include that our results are applicable to nonaffine parametric operator equations, dimensionally truncated conforming finite element discretized solutions of parametric PDEs, and even compositions of PDE solutions with smooth nonlinear quantities of interest. As a specific application of our method, we derive an improved dimension truncation error bound for elliptic PDEs with lognormally parameterized diffusion coefficients. Numerical examples support our theoretical findings. [ABSTRACT FROM AUTHOR]

The aim of this paper is to carry out convergence analysis and algorithm implementation of a novel sample-wise backpropagation method for training a class of stochastic neural networks (SNNs). The preliminary discussion on such an SNN framework was first introduced in [Archibald et al., Discrete Contin. Dyn. Syst. Ser. S, 15 (2022), pp. 2807--2835]. The structure of the SNN is formulated as a discretization of a stochastic differential equation (SDE). A stochastic optimal control framework is introduced to model the training procedure, and a sample-wise approximation scheme for the adjoint backward SDE is applied to improve the efficiency of the stochastic optimal control solver, which is equivalent to the backpropagation for training the SNN. The convergence analysis is derived by introducing a novel joint conditional expectation for the gradient process. Under the convexity assumption, our result indicates that the number of SNN training steps should be proportional to the square of the number of layers in the convex optimization case. In the implementation of the sample-based SNN algorithm with the benchmark MNIST dataset, we adopt the convolution neural network (CNN) architecture and demonstrate that our sample-based SNN algorithm is more robust than the conventional CNN. [ABSTRACT FROM AUTHOR]

Convection-diffusion-reaction equations model the conservation of scalar quantities. From the analytic point of view, solutions of these equations satisfy, under certain conditions, maximum principles, which represent physical bounds of the solution. That the same bounds are respected by numerical approximations of the solution is often of utmost importance in practice. The mathematical formulation of this property, which contributes to the physical consistency of a method, is called the discrete maximum principle (DMP). In many applications, convection dominates diffusion by several orders of magnitude. It is well known that standard discretizations typically do not satisfy the DMP in this convectiondominated regime. In fact, in this case it turns out to be a challenging problem to construct discretizations that, on the one hand, respect the DMP and, on the other hand, compute accurate solutions. This paper presents a survey on finite element methods, with the main focus on the convection-dominated regime, that satisfy a local or a global DMP. The concepts of the underlying numerical analysis are discussed. The survey reveals that for the steady-state problem there are only a few discretizations, all of them nonlinear, that at the same time both satisfy the DMP and compute reasonably accurate solutions, e.g., algebraically stabilized schemes. Moreover, most of these discretizations have been developed in recent years, showing the enormous progress that has been achieved lately. Similarly, methods based on algebraic stabilization, both nonlinear and linear, are currently the only finite element methods that combine the satisfaction of the global DMP and accurate numerical results for the evolutionary equations in the convection-dominated scenario. [ABSTRACT FROM AUTHOR]

This paper proposes a sketching strategy with spherical designs to equip the classical spherical basis function (SBF) approach for massive spherical data fitting. We conduct theoretical analysis and numerical verifications to demonstrate the feasibility of the proposed sketching strategy. From the theoretical side, we prove that sketching based on spherical designs can reduce the computational burden of the SBF approach without sacrificing its approximation capability. In particular, we provide upper and lower bounds for the proposed sketching strategy to fit noisy data on spheres. From the experimental side, we numerically illustrate the feasibility of the sketching strategy by showing its comparable fitting performance with the SBF approach. These interesting findings show that the proposed sketching strategy is effective in fitting massive and noisy data on spheres. [ABSTRACT FROM AUTHOR]

The noisy leaky integrate-and-fire (NLIF) model describes the voltage configurations of neuron networks with an interacting many-particles system at a microscopic level. When simulating neuron networks of large sizes, computing a coarse-grained mean-field Fokker--Planck equation solving the voltage densities of the networks at a macroscopic level practically serves as a feasible alternative in its high efficiency and credible accuracy when the interaction within the network remains relatively low. However, the macroscopic model fails to yield valid results of the networks when simulating considerably synchronous networks with active firing events. In this paper, we propose a multiscale solver for the NLIF networks, inheriting the macroscopic solver's low cost and the microscopic solver's high reliability. For each temporal step, the multiscale solver uses the macroscopic solver when the firing rate of the simulated network is low, while it switches to the microscopic solver when the firing rate tends to blow up. Moreover, the macroscopic and microscopic solvers are integrated with a high-precision switching algorithm to ensure the accuracy of the multiscale solver. The validity of the multiscale solver is analyzed from two perspectives: first, we provide practically sufficient conditions that guarantee the mean-field approximation of the macroscopic model and present rigorous numerical analysis on simulation errors when coupling the two solvers; second, the numerical performance of the multiscale solver is validated through simulating several large neuron networks, including networks with either instantaneous or periodic input currents which prompt active firing events over some time. [ABSTRACT FROM AUTHOR]

In this paper, we analyze the celebrated Haken--Kelso--Bunz model, describing the dynamics of bimanual coordination, in the presence of delay. We study the linear dynamics, stability, nonlinear behavior, and bifurcations of this model by both theoretical and numerical analysis. We calculate in-phase and antiphase limit cycles as well as quasi-periodic solutions via double Hopf bifurcation analysis and center manifold reduction. Moreover, we uncover further details on the global dynamic behavior by numerical continuation, including the occurrence of limit cycles in phase quadrature and 1-1 locking of quasi-periodic solutions. [ABSTRACT FROM AUTHOR]

Quasiperiodic systems are important space-filling ordered structures, without decay and translational invariance. How to solve quasiperiodic systems accurately and efficiently is a great challenge. A useful approach, the projection method (PM) [J. Comput. Phys., 256 (2014), pp. 428--440], has been proposed to compute quasiperiodic systems. Various studies have demonstrated that the PM is an accurate and efficient method to solve quasiperiodic systems. However, there is a lack of theoretical analysis of the PM. In this paper, we present a rigorous convergence analysis of the PM by establishing a mathematical framework of quasiperiodic functions and their high-dimensional periodic functions. We also give a theoretical analysis of the quasiperiodic spectral method (QSM) based on this framework. Results demonstrate that the PM and QSM both have exponential decay, and the QSM (PM) is a generalization of the periodic Fourier spectral (pseudospectral) method. Then, we analyze the computational complexity of the PM and QSM in calculating quasiperiodic systems. The PM can use a fast Fourier transform, while the QSM cannot. Moreover, we investigate the accuracy and efficiency of the PM, QSM, and periodic approximation method in solving the linear time-dependent quasiperiodic Schr\"odinger equation. [ABSTRACT FROM AUTHOR]

In this paper, a coupled electro-calcium model was developed and implemented to computationally explore the effects of neuronal synapse loss, in particular in the context of Alzheimer's disease. Established parameters affected by Alzheimer's disease, such as synapse loss, calcium leaks at deteriorating synaptic contacts, and downregulation of the calcium buffer calbindin, are subject to this study. Reconstructed neurons are used to define the computational domain for a system of PDEs and ODEs, discretized by finite differences and solved with a semi-implicit second-order time integrator. The results show neuronal resilience during synapse loss. When incorporating calcium leaks at affected synapses, neurons lose their ability to produce synapse-tonucleus calcium signals, necessary for learning, plasticity, and neuronal survival. Downregulation of calbindin concentrations partially recovers the signaling pathway to the cell nucleus. These results could define future research pathways toward stabilizing the calcium signaling pathways during Alzheimer's disease. The coupled electro-calcium model was implemented and solved using MATLAB https://github.com/NeuroBox3D/CalcSim. [ABSTRACT FROM AUTHOR]

In this paper we consider the distributed optimal control of the Navier—Stokes equations in the presence of pointwise mixed control-state constraints. After deriving a first order necessary condition, the regularity of the mixed constraint multiplier is investigated. Second order sufficient optimality conditions are studied as well. In the last part of the paper, a semismooth Newton method is applied for the numerical solution of the control problem. The convergence of the method is proved and numerical experiments are carried out. [ABSTRACT FROM AUTHOR]

Stochastic algebraic Riccati equations, also known as rational algebraic Riccati equations, arising in linear-quadratic optimal control for stochastic linear time-invariant systems, were considered to be not easy to solve. The state-of-the-art numerical methods mostly rely on differentiability or continuity, such as the Newton-type method, the LMI method, or the homotopy method. In this paper, we will build a novel theoretical framework and reveal the intrinsic algebraic structure appearing in this kind of algebraic Riccati equation. This structure guarantees that to solve them is almost as easy as to solve deterministic/classical ones, which will shed light on the theoretical analysis and numerical algorithm design for this topic. [ABSTRACT FROM AUTHOR]

In this paper, we investigate the energy minimization model arising in the ensemble Kohn--Sham density functional theory for metallic systems, in which a pseudo-eigenvalue matrix and a general smearing approach are involved. We study the invariance of the energy functional and the existence of the minimizer of the ensemble Kohn--Sham model. We propose an adaptive twoparameter step size strategy and the corresponding preconditioned conjugate gradient methods to solve the energy minimization model. Under some mild but reasonable assumptions, we prove the global convergence for the gradients of the energy functional produced by our algorithms. Numerical experiments show that our algorithms are efficient, especially for large scale metallic systems. In particular, our algorithms produce convergent numerical approximations for some metallic systems, for which the traditional self-consistent field iterations fail to converge. [ABSTRACT FROM AUTHOR]

Zeros of rational transfer function matrices R(λ) are the eigenvalues of associated polynomial system matrices P(λ) under minimality conditions. In this paper, we define a structured condition number for a simple eigenvalue λ0 of a (locally) minimal polynomial system matrix P(λ), which in turn is a simple zero λ0 of its transfer function matrix R(λ). Since any rational matrix can be written as the transfer function of a polynomial system matrix, our analysis yields a structured perturbation theory for simple zeros of rational matrices R(λ). To capture all the zeros of R(λ), regardless of whether they are poles, we consider the notion of root vectors. As corollaries of the main results, we pay particular attention to the special case of λ 0 being not a pole of R(λ) since in this case the results get simpler and can be useful in practice. We also compare our structured condition number with Tisseur's unstructured condition number for eigenvalues of matrix polynomials and show that the latter can be unboundedly larger. Finally, we corroborate our analysis by numerical experiments. [ABSTRACT FROM AUTHOR]

Over the last 10 years, results from [J. M. Melenk and S. Sauter, Math. Comp., 79 (2010), pp. 1871--1914], [J. M. Melenk and S. Sauter, SIAM J. Numer. Anal., 49 (2011), pp. 1210--1243], [S. Esterhazy and J. M. Melenk, Numerical Analysis of Multiscale Problems, Springer, New York, 2012, pp. 285--324] and [J. M. Melenk, A. Parsania, and S. Sauter, J. Sci. Comput., 57 (2013), pp. 536--581] decomposing high-frequency Helmholtz solutions into "low-" and "high-" frequency components have had a large impact in the numerical analysis of the Helmholtz equation. These results have been proved for the constant-coefficient Helmholtz equation in either the exterior of a Dirichlet obstacle or an interior domain with an impedance boundary condition. Using the Helffer--Sjöstrand functional calculus [B. Helffer and J. Sjöstrand, Schrödinger Operators, Springer, Berlin, 1989, pp. 118--197] this paper proves analogous decompositions for scattering problems fitting into the black-box scattering framework of Sjöstrand and Zworski [J. Amer. Math. Soc., 4 (1991), pp. 729--769] thus covering Helmholtz problems with variable coefficients, impenetrable obstacles, and penetrable obstacles all at once. These results allow us to prove new frequency-explicit convergence results for (i) the hp-finite-element method (hp-FEM) applied to the variable-coefficient Helmholtz equation in the exterior of an analytic Dirichlet obstacle, where the coefficients are analytic in a neighborhood of the obstacle, and (ii) the h-FEM applied to the Helmholtz penetrable-obstacle transmission problem. In particular, the result in (i) shows that the hp-FEM applied to this problem does not suffer from the pollution effect. [ABSTRACT FROM AUTHOR]

In this paper, we study the problem of initial data identification for weak-entropy solutions of the one-dimensional Burgers equation. This problem consists in identifying the set of initial data evolving to a given target at a final time. Due to the time-irreversibility of the Burgers equation, some target functions are unattainable from solutions of this equation, making the identification problem under consideration ill-posed. To get around this issue, we introduce a nonsmooth optimization problem, which consists in minimizing the difference between the predictions of the Burgers equation and the observations of the system at a final time in L² (ℝ) norm. Here, we characterize the set of minimizers of the aforementioned nonsmooth optimization problem. One of the minimizers is the backward entropy solution, constructed using a backward-forward method. Some simulations are given using a wave-front tracking algorithm. [ABSTRACT FROM AUTHOR]

The scattering of linear waves by periodic structures is a crucial phenomena in many branches of applied physics and engineering. In this paper we establish rigorous analytic results necessary for the proper numerical analysis of a class of high-order perturbation of surfaces/asymptotic waveform evaluation (HOPS/AWE) methods for numerically simulating scattering returns from periodic diffraction gratings. More specifically, we prove a theorem on existence and uniqueness of solutions to a system of partial differential equations which model the interaction of linear waves with a periodic two-layer structure. Furthermore, we establish joint analyticity of these solutions with respect to both geometry and frequency perturbations. This result provides hypotheses under which a rigorous numerical analysis could be conducted on our recently developed HOPS/AWE algorithm. [ABSTRACT FROM AUTHOR]

This paper proposes and examines a finite element method (FEM) for a Maxwell obstacle problem in electric shielding. The model is given by a coupled system comprising the Faraday equation and an evolutionary variational inequality (VI) of Ampère-Maxwell-type. Based on the leapfrog (Yee) time-stepping and the Nédélec edge elements, we set up a fully discrete FEM where the obstacle is discretized in such a way that no additional nonlinear solver is required for the computation of the discrete VI. While the L²-stability is achieved for the discrete solutions and the associated difference quotients, the scheme only guarantees the L¹-stability for the discrete magnetic curl field in the obstacle region. The lack of the global L²-stability for the magnetic curl field is justified by the low regularity issue in Maxwell obstacle problems and turns to be the main challenge in the convergence analysis. Our convergence proof consists of two main stages. First, exploiting the L¹-stability in the obstacle region, we derive a convergence result towards a weaker system involving smooth feasible test functions. In the second step, we recover the original system by enlarging the feasible test function set through a specific constraint preserving mollification process in the spirit of Ern and Guermond [Comput. Methods Appl. Math., 16 (2016), pp. 51-75]. This paper is closed by three-dimensional numerical results of the proposed FEM confirming the theoretical convergence result and, in particular, the Faraday shielding effect. [ABSTRACT FROM AUTHOR]

Saddle dynamics is a time continuous dynamics to efficiently compute the any-index saddle points and construct the solution landscape. In practice, the saddle dynamics needs to be discretized for numerical computations, while the corresponding numerical analyses are rarely studied in the literature, especially for the high-index cases. In this paper we propose the convergence analysis of discrete high-index saddle dynamics. To be specific, we prove the local linear convergence rates of numerical schemes of high-index saddle dynamics, which indicates that the local curvature in the neighborhood of the saddle point and the accuracy of computing the eigenfunctions are main factors that affect the convergence of discrete saddle dynamics. The proved results serve as compensations for the convergence analysis of high-index saddle dynamics and are substantiated by numerical experiments. [ABSTRACT FROM AUTHOR]

The goal of this paper is to relate the capacitance matrix formalism to the tightbinding approximation. By doing so, we open the way to the use of mathematical techniques and tools from condensed matter theory in the mathematical and numerical analysis of metamaterials, in particular for the understanding of their topological properties. We first study how the capacitance matrix formalism, both when the material parameters are static and modulated, can be posed in a Hamiltonian form. Then, we use this result to compare this formalism to the tight-binding approximation. We prove that the correspondence between the capacitance formulation and the tight-binding approximation holds only in the case of dilute resonators. On the other hand, the tight-binding model is often coupled with a nearest-neighbor approximation, whereby long-range interactions are neglected. Even in the dilute case, we show that long-range interactions between subwavelength resonators are relatively strong and nearest-neighbor approximations are not generally appropriate. [ABSTRACT FROM AUTHOR]

The aim of this paper is to discuss the mathematical modeling of Brownian active particle systems, a recently popular paradigmatic system for self-propelled particles. We present four microscopic models with different types of repulsive interactions between particles and their associated macroscopic models, which are formally obtained using different coarse-graining methods. The macroscopic limits are integro-differential equations for the density in phase space (positions and orientations) of the particles and may include nonlinearities in both the diffusive and advective components. In contrast to passive particles, systems of active particles can undergo phase separation without any attractive interactions, a mechanism known as motility-induced phase separation (MIPS). We explore the onset of such a transition for each model in the parameter space of occupied volume fraction and P'eclet number via a linear stability analysis and numerical simulations at both the microscopic and macroscopic levels. We establish that one of the models, namely, the mean-field model which assumes long-range repulsive interactions, cannot explain the emergence of MIPS. In contrast, MIPS is observed for the remaining three models that assume short-range interactions that localize the interaction terms in space. [ABSTRACT FROM AUTHOR]

In this paper, we study the propagation speeds of reaction-diffusion-advection fronts in time-periodic cellular and chaotic flows with Kolmogorov-Petrovsky-Piskunov (KPP) nonlinearity. We first apply the variational principle to reduce the computation of KPP front speeds to a principal eigenvalue problem of a linear advection-diffusion operator with space-time periodic coefficient on a periodic domain. To this end, we develop efficient Lagrangian particle methods to compute the principal eigenvalue through the Feynman-Kac formula. By estimating the convergence rate of Feynman-Kac semigroups and the operator splitting method for approximating the linear advection-diffusion solution operators, we obtain convergence analysis for the proposed numerical method. Finally, we present numerical results to demonstrate the accuracy and efficiency of the proposed method in computing KPP front speeds in time-periodic cellular and chaotic flows, especially the time-dependent Arnold-Beltrami-Childress flow and time-dependent Kolmogorov flow in three-dimensional space. [ABSTRACT FROM AUTHOR]

Motivated by a wide range of real-world problems whose solutions exhibit boundary and interior layers, the numerical analysis of discretizations of singularly perturbed differential equations is an established subdiscipline within the study of the numerical approximation of solutions to differential equations. Consequently, much is known about how to accurately and stably discretize such equations on a priori adapted meshes in order to properly resolve the layer structure present in their continuum solutions. However, despite being a key step in the numerical simulation process, much less is known about the efficient and accurate solution of the linear systems of equations corresponding to these discretizations. In this paper, we discuss problems associated with the application of direct solvers to these discretizations, and we propose a preconditioning strategy that is tuned to the matrix structure induced by using layer-adapted meshes for convection-diffusion equations, proving a strong condition-number bound on the preconditioned system in one spatial dimension and a weaker bound in two spatial dimensions. Numerical results confirm the efficiency of the resulting preconditioners in one and two dimensions, with time-to-solution of less than one second for representative problems on 1024 x 1024 meshes and up to 40 x speedup over standard sparse direct solvers. [ABSTRACT FROM AUTHOR]

This paper is concerned with resolvent estimates on the real axis for the Helmholtz equation posed in the exterior of a bounded obstacle with Dirichlet boundary conditions when the obstacle is trapping. There are two resolvent estimates for this situation currently in the literature: (i) in the case of elliptic trapping the general "worst case" bound of exponential growth applies, and examples show that this growth can be realized through some sequence of wavenumbers; (ii) in the prototypical case of hyperbolic trapping where the Helmholtz equation is posed in the exterior of two strictly convex obstacles (or several obstacles with additional constraints) the nontrapping resolvent estimate holds with a logarithmic loss. This paper proves the first resolvent estimate for parabolic trapping by obstacles, studying a class of obstacles the prototypical example of which is the exterior of two squares (in two dimensions) or two cubes (in three dimensions), whose sides are parallel. We show, via developments of the vector-field/multiplier argument of Morawetz and the first application of this methodology to trapping configurations, that a resolvent estimate holds with a polynomial loss over the nontrapping estimate. We use this bound, along with the other trapping resolvent estimates, to prove results about integral equation formulations of the boundary value problem in the case of trapping. Feeding these bounds into existing frameworks for analyzing finite and boundary element methods, we obtain the first wavenumber-explicit proofs of convergence for numerical methods for solving the Helmholtz equation in the exterior of a trapping obstacle. [ABSTRACT FROM AUTHOR]

This paper is mainly concerned with a two-species competition model in open advective environments, where individuals cannot pass through the upstream boundary and do not return to the habitat after leaving the downstream boundary. By the theory of principal eigenvalue, we first obtain two critical curves (Γ1 and Γ2) in the plane of bifurcation parameters that sharply determine the local stability of the two semitrivial steady states. Then under various conditions on given parameters, we discuss the global dynamics via different techniques, including the comparison principle for eigenvalues and perturbation and compactness arguments, and show that both competitive exclusion and coexistence are possible. For general values of parameters, we take both analytic and numerical approaches to further understand the potential behaviors of Γ1 and Γ2, and we numerically observe that in addition to the competitive exclusion and coexistence, the bistable phenomenon is also possible, which is different from the known results of competitive ODE and reaction-diffusion systems (where bistability is impossible). The implication of our numerical results on future work is also discussed. [ABSTRACT FROM AUTHOR]

Nonlinear multiscale problems are ubiquitous in materials science and biology. Complicated interactions between nonlinearities and (nonseparable) multiple scales pose a major challenge for analysis and simulation. In this paper, we study the numerical homogenization for multiscale elliptic PDEs with monotone nonlinearity, in particular the Leray--Lions problem (a prototypical example is the p-Laplacian equation), where the nonlinearity cannot be parameterized with low dimensional parameters, and the linearization error is nonnegligible. We develop the iterated numerical homogenization scheme by combining numerical homogenization methods for linear equations and the so-called quasi-norm based iterative approach for monotone nonlinear equation. We propose a residual regularized nonlinear iterative method and, in addition, develop the sparse updating method for the efficient update of coarse spaces. A number of numerical results are presented to complement the analysis and validate the numerical method. [ABSTRACT FROM AUTHOR]

In this paper, we consider the generalized eigenvalue problem (GEP) for bidiagonalproduct (BP) pairs with sign regularity (SR), which include structured pairs associated with illconditioned matrices, such as Cauchy and Vandermonde matrices, that arise in many applications. A sufficient and necessary condition is provided for an SR-BP pair to be regular. For regular SR-BP pairs having both matrices singular, we establish sharp relative perturbation bounds to show that all the generalized eigenvalues including infinite and zero ones can be accurately determined by their BP representations. By operating on the BP representations, a new method is developed to accurately compute generalized eigenvalues of such a regular SR-BP pair. Our method first transforms the pair into an equivalent pair with a certain sign regularity. Then a technique is proposed to implicitly deflate all the infinite eigenvalues. After deflating infinite eigenvalues, finite eigenvalues are computed by reducing the deflated GEP into a standard eigenvalue problem. An attractive property of our method is that all the generalized eigenvalues are returned to high relative accuracy especially, all the infinite and zero eigenvalues are computed exactly. Error analysis and numerical experiments are performed to confirm the claimed high relative accuracy. [ABSTRACT FROM AUTHOR]

This paper proposes a simple and rigorous Fourier-matching method to study transverse-magnetic-polarized electro-magnetic resonances by a perfectly conducting slab with a finite number of subwavelength slits of width h\ll 1. Since variable separation is applicable in the region outside the slits, by Fourier transforming its governing equation, a magnetic field can be represented in terms of its derivative on the aperture. Next, inside each slit where variable separation is still available, the field can be represented as a Fourier series in terms of a countable set of basis functions with unknown Fourier coefficients. Finally, by matching the two subdomain representations on the aperture, we establish a linear system of an infinite number of equations governing the countable Fourier coefficients; the unknowns are further rescaled to be in the standard l² space. By the asymptotic expansion of each entry of the coefficient matrix, we rigorously show that its certain principal submatrix is invertible so that the infinite-dimensional linear system can be reduced to a finite-dimensional linear system. Resonance frequencies are exactly those frequencies making the linear system rank-deficient. This in turn leads to an asymptotic formula of accuracy O (h³ log h) for computing the resonance frequencies. We emphasize that the new formula is more accurate than all existing results and is the first formula for slits of number more than two to the best of our knowledge. Numerical experiments are carried out finally to validate the proposed formula and demonstrate its accuracy. [ABSTRACT FROM AUTHOR]

In this paper, we study the problem of computing the effective diffusivity for particles moving in chaotic flows. Instead of solving a convection-diffusion type cell problem in the Eulerian formulation (arising from homogenization theory for parabolic equations), we compute the motion of particles in the Lagrangian formulation, which is modeled by stochastic differential equations (SDEs). A robust numerical integrator based on a splitting method was proposed to solve the SDEs and rigorous error analysis for the numerical integrator was provided using the backward error analysis technique in our previous work. However, the upper bound on the error estimate is not sharp. To improve our result, we propose a new and uniform in time error analysis for the numerical integrator that allows us to get rid of the exponential growth factor in our previous error estimate. Our new error analysis is based on a probabilistic approach, which interprets the solution process generated by our numerical integrator as a Markov process. By exploring the ergodicity of the solution process, we prove the convergence analysis of our method in computing effective diffusivity over infinite time. We present numerical results to verify the accuracy and efficiency of the proposed method in computing effective diffusivity for several chaotic flows, especially the Arnold--Beltrami--Childress flow and Kolmogorov flow in three-dimensional space. [ABSTRACT FROM AUTHOR]

We consider a class of stochastic path-dependent volatility models where the stochastic volatility, whose square follows the Cox-Ingersoll-Ross model, is multiplied by a (leverage) function of the spot process, its running maximum, and time. We propose a Monte Carlo simulation scheme which combines a log-Euler scheme for the spot process with the full truncation Euler scheme or the backward Euler-Maruyama scheme for the squared stochastic volatility component. Under some mild regularity assumptions and a condition on the Feller ratio, we establish the strong convergence with order 1/2 (up to a logarithmic factor) of the approximation process up to a critical time. The model studied in this paper contains as special cases Heston-type stochastic-local volatility models, the state of the art in derivative pricing, and a relatively new class of path-dependent volatility models. The present paper is the rst to prove the convergence of the popular Euler schemes with a positive rate, which is moreover consistent with that for Lipschitz coeffcients and hence optimal. [ABSTRACT FROM AUTHOR]

In this paper, we discuss the numerical solution of certain matrix-valued PDEs. Such PDEs arise, for example, when constructing a Riemannian contraction metric for a dynamical system given by an autonomous ODE. We develop and analyze a new meshfree discretization scheme using kernel-based approximation spaces. However, since these approximation spaces have now to be matrix-valued, the kernels we need to use are fourth-order tensors. We will review and extend recent results on even more general reproducing kernel Hilbert spaces. We will then apply this general theory to solve a matrix-valued PDE and derive error estimates for the approximate solution. The paper ends with applications to typical examples from dynamical systems. [ABSTRACT FROM AUTHOR]

In this paper we study the level set formulations of certain geometric evolution equations from a numerical point of view. Specifically, we consider the flow by powers greater than one of the mean curvature (PMCF) and the inverse mean curvature flow (IMCF). Since the corresponding equations in level set form are quasi-linear, degenerate, and especially possibly singular a regularization method is used in the literature to approximate these equations to overcome the singularities of the equations. The regularized equations depend both on an regularization parameter and in case of the ICMF additionally on further parameters. Motivated by Feng, Neilan, and Prohl [Numer. Math., 108 (2007), pp. 93-119], who study the finite element approximation of IMCF, we prove error estimates for the finite element approximation of the regularized equations for PMCF. We validate the rates with numerical examples. Additionally, the regularization error in the rotationally symmetric case for both flows is numerically analyzed in a two-dimensional setting. Therefore, in the case of IMCF we fix the additional parameters. Furthermore, having the goal to estimate the regularization error we derive barriers for the regularized level set IMCF respecting all parameters and specify them further in a rotationally symmetric simplified case. At the end of the paper we present simulations in the three-dimensional case. [ABSTRACT FROM AUTHOR]

The purpose of this paper is to investigate mathematically the fundamental mechanism of gas discharge. Townsend discovered that α- and γ-mechanisms are essential for a gas ionization process. Morrow, Degond, and Lucquin-Desreux derived a mathematical model taking these two mechanisms into account. This model consists of nonlinear hyperbolic, parabolic, and elliptic partial differential equations. In this paper, we establish a framework for analyzing this model rigorously, because no mathematical result has been announced. More precisely, we show the unique existence of a time-local solution to an initial-boundary value problem over an unbounded domain. The main difficulty of this problem lies in the fact that the direction of boundary characteristics of a hyperbolic equation degenerates at infinite distance. Generally speaking, the regularity of solutions to linearized hyperbolic equations with degenerate boundary characteristics may be lost near the boundary. This implies the possibility that loss on the derivatives arises at each step of the inductive scheme to solve the nonlinear problem. We first use the weighted Sobolev spaces in the construction of solutions to make clear the direction of characteristics. Furthermore, we reduce the initial-boundary value problem for the hyperbolic equation to an initial value problem by applying several extension operators. This reduction enables us to avoid the loss on the derivatives. [ABSTRACT FROM AUTHOR]

This paper is concerned with efficient representations and approximations of the solution to the scattering problem by a system of strongly coupled plasmonic particles. Three schemes are developed: the first is the resonant expansion which uses the resonant modes of the system of particles computed by a conformal transformation, the second is the hybridized resonant expansion which uses linear combinations of the resonant modes for each of the particles in the system as a basis to represent the solution, and the last one is the multipole expansion with respect to the origin. By considering a system formed by two plasmonic particles of circular shape, we demonstrate the relations between these expansion schemes and their advantages and disadvantages both analytically and numerically. In particular, we emphasize the efficiency of the resonant expansion scheme in approximating the near field of the system of particles. The difference between these plasmonic particle systems and the nonresonant dielectric particle system is also highlighted. The paper provides a guidance on the challenges for numerical simulations of strongly coupled plasmonic systems. [ABSTRACT FROM AUTHOR]

The paper studies the convergence properties of (continuous-time) best-response dynamics from game theory. Despite their fundamental role in game theory, best-response dynamics are poorly understood in many games of interest due to the discontinuous, set-valued nature of the best-response map. The paper focuses on elucidating several important properties of best-response dynamics in the class of multiagent games known as potential games--a class of games with fundamental importance in multiagent systems and distributed control. It is shown that in almost every potential game and for almost every initial condition, the best-response dynamics (i) have a unique solution, (ii) converge to pure-strategy Nash equilibria, and (iii) converge at an exponential rate. [ABSTRACT FROM AUTHOR]

In our recent work [Z. Peng et al., J. Comput. Phys., 415 (2020), 109485], a family of high-order asymptotic preserving (AP) methods, termed IMEX-LDG methods, are designed to solve some linear kinetic transport equations, including the one-group transport equation in slab geometry and the telegraph equation, in a diffusive scaling. As the Knudsen number ε goes to zero, the limiting schemes are implicit discretizations to the limiting diffusive equation. Both Fourier analysis and numerical experiments imply the methods are unconditionally stable in the diffusive regime when ε « 1. In this paper, we develop an energy approach to establish the numerical stability of the IMEX1-LDG method, the subfamily of the methods that is first-order accurate in time and arbitrary order in space, for the model with general material properties. Our analysis is the first to simultaneously confirm unconditional stability when ε « I and the uniform stability property with respect to e. To capture the unconditional stability, we propose a novel discrete energy and explore various stabilization mechanisms of the method and their relative contributions in different regimes. A general form of the weight function, introduced to obtain the unconditional stability for ε « 1, is also for the first time considered in such stability analysis. Based on uniform stability, a rigorous asymptotic analysis is then carried out to show the AP property. [ABSTRACT FROM AUTHOR]

This paper provides an a priori error analysis of a localized orthogonal decomposition method for the numerical stochastic homogenization of a model random diffusion problem. If the uniformly elliptic and bounded random coefficient field of the model problem is stationary and satisfies a quantitative decorrelation assumption in the form of the spectral gap inequality, then the expected L² error of the method can be estimated, up to logarithmic factors, by H + (ε/H)d/2, ε being the small correlation length of the random coefficient and H the width of the coarse finite element mesh that determines the spatial resolution. The proof bridges recent results of numerical homogenization and quantitative stochastic homogenization. [ABSTRACT FROM AUTHOR]

The normalized gradient ow, i.e., the gradient ow with discrete normalization (GFDN) introduced in [W. Bao and Q. Du, SIAM J. Sci. Comput., 25 (2004), pp. 1674{1697] or the imaginary time evolution method is one of the most popular techniques for computing the ground states of Bose{Einstein condensates (BECs). In this paper, we revisit the time discretizations for the GFDN and its generalization to the multicomponent BECs. Several widely used time discretizations are demonstrated to be not accurate for computing the ground state solution in the general case, especially for the multicomponent BECs with two or more constraints even for the most accepted linearized backward Euler schemes. More precisely, these schemes usually converge to a solution with an error depending on the time step size. To accurately and efficiently compute the ground state solution of BECs, we propose the gradient ow with Lagrange multiplier (GFLM) method which can be viewed as the modified GFDN by introducing the explicit Lagrange multiplier terms or an approximation of the continuous normalized gradient ow. Through analysis and numerical computation, we clarify that, in order to accurately compute the ground state solution, the GFDN method must be discretized in very special ways in time such as the linearized backward Euler scheme, while the GFLM can be discretized by various schemes and works for multicomponent BEC with multiple constraints. We propose a simple backward-forward Euler scheme to discretize the GFLM, where only a linear system with constant coefficients needs to be solved at each time step. Taking the spin-1 BEC as an example, we demonstrate how the traditional GFDN approach fails in computing the ground state solution of general multicomponent BECs accurately. The proposed GFLM with backward-forward Euler discretization is shown to be accurate and efficient. Extensive numerical results are provided to verify our claims. [ABSTRACT FROM AUTHOR]

The aim of this paper is to develop and analyze high-order time stepping schemes for approximately solving semilinear subdiffusion equations. We apply the convolution quadrature generated by k-step backward differentiation formula (BDFk) to discretize the time-fractional derivative with order α ∈ (0, 1) and modify the starting steps in order to achieve optimal convergence rate. This method has already been well-studied for the linear fractional evolution equations in Jin, Li, and Zhou [SIAM J. Sci. Comput., 39 (2017), pp. A3129-A3152], while the numerical analysis for the nonlinear problem is still missing in the literature. By splitting the nonlinear potential term into an irregular linear part and a smoother nonlinear part and using the generating function technique, we prove that the convergence order of the corrected BDFk scheme is O(τmin(k,1+2α-ε)) without imposing further assumption on the regularity of the solution. Numerical examples are provided to support our theoretical results. [ABSTRACT FROM AUTHOR]

In this paper, we will introduce a high order numerical method to solve the scattering problems with nonperiodic incident fields and (locally perturbed) periodic surfaces. For the problems we are considering, the classical methods to treat quasi-periodic scattering problems no longer work, while a Bloch transform based numerical method was proposed in [A. Lechleiter and R. Zhang, SIAM J. Sci. Comput., 39 (2017), pp. B819{B839]. This numerical method, on one hand, is able to solve this kind of problem convergently; on the other hand, it takes up a lot of time and memory during the computation. The motivation of this paper is to improve this numerical method from the regularity results of the Bloch transform of the total field, which have been studied in [R. Zhang, arXiv:1708.0756, 2017]. As the set of the singularities of the total field is discrete in R and finite in one periodic cell, we are able to improve the numerical method by designing a proper integration contour with special conditions at the singularities. With a good choice of the transformation, we can prove that the new numerical method could possess a super algebraic convergence rate. This new method improves the efficiency significantly. At the end of this paper, several numerical results will be provided to show the fast convergence of the new method. The method also provides a possibility to solve more complicated problems efficiently, e.g., three-dimensional problems or electromagnetic scattering problems. [ABSTRACT FROM AUTHOR]

The recent development of the quadrature-based moment equations (QBME) for the simulation of rarefied gases using the Boltzmann equation led to a promising hyperbolic moment model. This paper deals with the two-dimensional QBME model, presents its derivation based on different possible expansions, and gives explicit equations for the QBME model and related models before describing a numerical method to solve the nonconservative PDE system on two-dimensional, unstructured grids. The first simulations using the two-dimensional QBME are shown in this paper using the flow past a cylinder and a forward facing step test case. The results indicate the applicability of the QBME models for rarefied gas flows as well as convergence to the Euler solution in the case of vanishing Knudsen number. [ABSTRACT FROM AUTHOR]

This paper deals with a general form of variational problems in Banach spaces which encompasses variational inequalities as well as minimization problems. We prove a characterization of local error bounds for the distance to the (primal-dual) solution set and give a sufficient condition for such an error bound to hold. In the second part of the paper, we consider an algorithm of augmented Lagrangian type for the solution of such variational problems. We give some global convergence properties of the method and then use the error bound theory to provide estimates for the rate of convergence and to deduce boundedness of the sequence of penalty parameters. Finally, numerical results for optimal control, Nash equilibrium problems, and elliptic parameter estimation problems are presented. [ABSTRACT FROM AUTHOR]

In [SIAM J. Matrix Anal. Appl., 38 (2017), pp. 649–655], Crouzeix and Palencia show that the closure of the numerical range of a square matrix or linear operator A is a (1 +√2)spectral set for A; that is, for any function f analytic in the interior of the numerical range W(A) and continuous on its boundary, the inequality ‖f(A)‖ ≤ (1 +√2)‖f‖ W(A) holds, where the norm on the left is the operator 2-norm and ‖f‖ W(A) on the right denotes the supremum of ∣f(z)∣ over z ∈ W(A). In this paper, we show how the arguments in their paper can be extended to show that other regions in the complex plane that do not necessarily contain W(A) are K-spectral sets for a value of K that may be close to 1 +√2. We also find some special cases in which the constant (1 +√2) for W(A) can be replaced by 2, which is the value conjectured by Crouzeix. [ABSTRACT FROM AUTHOR]

Scattering problems for periodic structures have been studied a lot in the past few years. A main idea for numerical solution methods is to reduce such problems to one periodicity cell. In contrast to periodic settings, scattering from locally perturbed periodic surfaces is way more challenging. In this paper, we introduce and analyze a new numerical method to simulate scattering from locally perturbed periodic structures based on the Bloch transform. As this transform is applied only in periodic domains, we first rewrite the scattering problem artificially in a periodic domain. With the help of the Bloch transform, we next transform this problem into a coupled family of quasi-periodic problems posed in the periodicity cell. A numerical scheme then approximates the family of quasi-periodic solutions (we rely on the finite element method) and backtransformation provides the solution to the original scattering problem. In this paper, we give convergence analysis and error bounds for a Galerkin discretization in the spatial and the quasi periodicity's unit cells. We also provide a simple and efficient way of implementation that does not require numerical integration in the quasi periodicity, together with numerical examples for scattering from locally perturbed periodic surfaces computed by this scheme. [ABSTRACT FROM AUTHOR]

This paper provides a unified and detailed presentation of root-node-style algebraic multigrid (AMG). AMG is a popular and effective iterative method for solving large, sparse linear systems that arise from discretizing partial differential equations. However, while AMG is designed for symmetric positive definite (SPD) matrices, certain SPD problems, such as anisotropic diffusion, are still not adequately addressed by existing methods. Non-SPD problems pose an even greater challenge, and in practice AMG is often not considered as a solver for such problems. The focus of this paper is on so-called root-node AMG, which can be viewed as a combination of classical and aggregation-based multigrid. An algorithm for root-node AMG is outlined, and a filtering strategy is developed, which is able to control the cost of using root-node AMG, particularly on difficult problems. New theoretical motivation is provided for root-node and energy-minimization as applied to symmetric as well nonsymmetric systems. Numerical results are then presented demonstrating the robust ability of root-node AMG to solve nonsymmetric problems, systems-based problems, and difficult SPD problems, including strongly anisotropic diffusion, convection-diffusion, and upwind steady-state transport, in a scalable manner. New detailed estimates of the computational cost of the setup and solve phases are given for each example, providing additional support for root-node AMG over alternative methods. [ABSTRACT FROM AUTHOR]

A framework, PhaseLift, was recently proposed to solve the phase retrieval problem. In this framework, the problem is solved by optimizing a cost function over the set of complex Hermitian positive semidefinite matrices. This approach to phase retrieval motivates a more general consideration of optimizing cost functions on semidefinite Hermitian matrices where the desired minimizers are known to have low rank. This paper considers an approach based on an alternative cost function defined on a union of appropriate manifolds. It is related to the original cost function in a manner that preserves the ability to find a global minimizer and is significantly more efficient computationally. A rank-based optimality condition for stationary points is given and optimization algorithms based on state-of-the-art Riemannian optimization and dynamically reducing rank are proposed. Empirical evaluations are performed using the PhaseLift problem. The new approach is shown to be an effective method of phase retrieval with computational efficiency increased substantially compared to a state-of-the-art algorithm, the Wirtinger flow algorithm. A preliminary version of this paper can be found in [W. Huang, K. A. Gallivan, and X. Zhang, Procedia Comput. Sci., 80 (2016), pp. 1125-1134]. [ABSTRACT FROM AUTHOR]

This paper is focused on the residual distribution (RD) interpretation of the Dobrev, Kolev, and Rieben scheme [SIAM J. Sci. Comput., 34 (2012), pp. B606-B641] for the numerical solution of the Euler equations in Lagrangian form. The first ingredient of the original scheme is the staggered grid formulation which uses continuous node-based finite element approximations for the kinematic variables and cell-centered discontinuous finite elements for the thermodynamic parameters. The second ingredient of the Dobrev et al. scheme is an artificial viscosity technique applied in order to make possible the computation of strong discontinuities. The aim of this paper is to provide an efficient mass matrix diagonalization method in order to avoid the inversion of the global sparse mass matrix while keeping all the accuracy properties and to construct a parameter-free stabilization of the scheme to get rid of the artificial viscosity. In addition, we study the conservation and entropy properties of the constructed RD scheme. To demonstrate the robustness of the proposed RD scheme, we solve several one-dimensional shock tube problems from rather mild to very strong ones. This paper also illustrates a general technique that enables, from a nonconservative formulation of a system that has a conservative formulation, design of a numerical approximation that will provably give sequences of solution converging to a weak solution of the problem. This enables us to use directly variables that are more pertinent, from an engineering point of view, than the standard conserved variables: here the specific internal energy. [ABSTRACT FROM AUTHOR]

In this paper we are concerned with solvers for the systems arising from the plane wave discretizations of three dimensional Helmholtz equations with large wave numbers. For simplicity, we consider only the plane wave weighted least squares (PWLS) method for Helmholtz equations. The main goal of this paper is to construct efficient multilevel preconditioners for solving the resulting Helmholtz systems. To this end, we first build a multilevel space decomposition for the plane wave discretization space based on overlapping domain decompositions. Then, corresponding to the space decomposition, we construct two additive multilevel preconditioners with smoothers for the underlying Helmholtz systems. In these preconditioners, each subproblem to be solved has a very small number of degrees of freedom, which just equals the number of the plane wave basis functions on one element. Moreover, the preconditioners possess the optimal computational complexity per iteration. We apply the proposed multilevel preconditioners with a constant coarse mesh size to solve the systems generated by the PWLS discretization for three dimensional Helmholtz equations, and we find that the new preconditioners possess nearly stable convergence, i.e., the iteration counts of the preconditioned CG methods with the preconditioners increase very slowly when the wave number increases (and the fine mesh size decreases). [ABSTRACT FROM AUTHOR]

In this paper, we propose a new stabilized projection-based proper orthogonal decomposition reduced order model (POD-ROM) for the numerical simulation of incompressible flows. The new method draws inspiration from successful numerical stabilization techniques used in the context of finite element (FE) methods, such as local projection stabilization (LPS). In particular, the new LPS-ROM is a velocity-pressure ROM that uses pressure modes as well to compute the reduced order pressure needed for instance in the computation of relevant quantities, such as drag and lift forces on bodies in the flow. The new LPS-ROM circumvents the standard discrete inf-sup condition for the POD velocity-pressure spaces, whose fulfillment can be rather expensive in realistic applications in computational fluid dynamics (CFD). Also, the velocity modes do not have to be either strongly or weakly divergence-free, which allows the use of snapshots generated, for instance, with penalty or projection-based stabilized methods. The numerical analysis of the fully Navier-- Stokes discretization for the new LPS-ROM is presented by mainly deriving the corresponding error estimates. Numerical studies are performed to discuss the accuracy and performance of the new LPS-ROM on a two-dimensional laminar unsteady flow past a circular obstacle. [ABSTRACT FROM AUTHOR]

The current paper proposes a parametric eigenvalue problem for the Helmholtz equation on the plane based on two well-known physical models of emission from two-dimensional (2-D) microcavity lasers. The first model has been developed for passive, either lossy or lossless, open cavities and is usually referred to as the Complex-Frequency Eigenvalue Problem. The second model, named the Lasing Eigenvalue Problem (LEP), has been tailored to characterize emission from open cavities, filled in with gain material, on the threshold of nonattenuating in time radiation. Our generalized model contains both of them as special cases. We reduce the original problem to a nonlinear eigenvalue problem for a set of boundary integral equations with weakly singular kernels and formulate it as a parametric eigenvalue problem for a holomorphic Fredholm operator-valued function, which is convenient for mathematical and numerical analysis. Using this formulation and fundamental results of the theory of holomorphic operator-valued functions, we study the properties of the spectrum. For numerical solution of the problem, we propose a Nyström method, prove its convergence, and derive error estimates in the eigenvalue approximation. Combining this discretization with the residual inverse iteration technique, we compute approximate solutions of LEP and compare them with the exact ones and with the results obtained using other numerical methods. [ABSTRACT FROM AUTHOR]