A Normal Mode Stability Analysis of Numerical Interface Conditions for Fluid/Structure Interaction

In multi-physics computations where a compressible fluid is coupled with a linearly elastic solid, it is standard to enforce continuity of the normal velocities and of the normal stresses at the interface between the fluid and the solid. In a numerical scheme, there are many ways that velocity- and stress-continuity can be enforced in the discrete approximation. This paper performs a normal mode stability analysis of the linearized problem to investigate the stability of different numerical interface conditions for a model problem approximated by upwind type finite difference schemes. The analysis shows that depending on the ratio of densities between the solid and the fluid, some numerical interface conditions are stable up to the maximal CFL-limit, while other numerical interface conditions suffer from a severe reduction of the stable CFL-limit. The paper also presents a new interface condition, obtained as a simplified characteristic boundary condition, that is proved to not suffer from any reduction of the stable CFL-limit. Numerical experiments in one space dimension show that the new interface condition is stable also for computations with the non-linear Euler equations of compressible fluid flow coupled with a linearly elastic solid.     

An Edge-Based Smoothed Finite Element Method with TBC for the Elastic Wave Scattering by an Obstacle

Elastic wave scattering has received ever-increasing attention in military and medical fields due to its high-precision solution. In this paper, an edge-based smoothed finite element method (ES-FEM) combined with the transparent boundary condition (TBC) is proposed to solve the elastic wave scattering problem by a rigid obstacle with smooth surface, which is embedded in an isotropic and homogeneous elastic medium in two dimensions. The elastic wave scattering problem satisfies Helmholtz equations with coupled boundary conditions obtained by Helmholtz decomposition. Firstly, the TBC of the elastic wave scattering is constructed by using the analytical solution to Helmholtz equations, which can truncate the boundary value problem (BVP) in an unbounded domain into the BVP in a bounded domain. Then the formulations of ES-FEM with the TBC are derived for Helmholtz equations with coupled boundary conditions. Finally, several numerical examples illustrate that the proposed ES-FEM with the TBC (ES-FEM-TBC) can work effectively and obtain more stable and accurate solution than the standard FEM with the TBC (FEM-TBC) for the elastic wave scattering problem.     

Super-Grid Modeling of the Elastic Wave Equation in Semi-Bounded Domains

We develop a super-grid modeling technique for solving the elastic wave equation in semi-bounded two- and three-dimensional spatial domains. In this method, waves are slowed down and dissipated in sponge layers near the far-field boundaries. Mathematically, this is equivalent to a coordinate mapping that transforms a very large physical domain to a significantly smaller computational domain, where the elastic wave equation is solved numerically on a regular grid. To damp out waves that become poorly resolved because of the coordinate mapping, a high order artificial dissipation operator is added in layers near the boundaries of the computational domain. We prove by energy estimates that the super-grid modeling leads to a stable numerical method with decreasing energy, which is valid for heterogeneous material properties and a free surface boundary condition on one side of the domain. Our spatial discretization is based on a fourth order accurate finite difference method, which satisfies the principle of summation by parts. We show that the discrete energy estimate holds also when a centered finite difference stencil is combined with homogeneous Dirichlet conditions at several ghost points outside of the far-field boundaries. Therefore, the coefficients in the finite difference stencils need only be boundary modified near the free surface. This allows for improved computational efficiency and significant simplifications of the implementation of the proposed method in multi-dimensional domains. Numerical experiments in three space dimensions show that the modeling error from truncating the domain can be made very small by choosing a sufficiently wide super-grid damping layer. The numerical accuracy is first evaluated against analytical solutions of Lamb's problem, where fourth order accuracy is observed with a sixth order artificial dissipation. We then use successive grid refinements to study the numerical accuracy in the more complicated motion due to a point moment tensor source in a regularized layered material.     

Operator Factorization for Multiple-Scattering Problems and an Application to Periodic Media

This work concerns multiple-scattering problems for time-harmonic equations in a reference generic media. We consider scatterers that can be sources, obstacles or compact perturbations of the reference media. Our aim is to restrict the computational domain to small compact domains containing the scatterers. We use Robin-to-Robin (RtR) operators (in the most general case) to express boundary conditions for the interior problem. We show that one can always factorize the RtR map using only operators defined using single-scatterer problems. This factorization is based on a decomposition of the diffracted field, on the whole domain where it is defined. Assuming that there exists a good method for solving single-scatterer problems, it then gives a convenient way to compute RtR maps for a random number of scatterers.     

