Showing papers in "Communications on Applied Mathematics and Computation in 2022"

Dual Quaternions and Dual Quaternion Vectors

Abstract: We introduce a total order and an absolute value function for dual numbers. The absolute value function of dual numbers takes dual number values, and has properties similar to those of the absolute value function of real numbers. We define the magnitude of a dual quaternion, as a dual number. Based upon these, we extend 1-norm, $$\infty$$ -norm, and 2-norm to dual quaternion vectors.

Quinpi: Integrating Conservation Laws with CWENO Implicit Methods

Abstract: Abstract Many interesting applications of hyperbolic systems of equations are stiff, and require the time step to satisfy restrictive stability conditions. One way to avoid small time steps is to use implicit time integration. Implicit integration is quite straightforward for first-order schemes. High order schemes instead also need to control spurious oscillations, which requires limiting in space and time also in the linear case. We propose a framework to simplify considerably the application of high order non-oscillatory schemes through the introduction of a low order implicit predictor, which is used both to set up the nonlinear weights of a standard high order space reconstruction, and to achieve limiting in time. In this preliminary work, we concentrate on the case of a third-order scheme, based on diagonally implicit Runge Kutta ( $$\mathsf {DIRK}$$ DIRK ) integration in time and central weighted essentially non-oscillatory ( $$\mathsf {CWENO}$$ CWENO ) reconstruction in space. The numerical tests involve linear and nonlinear scalar conservation laws.

Standard Dual Quaternion Optimization and Its Applications in Hand-Eye Calibration and SLAM

Abstract: Several common dual quaternion functions, such as the power function, the magnitude function, the 2-norm function, and the kth largest eigenvalue of a dual quaternion Hermitian matrix, are standard dual quaternion functions, i.e., the standard parts of their function values depend upon only the standard parts of their dual quaternion variables. Furthermore, the sum, product, minimum, maximum, and composite functions of two standard dual functions, the logarithm and the exponential of standard unit dual quaternion functions, are still standard dual quaternion functions. On the other hand, the dual quaternion optimization problem, where objective and constraint function values are dual numbers but variables are dual quaternions, naturally arises from applications. We show that to solve an equality constrained dual quaternion optimization (EQDQO) problem, we only need to solve two quaternion optimization problems. If the involved dual quaternion functions are all standard, the optimization problem is called a standard dual quaternion optimization problem, and some better results hold. Then, we show that the dual quaternion optimization problems arising from the hand-eye calibration problem and the simultaneous localization and mapping (SLAM) problem are equality constrained standard dual quaternion optimization problems.

Global Solutions to Nonconvex Problems by Evolution of Hamilton-Jacobi PDEs

Abstract: Computing tasks may often be posed as optimization problems. The objective functions for real-world scenarios are often nonconvex and/or nondifferentiable. State-of-the-art methods for solving these problems typically only guarantee convergence to local minima. This work presents Hamilton-Jacobi-based Moreau adaptive descent (HJ-MAD), a zero-order algorithm with guaranteed convergence to global minima, assuming continuity of the objective function. The core idea is to compute gradients of the Moreau envelope of the objective (which is “piece-wise convex”) with adaptive smoothing parameters. Gradients of the Moreau envelope (i.e., proximal operators) are approximated via the Hopf-Lax formula for the viscous Hamilton-Jacobi equation. Our numerical examples illustrate global convergence.

A Fixed-Point Fast Sweeping WENO Method with Inverse Lax-Wendroff Boundary Treatment for Steady State of Hyperbolic Conservation Laws

