Search

This paper is devoted to the robust approximation with a variational phase field approach of multiphase mean curvature flows with possibly highly contrasted mobilities. The case of harmonically additive mobilities has been addressed recently using a suitable metric to define the gradient flow of the phase field approximate energy. We generalize this approach to arbitrary nonnegative mobilities using a decomposition as sums of harmonically additive mobilities. We establish the consistency of the resulting method by analyzing the sharp interface limit of the flow: a formal expansion of the phase field shows that the method is of second order. We propose a simple numerical scheme to approximate the solutions to our new model. Finally, we present some numerical experiments in dimensions 2 and 3 that illustrate the interest and effectiveness of our approach, in particular for approximating flows in which the mobility of some phases is zero., Comment: 22 pages, 8 figures

In this paper, we analyze an operator splitting scheme of the nonlinear heat equation in $\Omega\subset\mathbb{R}^d$ ($d\geq 1$): $\partial_t u = \Delta u + \lambda |u|^{p-1} u$ in $\Omega\times(0,\infty)$, $u=0$ in $\partial\Omega\times(0,\infty)$, $u ({\bf x},0) =\phi ({\bf x})$ in $\Omega$. where $\lambda\in\{-1,1\}$ and $\phi \in W^{1,q}(\Omega)\cap L^{\infty} (\Omega)$ with $2\leq p < \infty$ and $d(p-1)/20$. Finally, we give some numerical examples to confirm the reliability of the analyzed result.

In the context of mathematical modeling, it is sometimes convenient to integrate models of different nature. These types of combinations, however, might entail difficulties even when individual models are well-understood, particularly in relation to the well-posedness of the ensemble. In this article, we focus on combining two classes of dissimilar diffusive models: the first one defined over a continuum and the second one based on discrete equations that connect average values of the solution over disjoint subdomains. For stationary problems, we show unconditional stability of the linked problems and then the stability and convergence of its discretized counterpart when mixed finite elements are used to approximate the model on the continuum. The theoretical results are highlighted with numerical examples illustrating the effects of linking diffusive models. As a side result, we show that the methods introduced in this article can be used to infer the solution of diffusive problems with incomplete data., 22 pages, 7 figures, preprint. International Journal for Numerical Methods in Engineering (2023)

In this article, we investigate the estimation of the transient mold-slab heat flux in continuous casting molds given some thermocouples measurements in the mold plates. Mathematically, we can see this problem as the estimation of a Neumann boundary condition given pointwise state observations in the interior of the domain. We formulate it in a deterministic inverse problem setting. After introducing the industrial problem, we present the mold thermal model and related assumptions. Then, we formulate the boundary heat flux estimation problem in a deterministic inverse problem setting using a sequential approach according to the sequentiality of the temperature measurements. We consider different formulations of the inverse problem. For each one, we develop novel direct methodologies exploiting a space parameterization of the heat flux and the linearity of the mold model. We construct these methods to be divided into a computationally expensive offline phase that can be computed before the process starts, and a cheaper online phase to be performed during the casting process. To conclude, we test the performance of the proposed methods in two benchmark cases Agencia Estatal de Investigación, Grant/Award Number: PID2019-105615RBI00/AEI; European Research Council, Grant/Award Number:765374; H2020 Marie Skłodowska-Curie Actions, Grant/Award Number: 681447; Ministerio de Economía, Industria y Competitividad, Gobierno de España, Grant/Award Number: MTM2015-68275-R SI

The phase field approach is widely used to model fracture behaviors due to the absence of the need to track the crack topology and the ability to predict crack nucleation and branching. In this work, the asynchronous variational integrators (AVI) is adapted for the phase field approach of dynamic brittle fractures. The AVI is derived from Hamilton's principle and allows each element in the mesh to have its own local time step that may be different from others'. While the displacement field is explicitly updated, the phase field is implicitly solved, with upper and lower bounds strictly and conveniently enforced. In particular, the AT1 and AT2 variants are equally easily implemented. Three benchmark problems are used to study the performances of both AT1 and AT2 models and the results show that AVI for phase field approach significantly speeds up the computational efficiency and successfully captures the complicated dynamic fracture behavior.

The large-deformation mechanics and multiphysics of continuous or fracturing partially saturated porous media under static and dynamic loads are significant in engineering and science. This article is devoted to a computational coupled large-deformation periporomechanics paradigm assuming passive air pressure for modeling dynamic failure and fracturing in variably saturated porous media. The coupled governing equations for bulk and fracture material points are formulated in the current/deformed configuration through the updated Lagrangian-Eulerian framework. It is hypothesized that the horizon of a mixed material point remains spherical, and its neighbor points are determined in the current configuration. As a significant contribution, the mixed interface/phreatic material points near the phreatic line are explicitly considered for modeling the transition from partial to full saturation (vice versa) through the mixed peridynamic state concept. We have formulated the coupled constitutive correspondence principle and stabilization scheme in the updated Lagrangian-Eulerian framework for bulk and interface points. We numerically implement the coupled large deformation periporomechanics through a fully implicit fractional-step algorithm in time and a hybrid updated Lagrangian-Eulerian meshfree method in space. Numerical examples are presented to validate the implemented stabilized computational coupled large deformation periporomechanics and demonstrate its efficacy and robustness in modeling dynamic failure and fracturing in variably saturated porous media.

The viscous flow of two immiscible fluids in a porous medium on the Darcy scale is governed by a system of nonlinear parabolic equations. If infinite mobility of one phase can be assumed (e.g. in soil layers in contact with the atmosphere) the system can be substituted by the scalar Richards model. Thus, the domain of the porous medium may be partitioned into disjoint subdomains with either the full two-phase or the simplified Richards model dynamics. Extending the one-model approach from [1, 2] we suggest coupling conditions for this hybrid model approach. Based on an Euler implicit discretisation, a linear iterative (-type) domain decomposition scheme is proposed, and proven to be convergent. The theoretical findings are verified by a comparative numerical study that in particular confirms the efficiency of the hybrid ansatz as compared to full two-phase model computations., 26 pages

We present a probabilistic analysis of two Krylov subspace methods for solving linear systems. We prove a central limit theorem for norms of the residual vectors that are produced by the conjugate gradient and MINRES algorithms when applied to a wide class of sample covariance matrices satisfying some standard moment conditions. The proof involves establishing a four moment theorem for the so-called spectral measure, implying, in particular, universality for the matrix produced by the Lanczos iteration. The central limit theorem then implies an almost-deterministic iteration count for the iterative methods in question.

In this paper, we present and analyze a linear fully discrete second order scheme with variable time steps for the phase field crystal equation. More precisely, we construct a linear adaptive time stepping scheme based on the second order backward differentiation formulation (BDF2) and use the Fourier spectral method for the spatial discretization. The scalar auxiliary variable approach is employed to deal with the nonlinear term, in which we only adopt a first order method to approximate the auxiliary variable. This treatment is extremely important in the derivation of the unconditional energy stability of the proposed adaptive BDF2 scheme. However, we find for the first time that this strategy will not affect the second order accuracy of the unknown phase function $\phi^{n}$ by setting the positive constant $C_{0}$ large enough such that $C_{0}\geq 1/\Dt.$ The energy stability of the adaptive BDF2 scheme is established with a mild constraint on the adjacent time step radio $\gamma_{n+1}:=\Dt_{n+1}/\Dt_{n}\leq 4.8645$. Furthermore, a rigorous error estimate of the second order accuracy of $\phi^{n}$ is derived for the proposed scheme on the nonuniform mesh by using the uniform $H^{2}$ bound of the numerical solutions. Finally, some numerical experiments are carried out to validate the theoretical results and demonstrate the efficiency of the fully discrete adaptive BDF2 scheme., Comment: 21 pages, 5 figures

We prove several optimal-order error estimates for a finite-element method applied to an inhomogeneous Robin boundary value problem (BVP) for the Poisson equation defined in a smooth bounded domain in $\mathbb{R}^n$, $n=2,3$. The boundary condition is weakly imposed using Nitsche's method. The Robin BVP is interpreted as the classical penalty method with the penalty parameter $\varepsilon$. The optimal choice of the mesh size $h$ relative to $\varepsilon$ is a non-trivial issue. This paper carefully examines the dependence of $\varepsilon$ on error estimates. Our error estimates require no unessential regularity assumptions on the solution. Numerical examples are also reported to confirm our results., 17 pages, 3 figures

The standard in rod finite element formulations is the Bubnov-Galerkin projection method, where the test functions arise from a consistent variation of the ansatz functions. This approach becomes increasingly complex when highly nonlinear ansatz functions are chosen to approximate the rod's centerline and cross-section orientations. Using a Petrov-Galerkin projection method, we propose a whole family of rod finite element formulations where the nodal generalized virtual displacements and generalized velocities are interpolated instead of using the consistent variations and time derivatives of the ansatz functions. This approach leads to a significant simplification of the expressions in the discrete virtual work functionals. In addition, independent strategies can be chosen for interpolating the nodal centerline points and cross-section orientations. We discuss three objective interpolation strategies and give an in-depth analysis concerning locking and convergence behavior for the whole family of rod finite element formulations., arXiv admin note: text overlap with arXiv:2301.05595