Analysis of High-Order Absorbing Boundary Conditions for the Schrödinger Equation

The paper is concerned with the numerical solution of Schrödinger equations on an unbounded spatial domain. High-order absorbing boundary conditions for one-dimensional domain are derived, and the stability of the reduced initial boundary value problem in the computational interval is proved by energy estimate. Then a second order finite difference scheme is proposed, and the convergence of the scheme is established as well. Finally, numerical examples are reported to confirm our error estimates of the numerical methods.     

Application of the LS-STAG Immersed Boundary/ Cut-Cell Method to Viscoelastic Flow Computations

This paper presents the extension of a well-established Immersed Boundary (IB)/cut-cell method, the LS-STAG method (Y. Cheny & O. Botella, J. Comput. Phys. Vol. 229, 1043-1076, 2010), to viscoelastic flow computations in complex geometries. We recall that for Newtonian flows, the LS-STAG method is based on the finite-volume method on staggered grids, where the IB boundary is represented by its level-set function. The discretization in the cut-cells is achieved by requiring that global conservation properties equations be satisfied at the discrete level, resulting in a stable and accurate method and, thanks to the level-set representation of the IB boundary, at low computational costs.In the present work, we consider a general viscoelastic tensorial equation whose particular cases recover well-known constitutive laws such as the Oldroyd-B, White-Metzner and Giesekus models. Based on the LS-STAG discretization of the Newtonian stresses in the cut-cells, we have achieved a compatible velocity-pressure-stress discretization that prevents spurious oscillations of the stress tensor. Applications to popular benchmarks for viscoelastic fluids are presented: the four-to-one abrupt planar contraction flows with sharp and rounded re-entrant corners, for which experimental and numerical results are available. The results show that the LS-STAG method demonstrates an accuracy and robustness comparable to body-fitted methods.     

A Kernel-Free Boundary Integral Method for Variable Coefficients Elliptic PDEs

This work proposes a generalized boundary integral method for variable coefficients elliptic partial differential equations (PDEs), including both boundary value and interface problems. The method is kernel-free in the sense that there is no need to know analytical expressions for kernels of the boundary and volume integrals in the solution of boundary integral equations. Evaluation of a boundary or volume integral is replaced with interpolation of a Cartesian grid based solution, which satisfies an equivalent discrete interface problem, while the interface problem is solved by a fast solver in the Cartesian grid. The computational work involved with the generalized boundary integral method is essentially linearly proportional to the number of grid nodes in the domain. This paper gives implementation details for a second-order version of the kernel-free boundary integral method in two space dimensions and presents numerical experiments to demonstrate the efficiency and accuracy of the method for both boundary value and interface problems. The interface problems demonstrated include those with piecewise constant and large-ratio coefficients and the heterogeneous interface problem, where the elliptic PDEs on two sides of the interface are of different types.     

A Scalable Domain Decomposition Method for Ultra-Parallel Arterial Flow Simulations

Ultra-parallel flow simulations on hundreds of thousands of processors require new multi-level domain decomposition methods. Here we present such a new two-level method that has features both of discontinuous and continuous Galerkin formulations. Specifically, at the coarse level the domain is subdivided into several big patches and within each patch a spectral element discretization (fine level) is employed. New interface conditions for the Navier-Stokes equations are developed to connect the patches, relaxing the C⁰continuity and minimizing data transfer at the patch interface. We perform several 3D flow simulations of a benchmark problem and of arterial flows to evaluate the performance of the new method and investigate its accuracy.     

Quasi-Lagrangian Acceleration of Eulerian Methods

