Singularities of the Radon transform

Singularities of the Radon transform of a piecewise smooth function $f(x)$, $x\in R^n$, $n\geq 2$, are calculated. If the singularities of the Radon transform are known, then the equations of the surfaces of discontinuity of $f(x)$ are calculated by applying the Legendre transform to the functions, which appear in the equations of the discontinuity surfaces of the Radon transform of $f(x)$; examples are given. Numerical aspects of the problem of finding discontinuities of $f(x)$, given the discontinuities of its Radon transform, are discussed.Comment: 7 page

The Vortex Solution in the (2+1)-Dimensional Yang-Mills-Chern-Simons Theory at High Temperature

The vortex-like solution to the non-linear field equations in a two-dimensional SU(2) gauge theory with the Chern-Simons mass term is found at high temperature. It is derived from the effective Lagrangian including the leading order finite temperature corrections. The discovered field configuration possesses the finite energy and the quantized magnetic flux. At the centre of the vortex the point charge is located which is surrounded by the distributed charge of the opposite sign and the vortex is neutral as a whole. At high temperature the energy of the vortex is negative and it corresponds to the ground state. The derived solution is considered to be a result of heating the lattice vacuum structure formed at zero temperature.Comment: 9 pages, LaTeX, no figures, a4, cite.st

S-fraction multiscale finite-volume method for spectrally accurate wave propagation

We develop a method for numerical time-domain wave propagation based on the model order reduction approach. The method is built with high-performance computing (HPC) implementation in mind that implies a high level of parallelism and greatly reduced communication requirements compared to the traditional high-order finite-difference time-domain (FDTD) methods. The approach is inherently multiscale, with a reference fine grid model being split into subdomains. For each subdomain the coarse scale reduced order models (ROMs) are precomputed off-line in a parallel manner. The ROMs approximate the Neumann-to-Dirichlet (NtD) maps with high (spectral) accuracy and are used to couple the adjacent subdomains on the shared boundaries. The on-line part of the method is an explicit time stepping with the coupled ROMs. To lower the on-line computation cost the reduced order spatial operator is sparsified by transforming to a matrix Stieltjes continued fraction (S-fraction) form. The on-line communication costs are also reduced due to the ROM NtD map approximation properties. Another source of performance improvement is the time step length. Properly chosen ROMs substantially improve the Courant-Friedrichs-Lewy (CFL) condition. This allows the CFL time step to approach the Nyquist limit, which is typically unattainable with traditional schemes that have the CFL time step much smaller than the Nyquist sampling rate.Comment: 5 pages, 3 figure

A nonlinear method for imaging with acoustic waves via reduced order model backprojection

We introduce a novel nonlinear imaging method for the acoustic wave equation based on data-driven model order reduction. The objective is to image the discontinuities of the acoustic velocity, a coefficient of the scalar wave equation from the discretely sampled time domain data measured at an array of transducers that can act as both sources and receivers. We treat the wave equation along with transducer functionals as a dynamical system. A reduced order model (ROM) for the propagator of such system can be computed so that it interpolates exactly the measured time domain data. The resulting ROM is an orthogonal projection of the propagator on the subspace of the snapshots of solutions of the acoustic wave equation. While the wavefield snapshots are unknown, the projection ROM can be computed entirely from the measured data, thus we refer to such ROM as data-driven. The image is obtained by backprojecting the ROM. Since the basis functions for the projection subspace are not known, we replace them with the ones computed for a known smooth kinematic velocity model. A crucial step of ROM construction is an implicit orthogonalization of solution snapshots. It is a nonlinear procedure that differentiates our approach from the conventional linear imaging methods (Kirchhoff migration and reverse time migration - RTM). It resolves all dynamical behavior captured by the data, so the error from the imperfect knowledge of the velocity model is purely kinematic. This allows for almost complete removal of multiple reflection artifacts, while simultaneously improving the resolution in the range direction compared to conventional RTM.Comment: 33 pages, 8 figure