The classical Smagorinsky model's solution is an approximation to a (resolved) mean velocity. Since it is an eddy viscosity model, it cannot represent a flow of energy from unresolved fluctuations to the (resolved) mean velocity. This model has recently been modified to incorporate this flow and still be well-posed. Herein we first develop some basic properties of the modified model. Next, we perform a complete numerical analysis of two algorithms for its approximation. They are tested and proven to be effective., 26pages,4 figures, 2 table

Distributing points on a (possibly high-dimensional) sphere with minimal energy is a long-standing problem in and outside the field of mathematics. This paper considers a novel energy function that arises naturally from statistics and combinatorial optimization, and studies its theoretical properties. Our result solves both the exact optimal spherical point configurations in certain cases and the minimal energy asymptotics under general assumptions. Connections between our results and the L1-Principal Component analysis and Quasi-Monte Carlo methods are also discussed., 24 pages, 5 figures. To appear in Random Structures & Algorithms

This work deals with shape optimization for contact mechanics. More specifically, the linear elasticity model is considered under the small deformations hypothesis, and the elastic body is assumed to be in contact (sliding or with Tresca friction) with a rigid foundation. The mathematical formulations studied are two regularized versions of the original variational inequality: the penalty formulation and the augmented Lagrangian formulation. In order to get the shape derivatives associated to those two non-differentiable formulations, we follow an approach based on directional derivatives \cite{chaudet2020shape,chaudet2021shape}. This allows us to develop a gradient-based topology optimization algorithm, built on these derivatives and a level-set representation of shapes. The algorithm also benefits from a mesh-cutting technique, which gives an explicit representation of the shape at each iteration, and enables us to apply the boundary conditions strongly on the contact zone. The different steps of the method are detailed. Then, to validate the approach, some numerical results on two-dimensional and three-dimensional benchmarks are presented.

Gradient-type iterative methods for solving Hermitian eigenvalue problems can be accelerated by using preconditioning and deflation techniques. A preconditioned steepest descent iteration with implicit deflation (PSD-id) is one of such methods. The convergence behavior of the PSD-id is recently investigated based on the pioneering work of Samokish on the preconditioned steepest descent method (PSD). The resulting non-asymptotic estimates indicate a superlinear convergence of the PSD-id under strong assumptions on the initial guess. The present paper utilizes an alternative convergence analysis of the PSD by Neymeyr under much weaker assumptions. We embed Neymeyr's approach into the analysis of the PSD-id using a restricted formulation of the PSD-id. More importantly, we extend the new convergence analysis of the PSD-id to a practically preferred block version of the PSD-id, or BPSD-id, and show the cluster robustness of the BPSD-id. Numerical examples are provided to validate the theoretical estimates., Comment: 26 pages, 10 figures

We propose a mathematical and numerical model for the simulation of the heart function that couples cardiac electrophysiology, active and passive mechanics and hemodynamics, and includes reduced models for cardiac valves and the circulatory system. Our model accounts for the major feedback effects among the different processes that characterize the heart function, including electro-mechanical and mechano-electrical feedback as well as force-strain and force-velocity relationships. Moreover, it provides a three-dimensional representation of both the cardiac muscle and the hemodynamics, coupled in a fluid-structure interaction (FSI) model. By leveraging the multiphysics nature of the problem, we discretize it in time with a segregated electrophysiology-force generation-FSI approach, allowing for efficiency and flexibility in the numerical solution. We employ a monolithic approach for the numerical discretization of the FSI problem. We use finite elements for the spatial discretization of partial differential equations. We carry out a numerical simulation on a realistic human left heart model, obtaining results that are qualitatively and quantitatively in agreement with physiological ranges and medical images. This article is protected by copyright. All rights reserved.

Numerical stabilization is often used to eliminate (alleviate) the spurious oscillations generally produced by full order models (FOMs) in under-resolved or marginally-resolved simulations of convection-dominated flows. In this article, we investigate the role of numerical stabilization in reduced order models (ROMs) of marginally-resolved, convection-dominated incompressible flows. Specifically, we investigate the FOM-ROM consistency, that is, whether the numerical stabilization is beneficial both at the FOM and the ROM level. As a numerical stabilization strategy, we focus on the evolve-filter-relax (EFR) regularization algorithm, which centers around spatial filtering. To investigate the FOM-ROM consistency, we consider two ROM strategies: (i) the EFR-noEFR, in which the EFR stabilization is used at the FOM level, but not at the ROM level; and (ii) the EFR-EFR, in which the EFR stabilization is used both at the FOM and at the ROM level. We compare the EFR-noEFR with the EFR-EFR in the numerical simulation of a 2D incompressible flow past a circular cylinder in the convection-dominated, marginally-resolved regime. We also perform model reduction with respect to both time and Reynolds number. Our numerical investigation shows that the EFR-EFR is more accurate than the EFR-noEFR, which suggests that FOM-ROM consistency is beneficial in convection-dominated, marginally-resolved flows.

In this article we formulate and implement a computational multiphase periporomechanics model for unguided fracturing in unsaturated porous media. The same governing equation for the solid phase applies on and off cracks. Crack formation in this framework is autonomous, requiring no prior estimates of crack topology. As a new contribution, an energy-based criterion for arbitrary crack formation is formulated using the peridynamic effective force state for unsaturated porous media. Unsaturated fluid flow in the fracture space is modeled in a simplified way in line with the nonlocal formulation of unsaturated fluid flow in bulk. The formulated unsaturated fracturing periporomechanics is numerically implemented through a fractional step algorithm in time and a two-phase mixed meshless method in space. The two-stage operator split converts the coupled periporomechanics problem into an undrained deformation and fracture problem and an unsaturated fluid flow in the deformed skeleton configuration. Numerical simulations of in-plane open and shear cracking are conducted to validate the accuracy and robustness of the fracturing unsaturated periporomechanics model. Then numerical examples of wing cracking and non-planar cracking in unsaturated soil specimens are presented to demonstrate the efficacy of the proposed multiphase periporomechanics model for unguided cracking in unsaturated porous media.

Let $A$ be an $n\times n$ real matrix, and let $M$ be an $n\times n$ random matrix whose entries are i.i.d sub-Gaussian random variables with mean $0$ and variance $1$. We make two contributions to the study of $s_n(A+M)$, the smallest singular value of $A+M$. (1) We show that for all $\epsilon \geq 0$, $$\mathbb{P}[s_n(A + M) \leq \epsilon] = O(\epsilon \sqrt{n}) + 2e^{-\Omega(n)},$$ provided only that $A$ has $\Omega (n)$ singular values which are $O(\sqrt{n})$. This extends a well-known result of Rudelson and Vershynin, which requires all singular values of $A$ to be $O(\sqrt{n})$. (2) We show that any bound of the form $$\sup_{\|{A}\|\leq n^{C_1}}\mathbb{P}[s_n(A+M)\leq n^{-C_3}] \leq n^{-C_2}$$ must have $C_3 = \Omega (C_1 \sqrt{C_2})$. This complements a result of Tao and Vu, who proved such a bound with $C_3 = O(C_1C_2 + C_1 + 1)$, and counters their speculation of possibly taking $C_3 = O(C_1 + C_2)$., Comment: 15 pages; comments welcome!

In the wake of the 2020 COVID-19 epidemic, much work has been performed on the development of mathematical models for the simulation of the epidemic, and of disease models generally. Most works follow the susceptible-infected-removed (SIR) compartmental framework, modeling the epidemic with a system of ordinary differential equations. Alternative formulations using a partial differential equation (PDE) to incorporate both spatial and temporal resolution have also been introduced, with their numerical results showing potentially powerful descriptive and predictive capacity. In the present work, we introduce a new variation to such models by using delay differential equations (DDEs). The dynamics of many infectious diseases, including COVID-19, exhibit delays due to incubation periods and related phenomena. Accordingly, DDE models allow for a natural representation of the problem dynamics, in addition to offering advantages in terms of computational time and modeling, as they eliminate the need for additional, difficult-to-estimate, compartments (such as exposed individuals) to incorporate time delays. In the present work, we introduce a DDE epidemic model in both an ordinary- and partial differential equation framework. We present a series of mathematical results assessing the stability of the formulation. We then perform several numerical experiments, validating both the mathematical results and establishing model’s ability to reproduce measured data on realistic problems.