We present a simple and efficient strategy for the acceleration of explicit Eulerian methods for multidimensional hyperbolic systems of conservation laws. The strategy is based on the Galilean invariance of dynamic equations and optimization of the reference frame, in which the equations are numerically solved. The optimal reference frame moves (locally in time) with the average characteristic speed of the system, and, in this sense, the resulting method is quasi-Lagrangian. This leads to the acceleration of the numerical computations thanks to the optimal CFL condition and automatic adjustment of the computational domain to the evolving part of the solution. We show that our quasi-Lagrangian acceleration procedure may also reduce the numerical dissipation of the underlying Eulerian method. This leads to a significantly enhanced resolution, especially in the supersonic case. We demonstrate a great potential of the proposed method on a number of numerical examples.     

Numerical Simulation of Compressible Vortical Flows Using a Conservative Unstructured-Grid Adaptive Scheme

A two-dimensional numerical scheme for the compressible Euler equations is presented and applied here to the simulation of exemplary compressible vortical flows. The proposed approach allows to perform computations on unstructured moving grids with adaptation, which is required to capture complex features of the flow-field. Grid adaptation is driven by suitable error indicators based on the Mach number and by element-quality constraints as well. At the new time level, the computational grid is obtained by a suitable combination of grid smoothing, edge-swapping, grid refinement and de-refinement. The grid modifications—including topology modification due to edge-swapping or the insertion/deletion of a new grid node—are interpreted at the flow solver level as continuous (in time) deformations of suitably-defined node-centered finite volumes. The solution over the new grid is obtained without explicitly resorting to interpolation techniques, since the definition of suitable interface velocities allows one to determine the new solution by simple integration of the Arbitrary Lagrangian-Eulerian formulation of the flow equations. Numerical simulations of the steady oblique-shock problem, of the steady transonic flow and of the start-up unsteady flow around the NACA 0012 airfoil are presented to assess the scheme capabilities to describe these flows accurately.     

Benchmark Computations of the Phase Field Crystal and Functionalized Cahn-Hilliard Equations via Fully Implicit,Nesterov Accelerated Schemes

We introduce a fast solver for the phase field crystal (PFC) and functionalized Cahn-Hilliard (FCH) equations with periodic boundary conditions on a rectangular domain that features the preconditioned Nesterov’s accelerated gradient descent (PAGD) method. We discretize these problems with a Fourier collocation method in space, and employ various second-order schemes in time. We observe a significant speedup with this solver when compared to the preconditioned gradient descent (PGD) method. With the PAGD solver, fully implicit, second-order-in-time schemes are not only feasible to solve the PFC and FCH equations, but also do so more efficiently than some semi-implicit schemes in some cases where accuracy issues are taken into account. Benchmark computations of four different schemes for the PFC and FCH equations are conducted and the results indicate that, for the FCH experiments, the fully implicit schemes (midpoint rule and BDF2 equipped with the PAGD as a nonlinear time marching solver) perform better than their IMEX versions in terms of computational cost needed to achieve a certain precision. For the PFC, the results are not as conclusive as in the FCH experiments, which, we believe, is due to the fact that the nonlinearity in the PFC is milder nature compared to the FCH equation. We also discuss some practical matters in applying the PAGD. We introduce an averaged Newton preconditioner and a sweeping-friction strategy as heuristic ways to choose good preconditioner parameters. The sweeping-friction strategy exhibits almost as good a performance as the case of the best manually tuned parameters.     

A Finite-Volume Method for Nonlinear Nonlocal Equations with a Gradient Flow Structure

We propose a positivity preserving entropy decreasing finite volume scheme for nonlinear nonlocal equations with a gradient flow structure. These properties allow for accurate computations of stationary states and long-time asymptotics demonstrated by suitably chosen test cases in which these features of the scheme are essential. The proposed scheme is able to cope with non-smooth stationary states, different time scales including metastability, as well as concentrations and self-similar behavior induced by singular nonlocal kernels. We use the scheme to explore properties of these equations beyond their present theoretical knowledge.     

Source-Independent Amplitude-Semblance Full-Waveform Inversion Using a Hybrid Time- and Frequency-Domain Approach