Abstract: Abstract Fixed-point fast sweeping WENO methods are a class of efficient high-order numerical methods to solve steady-state solutions of hyperbolic partial differential equations (PDEs). The Gauss-Seidel iterations and alternating sweeping strategy are used to cover characteristics of hyperbolic PDEs in each sweeping order to achieve fast convergence rate to steady-state solutions. A nice property of fixed-point fast sweeping WENO methods which distinguishes them from other fast sweeping methods is that they are explicit and do not require inverse operation of nonlinear local systems. Hence, they are easy to be applied to a general hyperbolic system. To deal with the difficulties associated with numerical boundary treatment when high-order finite difference methods on a Cartesian mesh are used to solve hyperbolic PDEs on complex domains, inverse Lax-Wendroff (ILW) procedures were developed as a very effective approach in the literature. In this paper, we combine a fifth-order fixed-point fast sweeping WENO method with an ILW procedure to solve steady-state solution of hyperbolic conservation laws on complex computing regions. Numerical experiments are performed to test the method in solving various problems including the cases with the physical boundary not aligned with the grids. Numerical results show high-order accuracy and good performance of the method. Furthermore, the method is compared with the popular third-order total variation diminishing Runge-Kutta (TVD-RK3) time-marching method for steady-state computations. Numerical examples show that for most of examples, the fixed-point fast sweeping method saves more than half CPU time costs than TVD-RK3 to converge to steady-state solutions.

Efficient WENO-Based Prolongation Strategies for Divergence-Preserving Vector Fields

Abstract: Abstract Adaptive mesh refinement (AMR) is the art of solving PDEs on a mesh hierarchy with increasing mesh refinement at each level of the hierarchy. Accurate treatment on AMR hierarchies requires accurate prolongation of the solution from a coarse mesh to a newly defined finer mesh. For scalar variables, suitably high-order finite volume WENO methods can carry out such a prolongation. However, classes of PDEs, such as computational electrodynamics (CED) and magnetohydrodynamics (MHD), require that vector fields preserve a divergence constraint. The primal variables in such schemes consist of normal components of the vector field that are collocated at the faces of the mesh. As a result, the reconstruction and prolongation strategies for divergence constraint-preserving vector fields are necessarily more intricate. In this paper we present a fourth-order divergence constraint-preserving prolongation strategy that is analytically exact. Extension to higher orders using analytically exact methods is very challenging. To overcome that challenge, a novel WENO-like reconstruction strategy is invented that matches the moments of the vector field in the faces, where the vector field components are collocated. This approach is almost divergence constraint-preserving, therefore, we call it WENO-ADP. To make it exactly divergence constraint-preserving, a touch-up procedure is developed that is based on a constrained least squares (CLSQ) method for restoring the divergence constraint up to machine accuracy. With the touch-up, it is called WENO-ADPT. It is shown that refinement ratios of two and higher can be accommodated. An item of broader interest in this work is that we have also been able to invent very efficient finite volume WENO methods, where the coefficients are very easily obtained and the multidimensional smoothness indicators can be expressed as perfect squares. We demonstrate that the divergence constraint-preserving strategy works at several high orders for divergence-free vector fields as well as vector fields, where the divergence of the vector field has to match a charge density and its higher moments. We also show that our methods overcome the late time instability that has been known to plague adaptive computations in CED.

Space-Fractional Diffusion with Variable Order and Diffusivity: Discretization and Direct Solution Strategies

Abstract: We consider the multidimensional space-fractional diffusion equations with spatially varying diffusivity and fractional order. Significant computational challenges are encountered when solving these equations due to the kernel singularity in the fractional integral operator and the resulting dense discretized operators, which quickly become prohibitively expensive to handle because of their memory and arithmetic complexities. In this work, we present a singularity-aware discretization scheme that regularizes the singular integrals through a singularity subtraction technique adapted to the spatial variability of diffusivity and fractional order. This regularization strategy is conveniently formulated as a sparse matrix correction that is added to the dense operator, and is applicable to different formulations of fractional diffusion equations. We also present a block low rank representation to handle the dense matrix representations, by exploiting the ability to approximate blocks of the resulting formally dense matrix by low rank factorizations. A Cholesky factorization solver operates directly on this representation using the low rank blocks as its atomic computational tiles, and achieves high performance on multicore hardware. Numerical results show that the singularity treatment is robust, substantially reduces discretization errors, and attains the first-order convergence rate allowed by the regularity of the solutions. They also show that considerable savings are obtained in storage ( $$O(N^{1.5})$$ ) and computational cost ( $$O(N^2)$$ ) compared to dense factorizations. This translates to orders-of-magnitude savings in memory and time on multidimensional problems, and shows that the proposed methods offer practical tools for tackling large nonlocal fractional diffusion simulations.