Knitting is an effective technique for producing complex three-dimensional surfaces owing to the inherent flexibility of interlooped yarns and recent advances in manufacturing providing better control of local stitch patterns. Fully yarn-level modelling of large-scale knitted membranes is not feasible. Therefore, we use a two-scale homogenisation approach and model the membrane as a Kirchhoff-Love shell on the macroscale and as Euler-Bernoulli rods on the microscale. The governing equations for both the shell and the rod are discretised with cubic B-spline basis functions. For homogenisation we consider only the in-plane response of the membrane. The solution of the nonlinear microscale problem requires a significant amount of time due to the large deformations and the enforcement of contact constraints, rendering conventional online computational homogenisation approaches infeasible. To sidestep this problem, we use a pre-trained statistical Gaussian Process Regression (GPR) model to map the macroscale deformations to macroscale stresses. During the offline learning phase, the GPR model is trained by solving the microscale problem for a sufficiently rich set of deformation states obtained by either uniform or Sobol sampling. The trained GPR model encodes the nonlinearities and anisotropies present in the microscale and serves as a material model for the membrane response of the macroscale shell. The bending response can be chosen in dependence of the mesh size to penalise the fine out-of-plane wrinkling of the membrane. After verifying and validating the different components of the proposed approach, we introduce several examples involving membranes subjected to tension and shear to demonstrate its versatility and good performance., 24 pages, 14 figures

Floating offshore structures often exhibit low-frequency oscillatory motions in the horizontal plane, with amplitudes in the same order as their characteristic dimensions and larger than the corresponding wave-frequency responses, making the traditional formulations in an inertial coordinate system inconsistent and less applicable. To address this issue, we extend and explore an alternative formulation completely based on a non-inertial body-fixed coordinate system. Unlike the traditional seakeeping models, this formulation consistently allows for large-amplitude horizontal motions. A numerical model based on a higher-order boundary element is applied to solve the resulting boundary-value problems in the time domain. A recently developed new set of explicit time-integration methods, which do not necessitate the use of upwind schemes for spatial derivatives, are adopted to deal with the convective-type free-surface conditions. To suppress the weak saw-tooth instabilities on the free surface in time marching, we also present novel low-pass filters based on optimized weighted least squares, which are in principle applicable for both structured and unstructured meshes. The presented schemes for the convection equation and the low-pass filter are also relevant for other engineering fields dealing with similar mathematical problems. For ship seakeeping and added resistance analyses, we show that the present computational model does not need to use soft-springs for surge and sway, in contrast to the traditional models. For a floating monopile, the importance of consistently taking into account the effects of large horizontal motions is demonstrated considering the bichromatic incident waves. The present model is considered as a complete second-order wave-load model, as all the second-order wave loads, including the sum-frequency and difference-frequency components, are solved simultaneously.

This article develops a new algorithm named TTRISK to solve high-dimensional risk-averse optimization problems governed by differential equations (ODEs and/or PDEs) under uncertainty. As an example, we focus on the so-called Conditional Value at Risk (CVaR), but the approach is equally applicable to other coherent risk measures. Both the full and reduced space formulations are considered. The algorithm is based on low rank tensor approximations of random fields discretized using stochastic collocation. To avoid non-smoothness of the objective function underpinning the CVaR, we propose an adaptive strategy to select the width parameter of the smoothed CVaR to balance the smoothing and tensor approximation errors. Moreover, unbiased Monte Carlo CVaR estimate can be computed by using the smoothed CVaR as a control variate. To accelerate the computations, we introduce an efficient preconditioner for the KKT system in the full space formulation.The numerical experiments demonstrate that the proposed method enables accurate CVaR optimization constrained by large-scale discretized systems. In particular, the first example consists of an elliptic PDE with random coefficients as constraints. The second example is motivated by a realistic application to devise a lockdown plan for United Kingdom under COVID-19. The results indicate that the risk-averse framework is feasible with the tensor approximations under tens of random variables., 33 pages

Funder: Narodowe Centrum Nauki; Id: http://dx.doi.org/10.13039/501100004281, Funder: Simons Foundation; Id: http://dx.doi.org/10.13039/100000893, We introduce and develop a theory of orthogonality with respect to Sobolev inner products on the real line for sequences of functions with a tridiagonal, skew‐Hermitian differentiation matrix. While a theory of such L2 ‐orthogonal systems is well established, Sobolev orthogonality requires new concepts and their analysis. We characterize such systems completely as appropriately weighted Fourier transforms of orthogonal polynomials and present a number of illustrative examples, inclusive of a Sobolev‐orthogonal system whose leading N coefficients can be computed in O ( N log N ) $ \mathcal{O} (N\log N)$ operations.

There is growing awareness that errors in the model equations cannot be ignored in data assimilation methods such as four-dimensional variational assimilation (4D-Var). If allowed for, more information can be extracted from observations, longer time windows are possible, and the minimization process is easier, at least in principle. Weak constraint 4D-Var estimates the model error and minimizes a series of linear least-squares cost functions using the conjugate gradient (CG) method; minimising each cost function is called an inner loop. CG needs preconditioning to improve its performance. In previous work, limited memory preconditioners (LMPs) have been constructed using approximations of the eigenvalues and eigenvectors of the Hessian in the previous inner loop. If the Hessian changes signicantly in consecutive inner loops, the LMP may be of limited usefulness. To circumvent this, we propose using randomised methods for low rank eigenvalue decomposition and use these approximations to cheaply construct LMPs using information from the current inner loop. Three randomised methods are compared. Numerical experiments in idealized systems show that the resulting LMPs perform better than the existing LMPs. Using these methods may allow more efficient and robust implementations of incremental weak constraint 4D-Var.

In this work we introduce a novel semi-implicit structure-preserving finite-volume/finite-difference scheme for the viscous and resistive equations of magnetohydrodynamics (MHD) based on an appropriate 3-split of the governing PDE system, which is decomposed into a first convective subsystem, a second subsystem involving the coupling of the velocity field with the magnetic field and a third subsystem involving the pressure-velocity coupling. The nonlinear convective terms are discretized explicitly, while the remaining two subsystems accounting for the Alfven waves and the magneto-acoustic waves are treated implicitly. The final algorithm is at least formally constrained only by a mild CFL stability condition depending on the velocity field of the pure hydrodynamic convection. To preserve the divergence-free constraint of the magnetic field exactly at the discrete level, a proper set of overlapping dual meshes is employed. The resulting linear algebraic systems are shown to be symmetric and therefore can be solved by means of an efficient standard matrix-free conjugate gradient algorithm. One of the peculiarities of the presented algorithm is that the magnetic field is defined on the edges of the main grid, while the electric field is on the faces. The final scheme can be regarded as a novel shock-capturing, conservative and structure preserving semi-implicit scheme for the nonlinear viscous and resistive MHD equations. Several numerical tests are presented to show the main features of our novel solver: linear-stability in the sense of Lyapunov is verified at a prescribed constant equilibrium solution; a 2nd-order of convergence is numerically estimated; shock-capturing capabilities are proven against a standard set of stringent MHD shock-problems; accuracy and robustness are verified against a nontrivial set of 2- and 3-dimensional MHD problems., 44 pages, 22 figures, 2 tables

In this paper we design efficient quadrature rules for finite element discretizations of nonlocal diffusion problems with compactly supported kernel functions. Two of the main challenges in nonlocal modeling and simulations are the prohibitive computational cost and the nontrivial implementation of discretization schemes, especially in three-dimensional settings. In this work we circumvent both challenges by introducing a parametrized mollifying function that improves the regularity of the integrand, utilizing an adaptive integration technique, and exploiting parallelization. We first show that the "mollified" solution converges to the exact one as the mollifying parameter vanishes, then we illustrate the consistency and accuracy of the proposed method on several two- and three-dimensional test cases. Furthermore, we demonstrate the good scaling properties of the parallel implementation of the adaptive algorithm and we compare the proposed method with recently developed techniques for efficient finite element assembly.

A matrix $A\in\mathbb{C}^{n\times n}$ is diagonalizable if it has a basis of linearly independent eigenvectors. Since the set of nondiagonalizable matrices has measure zero, every $A\in \mathbb{C}^{n\times n}$ is the limit of diagonalizable matrices. We prove a quantitative version of this fact conjectured by E.B. Davies: for each $\delta\in (0,1)$, every matrix $A\in \mathbb{C}^{n\times n}$ is at least $\delta\|A\|$-close to one whose eigenvectors have condition number at worst $c_n/\delta$, for some constants $c_n$ dependent only on $n$. Our proof uses tools from random matrix theory to show that the pseudospectrum of $A$ can be regularized with the addition of a complex Gaussian perturbation. Along the way, we explain how a variant of a theorem of \'Sniady implies a conjecture of Sankar, Spielman and Teng on the optimal constant for smoothed analysis of condition numbers., Comment: 17pp. Fixed a mistake in the appendix, all results are unchanged. Accepted version, to appear in CPAM