Distance preserving model order reduction of graph-Laplacians and cluster analysis

Graph-Laplacians and their spectral embeddings play an important role in multiple areas of machine learning. This paper is focused on graph-Laplacian dimension reduction for the spectral clustering of data as a primary application. Spectral embedding provides a low-dimensional parametrization of the data manifold which makes the subsequent task (e.g., clustering) much easier. However, despite reducing the dimensionality of data, the overall computational cost may still be prohibitive for large data sets due to two factors. First, computing the partial eigendecomposition of the graph-Laplacian typically requires a large Krylov subspace. Second, after the spectral embedding is complete, one still has to operate with the same number of data points. For example, clustering of the embedded data is typically performed with various relaxations of k-means which computational cost scales poorly with respect to the size of data set. In this work, we switch the focus from the entire data set to a subset of graph vertices (target subset). We develop two novel algorithms for such low-dimensional representation of the original graph that preserves important global distances between the nodes of the target subset. In particular, it allows to ensure that target subset clustering is consistent with the spectral clustering of the full data set if one would perform such. That is achieved by a properly parametrized reduced-order model (ROM) of the graph-Laplacian that approximates accurately the diffusion transfer function of the original graph for inputs and outputs restricted to the target subset. Working with a small target subset reduces greatly the required dimension of Krylov subspace and allows to exploit the conventional algorithms (like approximations of k-means) in the regimes when they are most robust and efficient.Comment: 28 pages, 10 figure

Multi-scale S-fraction reduced-order models for massive wavefield simulations

We developed a novel reduced-order multi-scale method for solving large time-domain wavefield simulation problems. Our algorithm consists of two main stages. During the first "off-line" stage the fine-grid operator (of the graph Laplacian type} is partitioned on coarse cells (subdomains). Then projection-type multi-scale reduced order models (ROMs) are computed for the coarse cell operators. The off-line stage is embarrassingly parallel as ROM computations for the subdomains are independent of each other. It also does not depend on the number of simulated sources (inputs) and it is performed just once before the entire time-domain simulation. At the second "on-line" stage the time-domain simulation is performed within the obtained multi-scale ROM framework. The crucial feature of our formulation is the representation of the ROMs in terms of matrix Stieltjes continued fractions (S-fractions). The layered structure of the S-fraction introduces several hidden layers in the ROM representation, that results in the block-tridiagonal dynamic system within each coarse cell. This allows us to sparsify the obtained multi-scale subdomain operator ROMs and to reduce the communications between the adjacent subdomains which is highly beneficial for a parallel implementation of the on-line stage. Our approach suits perfectly the high performance computing architectures, however in this paper we present rather promising numerical results for a serial computing implementation only. These results include 3D acoustic and multi-phase anisotropic elastic problems.Comment: 31 pages, 11 figure

Robust nonlinear processing of active array data in inverse scattering via truncated reduced order models

We introduce a novel algorithm for nonlinear processing of data gathered by an active array of sensors which probes a medium with pulses and measures the resulting waves. The algorithm is motivated by the application of array imaging. We describe it for a generic hyperbolic system that applies to acoustic, electromagnetic or elastic waves in a scattering medium modeled by an unknown coefficient called the reflectivity. The goal of imaging is to invert the nonlinear mapping from the reflectivity to the array data. Many existing imaging methodologies ignore the nonlinearity i.e., operate under the assumption that the Born (single scattering) approximation is accurate. This leads to image artifacts when multiple scattering is significant. Our algorithm seeks to transform the array data to those corresponding to the Born approximation, so it can be used as a pre-processing step for any linear inversion method. The nonlinear data transformation algorithm is based on a reduced order model defined by a proxy wave propagator operator that has four important properties. First, it is data driven, meaning that it is constructed from the data alone, with no knowledge of the medium. Second, it can be factorized in two operators that have an approximately affine dependence on the unknown reflectivity. This allows the computation of the Fr\'{e}chet derivative of the reflectivity to the data mapping which gives the Born approximation. Third, the algorithm involves regularization which balances numerical stability and data fitting with accuracy of the order of the standard deviation of additive data noise. Fourth, the algebraic nature of the algorithm makes it applicable to scalar (acoustic) and vectorial (elastic, electromagnetic) wave data without any specific modifications.Comment: 26 pages, 6 figure