Singularity Formation for the General Poiseuille Flow of Nematic Liquid Crystals

Abstract: We consider the Poiseuille flow of nematic liquid crystals via the full Ericksen-Leslie model. The model is described by a coupled system consisting of a heat equation and a quasilinear wave equation. In this paper, we will construct an example with a finite time cusp singularity due to the quasilinearity of the wave equation, extended from an earlier result on a special case.

A Space-Time Interior Penalty Discontinuous Galerkin Method for the Wave Equation

Abstract: Abstract A new higher-order accurate space-time discontinuous Galerkin (DG) method using the interior penalty flux and discontinuous basis functions, both in space and in time, is presented and fully analyzed for the second-order scalar wave equation. Special attention is given to the definition of the numerical fluxes since they are crucial for the stability and accuracy of the space-time DG method. The theoretical analysis shows that the DG discretization is stable and converges in a DG-norm on general unstructured and locally refined meshes, including local refinement in time. The space-time interior penalty DG discretization does not have a CFL-type restriction for stability. Optimal order of accuracy is obtained in the DG-norm if the mesh size h and the time step $$\Delta t$$ Δ t satisfy $$h\cong C\Delta t$$ h ≅ C Δ t , with C a positive constant. The optimal order of accuracy of the space-time DG discretization in the DG-norm is confirmed by calculations on several model problems. These calculations also show that for p th-order tensor product basis functions the convergence rate in the $$L^\infty$$ L ∞ and $$L^2$$ L 2 -norms is order $$p+1$$ p + 1 for polynomial orders $$p=1$$ p = 1 and $$p=3$$ p = 3 and order p for polynomial order $$p=2$$ p = 2 .

A Space-Time Interior Penalty Discontinuous Galerkin Method for the Wave Equation

Abstract: Abstract A new higher-order accurate space-time discontinuous Galerkin (DG) method using the interior penalty flux and discontinuous basis functions, both in space and in time, is presented and fully analyzed for the second-order scalar wave equation. Special attention is given to the definition of the numerical fluxes since they are crucial for the stability and accuracy of the space-time DG method. The theoretical analysis shows that the DG discretization is stable and converges in a DG-norm on general unstructured and locally refined meshes, including local refinement in time. The space-time interior penalty DG discretization does not have a CFL-type restriction for stability. Optimal order of accuracy is obtained in the DG-norm if the mesh size h and the time step $$\Delta t$$ Δ t satisfy $$h\cong C\Delta t$$ h ≅ C Δ t , with C a positive constant. The optimal order of accuracy of the space-time DG discretization in the DG-norm is confirmed by calculations on several model problems. These calculations also show that for p th-order tensor product basis functions the convergence rate in the $$L^\infty$$ L ∞ and $$L^2$$ L 2 -norms is order $$p+1$$ p + 1 for polynomial orders $$p=1$$ p = 1 and $$p=3$$ p = 3 and order p for polynomial order $$p=2$$ p = 2 .

Numerical Approach of a Coupled Pressure-Saturation Model Describing Oil-Water Flow in Porous Media

Abstract: Abstract Two-phase flow in porous media is a very active field of research, due to its important applications in groundwater pollution, $${\text {CO}}_2$$ CO 2 sequestration, or oil and gas production from petroleum reservoirs, just to name a few of them. Fractional flow equations, which make use of Darcy’s law, for describing the movement of two immiscible fluids in a porous medium, are among the most relevant mathematical models in reservoir simulation. This work aims to solve a fractional flow model formed by an elliptic equation, representing the spatial distribution of the pressure, and a hyperbolic equation describing the space-time evolution of water saturation. The numerical solution of the elliptic part is obtained using a finite-element (FE) scheme, while the hyperbolic equation is solved by means of two different numerical approaches, both in the finite-volume (FV) framework. One is based on a monotonic upstream-centered scheme for conservation laws (MUSCL)-Hancock scheme, whereas the other makes use of a weighted essentially non-oscillatory (ENO) reconstruction. In both cases, a first-order centered (FORCE)- $$\alpha$$ α numerical scheme is applied for intercell flux reconstruction, which constitutes a new contribution in the field of fractional flow models describing oil-water movement. A relevant feature of this work is the study of the effect of the parameter $$\alpha$$ α on the numerical solution of the models considered. We also show that, in the FORCE- $$\alpha$$ α method, when the parameter $$\alpha$$ α increases, the errors diminish and the order of accuracy is more properly attained, as verified using a manufactured solution technique.