In this work, we focus on a fractional differential equation in Riesz form discretized by a polynomial B-spline collocation method. For an arbitrary polynomial degree $p$, we show that the resulting coefficient matrices possess a Toeplitz-like structure. We investigate their spectral properties via their symbol and we prove that, like for second order differential problems, also in this case the given matrices are ill-conditioned both in the low and high frequencies for large $p$. More precisely, in the fractional scenario the symbol has a single zero at $0$ of order $\alpha$, with $\alpha$ the fractional derivative order that ranges from $1$ to $2$, and it presents an exponential decay to zero at $\pi$ for increasing $p$ that becomes faster as $\alpha$ approaches $1$. This translates in a mitigated conditioning in the low frequencies and in a deterioration in the high frequencies when compared to second order problems. Furthermore, the derivation of the symbol reveals another similarity of our problem with a classical diffusion problem. Since the entries of the coefficient matrices are defined as evaluations of fractional derivatives of the B-spline basis at the collocation points, we are able to express the central entries of the coefficient matrix as inner products of two fractional derivatives of cardinal B-splines. Finally, we perform a numerical study of the approximation behavior of polynomial B-spline collocation. This study suggests that, in line with non-fractional diffusion problems, the approximation order for smooth solutions in the fractional case is $p+2-\alpha$ for even $p$, and $p+1-\alpha$ for odd $p$., Comment: 23 pages, 13 figures

In this paper, we are concerned with the inversion of circulant matrices and their quantized tensor-train (QTT) structure. In particular, we show that the inverse of a complex circulant matrix $A$, generated by the first column of the form $(a_0,\dots,a_{m-1},0,\dots,0,a_{-n},\dots, a_{-1})^\top$ admits a QTT representation with the QTT ranks bounded by $(m+n)$. Under certain assumptions on the entries of $A$, we also derive an explicit QTT representation of $A^{-1}$. The latter can be used, for instance, to overcome stability issues arising when numerically solving differential equations with periodic boundary conditions in the QTT format.

This note introduces a simple metric for benchmarking shock-capturing schemes. This metric is especially focused on the shock-capturing overshoots, which may undermine the robustness of numerical simulations, as well as the reliability of numerical results. The idea is to numerically solve the model linear advection equation with an initial condition of a square wave characterized with different wavenumbers. After one step of temporal evolution, the exact numerical overshoot error can be easily determined and shown as a function of the CFL number and the reduced wavenumber. With the overshoot error quantified by the present metric, a number of representative shock-capturing schemes are analyzed accordingly, and several findings including the amplitude of overshoots non-monotonously varying with the CFL number, and the amplitude of overshoots significantly depending on the reduced wavenumber of the square waves (discontinuities), are newly discovered, which are not before., 12 pages, 8 figures

The goal of the present paper is to understand the impact of numerical schemes for the reconstruction of data at cell faces in finite-volume methods, and to assess their interaction with the quadrature rule used to compute the average over the cell volume. Here, third-, fifth- and seventh-order WENO-Z schemes are investigated. On a problem with a smooth solution, the theoretical order of convergence rate for each method is retrieved, and changing the order of the reconstruction at cell faces does not impact the results, whereas for a shock-driven problem all the methods collapse to first-order. Study of the decay of compressible homogeneous isotropic turbulence reveals that using a high-order quadrature rule to compute the average over a finite volume cell does not improve the spectral accuracy and that all methods present a second-order convergence rate. However the choice of the numerical method to reconstruct data at cell faces is found to be critical to correctly capture turbulent spectra. In the context of simulations with finite-volume methods of practical flows encountered in engineering applications, it becomes apparent that an efficient strategy is to perform the average integration with a low-order quadrature rule on a fine mesh resolution, whereas high-order schemes should be used to reconstruct data at cell faces., arXiv admin note: text overlap with arXiv:1902.06665

We solve Maxwell's eigenvalue problem via isogeometric boundary elements and a contour integral method. We discuss the analytic properties of the discretisation, outline the implementation, and showcase numerical examples.

This is the peer reviewed version of the following article: Cavaliere, F. [et al.]. Nonintrusive reduced order model for parametric solutions of inertia relief problems. "International journal for numerical methods in engineering", 30 Agost 2021, vol. 122, núm. 16, p. 4270-4291., which has been published in final form at DOI:10.1002/nme.6702. This article may be used for non-commercial purposes in accordance with Wiley Terms and Conditions for Self-Archiving. The Inertia Relief (IR) technique is widely used by industry and produces equilibrated loads allowing to analyze unconstrained systems without resorting to the more expensive full dynamic analysis. The main goal of this work is to develop a computational framework for the solution of unconstrained parametric structural problems with IR and the Proper Generalized Decomposition (PGD) method. First, the IR method is formulated in a parametric setting for both material and geometric parameters. A reduced order model using the encapsulated PGD suite is then developed to solve the parametric IR problem, circumventing the so-called curse of dimensionality. With just one offline computation, the proposed PGD-IR scheme provides a computational vademecum that contains all the possible solutions for a predefined range of the parameters. The proposed approach is nonintrusive and it is therefore possible to be integrated with commercial finite element (FE) packages. The applicability and potential of the developed technique is shown using a three-dimensional test case and a more complex industrial test case. The first example is used to highlight the numerical properties of the scheme, whereas the second example demonstrates the potential in a more complex setting and it shows the possibility to integrate the proposed framework within a commercial FE package. In addition, the last example shows the possibility to use the generalized solution in a multi-objective optimization setting. This project is part of the Marie Sk lodowska-Curie ITN-EJD ProTechTion funded by the European Union Horizon 2020 research and innovation program with grant number 764636. The work of Fabiola Cavaliere, Sergio Zlotnik and Pedro D´ıez is partially supported by the Spanish Ministry of Economy and Competitiveness, Spain (Grant number: DPI2017-85139- C2-2-R) and by the Generalitat de Catalunya, Spain (Grant number: 2017-SGR-1278). Ruben Sevilla also acknowledges the support of the Engineering and Physical Sciences Research Council (Grant number: EP/P033997/1).

Using an autoencoder for dimensionality reduction, this paper presents a novel projection-based reduced-order model for eigenvalue problems. Reduced-order modelling relies on finding suitable basis functions which define a low-dimensional space in which a high-dimensional system is approximated. Proper orthogonal decomposition (POD) and singular value decomposition (SVD) are often used for this purpose and yield an optimal linear subspace. Autoencoders provide a nonlinear alternative to POD/SVD, that may capture, more efficiently, features or patterns in the high-fidelity model results. Reduced-order models based on an autoencoder and a novel hybrid SVD-autoencoder are developed. These methods are compared with the standard POD-Galerkin approach and are applied to two test cases taken from the field of nuclear reactor physics., Comment: 35 pages, 33 figures

In this paper, we present consistent and inconsistent discontinuous Galerkin methods for incompressible Euler and Navier-Stokes equations with the kinematic pressure, Bernoulli function and EMAC function. Semi- and fully discrete energy stability of the proposed dG methods are proved in a unified fashion. Conservation of total energy, linear and angular momentum is discussed with both central and upwind fluxes. Numerical experiments are presented to demonstrate our findings and compare our schemes with conventional schemes in the literature in both unsteady and steady problems. Numerical results show that global conservation of the physical quantities may not be enough to demonstrate the performance of the schemes, and our schemes are competitive and able to capture essential physical features in several benchmark problems., Comment: 32 pages, 24 figures

