emgr - The Empirical Gramian Framework

System Gramian matrices are a well-known encoding for properties of input-output systems such as controllability, observability or minimality. These so-called system Gramians were developed in linear system theory for applications such as model order reduction of control systems. Empirical Gramian are an extension to the system Gramians for parametric and nonlinear systems as well as a data-driven method of computation. The empirical Gramian framework - emgr - implements the empirical Gramians in a uniform and configurable manner, with applications such as Gramian-based (nonlinear) model reduction, decentralized control, sensitivity analysis, parameter identification and combined state and parameter reduction

An inexact Newton-Krylov algorithm for constrained diffeomorphic image registration

We propose numerical algorithms for solving large deformation diffeomorphic image registration problems. We formulate the nonrigid image registration problem as a problem of optimal control. This leads to an infinite-dimensional partial differential equation (PDE) constrained optimization problem. The PDE constraint consists, in its simplest form, of a hyperbolic transport equation for the evolution of the image intensity. The control variable is the velocity field. Tikhonov regularization on the control ensures well-posedness. We consider standard smoothness regularization based on $H^1$- or $H^2$-seminorms. We augment this regularization scheme with a constraint on the divergence of the velocity field rendering the deformation incompressible and thus ensuring that the determinant of the deformation gradient is equal to one, up to the numerical error. We use a Fourier pseudospectral discretization in space and a Chebyshev pseudospectral discretization in time. We use a preconditioned, globalized, matrix-free, inexact Newton-Krylov method for numerical optimization. A parameter continuation is designed to estimate an optimal regularization parameter. Regularity is ensured by controlling the geometric properties of the deformation field. Overall, we arrive at a black-box solver. We study spectral properties of the Hessian, grid convergence, numerical accuracy, computational efficiency, and deformation regularity of our scheme. We compare the designed Newton-Krylov methods with a globalized preconditioned gradient descent. We study the influence of a varying number of unknowns in time. The reported results demonstrate excellent numerical accuracy, guaranteed local deformation regularity, and computational efficiency with an optional control on local mass conservation. The Newton-Krylov methods clearly outperform the Picard method if high accuracy of the inversion is required.Comment: 32 pages; 10 figures; 9 table

Reconstruction of hidden 3D shapes using diffuse reflections

We analyze multi-bounce propagation of light in an unknown hidden volume and demonstrate that the reflected light contains sufficient information to recover the 3D structure of the hidden scene. We formulate the forward and inverse theory of secondary and tertiary scattering reflection using ideas from energy front propagation and tomography. We show that using careful choice of approximations, such as Fresnel approximation, greatly simplifies this problem and the inversion can be achieved via a backpropagation process. We provide a theoretical analysis of the invertibility, uniqueness and choices of space-time-angle dimensions using synthetic examples. We show that a 2D streak camera can be used to discover and reconstruct hidden geometry. Using a 1D high speed time of flight camera, we show that our method can be used recover 3D shapes of objects "around the corner"

Quantifying Nonlinearity Susceptibility via Site-Response Modeling Uncertainty at Three Sites in the Los Angeles Basin

