CAAM Technical Reports

Permanent URI for this collection

https://hdl.handle.net/1911/101542

Computational and Applied Mathematics (CAAM) technical reports 1981-present. Use the "Browse This Collection" sidebar to filter by date, author, etc.

Browse

Now showing 1 - 20 of 719

Discretization of Multipole Sources in a Finite Difference Setting for Wave Propagation Problems
(2018-06-20) Bencomo, Mario J.; Symes, William
Seismic sources are commonly idealized as point-sources due to their small spatial extent relative to seismic wavelengths. The acoustic isotropic point-radiator is inadequate as a model of seismic wave generation for seismic sources that are known to exhibit directivity. Therefore, accurate modeling of seismic wavefields must include source representations generating anisotropic radiation patterns. Such seismic sources can be modeled as linear combinations of multipole point-sources. In this paper we present a method for discretizing multipole sources in a finite difference setting, an extension of the moment matching conditions developed for the Dirac delta function in other applications. We also provide the necessary analysis and numerical evidence to demonstrate the accuracy of our singular source approximations. In particular, we develop a weak convergence theory for the discretization of a family of symmetric hyperbolic systems of first-order partial differential equations, with singular source terms, solved via staggered-grid finite difference methods. Numerical experiments demonstrate a stronger result than what is presented in our convergence theory, namely, optimal convergence rates of numerical solutions are achieved point-wise in space away from the source if an appropriate source discretization is used.
Nonlinear Waveform Inversion with Surface-Oriented Extended Modeling
(2017-03) Terentyev, Igor
This thesis investigates surface-oriented model extension approach to nonlinear full waveform inversion (FWI). Conventional least-squares (LS) approach is capable of reconstructing highly detailed models of subsurface. Resolution requirements of the realistic problems dictate the use of local descent methods to solve the LS optimization problem. However, in the setting of any characteristic seismic problem, LS objective functional has numerous local extrema, rendering descent methods unsuitable when initial estimate is not kinematically accurate. The aim of my work is to improve convexity properties of the objective functional. I use the extended modeling approach, and construct an extended optimization functional incorporating differential semblance condition. An important advantage of surface-oriented extensions is that they do not increase the computational complexity of the forward modeling. This approach blends FWI technique with migration velocity analysis (MVA) capability to recover long scale velocity model, producing optimization problems that combine global convergence properties of the MVA with data fitting approach of FWI. In particular, it takes into account nonlinear physical effects, such as multiple reflections. I employ variable projection approach to solve the extended optimization problem. I validate the method on synthetic models for the constant density acoustics problem.
Efficient estimation of coherent risk measures for risk-averse optimization problems governed by partial differential equations with random inputs
(2017-05) Takhtaganov, Timur
This thesis assesses and designs structure-exploiting methods for the efficient estimation of risk measures of quantities of interest in the context of optimization of partial differential equations (PDEs) with random inputs. Risk measures of the quantities of interest arise as objective functions or as constraints in the PDE-constrained optimization problems under uncertainty. A single evaluation of a risk measure requires numerical integration in a high-dimensional parameter space, which requires the solution of the PDE at many parameter samples. When the integrand is smooth in the random parameters, efficient methods, such as sparse grids, exist that substantially reduce the sample size. Unfortunately, many risk-averse formulations, such as semideviation and Conditional Value-at-Risk, introduce a non-smoothness in integrand. This work demonstrates that naive application of sparse grids and other smoothness-exploiting approaches is not beneficial in the risk-averse case. For the widely used class of coherent risk measures, this thesis proposes a new method for evaluating risk-averse objectives based on the biconjugate representation of coherent risk functions and importance sampling. The method is further enhanced by utilizing reduced order models of the PDEs under consideration. The proposed method leads to substantial reduction in the number of PDE solutions required to accurately estimate coherent risk measures. The performance of existing and of the new methods for the estimation of risk measures is demonstrated on examples of risk-averse PDE-constrained optimization problems. The resulting method can substantially reduce the number of PDE solutions required to solve optimization problems, and, therefore, enlarge the applicability of important risk measures for PDE-constrained optimization problems under uncertainty.
Analysis of inverse boundary value problems for elastic waves
(2018-04) Zhai, Jian
In seismology, people use waves generated by earthquakes, artificial explosions, or even “noises”, to detect the Earth’s interior structure. The waves traveling in rocks, which are the main components of Earth’s crust and mantle, are elastic waves. A key problem in seismology is whether the interior structures of interest could be uniquely determined by the measurements people can collect. Since in most situations, people can only record the vibrations at the surface, such problems usually are formulated as inverse boundary value problems for the elastic wave equation. We study several inverse boundary value problems for the elastic wave equation and give some uniqueness and stability results for them. In Chapter 2, we consider the inverse problem of determining the Lamé parameters and the density of a three-dimensional elastic body from the time-harmonic Dirichletto-Neumann map. We prove uniqueness and Lipschitz stability of this inverse problem when the Lamé parameters and the density are assumed to be piecewise constant on a given domain partition. In Chapter 3 and 4, we study the recovery of the density and stiffness tensors of a three-dimensional domain from the dynamical Dirichlet-to-Neumann map. We give explicit reconstruction schemes for the determination of stiffness tensors and the density at the boundary under certain assumptions. In Chapter 4, we prove uniqueness of piecewise analytic parameters in the interior for which we have boundary determination. In Chapter 5, we give a semiclassical description of surface waves or modes in an elastic medium that is stratified near its boundary at some scale comparable to the wave length. The analysis is based on the work of Colin de Verdière [48] on acoustic surface waves. In Chapter 6, we consider a Riemannian manifold in dimension n ≥ 3 with strictly convex boundary. We prove the local invertibility, up to potential fields, of the geodesic X-ray transform on tensor fields of order four near a boundary point. Under the condition that the manifold can be foliated with a continuous family of strictly convex surfaces, the local invertibility implies a global result.
Subset Selection and Feature Identification in the Electrocardiogram
(2018-04) Hendryx, Emily
Each feature in the electrocardiogram (ECG) corresponds to a different part of the cardiac cycle. Tracking changes in these features over long periods of time can offer insight regarding changes in a patient’s clinical status. However, the automated identification of features in some patient populations, such as the pediatric congenital heart disease population, remains a nontrivial task that has yet to be mastered. Working toward a solution to this problem, this thesis outlines an overall framework for the identification of individual features in the ECGs of different populations. With a goal of applying part of this framework retrospectively to large sets of patient data, we focus primarily on the selection of relevant subsets of ECG beats for subsequent interpretation by clinical experts. We demonstrate the viability of the discrete empirical interpolation method (DEIM) in identifying representative subsets of beat morphologies relevant for future classification models. The success of DEIM applied to data sets from a variety of contexts is compared to results from related approaches in numerical linear algebra, as well as some more common clustering algorithms. We also present a novel extension of DEIM, called E-DEIM, in which additional representative data points can be identified as important without being limited by the rank of the corresponding data matrix. This new algorithm is evaluated on two different data sets to demonstrate its use in multiple settings, even beyond medicine. With DEIM and its related methods identifying beat-class representatives, we then propose an approach to automatically extend physician expertise on the selected beat morphologies to new and unlabeled beats. Using a fuzzy classification scheme with dynamic time warping, we are able to provide preliminary results suggesting further pursuit of this framework in application to patient data.
Graph Coloring, Zero Forcing, and Related Problems
(2017-05) Brimkov, Boris
This thesis investigates several problems related to classical and dynamic coloring of graphs, and enumeration of graph attributes. In the first part of the thesis, I present new efficient methods to compute the chromatic and flow polynomials of specific families of graphs. The chromatic and ow polynomials count the number of ways to color and assign ow to the graph, and encode information relevant to various physical applications. The second part of the thesis focuses on zero forcing-- a dynamic graph coloring process whereby at each discrete time step, a colored vertex with a single uncolored neighbor forces that neighbor to become colored. Zero forcing has applications in linear algebra, quantum control, and power network monitoring. A connected forcing set is a connected set of initially colored vertices which forces the entire graph to become colored; the connected forcing number is the cardinality of the smallest connected forcing set. I present a variety of structural results about connected forcing, such as the effects of vertex and edge operations on the connected forcing number, the relations between the connected forcing number and other graph parameters, and the computational complexity of connected forcing. I also give efficient algorithms for computing the connected forcing numbers of different families of graphs, and characterize the graphs with extremal connected forcing numbers. Finally, I investigate several enumeration problems associated with zero forcing, such as the exponential growth of certain families of forcing sets, relations of families of forcing sets to matroids and greedoids, and polynomials which count the number of distinct forcing sets of a given size.
An Inverse Free Projected Gradient Descent Method for the Generalized Eigenvalue Problem
(2017-05) Camacho, Frankie
The generalized eigenvalue problem is a fundamental numerical linear algebra problem whose applications are wide ranging. For truly large-scale problems, matrices themselves are often not directly accessible, but their actions as linear operators can be probed through matrix-vector multiplications. To solve such problems, matrix-free algorithms are the only viable option. In addition, algorithms that do multiple matrix-vector multiplications simultaneously (instead of sequentially), or so-called block algorithms, generally have greater parallel scalability that can prove advantageous on highly parallel, modern computer architectures. In this work, we propose and study a new inverse-free, block algorithmic framework for generalized eigenvalue problems that is based on an extension of a recent framework called eigpen|an unconstrained optimization formulation utilizing the Courant Penalty function. We construct a method that borrows several key ideas, including projected gradient descent, back-tracking line search, and Rayleigh-Ritz (RR) projection. We establish a convergence theory for this framework. We conduct numerical experiments to assess the performance of the proposed method in comparison to two well-known existing matrix-free algorithms, as well as to the popular solver arpack as a benchmark (even though it is not matrix-free). Our numerical results suggest that the new method is highly promising and worthy of further study and development.
Hermite Methods for the Simulation of Wave Propagation
(2017-05) Vargas, Arturo
Simulations of wave propagation play a crucial role in science and engineering. In applications of geophysics, they are the engine of many seismic imaging algorithms. For electrical engineers, they can be a useful tool for the design of radars and antennas. In these applications achieving high fidelity, simulations are challenging due to the inherent issues in modeling highly oscillatory waves and the associated high computational cost of high-resolution simulations. Thus the ideal numerical method should be able to capture high-frequency waves and be suitable for parallel computing. In both seismic applications and computational electromagnetics the Yee scheme, a finite difference time domain (FDTD) method, is the method of choice for structured grids. The scheme has the benefit of being easy to implement but performs poorly in the presence of high-frequency waves. High order accurate FDTD methods may be derived but ultimately rely on neighboring grid points when approximating derivative. In contrast to FDTD methods, the Hermite methods of Goodrich and co-authors (2006) use Hermite interpolation and a staggered (dual) grid to construct high order accurate numerical methods for first order hyperbolic equations. These methods achieve high order approximations in both time and space by reconstructing local polynomials within cells of the computational domain and employing Hermite-Taylor time stepping. The resulting schemes are able to evolve the solution locally within a cell making them ideal for parallel computing. Building on the original Hermite methods this thesis focuses on two goals: (1) the development of new Hermite methods and (2) their implementation on modern computing architectures. To accomplish the first objective, this thesis presents two variations of Hermite methods which are designed to simplify the scheme while preserving the favorable features. The first variation is a family of Hermite methods which do not require a dual grid. These methods eliminate the need for storing dual coefficients while maintaining optimal convergence rates. The second type of variation are Hermite methods which use leapfrog time-stepping. These schemes propagate the solution with less computation than the original scheme and may be used for either first or second order equations. To address the second objective, this thesis presents algorithms which take advantage of the many-core architecture of graphics processing units (GPU). As threedimensional simulations can easily exceed the memory of a single GPU, techniques for partitioning the data across multiple GPUs are presented. Finally, this thesis presents numerical results and performance studies which confim the accuracy and efficiency of the proposed Hermite methods for linear and nonlinear wave equations.
Novel Techniques for the Zero-Forcing and p-Median Graph Location Problems
(2017-05) Fast, Caleb C.
This thesis presents new methods for solving two graph location problems, the p-Median problem and the zero-forcing problem. For the p-median problem, I present a branch decomposition based method that finds the best p-median solution that is limited to some input support graph. The algorithm can be used to either find an integral solution from a fractional linear programming solution, or it can be used to improve on the solutions given by a pool of heuristics. In either use, the algorithm compares favorably in running time or solution quality to state-of-the-art heuristics. For the zero-forcing problem, this thesis gives both theoretical and computational results. In the theoretical section, I show that the branchwidth of a graph is a lower bound on its zero-forcing number, and I introduce new bounds on the zero-forcing iteration index for cubic graphs. This thesis also introduces a special type of graph structure, a zero-forcing fort, that provides a powerful tool for the analysis and modeling of zero-forcing problems. In the computational section, I introduce multiple integer programming models for finding minimum zero-forcing sets and integer programming and combinatorial branch and bound methods for finding minimum connected zero-forcing sets. While the integer programming methods do not perform better than the best combinatorial method for the basic zero-forcing problem, they are easily adapted to the connected zero-forcing problem, and they are the best methods for the connected zero-forcing problem.
A Parallel-In-Time Gradient-Type Method For Optimal Control Problems
(2017-05) Deng, Xiaodi
This thesis proposes and analyzes a new parallel-in-time gradient-type method for time-dependent optimal control problems. When the classical gradient method is applied to such problems, each iteration requires the forward solution of the state equations followed by the backward solution of the adjoint equations before the gradient can be computed and control can be updated. The solution of the state/adjoint equations is computationally expensive and consumes most of the computation time. The proposed parallel-in-time gradient-type method introduces parallelism by splitting the time domain into N subdomains and executes the forward and backward computation in each time subdomain in parallel using state and adjoint variables at time subdomain boundaries from the last optimization iteration as initial values. The proposed method is generalized to allow different time domain partitions for forward/backward computations and overlapping time subdomains. Convergence analyses for the parallel-in-time gradient-type method applied to discrete-time optimal control problems are established. For linear-quadratic problems, the method is interpreted as a multiple-part splitting method and convergence is proven by analyzing the corresponding iteration matrix. In addition, the connection of the new method to the multiple shooting reformulation of the problem is revealed and an alternative convergence proof based on this connection is established. For general non-linear problems, the new method is combined with metric projection to handle bound constraints on the controls and convergence of the method is proven. Numerically, the parallel-in-time gradient-type method is applied to linear-quadratic optimal control problems and to a well-rate optimization problem governed by a system of highly non-linear partial differential equations. For linear-quadratic problems, the method exhibits strong scaling with up to 50 cores. The parallelization in time is on top of the existing parallelization in space to solve the state/adjoint equations. This is exploited in the well-rate optimization problem. If the existing parallelism in space scales well up to K processors, the addition of time domain decomposition by the proposed method scales well up to K X N processors for small to moderate number N of time subdomains.
Bilevel Clique Interdiction and Related Problems
(2017-05) Becker, Timothy Joseph
I introduce a formulation of the bilevel clique interdiction problem. Interdiction, a military term, describes the removal of enemy resources. The single level clique interdiction problem describes the attempt of an attacker to interdict a maximum number of cliques. The bilevel form of the problem introduces a defender who attempts to minimize the number of cliques interdicted by the attacker. An algorithm and formulation for the bilevel clique interdiction problem has not previously been investigated. I start by introducing a formulation and a column-generation algorithm to solve the problem of bilevel interdiction of a minimum clique transversal and move forward to the creation of a delayed row-and-column generation algorithm for bilevel clique interdiction. Next, I introduce a formulation and algorithm to solve the bilevel interdiction of a maximum stable set problem. Bilevel interdiction of a maximum stable set is choosing a maximum stable set, but with a defender who is attempting to minimize the maximum stable set that can be chosen by the interdictor. I introduce a deterministic formulation and a delayed column generation algorithm. Additionally, I introduce a stochastic formulation of the problem. I solve this problem using a cross-decomposition method that involves L-shaped cuts into a master problem as well as new "clique" cuts for the inner problem. Lastly, I define new classes of valid inequalities and facets for the clique transversal polytope. The valid inequalities come from two graph structures who have a closed form for their vertex cover number, which we use as a specic case for nding a minimum clique transversal. The first class of facets are just the maximal clique constraints of the clique transversal polytope. The next class contains an odd hole with distinct cliques on each edge of the hole. Another similar class contains an odd clique with distinct maximal cliques on the edges of one of its spanning cycles. The fourth class contains a clique with distinct maximal cliques on every edge of the initial clique, while the last class is a prism graph with distinct maximal cliques on every edge of the prism.
GPU-Accelerated Discontinuous Galerkin Methods on Hybrid Meshes: Applications in Seismic Imaging
(2017-05) Wang, Zheng
Seismic imaging is a geophysical technique assisting in the understanding of subsurface structure on a regional and global scale. With the development of computer technology, computationally intensive seismic algorithms have begun to gain attention in both academia and industry. These algorithms typically produce high-quality subsurface images or models, but require intensive computations for solving wave equations. Achieving high-fidelity wave simulations is challenging: first, numerical wave solutions may suffer from dispersion and dissipation errors in long-distance propagations; second, the efficiency of wave simulators is crucial for many seismic applications. High-order methods have advantages of decreasing numerical errors efficiently and hence are ideal for wave modelings in seismic problems. Various high order wave solvers have been studied for seismic imaging. One of the most popular solvers is the finite difference time domain (FDTD) methods. The strengths of finite difference methods are the computational efficiency and ease of implementation, but the drawback of FDTD is the lack of geometric flexibility. It has been shown that standard finite difference methods suffer from first order numerical errors at sharp media interfaces. In contrast to finite difference methods, discontinuous Galerkin (DG) methods, a class of high-order numerical methods built on unstructured meshes, enjoy geometric flexibility and smaller interface errors. Additionally, DG methods are highly parallelizable and have explicit semi-discrete form, which makes DG suitable for large-scale wave simulations. In this dissertation, the discontinuous Galerkin methods on hybrid meshes are developed and applied to two seismic algorithms---reverse time migration (RTM) and full waveform inversion (FWI). This thesis describes in depth the steps taken to develop a forward DG solver for the framework that efficiently exploits the element specific structure of hexahedral, tetrahedral, prismatic and pyramidal elements. In particular, we describe how to exploit the tensor-product property of hexahedral elements, and propose the use of hex-dominant meshes to speed up the computation. The computational efficiency is further realized through a combination of graphics processing unit (GPU) acceleration and multi-rate time stepping. As DG methods are highly parallelizable, we build the DG solver on multiple GPUs with element-specific kernels. Implementation details of memory loading, workload assignment and latency hiding are discussed in the thesis. In addition, we employ a multi-rate time stepping scheme which allows different elements to take different time steps. This thesis applies DG schemes to RTM and FWI to highlight the strengths of the DG methods. For DG-RTM, we adopt the boundary value saving strategy to avoid data movement on GPUs and utilize the memory load in the temporal updating procedure to produce images of higher qualities without a significant extra cost. For DG-FWI, a derivation of the DG-specific adjoint-state method is presented for the fully discretized DG system. Finally, sharp media interfaces are inverted by specifying perturbations of element faces, edges and vertices.
Numerical Methods and Applications for Reduced Models of Blood Flow
(2017-05) Puelz, Charles
The human cardiovascular system is a vastly complex collection of interacting components, including vessels, organ systems, valves, regulatory mechanisms, mi- crocirculations, remodeling tissue, and electrophysiological signals. Experimental, mathematical, and computational research efforts have explored various hemody- namic questions; the scope of this literature is a testament to the intricate nature of cardiovascular physiology. In this work, we focus on computational modeling of blood ow in the major vessels of the human body. We consider theoretical questions related to the numerical approximation of reduced models for blood ow, posed as nonlinear hyperbolic systems in one space dimension. Further, we apply this modeling framework to abnormal physiologies resulting from surgical intervention in patients with congenital heart defects. This thesis contains three main parts: (i) a discussion of the implementation and analysis for numerical discretizations of reduced models for blood ow, (ii) an investigation of solutions to different classes of models in the realm of smooth and discontinuous solutions, and (iii) an application of these mod- els within a multiscale framework for simulating ow in patients with hypoplastic left heart syndrome. The two numerical discretizations studied in this thesis are a characteristics{based method for approximating the Riemann{invariants of reduced blood ow models, and a discontinuous Galerkin scheme for approximating solutions to the reduced models directly. A priori error estimates are derived in particular cases for both methods. Further, two classes of hyperbolic systems for blood ow, namely the mass{momentum and the mass{velocity formulations, are systematically compared with each numerical method and physiologically relevant networks of ves- sels and boundary conditions. Lastly, closed loop vessel network models of various Fontan physiologies are constructed. Arterial and venous trees are built from net- works of one{dimensional vessels while the heart, valves, vessel junctions, and organ beds are modeled by systems of algebraic and ordinary differential equations.
Projection-Based Model Reduction in the Context of Optimization with Implicit PDE Constraints
(2017-05) Magruder, Caleb
I use reduced order models (ROMs) to substantially decrease the computational cost of Newton's method for large-scale time-dependent optimal control problems in settings where solving the implicit constraints and their adjoints can be prohibitively expensive. Reservoir management in particular can take weeks to solve the state and adjoint equations necessary to generate one gradient evaluation. Model order reduction has potential to reduce the cost; however, ROMs are valid only nearby their training and must be retrained as the optimization progresses. I will demonstrate that in the case of reservoir management, frequent retraining defeats of purpose the reduced order modeling the state equation altogether. Rather than generating a reduced order model to replace the original implicit constraint, I use the structure of Hessian and subspace-based ROMs to compute approximate reduced order Hessian information by recycling data from the full order state and adjoint solves in the gradient computation. The resulting method often enjoys convergence rates similar to those of Newton's Method, but at the computational cost of the gradient method. I demonstrate my approach on nonlinear parabolic optimal control problems on two cases: (1) a semilinear parabolic case with nonlinear reaction terms and (2) the well rate optimization problem in reservoir engineering. For semilinear parabolic case, I consider a cubic reaction term and an exponential reaction term modeling solid fuel injection. My results show a dramatic speedup in wall clock time for the entire optimization procedure. This overall acceleration includes the training and precomputation of the ROMs that is conventionally discounted as "offline" computational costs.
Representation and Estimation of Seismic Sources via Multipoles
(2017-05) Bencomo, Mario J.
Accurate representation and estimation of seismic sources are essential to the seismic inversion problem. General sources can be approximated by a truncated series of multipoles depending on the source anisotropy. Most research in the joint inversion of source and medium parameters assumes seismic sources can be modeled as isotropic point-sources resulting in an inability to fit the anisotropy observed in data, ultimately impacting the recovery of medium parameters. In this thesis I lay the groundwork for joint source-medium parameter inversion with potentially anisotropic seismic sources via full waveform inversion through three key contributions: a mathematical and computational framework for the modeling and inversion of sources via multipoles, construction and analysis of discretizations of multipole sources on regular grids, and preconditioners based on fractional time derivative/integral operators for the ill-conditioned source estimation subproblem. As an application of my multipole framework, I also study the efficacy of multipoles in modeling the airgun array source, the most common type of active source in marine seismic surveying. Inversion results recovered a dominating isotropic component of the multipole source model that accounted for 84% of the observed radiation pattern. An extra 10% of the observed output pressure field can be explained when incorporating dipole terms in the source representation, thus motivating the use of multipoles to capture source anisotropy.
Black Oil Simulation Utilizing a Central Finite Volume Scheme
(2016-05) Chinomona, Rujeko
Black-oil simulation is a valuable tool in predicting the multi-phase multi-component flow of fluids in reservoirs. This research validates the use of a central high resolution finite volume scheme developed by Kurganov and Tadmor (KT) to the black-oil model problem. The KT scheme is desirable because of its relative ease of implementation and its ability to generate high resolution solutions at low computational costs. In addition, convergence rates on simple hyperbolic conservation law problems provided in this thesis indicate that the KT scheme is second order. Results obtained from simulations are in alignment with published literature and simulations also accommo- date changing variables with predictable outcomes. The KT scheme can be applied to the black-oil model problem with increased confidence.
Born Waveform Inversion in Shot Coordinate Domain
(2016-05) Huang, Yin
The goal of this thesis is to integrate Born waveform inversion, variable projection algorithm and model extension concept to get a method that can improve the long scale background model updates reliably and efficiently from seismic data. Born waveform inversion is a partially linearized version of full waveform inversion based on Born (linearized) modeling, in which the earth model is separated into a smooth long scale background model and an oscillatory short scale reflectivity and both are updated to fit observed trace data. Because kinematic variables (background model) are updated, Born waveform inversion has the same feature as full waveform inversion: very sensitive to initial model when solved by gradient based optimization method (almost the only possible method because of the problem scale). Extended Born waveform inversion allows reflectivity to depend on additional parameters, potentially enlarging the convexity domain by enlarging the searching model space and permitting data fit throughout the inversion process and in turn reducing the sensitivity to initial model. Extended or not, the Born waveform inversion objective function is quadratic in the reflectivity, so that a nested optimization approach is available: minimize over reflectivity in an inner stage, then minimize the background-dependent result in a second, outer stage, which results in a reduced objective function of the background model only (VPE method). This thesis integrates the nested optimization approach into an inversion scheme and analyzes that the accuracy of the solution to the inner optimization is crucial for a robust outer optimization and both model extension and the nested optimization are necessary for a successful Born waveform inversion. And then we propose a flexibly preconditioned least squares migration scheme (FPCG) that significantly improves the convergence of iterative least squares migration and produces high resolution images with balanced amplitude. The proposed scheme also improves the efficiency of the solution to the inner stage of the nested optimization scheme and the accuracy of the gradient, and thus potentially improves the efficiency of the VPE method. However, a theoretical error estimate in the gradient computation of the VPE method is still hard to obtain, and we explain the reason and illustrate with numerical examples.
Nonnormality in Lyapunov Equations
(2016-05) Baker, Jonathan
The singular values of the solution to a Lyapunov equation determine the potential accuracy of the low-rank approximations constructed by iterative methods. Low-rank solutions are more accurate if most of the singular values are small, so a-priori bounds that describe coefficient matrix properties that correspond to rapid singular value decay are valuable. Previous bounds take similar forms, all of which weaken (quadratically) as the coefficient matrix departs from normality. Such bounds suggest that the farther from normal the coefficient matrix is, the slower the singular values of the solution will decay. This predicted slowing of decay is manifest in the ADI algorithm for solving Lyapunov equations, which converges more slowly when the coefficient is far from normal. However, simple examples typically exhibit an eventual acceleration of decay if the coefficient departs sufficiently from normality. This thesis shows that the decay acceleration principle is universal: decay always improves as departure from normality increases beyond a given threshold, specifically, as the numerical range of the coefficient matrix extends farther into the right half-plane. The concluding chapter gives examples showing that similar behavior can occur for general Sylvester equations, though the right-hand side plays a more important role.
Clique Generalizations and Related Problems
(2016-05) Wood, Cynthia Ivette
A large number of real-world problems can be model as optimization problems in graphs. The clique model was introduced to aid the study of network structure for social interaction. Each vertex represented an actor and the edges represented the relations between them. Nevertheless, the model has been shown to be restrictive for modeling real-world problems, since it leaves out subgraphs that do not have all possible edges. As a consequence, clique generalizations were introduced to overcome the disadvantages of the clique model. In this thesis, I present three computationally difficult combinatorial optimization problems related to clique generalization problems: co-2-plexes and k-cores. A k-core is a subgraph with minimum degree greater than or equal to k. In this work, I discuss the minimal k-core problem and the minimum k-core problem. I present a backtracking algorithm to find all minimal k-cores of a given undirected graph and its applications to the study of associative memory. The proposed method is a modification of the Bron and Kerbosch algorithm for finding all cliques of an undirected graph. In addition, I study the polyhedral structure of the k-core polytope. The minimum k-core problem is modeled as a binary integer program and relaxed as a linear program. Since the relaxation yields to a non-integral solution, cuts must be added in order to improve the solution. I show that edge and cycle transversals of the graph give valid inequalities for the convex hull of k-cores. A set of pairwise non-adjacent vertices defines a stable set. A stable set is the complement of a clique. A co-2-plex is a subgraph with degree less than or equal to one, and it is a stable set relaxation. I introduce a study of the maximum weighted co-2-plex (MWC2P) problem for {claw, bull}-free graphs and present two polynomial time algorithms to solve it. One of the algorithms transforms the original graph to solve an instance of the maximum weighted stable set problem utilizing Minty’s algorithm. The second algorithm is an extension of Minty’s algorithm and solves the problem in the original graph. All the algorithms discussed in this thesis were implemented and tested. Numerical results are provided for each one of them.
Accelerating Convergence by Augmented Rayleigh-Ritz Projections For Large-Scale Eigenpair Computation
(2016-01) Wen, Zaiwen; Zhang, Yin
Iterative algorithms for large-scale eigenpair computation are mostly based subspace projections consisting of two main steps: a subspace update (SU) step that generates bases for approximate eigenspaces, followed by a Rayleigh-Ritz (RR) projection step that extracts approximate eigenpairs. A predominant methodology for the SU step makes use of Krylov subspaces that builds orthonormal bases piece by piece in a sequential manner. On the other hand, block methods such as the classic (simultaneous) subspace iteration, allow higher levels of concurrency than what is reachable by Krylov subspace methods, but may suffer from slow convergence. In this work, we analyze the rate of convergence for a simple block algorithmic framework that combines an augmented Rayleigh-Ritz (ARR) procedure with the subspace iteration. Our main results are Theorem 4.5 and its corollaries which show that the ARR procedure can provide significant accelerations to convergence speed. Our analysis will offer useful guidelines for designing and implementing practical algorithms from this framework.

Browse

Recent Submissions