Global seismicity on all three solar system’s bodies with in situ measurements — Earth, Moon, and Mars — is due mainly to mechanical Rieger resonance (RR) of solar wind’s macroscopic flapping, driven by the well-known PRg=~154-day Rieger period and detected commonly in most heliophysical data types and the interplanetary magnetic field (IMF). Thus, InSight mission marsquakes rates are periodic with PRg as characterized by a very high (>>12) fidelity Φ=2.8·106 and by being the only 99%-significant spectral peak in the 385.8–64.3-nHz (1–180-day) band of highest planetary energies; the longest-span (v.9) release of raw data revealed the entire RR, excluding a tectonically active Mars. For check, I analyze rates of Oct 2015–Feb 2019, Mw5.6+ earthquakes, and all (1969–1977) Apollo mission moonquakes. To decouple magnetosphere and IMF effects, I study Earth and Moon seismicity during traversals of the Earth magnetotail vs. IMF. The analysis showed with 99–67% confidence and Φ>>12 fidelity that (an unspecified majority of) moonquakes and Mw5.6+ earthquakes also recur at Rieger periods. About half of the spectral peaks split but also into clusters that average to the usual Rieger periodicities, where magnetotail reconnecting clears the signal. Earlier claims that solar plasma dynamics could be seismogenic due to electrical surging or magnetohydrodynamic interactions between magnetically trapped plasma and water molecules embedded within solid matter are confirmed. This result calls for reinterpreting the seismicity phenomenon and for reliance on global magnitude scales. The predictability of solar-wind macroscopic dynamics is now within reach for the first time, which will benefit seismic and weather prediction and the safety of space missions., HIGHLIGHTS: • Seismicity on Mars, Moon, & Earth is caused mostly externally due to macroscopic dynamics of the solar wind • Seismicity phenomenon is also astrophysical in origin instead of being exclusively geophysical as believed by some • Marsquakes are almost exclusively forced by solar wind, confirming a general view that Mars is tectonically inactive • Moonquakes are mostly forced at solar-wind resonance times and not just during tides, which was a simplistic view • Connection between geomagnetism and strong seismicity, known for over a century, is explained for the first time • The first application of the Gauss-Vaníček method for rigorous spectral analysis in planetary seismology. New version: https://doi.org/10.48550/arXiv.2301.10800, {"references":["Abreu, J.A., Beer, J., Ferriz-Mas, A., McCracken, K.G., Steinhilber, F. (2012) Is there a planetary influence on solar activity? Astron. Astrophys. 548:A88. https://doi.org/10.1051/0004-6361/201219997","Alfvén, H. (1942) Existence of electromagnetic-hydrodynamic waves. Nature 150(3805):405–406. https://doi.org/10.1038%2F150405d0","Bai, T. (2003) Hot Spots for Solar Flares Persisting for Decades: Longitude Distributions of Flares of Cycles 19–23. Astrophys. J. 585:1114-1123. https://dx.doi.org/10.1086/346152","Bai T. and Cliver E. W. (1990) A 154 day periodicity in the occurrence rate of proton flares. Astrophys. J. 363:299–309. https://doi.org/10.1086/169342","Banerdt, W.B., Smrekar, S.E., Banfield, D. et al. (2020) Initial results from the InSight mission on Mars. Nat. Geosci. 13:183–189. https://doi.org/10.1038/s41561-020-0544-y","Bulow, R.C., Johnson, C.L., Bills, B.G., Shearer, P.M. (2007) Temporal and spatial properties of some deep moonquake clusters. J. Geophys. Res. 112:E09003. https://doi.org/10.1029/2006JE002847","Burrell, A.G., Sánchez-Cano, B., Witasse, O., Lester, M., Cartacci, M. (2020) Comparison of terrestrial and Martian TEC at dawn and dusk during solstices. Earth Planets Space 72:140 (2020). https://doi.org/10.1186/s40623-020-01258-3","Cane, H.V., Richardson, I.G., von Rosenvinge, T.T. (1998) Interplanetary magnetic field periodicity of ∼153 days. Geophys. Res. Lett. 25(24):4437–4440. https://doi.org/10.1029/1998GL900208","Ceylan, S., Clinton, J.F., Giardini, D. et al. (2021a) Companion guide to the marsquake catalog from InSight, Sols 0–478: Data content and non-seismic events. Phys. Earth Planet. Inter. 310:106597. https://doi.org/10.1016/j.pepi.2020.106597","Ceylan, S., Clinton, J.F., Horleston, A., et al. (2021b) The Seismicity on Mars as Recorded by Insight's Marsquake Service. Presentation at the Annual Meeting of the Seismological Society of America, 23 Apr. https://seismosoc.secure-platform.com/a/solicitations/24/sessiongallery/380/application/6767, https://www.seismosoc.org/news/seismicity-on-mars-full-of-surprises-in-first-continuous-year-of-data-collection","Chowdhury, P., Kudela, K. Moon, Y.J. (2016) A Study of heliospheric modulation and periodicities of Galactic cosmic rays during cycle 24. Sol. Phys. 291:581–602. https://doi.org/10.1007/s11207-015-0832-7","Chowdhury, P., Khan, M., Ray, P.C. (2009) Intermediate-term periodicities in sunspot areas during solar cycles 22 and 23. Mon. Not. R. Astron. Soc. 392(1):1159–1180. https://doi.org/10.1111%2Fj.1365-2966.2008.14117.x","Chowdhury, P., Ray, P.C., Ray, S. (2008) Periodicity of ~155 days in solar electron fluence. Indian J. Phys. 82:95–104.","Clinton, J.F., Ceylan, S., van Driel, M., Giardini, D. et al., (2021) The Marsquake catalogue from InSight, sols 0–478. Phys. Earth Planet. Inter. 310:106595 (corrected Supplementary Table 2). https://doi.org/10.1016/j.pepi.2020.106595","Connerney, J.E.P., Acuna, M.H., Wasilewski, P.J., Kletetschka, G., Ness, N.F., Reme, H., Lin, R.P., Mitchell, D.L. (2001) The global magnetic field of Mars and implications for crustal evolution. Geophys. Res. Lett. 28(21):4015–4018. https://doi.org/10.1029/2001GL013619","Craymer, M.R. (1998) The Least Squares Spectrum, Its Inverse Transform and Autocorrelation Function: Theory and Some Applications in Geodesy. Ph.D. Dissertation, University of Toronto, Canada. URI: https://hdl.handle.net/1807/12263","Das, R., Meneses, C. (2021) Scaling relations for energy magnitudes. J. Geophys. 64(1):1–11. https://n2t.net/ark:/88439/x063005","Das, R., Sharma, M.L., Wason, H.R., Choudhury, D., Gonzalez, G. (2019) A seismic moment magnitude scale. Bull. Seismol. Soc. Am. 109(4):1542–1555. https://doi.org/10.1785/0120180338","Den Hartog, J. (1985) Mechanical Vibrations (4th ed.) Dover Publications, USA. ISBN 9780486647852","Denskat, K.U., Beinroth, H.J., Neubauer, F.M. (1983) Interplanetary Magnetic Field Power Spectra With Frequencies from 2.4 X 10 to the-5th Hz to 470 Hz from HELIOS-Observations During Solar Minimum Conditions. J. Geophys. 54(1):60–67. https://n2t.net/ark:/88439/y063880","Dimitropoulou, M., Moussas, X., Strintzi, D. (2008) Enhanced Rieger type periodicities' detection in X-ray solar flares and statistical validation of Rossby waves' existence. Proc. Int. Astron. Union 4(S257):159–163. https://doi.org/10.1017/S1743921309029226","van Driel, M., Ceylan, S., Clinton, J. F., Giardini, D., Horleston, A., Margerin, L., et al. (2021) High‐frequency seismic events on Mars observed by InSight. J. Geophys. Res. Planets 126:e2020JE006670. https://doi.org/10.1029/2020JE006670","Dziewonski, A.M., Chou, T.-A., Woodhouse, J.H. (1981) Determination of earthquake source parameters from wave-form data for studies of global and regional seismicity. J. Geophys. Res. Solid Earth 86(B4):2825–2852. https://doi.org/10.1029/JB086iB04p02825","Ferrazzini, V., Aki, K. (1987) Slow waves trapped in a fluid‐filled infinite crack: Implication for volcanic tremor. J. Geophys. Res. 92(B9):9215–9223. https://doi.org/10.1029/JB092iB09p09215","Forgács–Dajka, E., Borkovits, T. (2007) Searching for mid-term variations in different aspects of solar activity – looking for probable common origins and studying temporal variations of magnetic polarities. Mon. Not. R. Astron. Soc. 374:282–291. https://doi.org/doi:10.1111/j.1365-2966.2006.11167.x","Frohlich, C., Nakamura, Y. (2009) The physical mechanisms of deep moonquakes and intermediate-depth earthquakes: How similar and how different? Phys. Earth Planet. Inter. 173(3–4):365–374 https://doi.org/10.1016/j.pepi.2009.02.004","Gupta, H.K. (Ed.) (2011) Encyclopedia of Solid Earth Geophysics. Springer. ISBN 9789048187010. https://doi.org/10.1007/978-90-481-8702-7","Gurgenashvili, E., Zaqarashvili, T.V., Kukhianidze, V., Oliver, R., Ballester, J.L., Dikpati, M., McIntosh, S.W. (2017) North–South Asymmetry in Rieger-type Periodicity during Solar Cycles 19–23. Astrophys. J. 845(2):137–148. https://dx.doi.org/10.3847/1538-4357/aa830a","Halekas, J.S., Poppe, A.R., Harada, Y., Bonnell, J.W., Ergun, R.E., McFadden, J.P. (2018) A tenuous lunar ionosphere in the geomagnetic tail. Geophys. Res. Lett. 45:9450–9459. https://doi.org/10.1029/2018GL079936","Harada, Y., Machida, S., Saito, Y., Yokota, S., Asamura, K. et al. (2010) Interaction between terrestrial plasma sheet electrons and the lunar surface: SELENE (Kaguya) observations, Geophys. Res. Lett. 37:L19202. https://doi.org/10.1029/2010GL044574","Jacoby, W.R. (2001) Successes and failures in geodynamics: from past to future. J. Geodyn. 32:3–27. https://doi.org/10.1016/S0264-3707(01)00026-6","Kanamori, H. (1977) The energy release in great earthquakes. J. Geophys. Res. Solid Earth Planets 82(29):2981–2987. https://doi.org/10.1029/JB082i020p02981","Kinkhabwala, A. (2013) Maximum Fidelity. Max Planck Institute of Molecular Physiology report. arXiv:1301.5186, Subject Field(s): Statistics Theory (math.ST); Instrumentation and Methods for Astrophysics (astro-ph.IM); High Energy Physics - Experiment (hep-ex). https://doi.org/10.48550/arXiv.1301.5186","Knapmeyer, M., Stähler, S.C., Daubar, I., Forget, F., Spiga, A., et al. (2021) Seasonal seismic activity on Mars. Earth Planet. Sci. Lett. 576:117171. https://doi.org/10.1016/j.epsl.2021.117171","Kundt, W., Jessner, A. (1986) Volcanoes, fountains, earthquakes, and continental motion — What causes them? J. Geophys. 60(1):33-40. https://n2t.net/ark:/88439/y049522","Kurochkin, N.E. (1998) Transient periodicity in solar activity. Astron. Astrophys. Trans. 15(1-4):277–279. https://doi.org/10.1080/10556799808201781","LaViolette, P.A. (1983) Galactic explosions, cosmic dust invasions, and climatic change. Ph.D. dissertation, Portland State University, pp.737. https://search.library.pdx.edu/permalink/f/p82vj0/CP71128014660001451","Li, W., Raeder, J., Thomsen, M.F., Lavraud, B. (2008) Solar wind plasma entry into the magnetosphere under northward IMF conditions. J. Geophys. Res. 113:A04204. https://doi.org/10.1029/2007JA012604","Lorenz, R.D., Panning, M. (2018) Empirical recurrence rates for ground motion signals on planetary surfaces. Icarus 303:273–279. https://doi.org/10.1016/j.icarus.2017.10.008","Lorenz, R.D., Nakamura, Y., Murphy, J.R. (2017) Viking‐2 seismometer measurements on Mars: PDS data archive and meteorological applications. Earth Space Sci. 4:681–688. https://doi.org/10.1002/2017EA000306","Love, J.J., Thomas, J.N. (2013) Insignificant solar-terrestrial triggering of earthquakes. Geophys. Res. Lett. 40:1165–1170. https://doi.org/10.1002/grl.50211","Nakamura, Y., Latham, G.V., Dorman, H.J., Harris, J.E. (1981) Passive Seismic Experiment, Long Period Event Catalog, Final Version. Institute for Geophysics, University of Texas at Austin. https://hdl.handle.net/2152/65671","Odintsov, S., Boyarchuk, K., Georgieva, K., Kirov, B., Atanasov, D. (2006) Long-period trends in global seismic and geomagnetic activity and their relation to solar activity. Phys. Chem. Earth 31(1–3):88–93. https://doi.org/10.1016/j.pce.2005.03.004","Ogawa, M. (2007) Mantle convection: A review. Fluid Dyn. Res. 40(6):379–398. https://doi.org/10.1016/j.fluiddyn.2007.09.001","Omerbashich, M. (2023) Earth as a time crystal: macroscopic nature of a quantum-scale phenomenon from transformative moderation of geomagnetic polarity, topography, and climate by precession resonance due to many-body entrainment. arXiv:2301.02578, Subject Field(s): Geophysics (physics.geo-ph); Earth & Planetary Astrophysics (astro-ph.EP); Chaotic Dynamics (nlin.CD); Quantum Physics (quant-ph). https://doi.org/10.48550/arXiv.2301.02578","Omerbashich, M. (2022) Detection and mapping of Earth body resonances with continuous GPS. J. Geophys. 64:12-33. https://n2t.net/ark:/88439/x073994","Omerbashich, M. (2021) Non-marine tetrapod extinctions solve extinction periodicity mystery. Hist. Biol. 34(1):188-191. https://doi.org/10.1080/08912963.2021.1907367","Omerbashich, M. (2020a) Earth body resonance. J. Geophys. 63:15-29. https://n2t.net/ark:/88439/x020219","Omerbashich, M. (2020b) Moon body resonance. J. Geophys. 63:30-42. https://n2t.net/ark:/88439/x034508","Omerbashich M. (2007) Magnification of mantle resonance as a cause of tectonics. Geod. Acta. 20(6):369–383. https://doi.org/10.3166/ga.20.369-383","Omerbashich, M. (2006a) Springtide-induced magnification of Earth mantle resonance causes tectonics and conceals universality of physics at all scales. arXiv. https://doi.org/10.48550/arXiv.physics/0608026","Omerbashich, M. (2006b) Gauss–Vaníček Spectral Analysis of the Sepkoski Compendium: No New Life Cycles. Comp. Sci. Eng. 8(4):26–30. https://doi.org/10.1109/MCSE.2006.68 (Erratum due to journal error, Comp. Sci. Eng. 9(4):5–6. https://doi.org/10.1109/MCSE.2007.79; full text: https://arxiv.org/abs/math-ph/0608014)","Omerbashich, M. (2004) Earth-model Discrimination Method. Ph.D. Dissertation, pp.129. ProQuest, USA. DOI: https://doi.org/10.6084/m9.figshare.12847304","Özgüç, A., Ataç, T. (1994) The 73-day periodicity of the flare index during the current solar cycle 22. Sol. Phys. 150:339–346. https://doi.org/10.1007/BF00712895","Pagiatakis, S. (1999) Stochastic significance of peaks in the least-squares spectrum. J. Geod. 73:67–78. https://doi.org/10.1007/s001900050220","Phillips, T. (2008) The Moon and the Magnetotail. NASA's Goddard Space Flight Center, NASA Space News, 16 Apr. http://web.archive.org/web/20210225013348/https://www.nasa.gov/topics/moonmars/features/magnetotail_080416.html","Press, W.H., Teukolsky, S.A., Vetterling, W.T., Flannery, B.P. (2007) Numerical Recipes: The Art of Scientific Computing (3rd Ed.). Cambridge University Press. ISBN 9780521880688","Ramstad, R., Brain, D.A., Dong, Y., Espley, J., Halekas, J., Jakosky, B. (2020) The global current systems of the Martian induced magnetosphere. Nat. Astron. 4:979–985. https://doi.org/10.1038/s41550-020-1099-y","Richter, F., McKenzie, D. (1977) Simple plate models of mantle convection. J. Geophys. 44(1):441-471. https://n2t.net/ark:/88439/y001916","Rieger, E., Share, G.H., Forrest, D.J., Kanbach, G., Reppin, C., Chupp, E.L. (1984) A 154-day periodicity in the occurrence of hard solar flares? Nature 312:623–625. https://doi.org/10.1038/312623a0","Ringler, A.T., Anthony, R.E., Wilson, D.C., Claycomb, A.C., Spritzer, J. (2020) Magnetic Field Variations in Alaska: Recording Space Weather Events on Seismic Stations in Alaska. Bull. Seismol. Soc. Am. 110(5):2530–2540. https://doi.org/10.1785/0120200019","Russell, C.T. (2001) The dynamics of planetary magnetospheres. Planet. Space Sci. 49(10–11):1005–1030. https://doi.org/10.1016/S0032-0633(01)00017-4","Shang, W.S., Tang, B.B., Shi, Q.Q., Tian, A.M., Zhou, X.-Y., Yao, Z.H., et al. (2020) Unusual location of the geotail magnetopause near lunar orbit: A case study. J. Geophys. Res. Space Phys. 125:e2019JA027401. https://doi.org/10.1029/2019JA027401","Shannon, C.E. (1948) A Mathematical Theory of Communication. Bell System Tech. J. 27:379–423, 623–656. https://doi.org/10.1002/j.1538-7305.1948.tb01338.x","Simpson, J.F. (1968) Solar activity as a triggering mechanism for earthquakes. Earth Planet. Sci. Lett. 3:417–425. https://doi.org/10.1016/0012-821X(67)90071-4","Singh, Y.P., Badruddin (2019) Study of the solar rotational period and its harmonics in solar activity, interplanetary, geomagnetic, and cosmic ray intensity indicators during solar polarity reversal periods. Sol. Phys. 294:27. https://doi.org/10.1007/s11207-019-1413-y","Stähler, S.C., Mittelholz, A., Perrin, C. et al. (2022) Tectonics of Cerberus Fossae unveiled by marsquakes. Nat. Astron. https://doi.org/10.1038/s41550-022-01803-y","Steeves, R.R. (1981). A statistical test for significance of peaks in the least squares spectrum. Collected Papers, Geodetic Sur-vey, Department of Energy, Mines and Resources. Surveys and Mapping Branch, Ottawa Canada, pp. 149–166. http://www2.unb.ca/gge/Research/GRL/LSSA/Literature/Steeves1981.pdf","Stevenson, D.J. (2008) A planetary perspective on the deep Earth. Nature 451:261‐265. https://doi.org/10.1038/nature06582","Suemoto, Y., Ikeda, T., Tsuji, T. (2020) Temporal variation and frequency dependence of seismic ambient noise on Mars from polarization analysis. Geophys. Res. Lett. 47:e2020GL087123. https://doi.org/10.1029/2020GL087123","Tape, C., Ringler; A.T., Hampton, D.L. (2020) Recording the Aurora at Seismometers across Alaska. Seismol. Res. Lett. (2020) 91(6):3039–3053. https://doi.org/10.1785/0220200161","Taylor, J., Hamilton, S. (1972) Some tests of the Vaníček Method of spectral analysis. Astrophys. Space Sci. 17:357–367. https://doi.org/10.1007/BF00642907","Trotignon, J.G, Parrot, M., Cerisier, J.C., Menvielle, M., Axford, W.I., Paetzold, M., Warnant, R., Wernik, A.W. (2000) The plasma environment of Mars: from the shocked solar wind down to the ionosphere. Planet. Space Sci. 48(12–14):1181–1191. https://doi.org/10.1016/S0032-0633(00)00102-1","Tucker, O.J., Farrell, W.M., Poppe, A.R. (2021) On the effect of magnetospheric shielding on the lunar hydrogen cycle. J. Geophys. Res. Planets 126:e2020JE006552. https://doi.org/10.1029/2020JE006552","Tzamali, M., Pagiatakis, S. (2021) GRACE-FO accelerometer data: An alternative approach using Least Squares Spectral Analysis. EGU General Assembly, 19–30 Apr. https://doi.org/10.5194/egusphere-egu21-6906","Vaníček, P. (1969) Approximate Spectral Analysis by Least-Squares Fit. Astrophys. Space Sci. 4(4):387–391. https://doi.org/10.1007/BF00651344","Vaníček, P. (1971) Further Development and Properties of the Spectral Analysis by Least-Squares Fit. Astrophys. Space Sci. 12(1):10–33. https://doi.org/10.1007/BF00656134","Walsh, A.P., Haaland, S., Forsyth, C., Keesee A.M., Kissinger, J., Li, K., Runov, A., Soucek, J., Walsh, B.M., Wing, S., Taylor, M.G.G.T. (2014) Dawn–dusk asymmetries in the coupled solar wind–magnetosphere–ionosphere system: a review. Ann. Geophys. 32:705–737. https://doi.org/10.5194/angeo-32-705-2014","Wang, H.Z., Zhang, J., Shi, Q.Q., Saito, Y., Degeling, A.W., Rae, I.J., Zong, Q.G., Wei, Y., Liu, J., Guo, R.L. (2021) Earth Wind as a Possible Exogenous Source of Lunar Surface Hydration. Astrophys. J. Lett. 907(2):L32. https://doi.org/10.3847/2041-8213/abd559","Wells, D.E., Vaníček, P., Pagiatakis, S. (1985) Least squares spectral analysis revisited. Department of Geodesy & Geomatics Engineering Technical Report 84, University of New Brunswick, Canada. http://www2.unb.ca/gge/Pubs/TR84.pdf","Yegorkin, A.V., Chernyshov, N.M. (1983) Peculiarities of Mantle Waves from Long-Range Profiles. J. Geophys. 54(1):30-34. https://n2t.net/ark:/88439/y090547"]}