The effects of near-surface soil stratigraphy on the amplitude and frequency content of ground motion are accounted for in most modern U.S. seismic design codes for building structures as a function of the soil conditions prevailing in the area of interest. Nonetheless, currently employed site-classification criteria do not adequately describe the nonlinearity susceptibility of soil formations, which prohibits the development of standardized procedures for the computationally efficient integration of nonlinear ground response analyses in broadband ground-motion simulations. In turn, the lack of a unified methodology for nonlinear site-response analyses affects both the prediction accuracy of site-specific ground-motion intensity measures and the evaluation of site-amplification factors when broadband simulations are used for the development of hybrid attenuation relations. In this article, we introduce a set of criteria for quantification of the nonlinearity susceptibility of soil profiles based on the site conditions and incident ground-motion characteristics, and we implement them to identify the least complex ground response prediction methodology required for the simulation of nonlinear site effects at three sites in the Los Angeles basin. The criteria are developed on the basis of a comprehensive nonlinear site-response modeling uncertainty analysis, which includes both detailed soil profile descriptions and statistical adequacy of ground-motion time histories. Approximate and incremental nonlinear models are implemented, and the limited site-response observations are initially compared to the ensemble site-response estimates. A suite of synthetic ground motions for rupture scenarios of weak, medium, and large magnitude events (M 3.5–7.5) is next generated, parametric studies are conducted for each fixed magnitude scenario by varying the source-to-site distance, and the variability introduced in ground-motion predictions is quantified for each nonlinear site-response methodology. A frequency index is developed to describe the frequency content of incident ground motion relative to the resonant frequencies of the soil profile, and this index is used in conjunction with the rock-outcrop acceleration peak amplitude (PGA_(RO)) to identify the site conditions and ground-motion characteristics where incremental nonlinear analyses should be employed in lieu of approximate methodologies. We show that the proposed intensity-frequency representation of ground motion may be implemented to describe the nonlinearity susceptibility of soil formations in broadband simulations by accounting both for the magnitude-distance-orientation characteristics of seismic motion and the profile stiffness characteristics. The synthetic ground-motion predictions are next used for the development of site-amplification factors for the alternative site-response methodologies, and the results are compared to published site factors of attenuation relations. For the site conditions investigated, currently established amplification factors compare well with synthetic simulations for class C and D site conditions, while long-period amplification factors are overestimated by a factor of 1.5 at the class E site, where site-specific nonlinear analyses should be employed for levels of PGA_(RO)>0.2g

Convex inner approximations of nonconvex semialgebraic sets applied to fixed-order controller design

We describe an elementary algorithm to build convex inner approximations of nonconvex sets. Both input and output sets are basic semialgebraic sets given as lists of defining multivariate polynomials. Even though no optimality guarantees can be given (e.g. in terms of volume maximization for bounded sets), the algorithm is designed to preserve convex boundaries as much as possible, while removing regions with concave boundaries. In particular, the algorithm leaves invariant a given convex set. The algorithm is based on Gloptipoly 3, a public-domain Matlab package solving nonconvex polynomial optimization problems with the help of convex semidefinite programming (optimization over linear matrix inequalities, or LMIs). We illustrate how the algorithm can be used to design fixed-order controllers for linear systems, following a polynomial approach

Space-time adaptive solution of inverse problems with the discrete adjoint method

Adaptivity in both space and time has become the norm for solving problems modeled by partial differential equations. The size of the discretized problem makes uniformly refined grids computationally prohibitive. Adaptive refinement of meshes and time steps allows to capture the phenomena of interest while keeping the cost of a simulation tractable on the current hardware. Many fields in science and engineering require the solution of inverse problems where parameters for a given model are estimated based on available measurement information. In contrast to forward (regular) simulations, inverse problems have not extensively benefited from the adaptive solver technology. Previous research in inverse problems has focused mainly on the continuous approach to calculate sensitivities, and has typically employed fixed time and space meshes in the solution process. Inverse problem solvers that make exclusive use of uniform or static meshes avoid complications such as the differentiation of mesh motion equations, or inconsistencies in the sensitivity equations between subdomains with different refinement levels. However, this comes at the cost of low computational efficiency. More efficient computations are possible through judicious use of adaptive mesh refinement, adaptive time steps, and the discrete adjoint method. This paper develops a framework for the construction and analysis of discrete adjoint sensitivities in the context of time dependent, adaptive grid, adaptive step models. Discrete adjoints are attractive in practice since they can be generated with low effort using automatic differentiation. However, this approach brings several important challenges. The adjoint of the forward numerical scheme may be inconsistent with the continuous adjoint equations. A reduction in accuracy of the discrete adjoint sensitivities may appear due to the intergrid transfer operators. Moreover, the optimization algorithm may need to accommodate state and gradient vectors whose dimensions change between iterations. This work shows that several of these potential issues can be avoided for the discontinuous Galerkin (DG) method. The adjoint model development is considerably simplified by decoupling the adaptive mesh refinement mechanism from the forward model solver, and by selectively applying automatic differentiation on individual algorithms. In forward models discontinuous Galerkin discretizations can efficiently handle high orders of accuracy, $h/p$-refinement, and parallel computation. The analysis reveals that this approach, paired with Runge Kutta time stepping, is well suited for the adaptive solutions of inverse problems. The usefulness of discrete discontinuous Galerkin adjoints is illustrated on a two-dimensional adaptive data assimilation problem