Showing papers in "SIAM Journal on Scientific Computing in 2015"
TL;DR: In this article, a general numerical framework to approximate so-lutions to linear programs related to optimal transport is presented, where the set of linear constraints can be split in an intersection of a few simple constraints, for which the projections can be computed in closed form.
Abstract: This article details a general numerical framework to approximate so-lutions to linear programs related to optimal transport. The general idea is to introduce an entropic regularization of the initial linear program. This regularized problem corresponds to a Kullback-Leibler Bregman di-vergence projection of a vector (representing some initial joint distribu-tion) on the polytope of constraints. We show that for many problems related to optimal transport, the set of linear constraints can be split in an intersection of a few simple constraints, for which the projections can be computed in closed form. This allows us to make use of iterative Bregman projections (when there are only equality constraints) or more generally Bregman-Dykstra iterations (when inequality constraints are in-volved). We illustrate the usefulness of this approach to several variational problems related to optimal transport: barycenters for the optimal trans-port metric, tomographic reconstruction, multi-marginal optimal trans-port and in particular its application to Brenier's relaxed solutions of in-compressible Euler equations, partial un-balanced optimal transport and optimal transport with capacity constraints.
567 citations
TL;DR: A sparsity oriented simulated annealing procedure with non-Gaussian random perturbation is proposed and the almost sure convergence of the combined algorithm (DCASA) to a global minimum is proved.
Abstract: We study minimization of the difference of $\ell_1$ and $\ell_2$ norms as a nonconvex and Lipschitz continuous metric for solving constrained and unconstrained compressed sensing problems. We establish exact (stable) sparse recovery results under a restricted isometry property (RIP) condition for the constrained problem, and a full-rank theorem of the sensing matrix restricted to the support of the sparse solution. We present an iterative method for $\ell_{1-2}$ minimization based on the difference of convex functions algorithm and prove that it converges to a stationary point satisfying the first-order optimality condition. We propose a sparsity oriented simulated annealing procedure with non-Gaussian random perturbation and prove the almost sure convergence of the combined algorithm (DCASA) to a global minimum. Computation examples on success rates of sparse solution recovery show that if the sensing matrix is ill-conditioned (non RIP satisfying), then our method is better than existing nonconvex compre...
349 citations
TL;DR: This research presents a novel, scalable, scalable and scalable approaches to solve the challenge of integrating NoSQL data stores to manage and manage distributed systems.
Abstract: United States. Air Force Office of Scientific Research. Dynamic Data-Driven Application Systems Program (Grant FA9550-11-1-0339)
179 citations
TL;DR: In this article, the authors present a novel error analysis that considers randomized algorithms within the subspace iteration framework and show with very high probability that highly accurate low-rank approximations as well as singular values ca...
Abstract: A classical problem in matrix computations is the efficient and reliable approximation of a given matrix by a matrix of lower rank. The truncated singular value decomposition (SVD) is known to provide the best such approximation for any given fixed rank. However, the SVD is also known to be very costly to compute. Among the different approaches in the literature for computing low-rank approximations, randomized algorithms have attracted researchers' attention recently due to their surprising reliability and computational efficiency in different application areas. Typically, such algorithms are shown to compute with very high probability low-rank approximations that are within a constant factor from optimal, and are known to perform even better in many practical situations. In this paper, we present a novel error analysis that considers randomized algorithms within the subspace iteration framework and show with very high probability that highly accurate low-rank approximations as well as singular values ca...
178 citations
TL;DR: A low-rank format called Block Low-Rank (BLR) is proposed, and it is explained how it can be used to reduce the memory footprint and the complexity of direct solvers for sparse matrices based on the multifrontal method.
Abstract: Matrices coming from elliptic Partial Differential Equations (PDEs) have been shown to have a low-rank property: well defined off-diagonal blocks of their Schur complements can be approximated by low-rank products. Given a suitable ordering of the matrix which gives to the blocks a geometrical meaning, such approximations can be computed using an SVD or a rank-revealing QR factorization. The resulting representation offers a substantial reduction of the memory requirement and gives efficient ways to perform many of the basic dense algebra operations. Several strategies have been proposed to exploit this property. We propose a low-rank format called Block Low-Rank (BLR), and explain how it can be used to reduce the memory footprint and the complexity of direct solvers for sparse matrices based on the multifrontal method. We present experimental results that show how the BLR format delivers gains that are comparable to those obtained with hierarchical formats such as Hierarchical matrices (H matrices) and Hierarchically Semi-Separable (HSS matrices) but provides much greater flexibility and ease of use which are essential in the context of a general purpose, algebraic solver.
170 citations
TL;DR: In this article, a space-time spectral method is presented for the numerical solution of the time fractional Fokker-planck initial-boundary value problem, which employs the Jacobi polynomials for the temporal discretization and Fourier-like basis functions for the spatial discretisation.
Abstract: The fractional Fokker--Planck equation is an important physical model for simulating anomalous diffusions with external forces. Because of the nonlocal property of the fractional derivative an interesting problem is to explore high accuracy numerical methods for fractional differential equations. In this paper, a space-time spectral method is presented for the numerical solution of the time fractional Fokker--Planck initial-boundary value problem. The proposed method employs the Jacobi polynomials for the temporal discretization and Fourier-like basis functions for the spatial discretization. Due to the diagonalizable trait of the Fourier-like basis functions, this leads to a reduced representation of the inner product in the Galerkin analysis. We prove that the time fractional Fokker--Planck equation attains the same approximation order as the time fractional diffusion equation developed in [X. Li and C. Xu, SIAM J. Numer. Anal., 47 (2009), pp. 2108--2131] by using the present method. That indicates an e...
165 citations
TL;DR: Numerical tests show that very few sweeps are needed to construct a factorization that is an effective preconditioner, and the amount of parallelism is large irrespective of the ordering of the matrix, and matrix ordering can be used to enhance the accuracy of the factorization rather than to increase parallelism.
Abstract: This paper presents a new fine-grained parallel algorithm for computing an incomplete LU factorization. All nonzeros in the incomplete factors can be computed in parallel and asynchronously, using one or more sweeps that iteratively improve the accuracy of the factorization. Unlike existing parallel algorithms, the amount of parallelism is large irrespective of the ordering of the matrix, and matrix ordering can be used to enhance the accuracy of the factorization rather than to increase parallelism. Numerical tests show that very few sweeps are needed to construct a factorization that is an effective preconditioner.
162 citations
TL;DR: This article aims to fill in the gap of the second-order accurate schemes for the time-fractional subdiffusion equation with unconditional stability by constructing two improved algorithms that exhibit better performances than the known ones.
Abstract: This article aims to fill in the gap of the second-order accurate schemes for the time-fractional subdiffusion equation with unconditional stability. Two fully discrete schemes are first proposed f...
147 citations
TL;DR: A high-order accurate numerical quadrature algorithm is presented for the evaluation of integrals over curved surfaces and volumes which are defined implicitly via a fixed isosurface of a given function restricted to a given hyperrectangle.
Abstract: A high-order accurate numerical quadrature algorithm is presented for the evaluation of integrals over curved surfaces and volumes which are defined implicitly via a fixed isosurface of a given function restricted to a given hyperrectangle. By converting the implicitly defined geometry into the graph of an implicitly defined height function, the approach leads to a recursive algorithm on the number of spatial dimensions which requires only one-dimensional root finding and one-dimensional Gaussian quadrature. The computed quadrature scheme yields strictly positive quadrature weights and inherits the high-order accuracy of Gaussian quadrature: a range of different convergence tests demonstrate orders of accuracy up to 20th order. Also presented is an application of the quadrature algorithm to a high-order embedded boundary discontinuous Galerkin method for solving partial differential equations on curved domains.
131 citations
TL;DR: A recently introduced approach for nonlinear model order reduction based on generalized moment matching using basic tensor calculus and the idea of two-sided interpolation methods is extended to this more general setting by employing the tensor structure of the Hessian.
Abstract: In this paper, we investigate a recently introduced approach for nonlinear model order reduction based on generalized moment matching. Using basic tensor calculus, we propose a computationally efficient way of computing reduced-order models. We further extend the idea of two-sided interpolation methods to this more general setting by employing the tensor structure of the Hessian. We investigate the use of oblique projections in order to preserve important system properties such as stability. We test one-sided and two-sided projection methods for different semi-discretized nonlinear partial differential equations and show their competitiveness when compared to proper orthogonal decomposition (POD).
124 citations
TL;DR: In this paper, a low-rank update of the prior covariance matrix is proposed to characterize and approximate the posterior distribution of the parameters in inverse problems, based on the leading eigendirections of the matrix pencil defined by the Hessian of the negative log-likelihood and the prior precision.
Abstract: In the Bayesian approach to inverse problems, data are often informative, relative to the prior, only on a low-dimensional subspace of the parameter space. Significant computational savings can be achieved by using this subspace to characterize and approximate the posterior distribution of the parameters. We first investigate approximation of the posterior covariance matrix as a low-rank update of the prior covariance matrix. We prove optimality of a particular update, based on the leading eigendirections of the matrix pencil defined by the Hessian of the negative log-likelihood and the prior precision, for a broad class of loss functions. This class includes the Forstner metric for symmetric positive definite matrices, as well as the Kullback--Leibler divergence and the Hellinger distance between the associated distributions. We also propose two fast approximations of the posterior mean and prove their optimality with respect to a weighted Bayes risk under squared-error loss. These approximations are dep...
TL;DR: In this paper, an infinite-dimensional deflation algorithm for systematically modifying the residual of a nonlinear PDE problem to eliminate known solutions from consideration is presented, which enables the Newton-Kantorovitch iteration to converge to several different solutions, even starting from the same initial guess.
Abstract: Nonlinear systems of partial differential equations (PDEs) may permit several distinct solutions. The typical current approach to finding distinct solutions is to start Newton's method with many different initial guesses, hoping to find starting points that lie in different basins of attraction. In this paper, we present an infinite-dimensional deflation algorithm for systematically modifying the residual of a nonlinear PDE problem to eliminate known solutions from consideration. This enables the Newton--Kantorovitch iteration to converge to several different solutions, even starting from the same initial guess. The deflated Jacobian is dense, but an efficient preconditioning strategy is devised, and the number of Krylov iterations is observed not to grow as solutions are deflated. The power of the approach is demonstrated on several problems from special functions, phase separation, differential geometry, and fluid mechanics that permit distinct solutions.
TL;DR: In this paper, the authors generalize existing Jacobi-Gauss-Lobatto collocation methods for variable-order fractional differential equations using a singular approximation basis in terms of weighted Jacobi polynomials.
Abstract: We generalize existing Jacobi--Gauss--Lobatto collocation methods for variable-order fractional differential equations using a singular approximation basis in terms of weighted Jacobi polynomials of the form $(1 \pm x)^\mu P_j^{a,b}(x)$, where ${\mu>-1}$. In order to derive the differentiation matrices of the variable-order fractional derivatives, we develop a three-term recurrence relation for both integrals and derivatives of these weighted Jacobi polynomials, hence extending the three-term recurrence relationship of Jacobi polynomials. The new spectral collocation method is applied to solve fractional ordinary and partial differential equations with endpoint singularities. We demonstrate that the singular basis enhances greatly the accuracy of the numerical solution by properly tuning the parameter $\mu$, even for cases where we do not know explicitly the form of singularity in the solution at the boundaries.
TL;DR: This work proposes improved rational approximants leading to FEAST variants with faster convergence, in particular, when using rational approximation based on the work of Zolotarev, and improves both convergence robustness and load balancing when FEAST runs on multiple search intervals in parallel.
Abstract: The FEAST method for solving large sparse eigenproblems is equivalent to subspace iteration with an approximate spectral projector and implicit orthogonalization. This relation allows to characterize the convergence of this method in terms of the error of a certain rational approximant to an indicator function. We propose improved rational approximants leading to FEAST variants with faster convergence, in particular, when using rational approximants based on the work of Zolotarev. Numerical experiments demonstrate the possible computational savings especially for pencils whose eigenvalues are not well separated and when the dimension of the search space is only slightly larger than the number of wanted eigenvalues. The new approach improves both convergence robustness and load balancing when FEAST runs on multiple search intervals in parallel.
TL;DR: The design and implementation of the AmgX library, which provides drop-in GPU acceleration of distributed algebraic multigrid (AMG) and preconditioned iterative methods, is discussed.
Abstract: The solution of large sparse linear systems arises in many applications, such as computational fluid dynamics and oil reservoir simulation. In realistic cases the matrices are often so large that they require large scale distributed parallel computing to obtain the solution of interest in a reasonable time. In this paper we discuss the design and implementation of the AmgX library, which provides drop-in GPU acceleration of distributed algebraic multigrid (AMG) and preconditioned iterative methods. The AmgX library implements both classical and aggregation-based AMG methods with different selector and interpolation strategies, along with a variety of smoothers and preconditioners, including block-Jacobi, Gauss--Seidel, and incomplete-LU factorization. The library contains many of the standard and flexible preconditioned Krylov subspace iterative methods, which can be combined with any of the available multigrid methods or simpler preconditioners. The parallelism in the aggregation scheme exploits parallel...
TL;DR: In this paper, an efficient and scalable low rank matrix completion algorithm was proposed, where the orthogonal matching pursuit method was extended from the vector case to the matrix case and a weight updating rule was introduced to reduce the time and storage complexity.
Abstract: In this paper, we propose an efficient and scalable low rank matrix completion algorithm. The key idea is to extend the orthogonal matching pursuit method from the vector case to the matrix case. We further propose an economic version of our algorithm by introducing a novel weight updating rule to reduce the time and storage complexity. Both versions are computationally inexpensive for each matrix pursuit iteration and find satisfactory results in a few iterations. Another advantage of our proposed algorithm is that it has only one tunable parameter, which is the rank. It is easy to understand and to use by the user. This becomes especially important in large-scale learning problems. In addition, we rigorously show that both versions achieve a linear convergence rate, which is significantly better than the previous known results. We also empirically compare the proposed algorithms with several state-of-the-art matrix completion algorithms on many real-world datasets, including the large-scale recommendati...
TL;DR: A Robin-type multiphysics domain decomposition method (DDM) for the steady-state Navier--Stokes--Darcy model with three interface conditions with well-posedness is proposed and analyzed by using a branch of nonsingular solutions.
Abstract: This paper proposes and analyzes a Robin-type multiphysics domain decomposition method (DDM) for the steady-state Navier--Stokes--Darcy model with three interface conditions. In addition to the two regular interface conditions for the mass conservation and the force balance, the Beavers--Joseph condition is used as the interface condition in the tangential direction. The major mathematical difficulty in adopting the Beavers--Joseph condition is that it creates an indefinite leading order contribution to the total energy budget of the system [Y. Cao et al., Comm. Math. Sci., 8 (2010), pp. 1--25; Y. Cao et al., SIAM J. Numer. Anal., 47 (2010), pp. 4239--4256]. In this paper, the well-posedness of the Navier--Stokes--Darcy model with Beavers--Joseph condition is analyzed by using a branch of nonsingular solutions. By following the idea in [Y. Cao et al., Numer. Math., 117 (2011), pp. 601--629], the three physical interface conditions are utilized together to construct the Robin-type boundary conditions on th...
TL;DR: In this article, a model-reduction methodology that preserves Lagrangian structure and achieves computational efficiency in the presence of high-order nonlinearities and arbitrary parameter dependence is proposed.
Abstract: This work proposes a model-reduction methodology that preserves Lagrangian structure and achieves computational efficiency in the presence of high-order nonlinearities and arbitrary parameter dependence. As such, the resulting reduced-order model retains key properties such as energy conservation and symplectic time-evolution maps. We focus on parameterized simple mechanical systems subjected to Rayleigh damping and external forces, and consider an application to nonlinear structural dynamics. To preserve structure, the method first approximates the system's “Lagrangian ingredients''---the Riemannian metric, the potential-energy function, the dissipation function, and the external force---and subsequently derives reduced-order equations of motion by applying the (forced) Euler--Lagrange equation with these quantities. From the algebraic perspective, key contributions include two efficient techniques for approximating parameterized reduced matrices while preserving symmetry and positive definiteness: matri...
TL;DR: A modification of this variational model includes additional constraints via penalty terms to enforce the irreversibility of the fracture as well as the applied displacement field to numerically compute the time-evolving minimizing solution.
Abstract: The quasi-static brittle fracture model proposed by G. Francfort and J.-J. Marigo can be $\Gamma$-approximated at each time evolution step by the Ambrosio--Tortorelli functional. In this paper, we focus on a modification of this functional which includes additional constraints via penalty terms to enforce the irreversibility of the fracture as well as the applied displacement field. Second, we build on this variational model an adapted discretization to numerically compute the time-evolving minimizing solution. We present the derivation of a novel a posteriori error estimator driving the anisotropic adaptive procedure. The main properties of these automatically generated meshes are to be very fine and strongly anisotropic in a very thin neighborhood of the crack, but only far away from the crack tip, while they show a highly isotropic behavior in a neighborhood of the crack tip instead. As a consequence of these properties, the resulting discretizations follow very closely the propagation of the fracture,...
TL;DR: In this paper, the authors combine hierarchical and topological relationships between octants to design efficient recursive algorithms that operate on distributed forests of octrees for adaptive mesh refinement and coarsening.
Abstract: The forest-of-octrees approach to parallel adaptive mesh refinement and coarsening has recently been demonstrated in the context of a number of large-scale PDE-based applications. Efficient reference software has been made freely available to the public both in the form of the standalone \tt p4est library and more indirectly by the general-purpose finite element library \tt deal.II, which has been equipped with a \tt p4est backend. Although linear octrees, which store only leaf octants, have an underlying tree structure by definition, it is not fully exploited in previously published mesh-related algorithms. This is because tree branches are not explicitly stored and because the topological relationships in meshes, such as the adjacency between cells, introduce dependencies that do not respect the octree hierarchy. In this work, we combine hierarchical and topological relationships between octants to design efficient recursive algorithms that operate on distributed forests of octrees. We present three imp...
TL;DR: A Fourier method to solve backward stochastic differential equations (BSDEs) using a general theta-discretization of the time-integrands leading to an induction scheme with conditional expectations is developed.
Abstract: We develop a Fourier method to solve backward stochastic differential equations (BSDEs). A general theta-discretization of the time-integrands leads to an induction scheme with conditional expectations. These are approximated by using Fourier cosine series expansions, relying on the availability of a characteristic function. The method is applied to BSDEs with jumps. Numerical experiments demonstrate the applicability of BSDEs in financial and economic problems and show fast convergence of our efficient probabilistic numerical method.
TL;DR: A high-order and energy stable scheme is developed to simulate phase-field models by combining the semi-implicit spectral deferred correction (SDC) method and the energy stable convex splitting technique, which is found very useful for producing accurate numerical solutions at small time (dynamics) as well as long time (steady state) with reasonably large time stepsizes.
Abstract: A high-order and energy stable scheme is developed to simulate phase-field models by combining the semi-implicit spectral deferred correction (SDC) method and the energy stable convex splitting technique. The convex splitting scheme we use here is a linear unconditionally stable method but is only of first-order accuracy, so the SDC method can be used to iteratively improve the rate of convergence. However, it is found that the accuracy improvement may affect the overall energy stability which is intrinsic to the phase-field models. To compromise the accuracy and stability, a local $p$-adaptive strategy is proposed to enhance the accuracy by sacrificing some local energy stability in an acceptable level. The proposed strategy is found very useful for producing accurate numerical solutions at small time (dynamics) as well as long time (steady state) with reasonably large time stepsizes. Numerical experiments are carried out to demonstrate the high effectiveness of the proposed numerical strategy.
TL;DR: In this article, a well-balanced second order Godunov-type finite volume scheme for compressible Euler equations with gravity is presented, achieved by a specific combination of source term discretization, hydrostatic reconstruction, and numerical flux that exactly resolves stationary contacts.
Abstract: We present a novel well-balanced second order Godunov-type finite volume scheme for compressible Euler equations with gravity. The well-balanced property is achieved by a specific combination of source term discretization, hydrostatic reconstruction, and numerical flux that exactly resolves stationary contacts. The scheme is able to preserve isothermal and polytropic stationary solutions up to machine precision. It is applied on several examples using the numerical flux of Roe to demonstrate its well-balanced property and the improved resolution of small perturbations around the stationary solution.
TL;DR: A fast, parallel maximum clique algorithm for large sparse graphs that is designed to exploit characteristics of social and information networks and employs a branch-and-bound strategy with novel and aggressive pruning techniques.
Abstract: We present a fast, parallel maximum clique algorithm for large sparse graphs that is designed to exploit characteristics of social and information networks. The method exhibits a roughly linear runtime scaling over real-world networks ranging from a thousand to a hundred million nodes. In a test on a social network with 1.8 billion edges, the algorithm finds the largest clique in about 20 minutes. At its heart the algorithm employs a branch-and-bound strategy with novel and aggressive pruning techniques. The pruning techniques include the combined use of core numbers of vertices along with a good initial heuristic solution to remove the vast majority of the search space. In addition, the exploration of the search tree is parallelized. During the search, processes immediately communicate changes to upper and lower bounds on the size of the maximum clique. This exchange of information occasionally results in a superlinear speedup because tasks with large search spaces can be pruned by other processes. We de...
TL;DR: In this paper, a CUR approximate matrix factorization based on the discrete empirical interpolation method (DEIM) is proposed, and the accuracy tracks the optimal approximation error within a factor that depends on the conditioning of submatrices of the matrix.
Abstract: We derive a CUR approximate matrix factorization based on the discrete empirical interpolation method (DEIM). For a given matrix ${\bf A}$, such a factorization provides a low-rank approximate decomposition of the form ${\bf A} \approx \bf C \bf U \bf R$, where ${\bf C}$ and ${\bf R}$ are subsets of the columns and rows of ${\bf A}$, and ${\bf U}$ is constructed to make $\bf C\bf U \bf R $ a good approximation. Given a low-rank singular value decomposition ${\bf A} \approx \bf V \bf S \bf W^T$, the DEIM procedure uses ${\bf V}$ and ${\bf W}$ to select the columns and rows of ${\bf A}$ that form ${\bf C}$ and ${\bf R}$. Through an error analysis applicable to a general class of CUR factorizations, we show that the accuracy tracks the optimal approximation error within a factor that depends on the conditioning of submatrices of ${\bf V}$ and ${\bf W}$. For very large problems, ${\bf V}$ and ${\bf W}$ can be approximated well using an incremental QR algorithm that makes only one pass through ${\bf A}$. Numer...
TL;DR: This work introduces a low-rank in time technique that exploits the low-Rank nature of the solution of time-dependent PDE-constrained optimization problems and illustrates how three different problems can be rewritten and used within aLow-rank Krylov subspace solver with appropriate preconditioning.
Abstract: The solution of time-dependent PDE-constrained optimization problems is a chal- lenging task in numerical analysis and applied mathematics. All-at-once discretizations and corre- sponding solvers provide efficient methods to robustly solve the arising discretized equations. One of the drawbacks of this approach is the high storage demand for the vectors representing the discrete space-time cylinder. Here we introduce a low-rank in time technique that exploits the low-rank nature of the solution. The theoretical foundations for this approach originate in the numerical treatment of matrix equations and can be carried over to PDE-constrained optimization. We illustrate how three different problems can be rewritten and used within a low-rank Krylov subspace solver with appropriate preconditioning.
TL;DR: This study introduces two classes of regular and singular tempered fractional Sturm--Liouville problems of two kinds (TFSLP-I and TF SLP-II) of order $
u \in (0,2)$.
Abstract: Continuum-time random walk is a general model for particle kinetics, which allows for incorporating waiting times and/or non-Gaussian jump distributions with divergent second moments to account for Levy flights. Exponentially tempering the probability distribution of the waiting times and the anomalously large displacements results in tempered-stable Levy processes with finite moments, where the fluid (continuous) limit leads to the tempered fractional diffusion equation. The development of fast and accurate numerical schemes for such nonlocal problems requires a new spectral theory and suitable choice of basis functions. In this study, we introduce two classes of regular and singular tempered fractional Sturm--Liouville problems of two kinds (TFSLP-I and TFSLP-II) of order $
u \in (0,2)$. In the regular case, the corresponding tempered differential operators are associated with tempering functions $p_I(x) = \exp(2\tau) $ and $p_{II}(x) = \exp(-2\tau)$, $\tau \geq 0$, respectively, in the regular TFSLP-I...
TL;DR: This contribution is considering error estimators that are based on conservative flux reconstruction and provide an efficient and rigorous bound on the full error with respect to the weak solut...
Abstract: In this contribution we consider localized, robust, and efficient a posteriori error estimation of the localized reduced basis multiscale (LRBMS) method for parametric elliptic problems with possibly heterogeneous diffusion coefficient. The numerical treatment of such parametric multiscale problems is characterized by a high computational complexity, arising from the multiscale character of the underlying differential equation and the additional parameter dependence. The LRBMS method can be seen as a combination of numerical multiscale methods and model reduction using reduced basis (RB) methods to efficiently reduce the computational complexity with respect to the multiscale as well as the parametric aspect of the problem, simultaneously. In contrast to the classical residual based error estimators currently used in RB methods, we are considering error estimators that are based on conservative flux reconstruction and provide an efficient and rigorous bound on the full error with respect to the weak solut...
TL;DR: A alternating least squares (ALS) fit to an overrelaxation scheme inspired by the LMaFit method for matrix completion and both approaches aim at finding a tensor $A$ that fulfills the first order optimality conditions by a nonlinear Gauss--Seidel-type solver.
Abstract: We consider the problem of fitting a low rank tensor $A\in\mathbb{R}^{{\mathcal I}}$, ${\mathcal I} = \{1,\ldots,n\}^{d}$, to a given set of data points $\{M_i\in\mathbb{R} \mid i\in P\}$, $P\subset{\mathcal I}$. The low rank format under consideration is the hierarchical or tensor train or matrix product states format. It is characterized by rank bounds $r$ on certain matricizations of the tensor. The number of degrees of freedom is in ${\cal O}(r^2dn)$. For a fixed rank and mode size $n$ we observe that it is possible to reconstruct random (but rank structured) tensors as well as certain discretized multivariate (but rank structured) functions from a number of samples that is in ${\cal O}(\log N)$ for a tensor having $N=n^d$ entries. We compare an alternating least squares (ALS) fit to an overrelaxation scheme inspired by the LMaFit method for matrix completion. Both approaches aim at finding a tensor $A$ that fulfills the first order optimality conditions by a nonlinear Gauss--Seidel-type solver that c...
TL;DR: A novel numerical method to solve the problem of optimal transport and the related elliptic Monge--Ampere equation is introduced, one of the few numerical algorithms capable of solving this problem efficiently with the proper transport boundary condition.
Abstract: In this article we introduce a novel numerical method to solve the problem of optimal transport and the related elliptic Monge--Ampere equation. It is one of the few numerical algorithms capable of solving this problem efficiently with the proper transport boundary condition. The computation time scales well with the grid size and has the additional advantage that the target domain may be nonconvex. We present the method and several numerical experiments.