Semi-implicit time-stepping schemes for atmosphere and ocean models require elliptic solvers that work efficiently on modern supercomputers. This paper reports our study of the potential computational savings when using mixed precision arithmetic in the elliptic solvers. The essential components of a representative elliptic solver are run at precision levels as low as half (16 bits), and accompanied with a detailed evaluation of the impact of reduced precision on the solver convergence and the solution quality. A detailed inquiry into reduced precision requires a model configuration that is meaningful but cheaper to run and easier to evaluate than full atmosphere/ocean models. This study is therefore conducted in the context of a novel semi-implicit shallow-water model on the sphere, purposely designed to mimic numerical intricacies of modern all-scale weather and climate (W&C) models with the numerical stability independent on celerity of all wave motions. The governing algorithm of the shallow-water model is based on the non-oscillatory MPDATA methods for geophysical flows, whereas the resulting elliptic problem employs a strongly preconditioned non-symmetric Krylov-subspace solver GCR, proven in advanced atmospheric applications. The classical longitude/latitude grid is deliberately chosen to retain the stiffness of global W&C models posed in thin spherical shells as well as to better understand the performance of reduced-precision arithmetic in the vicinity of grid singularities. Precision reduction is done on a software level, using an emulator. The reduced-precision experiments are conducted for established dynamical-core test-cases, like the Rossby-Haurwitz wave number 4 and a zonal orographic flow. The study shows that selected key components of the elliptic solver, most prominently the preconditioning, can be performed at the level of half precision., 38 pages, 12 figures; to be submitted to the Journal of Computational Physics (JCP)