Untangling the nonlinearity in inverse scattering with data-driven reduced order models

The motivation of this work is an inverse problem for the acoustic wave equation, where an array of sensors probes an unknown medium with pulses and measures the scattered waves. The goal of the inversion is to determine from these measurements the structure of the scattering medium, modeled by a spatially varying acoustic impedance function. Many inversion algorithms assume that the mapping from the unknown impedance to the scattered waves is approximately linear. The linearization, known as the Born approximation, is not accurate in strongly scattering media, where the waves undergo multiple reflections before they reach the sensors in the array. Thus, the reconstructions of the impedance have numerous artifacts. The main result of the paper is a novel, linear-algebraic algorithm that uses a reduced order model (ROM) to map the data to those corresponding to the single scattering (Born) model. The ROM construction is based only on the measurements at the sensors in the array. The ROM is a proxy for the wave propagator operator, that propagates the wave in the unknown medium over the duration of the time sampling interval. The output of the algorithm can be input into any off-the-shelf inversion software that incorporates state of the art linear inversion algorithms to reconstruct the unknown acoustic impedance.Comment: 27 pages, 9 figure

Direct, nonlinear inversion algorithm for hyperbolic problems via projection-based model reduction

We estimate the wave speed in the acoustic wave equation from boundary measurements by constructing a reduced-order model (ROM) matching discrete time-domain data. The state-variable representation of the ROM can be equivalently viewed as a Galerkin projection onto the Krylov subspace spanned by the snapshots of the time-domain solution. The success of our algorithm hinges on the data-driven Gram--Schmidt orthogonalization of the snapshots that suppresses multiple reflections and can be viewed as a discrete form of the Marchenko--Gel'fand--Levitan--Krein algorithm. In particular, the orthogonalized snapshots are localized functions, the (squared) norms of which are essentially weighted averages of the wave speed. The centers of mass of the squared orthogonalized snapshots provide us with the grid on which we reconstruct the velocity. This grid is weakly dependent on the wave speed in traveltime coordinates, so the grid points may be approximated by the centers of mass of the analogous set of squared orthogonalized snapshots generated by a known reference velocity. We present results of inversion experiments for one- and two-dimensional synthetic models.Comment: 54 pages, 6 figures fixed typos and small errors expanded several sections to aid in understandin

A model reduction approach to numerical inversion for a parabolic partial differential equation

We propose a novel numerical inversion algorithm for the coefficients of parabolic partial differential equations, based on model reduction. The study is motivated by the application of controlled source electromagnetic exploration, where the unknown is the subsurface electrical resistivity and the data are time resolved surface measurements of the magnetic field. The algorithm presented in this paper considers inversion in one and two dimensions. The reduced model is obtained with rational interpolation in the frequency (Laplace) domain and a rational Krylov subspace projection method. It amounts to a nonlinear mapping from the function space of the unknown resistivity to the small dimensional space of the parameters of the reduced model. We use this mapping as a nonlinear preconditioner for the Gauss-Newton iterative solution of the inverse problem. The advantage of the inversion algorithm is twofold. First, the nonlinear preconditioner resolves most of the nonlinearity of the problem. Thus the iterations are less likely to get stuck in local minima and the convergence is fast. Second, the inversion is computationally efficient because it avoids repeated accurate simulations of the time-domain response. We study the stability of the inversion algorithm for various rational Krylov subspaces, and assess its performance with numerical experiments.Comment: 31 pages, 9 figures, 2 table