A Combination of Residual Distribution and the Active Flux Formulations or a New Class of Schemes That Can Combine Several Writings of the Same Hyperbolic Problem: Application to the 1D Euler Equations

Abstract: We show how to combine in a natural way (i.e., without any test nor switch) the conservative and non-conservative formulations of an hyperbolic system that has a conservative form. This is inspired from two different classes of schemes: the residual distribution one (Abgrall in Commun Appl Math Comput 2(3): 341-368, 2020), and the active flux formulations (Eyman and Roe in 49th AIAA Aerospace Science Meeting, 2011; Eyman in active flux. PhD thesis, University of Michigan, 2013; Helzel et al. in J Sci Comput 80(3): 35-61, 2019; Barsukow in J Sci Comput 86(1): paper No. 3, 34, 2021; Roe in J Sci Comput 73: 1094-1114, 2017). The solution is globally continuous, and as in the active flux method, described by a combination of point values and average values. Unlike the "classical" active flux methods, the meaning of the point-wise and cell average degrees of freedom is different, and hence follow different forms of PDEs; it is a conservative version of the cell average, and a possibly non-conservative one for the points. This new class of scheme is proved to satisfy a Lax-Wendroff-like theorem. We also develop a method to perform non-linear stability. We illustrate the behaviour on several benchmarks, some quite challenging.

Radon Measure Solutions to Riemann Problems for Isentropic Compressible Euler Equations of Polytropic Gases

Abstract: We solve the Riemann problems for isentropic compressible Euler equations of polytropic gases in the class of Radon measures, and the solutions admit the concentration of mass. It is found that under the requirement of satisfying the over-compressing entropy condition: (i) there is a unique delta shock solution, corresponding to the case that has two strong classical Lax shocks; (ii) for the initial data that the classical Riemann solution contains a shock wave and a rarefaction wave, or two shocks with one being weak, there are infinitely many solutions, each consists of a delta shock and a rarefaction wave; (iii) there are no delta shocks for the case that the classical entropy weak solutions consist only of rarefaction waves. These solutions are self-similar. Furthermore, for the generalized Riemann problem with mass concentrated initially at the discontinuous point of initial data, there always exists a unique delta shock for at least a short time. It could be prolonged to a global solution. Not all the solutions are self-similar due to the initial velocity of the concentrated point-mass (particle). Whether the delta shock solutions constructed satisfy the over-compressing entropy condition is clarified. This is the first result on the construction of singular measure solutions to the compressible Euler system of polytropic gases, that is strictly hyperbolic, and whose characteristics are both genuinely nonlinear. We also discuss possible physical interpretations and applications of these new solutions.

Modeling Fast Diffusion Processes in Time Integration of Stiff Stochastic Differential Equations

Abstract: Numerical algorithms for stiff stochastic differential equations are developed using linear approximations of the fast diffusion processes, under the assumption of decoupling between fast and slow processes. Three numerical schemes are proposed, all of which are based on the linearized formulation albeit with different degrees of approximation. The schemes are of comparable complexity to the classical explicit Euler-Maruyama scheme but can achieve better accuracy at larger time steps in stiff systems. Convergence analysis is conducted for one of the schemes, that shows it to have a strong convergence order of 1/2 and a weak convergence order of 1. Approximations arriving at the other two schemes are discussed. Numerical experiments are carried out to examine the convergence of the schemes proposed on model problems.

Perturbation Analysis for t-Product-Based Tensor Inverse, Moore-Penrose Inverse and Tensor System

Abstract: This paper establishes some perturbation analysis for the tensor inverse, the tensor Moore-Penrose inverse, and the tensor system based on the t-product. In the settings of structured perturbations, we generalize the Sherman-Morrison-Woodbury (SMW) formula to the t-product tensor scenarios. The SMW formula can be used to perform the sensitivity analysis for a multilinear system of equations.