The moving discontinuous Galerkin finite element method with interface condition enforcement (MDG-ICE) is applied to the case of viscous flows. This method uses a weak formulation that separately enforces the conservation law, constitutive law, and the corresponding interface conditions in order to provide the means to detect interfaces or under-resolved flow features. To satisfy the resulting overdetermined weak formulation, the discrete domain geometry is introduced as a variable, so that the method implicitly fits a priori unknown interfaces and moves the grid to resolve sharp, but smooth, gradients, achieving a form of anisotropic curvilinear $r$-adaptivity. This approach avoids introducing low-order errors that arise using shock capturing, artificial dissipation, or limiting. The utility of this approach is demonstrated with its application to a series of test problems culminating with the compressible Navier-Stokes solution to a Mach 5 viscous bow shock for a Reynolds number of $10^{5}$ in two-dimensional space. Time accurate solutions of unsteady problems are obtained via a space-time formulation, in which the unsteady problem is formulated as a higher dimensional steady space-time problem. The method is shown to accurately resolve and transport viscous structures without relying on numerical dissipation for stabilization., Comment: 29 pages, 11 figures

A stress equilibration procedure for linear elasticity is proposed and analyzed in this paper with emphasis on the behavior for (nearly) incompressible materials. Based on the displacement-pressure approximation computed with a stable finite element pair, it constructs an $H (\text{div})$-conforming, weakly symmetric stress reconstruction. Our focus is on the Taylor-Hood combination of continuous finite element spaces of polynomial degrees $k+1$ and $k$ for the displacement and the pressure, respectively. Our construction leads then to reconstructed stresses by Raviart-Thomas elements of degree $k$ which are weakly symmetric in the sense that its anti-symmetric part is zero tested against continuous piecewise polynomial functions of degree $k$. The computation is performed locally on a set of vertex patches covering the computational domain in the spirit of equilibration \cite{BraSch:08}. Due to the weak symmetry constraint, the local problems need to satisfy consistency conditions associated with all rigid body modes, in contrast to the case of Poisson's equation where only the constant modes are involved. The resulting error estimator is shown to constitute a guaranteed upper bound for the error with a constant that depends only on the shape regularity of the triangulation. Local efficiency, uniformly in the incompressible limit, is deduced from the upper bound by the residual error estimator.

We present a class of new explicit and stable numerical algorithms to solve the spatially discretized linear heat or diffusion equation. After discretizing the space and the time variables like conventional finite difference methods, we do not approximate the time derivatives by finite differences, but use constant neighbor and linear neighbour approximations to decouple the ordinary differential equations and solve them analytically. During this process, the timestep-size appears not in polynomial, but in exponential form with negative exponents, which guarantees stability. We compare the performance of the new methods with analytical and numerical solutions. According to our results, the methods are first and second order in time and can be much faster than the commonly used explicit or implicit methods, especially in the case of extremely large stiff systems., 21 pages

In fractured natural formations, the equations governing fluid flow and geomechanics are strongly coupled. Hydrodynamical properties depend on the mechanical configuration, and they are therefore difficult to accurately resolve using uncoupled methods. In recent years, significant research has focused on discretization strategies for these coupled systems, particularly in the presence of complicated fracture network geometries. In this work, we explore a finite-volume discretization for the multiphase flow equations coupled with a finite-element scheme for the mechanical equations. Fractures are treated as lower dimensional surfaces embedded in a background grid. Interactions are captured using the Embedded Discrete Fracture Model (EDFM) and the Embedded Finite Element Method (EFEM) for the flow and the mechanics, respectively. This non-conforming approach significantly alleviates meshing challenges. EDFM considers fractures as lower dimension finiten volumes which exchange fluxes with the rock matrix cells. The EFEM method provides, instead, a local enrichment of the finite-element space inside each matrix cell cut by a fracture element. Both the use of piecewise constant and piecewise linear enrichments are investigated. They are also compared to an Extended Finite Element (XFEM) approach. One key advantage of EFEM is the element-based nature of the enrichment, which reduces the geometric complexity of the implementation and leads to linear systems with advantageous properties. Synthetic numerical tests are presented to study the convergence and accuracy of the proposed method. It is also applied to a realistic scenario, involving a heterogeneous reservoir with a complex fracture distribution, to demonstrate its relevance for field applications.