Full-waveform inversion is a promising tool to produce accurate and high-resolution subsurface models. Conventional full-waveform inversion requires an accurate estimation of the source wavelet, and its computational cost is high. We develop a novel source-independent full-waveform inversion method using a hybrid time- and frequency-domain scheme to avoid the requirement of source wavelet estimation and to reduce the computational cost. We employ an amplitude-semblance objective function to not only effectively remove the source wavelet effect on full-waveform inversion, but also to eliminate the impact of the inconsistency of source wavelets among different shots gathers on full-waveform inversion. To reduce the high computational cost of full-waveform inversion in the time domain, we implement our new algorithm using a hybrid time- and frequency-domain approach. The forward and backward wave propagation operations are conducted in the time domain, while the frequency-domain wavefields are obtained during modeling using the discrete-time Fourier transform. The inversion process is conducted in the frequency domain for selected frequencies. We verify our method using synthetic seismic data for the Marmousi model. The results demonstrate that our novel source-independent full-waveform inversion produces accurate velocity models even if the source signature is incorrect. In addition, our method can significantly reduce the computational time using the hybrid time- and frequency-domain approach compared to the conventional full-waveform inversion in the time domain.     

Three-Dimensional Lattice Boltzmann Flux Solver and Its Applications to Incompressible Isothermal and Thermal Flows

A three-dimensional (3D) lattice Boltzmann flux solver (LBFS) is presented in this paper for the simulation of both isothermal and thermal flows. The present solver combines the advantages of conventional Navier-Stokes (N-S) solvers and lattice Boltzmann equation (LBE) solvers. It applies the finite volume method (FVM) to solve the N-S equations. Different from the conventional N-S solvers, its viscous and inviscid fluxes at the cell interface are evaluated simultaneously by local reconstruction of LBE solution. As compared to the conventional LBE solvers, which apply the lattice Boltzmann method (LBM) globally in the whole computational domain, it only applies LBM locally at each cell interface, and flow variables at cell centers are given from the solution of N-S equations. Since LBM is only applied locally in the 3D LBFS, the drawbacks of the conventional LBM, such as limitation to uniform mesh, tie-up of mesh spacing and time step, tedious implementation of boundary conditions, are completely removed. The accuracy, efficiency and stability of the proposed solver are examined in detail by simulating plane Poiseuille flow, lid-driven cavity flow and natural convection. Numerical results show that the LBFS has a second order of accuracy in space. The efficiency of the LBFS is lower than LBM on the same grids. However, the LBFS needs very less non-uniform grids to get grid-independence results and its efficiency can be greatly improved and even much higher than LBM. In addition, the LBFS is more stable and robust.     

Large-Scale Electromagnetic Modeling for Multiple Inhomogeneous Domains

We develop a new formulation of the integral equation (IE) method for three-dimensional (3D) electromagnetic (EM) field computation in large-scale models with multiple inhomogeneous domains. This problem arises in many practical applications including modeling the EM fields within the complex geoelectrical structures in geophysical exploration. In geophysical applications, it is difficult to describe an earth structure using the horizontally layered background conductivity model, which is required for the efficient implementation of the conventional IE approach. As a result, a large domain of interest with anomalous conductivity distribution needs to be discretized, which complicates the computations. The new method allows us to consider multiple inhomogeneous domains, where the conductivity distribution is different from that of the background, and to use independent discretizations for different domains. This reduces dramatically the computational resources required for large-scale modeling. In addition, using this method, we can analyze the response of each domain separately without an inappropriate use of the superposition principle for the EM field calculations. The method was carefully tested for the modeling the marine controlled-source electromagnetic (MCSEM) fields for complex geoelectric structures with multiple inhomogeneous domains, such as a seafloor with the rough bathymetry, salt domes, and reservoirs. We have also used this technique to investigate the return induction effects from regional geoelectrical structures, e.g., seafloor bathymetry and salt domes, which can distort the EM response from the geophysical exploration target.     

An Efficient Neural-Network and Finite-Difference Hybrid Method for Elliptic Interface Problems with Applications

A new and efficient neural-network and finite-difference hybrid method is developed for solving Poisson equation in a regular domain with jump discontinuities on embedded irregular interfaces. Since the solution has low regularity across the interface, when applying finite difference discretization to this problem, an additional treatment accounting for the jump discontinuities must be employed. Here, we aim to elevate such an extra effort to ease our implementation by machine learning methodology. The key idea is to decompose the solution into singular and regular parts. The neural network learning machinery incorporating the given jump conditions finds the singular solution, while the standard five-point Laplacian discretization is used to obtain the regular solution with associated boundary conditions. Regardless of the interface geometry, these two tasks only require supervised learning for function approximation and a fast direct solver for Poisson equation, making the hybrid method easy to implement and efficient. The two- and three-dimensional numerical results show that the present hybrid method preserves second-order accuracy for the solution and its derivatives, and it is comparable with the traditional immersed interface method in the literature. As an application, we solve the Stokes equations with singular forces to demonstrate the robustness of the present method.     

