Topic overview

math.NA

6807 works11497 researchers

Open map Browse papers

Map preview

Start with the graph, then narrow the list

6807works

11497researchers

Next steps

Use the topic as a working map

Open the full map for clusters, then return here to scan ranked papers and people.

Inspect nearby papers, researchers, institutions and communities without opening a separate graph page.

Building this map preview

BZPEER is loading the nearby papers, people, topics and institutions for this page.

preprint2016arXiv

Stable computations with flat radial basis functions using vector-valued rational approximations

One commonly finds in applications of smooth radial basis functions (RBFs) that scaling the kernels so they are `flat' leads to smaller discretization errors. However, the direct numerical approach for computing with flat RBFs (RBF-Direct) is severely ill-conditioned. We present an algorithm for bypassing this ill-conditioning that is based on a new method for rational approximation (RA) of vector-valued analytic functions with the property that all components of the vector share the same singularities. This new algorithm (RBF-RA) is more accurate, robust, and easier to implement than the Contour-Padé method, which is similarly based on vector-valued rational approximation. In contrast to the stable RBF-QR and RBF-GA algorithms, which are based on finding a better conditioned base in the same RBF-space, the new algorithm can be used with any type of smooth radial kernel, and it is also applicable to a wider range of tasks (including calculating Hermite type implicit RBF-FD stencils). We present a series of numerical experiments demonstrating the effectiveness of this new method for computing RBF interpolants in the flat regime. We also demonstrate the flexibility of the method

preprint2016arXiv

A fast platform for simulating flexible fiber suspensions applied to cell mechanics

We present a novel platform for the large-scale simulation of fibrous structures immersed in a Stokesian fluid and evolving under confinement or in free-space. One of the main motivations for this work is to study the dynamics of fiber assemblies within biological cells. For this, we also incorporate the key biophysical elements that determine the dynamics of these assemblies, which include the polymerization and depolymerization kinetics of fibers, their interactions with molecular motors and other objects, their flexibility, and hydrodynamic coupling. This work, to our knowledge, is the first technique to include many-body hydrodynamic interactions (HIs), and the resulting fluid flows, in cellular fiber assemblies. We use the non-local slender body theory to compute the fluid-structure interactions of the fibers and a second-kind boundary integral formulation for other rigid bodies and the confining boundary. A kernel-independent implementation of the fast multiple method is utilized for efficient evaluation of HIs. The deformation of the fibers is described by the nonlinear Euler--Bernoulli beam theory and their polymerization is modeled by the reparametrization of the dynamic e

preprint2016arXiv

On the equivalence between the Scheduled Relaxation Jacobi method and Richardson's non-stationary method

The Scheduled Relaxation Jacobi (SRJ) method is an extension of the classical Jacobi iterative method to solve linear systems of equations ($Au=b$) associated with elliptic problems. It inherits its robustness and accelerates its convergence rate computing a set of $P$ relaxation factors that result from a minimization problem. In a typical SRJ scheme, the former set of factors is employed in cycles of $M$ consecutive iterations until a prescribed tolerance is reached. We present the analytic form for the optimal set of relaxation factors for the case in which all of them are different, and find that the resulting algorithm is equivalent to a non-stationary generalized Richardson's method. Our method to estimate the weights has the advantage that the explicit computation of the maximum and minimum eigenvalues of the matrix $A$ is replaced by the (much easier) calculation of the maximum and minimum frequencies derived from a von Neumann analysis. This set of weights is also optimal for the general problem, resulting in the fastest convergence of all possible SRJ schemes for a given grid structure. We also show that with the set of weights computed for the optimal SRJ scheme for

preprint2016arXiv

The Discrete Stochastic Galerkin Method for Hyperbolic Equations with Non-smooth and Random Coefficients

We develop a general polynomial chaos (gPC) based stochastic Galerkin (SG) for hyperbolic equations with random and singular coefficients. Due to the singu- lar nature of the solution, the standard gPC-SG methods may suffer from a poor or even non convergence. Taking advantage of the fact that the discrete solution, by the central type finite difference or finite volume approximations in space and time for example, is smoother, we first discretize the equation by a smooth finite difference or finite volume scheme, and then use the gPC-SG approximation to the discrete system. The jump condition at the interface is treated using the immersed upwind methods introduced in [8, 12]. This yields a method that converges with the spectral accuracy for finite mesh size and time step. We use a linear hyperbolic equation with discontinuous and random coefficient, and the Liouville equation with discontinuous and random potential, to illustrate our idea, with both one and second order spatial discretizations. Spectral convergence is established for the first equation, and numerical examples for both equations show the desired accu- racy of the method.

preprint2016arXiv

An Asymptotic Preserving method for strongly anisotropic diffusion equations based on field line integration

In magnetized plasma, the magnetic field confines the particles around the field lines. The anisotropy intensity in the viscosity and heat conduction may reach the order of $10^{12}$. When the boundary conditions are periodic or Neumann, the strong diffusion leads to an ill-posed limiting problem. To remove the ill-conditionedness in the highly anisotropic diffusion equations, we introduce a simple but very efficient asymptotic preserving reformulation in this paper. The key idea is that, instead of discretizing the Neumann boundary conditions locally, we replace one of the Neumann boundary condition by the integration of the original problem along the field line, the singular $1/ε$ terms can be replaced by $O(1)$ terms after the integration, so that yields a well-posed problem. Small modifications to the original code are required and no change of coordinates nor mesh adaptation are needed. Uniform convergence with respect to the anisotropy strength $1/ε$ can be observed numerically and the condition number does not scale with the anisotropy.

preprint2016arXiv

Nonlinear Flows for Displacement Correction and Applications in Tomography

In this paper we derive nonlinear evolution equations associated with a class of non-convex energy functionals which can be used for correcting displacement errors in imaging data. We study properties of these filtering flows and provide experiments for correcting angular perturbations in tomographical data.

preprint2017arXiv

Weight-adjusted discontinuous Galerkin methods: wave propagation in heterogeneous media

Time-domain discontinuous Galerkin (DG) methods for wave propagation require accounting for the inversion of dense elemental mass matrices, where each mass matrix is computed with respect to a parameter-weighted L2 inner product. In applications where the wavespeed varies spatially at a sub-element scale, these matrices are distinct over each element, necessitating additional storage. In this work, we propose a weight-adjusted DG (WADG) method which reduces storage costs by replacing the weighted L2 inner product with a weight-adjusted inner product. This equivalent inner product results in an energy stable method, but does not increase storage costs for locally varying weights. A-priori error estimates are derived, and numerical examples are given illustrating the application of this method to the acoustic wave equation with heterogeneous wavespeed.

preprint2016arXiv

Adaptive enriched Galerkin methods for miscible displacement problems with entropy residual stabilization

We present a novel approach to the simulation of miscible displacement by employing adaptive enriched Galerkin finite element methods (EG) coupled with entropy residual stabilization for transport. In particular, numerical simulations of viscous fingering instabilities in heterogeneous porous media and Hele-Shaw cells are illustrated. EG is formulated by enriching the conforming continuous Galerkin finite element method (CG) with piecewise constant functions. The method provides locally and globally conservative fluxes, which is crucial for coupled flow and transport problems. Moreover, EG has fewer degrees of freedom in comparison with discontinuous Galerkin (DG) and an efficient flow solver has been derived which allows for higher order schemes. Dynamic adaptive mesh refinement is applied in order to save computational cost for large-scale three dimensional applications. In addition, entropy residual based stabilization for high order EG transport systems prevents any spurious oscillations. Numerical tests are presented to show the capabilities of EG applied to flow and transport.

preprint2016arXiv

Elastic Splines II: unicity of optimal s-curves and $G^2$ regularity of splines

Given points $P_1,P_2,\ldots,P_m$ in the complex plane, we are concerned with the problem of finding an interpolating curve with minimal bending energy (i.e., an optimal interpolating curve). It was shown previously that existence is assured if one requires that the pieces of the interpolating curve be s-curves. In the present article we also impose the restriction that these s-curves have chord angles not exceeding $π/2$ in magnitude. With this setup, we have identified a sufficient condition for the $G^2$ regularity of optimal interpolating curves. This sufficient condition relates to the stencil angles $\{ψ_j\}$, where $ψ_j$ is defined as the angular change in direction from segment $[P_{j-1},P_j]$ to segment $[P_j,P_{j+1}]$. A distinguished angle $Ψ$ ($\approx 37^\circ$) is identified, and we show that if the stencil angles satisfy $|ψ_j|<Ψ$, then optimal interpolating curves are globally $G^2$. As with the previous article, most of our effort is concerned with the geometric Hermite interpolation problem of finding an optimal s-curve which connects $P_1$ to $P_2$ with prescribed chord angles $(α,β)$. Whereas existence was previously shown, and sometimes uniqueness, the present

preprint2016arXiv

Rank structured approximation method for quasi--periodic elliptic problems

We consider an iteration method for solving an elliptic type boundary value problem $\mathcal{A} u=f$, where a positive definite operator $\mathcal{A}$ is generated by a quasi--periodic structure with rapidly changing coefficients (typical period is characterized by a small parameter $ε$) . The method is based on using a simpler operator $\mathcal{A}_0$ (inversion of $\mathcal{A}_0$ is much simpler than inversion of $\mathcal{A}$), which can be viewed as a preconditioner for $\mathcal{A}$. We prove contraction of the iteration method and establish explicit estimates of the contraction factor $q$. Certainly the value of $q$ depends on the difference between $\mathcal{A}$ and $\mathcal{A}_0$. For typical quasi--periodic structures, we establish simple relations that suggest an optimal $\mathcal{A}_0$ (in a selected set of "simple" structures) and compute the corresponding contraction factor. Further, this allows us to deduce fully computable two--sided a posteriori estimates able to control numerical solutions on any iteration. The method is especially efficient if the coefficients of $\mathcal{A}$ admit low rank representations and algebraic operations are performed in tenso

preprint2017arXiv

Mixed one-bit compressive sensing with applications to overexposure correction for CT reconstruction

When a measurement falls outside the quantization or measurable range, it becomes saturated and cannot be used in classical reconstruction methods. For example, in C-arm angiography systems, which provide projection radiography, fluoroscopy, digital subtraction angiography, and are widely used for medical diagnoses and interventions, the limited dynamic range of C-arm flat detectors leads to overexposure in some projections during an acquisition, such as imaging relatively thin body parts (e.g., the knee). Aiming at overexposure correction for computed tomography (CT) reconstruction, we in this paper propose a mixed one-bit compressive sensing (M1bit-CS) to acquire information from both regular and saturated measurements. This method is inspired by the recent progress on one-bit compressive sensing, which deals with only sign observations. Its successful applications imply that information carried by saturated measurements is useful to improve recovery quality. For the proposed M1bit-CS model, alternating direction methods of multipliers is developed and an iterative saturation detection scheme is established. Then we evaluate M1bit-CS on one-dimensional signal recovery tasks. In s

preprint2017arXiv

Solving delay differential equations through RBF collocation

A general and easy-to-code numerical method based on radial basis functions (RBFs) collocation is proposed for the solution of delay differential equations (DDEs). It relies on the interpolation properties of infinitely smooth RBFs, which allow for a large accuracy over a scattered and relatively small discretization support. Hardy's multiquadric is chosen as RBF and combined with the Residual Subsampling Algorithm of Driscoll and Heryudono for support adaptivity. The performance of the method is very satisfactory, as demonstrated over a cross-section of benchmark DDEs, and by comparison with existing general-purpose and specialized numerical schemes for DDEs.

preprint2016arXiv

Fast symmetric factorization of hierarchical matrices with applications

We present a fast direct algorithm for computing symmetric factorizations, i.e. $A = WW^T$, of symmetric positive-definite hierarchical matrices with weak-admissibility conditions. The computational cost for the symmetric factorization scales as $\mathcal{O}(n \log^2 n)$ for hierarchically off-diagonal low-rank matrices. Once this factorization is obtained, the cost for inversion, application, and determinant computation scales as $\mathcal{O}(n \log n)$. In particular, this allows for the near optimal generation of correlated random variates in the case where $A$ is a covariance matrix. This symmetric factorization algorithm depends on two key ingredients. First, we present a novel symmetric factorization formula for low-rank updates to the identity of the form $I+UKU^T$. This factorization can be computed in $\mathcal{O}(n)$ time if the rank of the perturbation is sufficiently small. Second, combining this formula with a recursive divide-and-conquer strategy, near linear complexity symmetric factorizations for hierarchically structured matrices can be obtained. We present numerical results for matrices relevant to problems in probability \& statistics (Gaussian processes), interp

preprint2017arXiv

Integrating Lipschitzian Dynamical Systems using Piecewise Algorithmic Differentiation

In this article we analyze a generalized trapezoidal rule for initial value problems with piecewise smooth right hand side $F:\R^n\to\R^n$. When applied to such a problem the classical trapezoidal rule suffers from a loss of accuracy if the solution trajectory intersects a nondifferentiability of $F$. The advantage of the proposed generalized trapezoidal rule is threefold: Firstly we can achieve a higher convergence order than with the classical method. Moreover, the method is energy preserving for piecewise linear Hamiltonian systems. Finally, in analogy to the classical case we derive a third order interpolation polynomial for the numerical trajectory. In the smooth case the generalized rule reduces to the classical one. Hence, it is a proper extension of the classical theory. An error estimator is given and numerical results are presented.

preprint2016arXiv

Nested domain decomposition with polarized traces for the 2D Helmholtz equation

We present a solver for the 2D high-frequency Helmholtz equation in heterogeneous, constant density, acoustic media, with online parallel complexity that scales empirically as $\mathcal{O}(\frac{N}{P})$, where $N$ is the number of volume unknowns, and $P$ is the number of processors, as long as $P = \mathcal{O}(N^{1/5})$. This sublinear scaling is achieved by domain decomposition, not distributed linear algebra, and improves on the $P =\mathcal{O}(N^{1/8})$ scaling reported earlier in [L. Zepeda-Núñez and L. Demanet, J. Comput. Phys., 308 (2016), pp. 347-388 ]. The solver relies on a two-level nested domain decomposition: a layered partition on the outer level, and a further decomposition of each layer in cells at the inner level. The Helmholtz equation is reduced to a surface integral equation (SIE) posed at the interfaces between layers, efficiently solved via a nested version of the polarized traces preconditioner [L. Zepeda-Núñez and L. Demanet, J. Comput. Phys., 308 (2016), pp. 347-388.]. The favorable complexity is achieved via an efficient application of the integral operators involved in the SIE.

preprint2016arXiv

An efficient threshold dynamics method for wetting on rough surfaces

The threshold dynamics method developed by Merriman, Bence and Osher (MBO) is an efficient method for simulating the motion by mean curvature flow when the interface is away from the solid boundary. Direct generalization of MBO-type methods to the wetting problem with interfaces intersecting the solid boundary is not easy because solving the heat equation in a general domain with a wetting boundary condition is not as efficient as it is with the original MBO method. The dynamics of the contact point also follows a different law compared with the dynamics of the interface away from the boundary. In this paper, we develop an efficient volume preserving threshold dynamics method for simulating wetting on rough surfaces. This method is based on minimization of the weighted surface area functional over an extended domain that includes the solid phase. The method is simple, stable with $O(N \log N)$ complexity per time step and is not sensitive to the inhomogeneity or roughness of the solid boundary.

preprint2017arXiv

Trust-Region Methods for Nonlinear Elliptic Equations with Radial Basis Functions

We consider the numerical solution of nonlinear elliptic boundary value problems with Kansa's method. We derive analytic formulas for the Jacobian and Hessian of the resulting nonlinear collocation system and exploit them within the framework of the trust-region algorithm. This ansatz is tested on semilinear, quasilinear and fully nonlinear elliptic PDEs (including Plateau's problem, Hele-Shaw flow and the Monge-Ampère equation) with excellent results. The new approach distinctly outperforms previous ones based on linearization or finite-difference Jacobians.

preprint2016arXiv

Palindromic discontinuous Galerkin method for kinetic equations with stiff relaxation

We present a high order scheme for approximating kinetic equations with stiff relaxation. The objective is to provide efficient methods for solving the underlying system of conservation laws. The construction is based on several ingredients: (i) a high order implicit upwind Discontinuous Galerkin approximation of the kinetic equations with easy-to-solve triangular linear systems; (ii) a second order asymptotic-preserving time integration based on symmetry arguments; (iii) a palindromic composition of the second order method for achieving higher orders in time. The method is then tested at orders 2, 4 and 6. It is asymptotic-preserving with respect to the stiff relaxation and accepts high CFL numbers.

preprint2017arXiv

Elastic Splines I: Existence

Given interpolation points $P_1,P_2,\ldots,P_n$ in the plane, it is known that there does not exist an interpolating curve with minimal bending energy, unless the given points lie sequentially along a line. We say than an interpolating curve is {\it admissable} if each piece, connecting two consecutive points $P_i$ and $P_{i+1}$, is an s-curve, where an {\it s-curve} is a planar curve which first turns at most $180^\circ$ in one direction and then turns at most $180^\circ$ in the opposite direction. Our main result is that among all admissable interpolating curves there exists a curve with minimal bending energy. We also prove, in a very constructive manner, the existence of an s-curve, with minimal bending energy, which connects two given unit tangent vectors.

preprint2015arXiv

Convergence Analysis for a Finite Element Approximation of a Steady Model for Electrorheological Fluids

In this paper we study the finite element approximation of systems of $p(\cdot)$-Stokes type, where $p(\cdot)$ is a (non constant) given function of the space variables. We derive --in some cases optimal-- error estimates for finite element approximation of the velocity and of the pressure, in a suitable functional setting.

preprint2016arXiv

A Study on Moving Mesh Finite Element Solution of the Porous Medium Equation

An adaptive moving mesh finite element method is studied for the numerical solution of the porous medium equation with and without variable exponents and absorption. The method is based on the so-called moving mesh partial differential equation approach and employs its newly developed implementation. Three types of mesh are considered, uniform and arclength-based and Hessian-based adaptive meshes. The method shows a first order convergence for uniform and arclength-based adaptive meshes and a second-order convergence for Hessian-based adaptive meshes. It is also shown that the method can be used for situations with complex free boundaries, emerging and splitting of free boundaries, and the porous medium equation with variable exponents and absorption. Two dimensional numerical results are presented.

preprint2016arXiv

A Positive and Entropy-Satisfying Finite Volume Scheme for the Baer-Nunziato Model

We present a relaxation scheme for approximating the entropy dissipating weak solutions of the Baer-Nunziato two-phase flow model. This relaxation scheme is straightforwardly obtained as an extension of the relaxation scheme designed in [16] for the isentropic Baer-Nunziato model and consequently inherits its main properties. To our knowledge, this is the only existing scheme for which the approximated phase fractions, phase densities and phase internal energies are proven to remain positive without any restrictive condition other than a classical fully computable CFL condition. For ideal gas and stiffened gas equations of state, real values of the phasic speeds of sound are also proven to be maintained by the numerical scheme. It is also the only scheme for which a discrete entropy inequality is proven, under a CFL condition derived from the natural sub-characteristic condition associated with the relaxation approximation. This last property, which ensures the non-linear stability of the numerical method, is satisfied for any admissible equation of state. We provide a numerical study for the convergence of the approximate solutions towards some exact Riemann solutions. The numeric

preprint2014arXiv

Finite element approximation of the $p(\cdot)$-Laplacian

We study a~priori estimates for the Dirichlet problem of the $p(\cdot)$-Laplacian, \[-\mathrm{div}(|\nabla v|^{p(\cdot)-2} \nabla v) = f. \] We show that the gradients of the finite element approximation with zero boundary data converges with rate $O(h^α)$ if the exponent $p$ is $α$-Hölder continuous. The error of the gradients is measured in the so-called quasi-norm, i.e. we measure the $L^2$-error of $|\nabla v|^{\frac{p-2}{2}} \nabla v$.

preprint2015arXiv

On harmonic analysis of vector-valued signals

A vector-valued signal in N dimensions is a signal whose value at any time instant is an N-dimensional vector, that is, an element of $\mathbb{R}^N$. The sum of an arbitrary number of such signals of the same frequency is shown to trace an ellipse in N-dimensional space, that is, to be confined to a plane. The parameters of the ellipse (major and minor axes, represented by N-dimensional vectors; and phase) are obtained algebraically in terms of the directions of oscillation of the constituent signals, and their phases. It is shown that the major axis of the ellipse can always be determined algebraically. That is, a vector, whose value can be computed algebraically (without decisions or comparisons of magnitude) from parameters of the constituent signals, always represents the major axis of the ellipse. The ramifications of this result for the processing and Fourier analysis of signals with vector values or samples are discussed, with reference to the definition of Fourier transforms, particularly discrete Fourier transforms, such as have been defined in several hypercomplex algebras, including Clifford algebras. The treatment in the paper, however, is entirely based on signals with

702 works