For most finite element simulations, boundary-conforming meshes have significant advantages in terms of accuracy or efficiency. This is particularly true for complex domains. However, with increased complexity of the domain, generating a boundary-conforming mesh becomes more difficult and time consuming. One might therefore decide to resort to an approach where individual boundary-conforming meshes are pieced together in a modular fashion to form a larger domain. This paper presents a stabilized finite element formulation for fluid and temperature equations on sliding meshes. It couples the solution fields of multiple subdomains whose boundaries slide along each other on common interfaces. Thus, the method allows to use highly tuned boundary-conforming meshes for each subdomain that are only coupled at the overlapping boundary interfaces. In contrast to standard overlapping or fictitious domain methods the coupling is broken down to few interfaces with reduced geometric dimension. The formulation consists of the following key ingredients: the coupling of the solution fields on the overlapping surfaces is imposed weakly using a stabilized version of Nitsche's method. It ensures mass and energy conservation at the common interfaces. Additionally, we allow to impose weak Dirichlet boundary conditions at the non-overlapping parts of the interfaces. We present a detailed numerical study for the resulting stabilized formulation. It shows optimal convergence behavior for both Newtonian and generalized Newtonian material models. Simulations of flow of plastic melt inside single-screw as well as twin-screw extruders demonstrate the applicability of the method to complex and relevant industrial applications.

We investigate the temporal accuracy of two generalized-$\alpha$ schemes for the incompressible Navier-Stokes equations. The conventional approach treats the pressure with the backward Euler method while discretizing the remainder of the Navier-Stokes equations with the generalized-$\alpha$ method. We developed a suite of numerical codes using inf-sup stable higher-order non-uniform rational B-spline (NURBS) elements for spatial discretization. In doing so, we are able to achieve very high spatial accuracy and, furthermore, to perform temporal refinement without consideration of the stabilization terms, which can degenerate for small time steps. Numerical experiments suggest that only first-order accuracy is achieved, at least for the pressure, in this aforesaid approach. Evaluating the pressure at the intermediate time step recovers second-order accuracy, and the numerical implementation, in fact, becomes simpler. Therefore, although the pressure can be viewed as a Lagrange multiplier enforcing the incompressibility constraint, its temporal discretization is not independent and should be subject to the generalized-$\alpha$ method in order to maintain second-order accuracy of the overall algorithm.

This paper is a rebuttal to the claim found in the literature that the MUSCL scheme cannot be third-order accurate for nonlinear conservation laws. We provide a rigorous proof for third-order accuracy of the MUSCL scheme based on a careful and detailed truncation error analysis. Throughout the analysis, the distinction between the cell average and the point value will be strictly made for the numerical solution as well as for the target operator. It is shown that the average of the solutions reconstructed at a face by Van Leer's kappa-scheme recovers a cubic solution exactly with kappa = 1/3, the same is true for the average of the nonlinear fluxes evaluated by the reconstructed solutions, and a dissipation term is already sufficiently small with a third-order truncation error. Finally, noting that the target spatial operator is a cell-averaged flux derivative, we prove that the leading truncation error of the MUSCL finite-volume scheme is third-order with kappa = 1/3. The importance of the diffusion scheme is also discussed: third-order accuracy will be lost when the third-order MUSLC scheme is used with a wrong fourth-order diffusion scheme for convection-diffusion problems. Third-order accuracy is verified by thorough numerical experiments for both steady and unsteady problems. This paper is intended to serve as a reference to clarify confusions about third-order accuracy of the MUSCL scheme, as a guide to correctly analyze and verify the MUSCL scheme for nonlinear equations, and eventually as the basis for clarifying third-order unstructured-grid schemes in a subsequent paper.

The Fragile Points Method (FPM) is an elementarily simple Galerkin meshless method, employing Point-based discontinuous trial and test functions only, without using element-based trial and test functions. In this study, the algorithmic formulations of FPM for linear elasticity are given in detail, by exploring the concepts of point stiffness matrices and numerical flux corrections. Advantages of FPM for simulating the deformations of complex structures, and for simulating complex crack propagations and rupture developments, are also thoroughly discussed. Numerical examples of deformation and stress analyses of benchmark problems, as well as of realistic structures with complex geometries, demonstrate the accuracy, efficiency and robustness of the proposed FPM. Simulations of crack initiation and propagations are also given in this study, demonstrating the advantages of the present FPM in modeling complex rupture and fracture phenomena. The crack and rupture propagation modeling in FPM is achieved without remeshing or augmenting the trial functions as in standard, extended or generalized FEM. The simulation of impact, penetration and other extreme problems by FPM will be discussed in our future papers.

In this paper we solve $m$-parameter eigenvalue problems ($m$EPs), with $m$ any natural number by representing the problem using Tensor-Trains (TT) and designing a method based on this format. $m$EPs typically arise when separation of variables is applied to separable boundary value problems. Often, methods for solving $m$EP are restricted to $m = 3$, due to the fact that, to the best of our knowledge, no available solvers exist for $m>3$ and reasonable size of the involved matrices. In this paper, we prove that computing the eigenvalues of a $m$EP can be recast into computing the eigenvalues of TT-operators. We adapted the algorithm in \cite{Dolgov2014a} for symmetric eigenvalue problems in TT-format to an algorithm for solving generic $m$EPs. This leads to a subspace method whose subspace dimension does not depend on $m$, in contrast to other subspace methods for $m$EPS. This allows us to tackle $m$EPs with $m > 3$ and reasonable size of the matrices. We provide theoretical results and report numerical experiments. The MATLAB code is publicly available., 17 pages, 14 figures

Low precision arithmetic, in particular half precision floating point arithmetic, is now available in commercial hardware. Using lower precision can offer significant savings in computation and communication costs with proportional savings in energy. Motivated by this, there has been a renewed interest in mixed precision iterative refinement for solving linear systems $Ax=b$, and new variants of GMRES-based iterative refinement have been developed. Each particular variant with a given combination of precisions leads to different condition number-based constraints for convergence of the backward and forward errors, and each has different performance costs. The constraints for convergence given in the literature are, as an artifact of the analyses, often overly strict in practice, and thus could lead a user to select a more expensive variant when a less expensive one would have sufficed. In this work, we develop a multistage mixed precision iterative refinement solver which aims to combine existing mixed precision approaches to balance performance and accuracy and improve usability. For an initial combination of precisions, the algorithm begins with the least expensive approach and convergence is monitored via inexpensive computations with quantities produced during the iteration. If slow convergence or divergence is detected using particular stopping criteria, the algorithm switches to use a more expensive, but more reliable variant. A novel aspect of our approach is that, unlike existing implementations, our algorithm first attempts to use ``stronger'' solvers for the solution update before resorting to increasing the precision(s). In some scenarios, this can avoid the need to refactorize the matrix in higher precision. We perform extensive numerical experiments on random dense problems and problems from real applications which confirm the benefits of the multistage approach., Comment: 30 pages

The multilinear Pagerank model [Gleich, Lim and Yu, 2015] is a tensor-based generalization of the Pagerank model. Its computation requires solving a system of polynomial equations that contains a parameter $\alpha \in [0,1)$. For $\alpha \approx 1$, this computation remains a challenging problem, especially since the solution may be non-unique. Extrapolation strategies that start from smaller values of $\alpha$ and `follow' the solution by slowly increasing this parameter have been suggested; however, there are known cases where these strategies fail, because a globally continuous solution curve cannot be defined as a function of $\alpha$. In this paper, we improve on this idea, by employing a predictor-corrector continuation algorithm based on a more general representation of the solutions as a curve in $\mathbb{R}^{n+1}$. We prove several global properties of this curve that ensure the good behavior of the algorithm, and we show in our numerical experiments that this method is significantly more reliable than the existing alternatives.

The early detection of terrorist threat objects, such as guns and knives, through improved metal detection, has the potential to reduce the number of attacks and improve public safety and security. To achieve this, there is considerable potential to use the fields applied and measured by a metal detector to discriminate between different shapes and different metals since, hidden within the field perturbation, is object characterisation information. The magnetic polarizability tensor (MPT) offers an economical characterisation of metallic objects and its spectral signature provides additional object characterisation information. The MPT spectral signature can be determined from measurements of the induced voltage over a range frequencies in a metal signature for a hidden object. With classification in mind, it can also be computed in advance for different threat and non-threat objects. In the article, we evaluate the performance of probabilistic and non-probabilistic machine learning algorithms, trained using a dictionary of computed MPT spectral signatures, to classify objects for metal detection. We discuss the importances of using appropriate features and selecting an appropriate algorithm depending on the classification problem being solved and we present numerical results for a range of practically motivated metal detection classification problems.