A Computational Model for Simulation of Shallow Water Waves by Elastic Deformations in the Topography

We propose a coupled model to simulate shallow water waves induced by elastic deformations in the bed topography. The governing equations consist of the depth-averaged shallow water equations including friction terms for the water free-surface and the well-known second-order elastostatics formulation for the bed deformation. The perturbation on the free-surface is assumed to be caused by a sudden change in the bottom beds. At the interface between the water flow and the bed topography, transfer conditions are implemented. Here, the hydrostatic pressure and friction forces are considered for the elastostatic equations whereas bathymetric forces are accounted for in the shallow water equations. The focus in the present study is on the development of a simple and accurate representation of the interaction between water waves and bed deformations in order to simulate practical shallow water flows without relying on complex partial differential equations with free boundary conditions. The effects of location and magnitude of the deformation on the flow fields and free-surface waves are investigated in details. Numerical simulations are carried out for several test examples on shallow water waves induced by sudden changes in the bed. The proposed computational model has been found to be feasible and satisfactory.     

Numerical Simulation of Unidirectional Stratified Flow by Moving Particle Semi-Implicit Method

Numerical simulation of stratified flow of two fluids between two infinite parallel plates using the Moving Particle Semi-implicit (MPS) method is presented. The developing process from entrance to fully development flow is captured. In the simulation, the computational domain is represented by various types of particles. Governing equations are described based on particles and their interactions. Grids are not necessary in any calculation steps of the simulation. The particle number density is implicitly required to be constant to satisfy incompressibility. The weight function is used to describe the interaction between different particles. The particle is considered to constitute the free interface if the particle number density is below a set point. Results for various combinations of density, viscosity, mass flow rates, and distance between the two parallel plates are presented. The proposed procedure is validated using the derived exact solution and the earlier numerical results from the Level-Set method. Furthermore, the evolution of the interface in the developing region is captured and compares well with the derived exact solutions in the developed region.     

A Riccati-equation-based algorithm for continuous-time optimal control problems

In this paper we consider continuous-time unconstrained optimal control problems. We propose a computational method which is essentially based on the closed-loop solutions of the linear quadratic optimal control problems. In the proposed algorithm, Riccati differential equations play an important role. We prove that accumulation points generated by the present algorithm, if they exist, satisfy the weak necessary conditions for optimality, under some assumptions including Kalman's sufficient conditions for the bounded Riccati solutions. In addition, we also propose the simple but effective technique to guarantee the boundedness of the solutions of Riccati equations. Lastly, we illustrate the usefulness of the present algorithm through simulation experiences. Copyright © 1998 John Wiley & Sons, Ltd.     

An alternating direction method of multipliers algorithm for symmetric model predictive control

This article presents an alternating direction method of multipliers (ADMM) algorithm for solving large‐scale model predictive control (MPC) problems that are invariant under the symmetric‐group. Symmetry was used to find transformations of the inputs, states, and constraints of the MPC problem that decompose the dynamics and cost. We prove an important property of the symmetric decomposition for the symmetric‐group that allows us to efficiently transform between the original and decomposed symmetric domains. This allows us to solve different subproblems of a baseline ADMM algorithm in different domains where the computations are less expensive. This reduces the computational cost of each iteration from quadratic to linear in the number of repetitions in the system. In addition, we show that the memory complexity for our ADMM algorithm is also linear in number of repetitions in the system, rather than the typical quadratic complexity. We demonstrate our algorithm for two case studies; battery balancing and heating, ventilation, and air conditioning. In both case studies, the symmetric algorithm reduced the computation‐time from minutes to seconds and memory usage from tens of megabytes to tens or hundreds of kilobytes, allowing the previously nonviable MPCs to be implemented in real time on embedded computers with limited computational and memory resources.