# Cornell Theory Center Technical Reports

## Permanent URI for this collection

This is a collection of technical reports from the Cornell Theory Center from the time period of 1990-2005. These reports are part of the NCSTRL collection of Computer Science Technical Reports.

## Browse

### Recent Submissions

Item Accurate Solution of Weighted Least Squares by Iterative MethodsBobrovnikova, Elena Y.; Vavasis, Stephen A. (Cornell University, 1997-02-06)We consider the weighted least-squares (WLS) problem with a very ill-conditioned weight matrix. Weighted least-squares problems arise in many applications including linear programming, electrical networks, boundary value problems, and structures. Because of roundoff errors, standard iterative methods for solving a WLS problem with ill-conditioned weights may not give the correct answer. Indeed, the difference between the true and computed solution (forward error) may be large. We propose an iterative algorithm, called MINRES-L, for solving WLS problems. The MINRES-L method is the application of MINRES, a Krylov-space method due to Paige and Saunders, to a certain layered linear system. Using a simplified model of the effects of round off error, we prove that MINRES-L gives answers with small forward error. We present computational experiments for some applications.Item Local correlation energies of two-electron atoms and model systemsHuang, Chien-Jung; Umrigar, C.J. (Cornell University, 1997-01)We present nearly-local definitions of the correlation energy density, and its potential and kinetic components, and evaluate them for several two-electron systems. This information should provide valuable guidance in constructing better correlation functionals than those in common use. In addition, we demonstrate that the quantum chemistry and the density functional definitions of the correlation energy rapidly approach one another with increasing atomic number.Item Quality Mesh Generation in Higher DimensionsMitchell, Scott A.; Vavasis, Stephen A. (Cornell University, 1996-12)We consider the problem of triangulating a d-dimensional region. Our mesh generation algorithm, called QMG, is a qradtree-based algorithm that can triangulate any polyhedral region including nonconvexregions with holes. Furthermore, our algorithm guarantees a bounded aspect ratio triangulation provided that the input domain itself has no sharp angles. Finally, our algorithm is guaranteed never to over refine the domain in the sense that the number of simplices produced by QMG is bounded above by a factor times the number produced by any competing algorithm, where the factor depends on the aspect ratio bound satisfied by the competing algorithm. The QMG algorithm has been implemented in C++ and is used as a mesh generator for the finite element method.Item Model fermion Monte Carlo with correlated pairs IIKalos, M.H.; Schmidt, K.E. (Cornell University, 1996-12)Correlated dynamics can produce stable algorithms for excited states of quantum many-body problems. We study a variety of harmonic oscillator problems to demonstrate the kinds of correlations needed. We show that marginally correct dynamics that produce a stable overlap with an antisymmetrictrial function give the correct fermion ground state.Item Symmetry, Nonlinear Bifurcation Analysis, and Parallel ComputationWohlever, J.C. (Cornell University, 1996-10)In the natural and engineering sciences the equations which model physical systems with symmetry often exhibit an invariance with respect to a particular group "G" of linear transformations. "G" is typically a linear representation of a symmetry group "g" which characterizes the symmetry of the physical system. In this work, we will discuss the natural parallelism which arises while seeking families of solutions to a specific class of nonlinear vector equations which display a special type of group invariance, referred to as equivariance. The inherent parallelism stems for a global de-coupling, due to symmetry, of the full nonlinear equations which effectively splits the original problem into a set of smaller problems. Numerical results from asymmetry-adapted numerical procedure, (MMcontcm.m), written in MultiMATLAB are discussed.Item Monte Carlo Optimization of Trial Wave Functions in Quantum Mechanics and Statistical MechanicsNightingale, M.P.; Umrigar, C.J. (Cornell University, 1996-10)This review covers applications of quantum Monte Carlo methods to quantum mechanical problems in the study of electronic and atomic structure, as well as applications to statistical mechanical problems both of static and dynamic nature. The common thread in all these applications is optimization of many-parameter trial states, which is done by minimization of the variance of the local energy or, more generally for arbitrary eigenvalue problems, minimization of the variance of the configurational eigenvalue.Item A critical assessment of the Self-Interaction Corrected Local Density Functional method and its algorithmic implementationGoedecker, S.; Umrigar, C.J. (Cornell University, 1996-10)We calculate the electronic structure of several atoms and small molecules by direct minimization of the Self-Interaction Corrected Local Density Approximation (SIC-LDA) functional. To do this we first derive an expression for the gradient of this functional under the constraint that the orbitals be orthogonal and show that previously given expressions do not correctly incorporate this constraint. In our atomic calculations the SIC-LDA yields total energies, ionization energies and charge densities that are superior to results obtained with the Local Density Approximation (LDA). However, for molecules SIC-LDA gives bond lengths and reaction energies that are inferior to those obtained from LDA. The nonlocal BLYP functional, which we include as a representative GGA functional, out performs both LDA and SIC-LDA forall ground state properties we considered.Item European Option Pricing with Fixed Transaction CostsAiyer, Ajay Subramanian (Cornell University, 1996-09)In this paper, we study the problem of European option pricing in the presence of fixed transaction costs. The problems of optimal portfolio selection and option pricing in the presence of proportional transaction costs has been extensively studied in the mathematical finance literature. However, much less is known when we have fixed transaction costs. In this paper, we show that calculating the price of an European optioninvolves calculating the value functions of two stochastic impulse control problems and we obtain the explicit expressions for the resultant quasi-variational ine qualities satisfied by the value functions and then carry out a numerical calculation of the option price.Item Optimal Portfolio Selection with Fixed Transactions Costs in the presence of Jumps and Random DriftAiyer, Ajay Subramanian (Cornell University, 1996-09)In this paper, we study the general problem of optimal portfolio selection with fixed transactions costs in the presence of jumps. We extend the analysis of Morton and Pliska to this setting by modeling the return processes of the risky assets in the investor's portfolio as jump-diffusion processes and derive the expression for the related optimal stopping time problem of a Markov process with jumps and explicitly solve it in the situation when the portfolio consists only of one risky asset. We also provide an asymptotic analysis of our model with one risky asset following the ideas of Wilmott and Atkinson. In the process, we also obtain a solution for the "Merton problem" generalized to the situation when there is credit risk. Finally, we consider the case where the drift of the stockprice process is random and unobservable and obtain expressions for the optimal trading policies.Item Non-normal Dynamics and Hydrodynamic StabilityBaggett, Jeffrey Scott (Cornell University, 1996-08)This thesis explores the interaction of non-normality and nonlinearity incontinuous dynamical systems. A solution beginning near a linearly stable fixed point may grow large by a linear mechanism, if the linearization is non-normal, until it is swept away by nonlinearities resulting in a much smaller basin of attraction than could possibly be predicted by the spectrum of the linearization. Exactly this situation occurs in certain linearly stable shear flows, where the linearization about the laminar flow may be highly non-normal leading to the transient growth of certain small disturbances by factors which scale with the Reynolds number. These issues are brought into focus in Chapter 1 through the study of atwo-dimensional model system of ordinary differential equations proposed by Trefethen, et al. [Science, 261, 1993]. In Chapter 2, two theorems are proved which show that the basin of attraction of a stable fixed point, in systems of differential equations combining a non-normal linear term with quadratic nonlinearities, can decrease rapidly as the degree of non-normality is increased, often faster than inverse linearly. Several different low-dimensional models of transition to turbulence are examined in Chapter 3. These models were proposed by more than a dozen authors for a wide variety of reasons, but they all incorporate non-normal linear terms and quadratic nonlinearities. Surprisingly, in most cases, the basin of attraction of the "laminar flow" shrinks much faster than the inverse Reynolds number. Transition to turbulence from optimally growing linear disturbances, streamwise vortices, is investigated in plane Poiseuille and plane Couette flows in Chapter4. An explanation is given for why smaller streamwise vortices can lead to turbulence in plane Poiseuille flow. In plane Poiseuille flow, the transient linear growth of streamwise streaks caused by non-normality leads directly to a secondary instability. Certain unbounded operators are so non-normal that the evolution of infinitesimal perturbations to the fixed point is entirely unrelated to the spectrum, even as i to infinity. Two examples of this phenomenonare presented in Chapter 5.Item Structure and Efficient Hessian CalculationColeman, Thomas F.; Verma, Arun (Cornell University, 1996-08)Modern methods for numerical optimization calculate (or approximate) the matrix of second derivatives, the Hessian matrix, at each iteration. The recent arrival of robust software for automatic differentiation allows for the possibility of automatically computing the Hessian matrix, and the gradient, given a code to evaluate the objective function itself. However, for large-scale problems direct application of automatic differentiation may be unacceptably expensive. Recent work has shown that this cost can be dramatically reduced in the presence of sparsity. In this paper we show that for structured problems it is possible to apply automatic differentiation tools in an economical way - even in the absence of sparsity in the Hessian.Item Stable and Efficient Solution of Weighted Least-Squares Problems with Applications in Interior Point MethodsHough, Patricia D. (Cornell University, 1996-08)In this thesis, we consider two closely related problems. The first is a full-rank weighted least-squares problem with a weight matrix that is positive definite, diagonal, and extremely ill conditioned. The ill-conditioning can cause standard algorithms to compute solutions with, in some cases, no digits of accuracy. Theory suggests the existence of an algorithm that will compute an accurate solution despite the ill-conditioning in the weight matrix. We describe a new algorithm, the Complete Orthogonal Decomposition (COD) Algorithm, for solving the weighted least-squares problem and show that it has this desirable property. In addition, the COD Algorithm is based on standard, well-understood techniques and is straightforward to implement. A natural application for the weighted least-squares problem described in the previous paragraph is interior point methods for linear programming. We discuss the problem in this context, and describe how the COD algorithm can be extended and used in this setting. Unlike other algorithms, this one is stable for interior point methods without assuming nondegeneracy in the linear programming instance. Computational experiments indicate that it is more reliable than other algorithms when the problem is near degenerate. The second problem involves a particular interior point algorithm. In 1994, Vavasis and Ye proposed a new primal-dual path-following interior point method, the layered-step interior point (LIP) method. This algorithm interleaves traditional steps with longer, layered least-squares (LLS) steps. Computation of a LLS step requires solving a weighted least-squares problem similar to the one described above, but the weight matrix also has the property that the weights fall into well-separated groups. This additional structure allows the problem to be broken down into smaller, constrained problems with well-conditioned weight matrices. The smaller problems can then be solved stably with standard algorithms, and the LLS step can be computed. Vavasis and Ye did not propose a particular algorithm for solving the LLS problem. In this thesis, we present an algorithm based on Cholesky factorization. The algorithm is such that a modified version of the sparse Cholesky code of Ng and Peyton of Oak Ridge National Laboratories can be used. Thus, the theoretical results are straight forward, and this algorithm proves to be accurate and efficient in practice.Item The Procrustes Problem for Orthogonal Stiefel MatricesBojanczyk, A.W.; Lutoborski, A. (Cornell University, 1996-08)(This abstract contains mathematical symbols that may not reproduce wellin ASCII text.) In this paper we consider the Procrustes problem on the manifold of orthogonal Stiefel matrices. That is, given matrices A epsilon R(m*k), B epsilon R(m*p), m greater tahn or equal to p greater than or equal to k, we seek the minimum of ||A - BQ||2 for all matrices Q epsilon R(p*k), QTQ = I(k*k). We introduce a class of relaxation methods for generating minimizing sequences and offer a geometric interpretation of these methods. Results of numerical experiments illustrating the convergence of the methods are given.Item Scalable Parallel Electronic Structure Calculations on the IBM SP2Goedecker, Stefan; Hoisie, Adolfy (Cornell University, 1996-08)We have developed a highly efficient and scalable electronic structure code for parallel computers using message passing. The algorithm takes advantages of the natural parallelism in quantum chemistry problems to obtain very high performance even on a large number of processors. Most of the terms which scale cubically with respect to the number of atoms have been eliminated allowing the treatment of very large systems. It uses one of the most precise versions of Density Functional Theory, namely Self-Interaction Corrected Density Functional Theory.Item On the Efficient Methods to Solve ODEs and BVPs Using Automatic DifferentiationVerma, Arun (Cornell University, 1996-08)A large number of physical phenomena are modeled by a system of ODEs or a system of implicit ODEs. We demonstrate applicability of automatic differentiation (AD) for solving: (1) Boundary value problems in ODEs and implicit ODEs. (2) Initial state and parameter estimation problems. The impact of using AD is two fold. Firstly, efficient methods for computing the gradient vectors and Jacobian matrices have been developed using AD. Secondly the process of getting derivatives via AD is robust, more user friendly, and provides error free derivatives. Furthermore, techniques using AD have been developed which exploit structure in the user's computation, and particularly the structure we observe in boundary value problems or state/parameter estimation problems. We demonstrate by a few experiments the efficiency gained by the usage of AD in solving these problems.Item Piecewise Differentiable Minimization for Ill-posed Inverse ProblemsLi, Yuying (Cornell University, 1996-08)Based on minimizing a piece wise differentiable lp function subject to a single inequality constraint, this paper discusses algorithms for a discretized regularization problem for ill-posed inverse problems. We examine computational challenges of solving this regularization problem. Possible minimization algorithms such as the steepest descent method, iteratively weighted least squares (IRLS) method and a recent globally convergent affine scaling Newton approach are considered. Limitations and efficiency of these algorithms are demonstrated using the geophysical travel time tomographic inversion and image restoration applications.Item Transformation Techniques for Toeplitz and Toeplitz-plus-Hankel Matrices Part II. AlgorithmsBojanczyk, A.W.; Heinig, George (Cornell University, 1996-08)In the first part of the paper transformations mapping Toeplitz and Toeplitz-plus-Hankel matrices into generalized Cauchy matrices were studied. In this second part fast algorithms for LU-factorization and inversion of generalized Cauchy matrices are discussed. It is shown that the combination of transformation pivoting techniques leads to algorithms for indefinite Toeplitz and Toeplitz-plus-Hankel matrices that are more stable than the classical ones Special attention is paid to the symmetric and hermitian cases.Item Transformation Techniques for Toeplitz and Toeplitz-plus-Hankel Matrices Part I. TransformationsBojanczyk, A.W.; Heinig, George (Cornell University, 1996-08)Transformations of the form A to C1AC2 are investigated that transform Toeplitz and Toeplitz-plus-Hankel matrices into generalized Cauchy matrices. C1 and C2 are matrices related to the discrete Fourier transformation or to various real trigonometric transformations. Combining these results with pivoting techniques, in part II algorithms for Toeplitz and Toeplitz-plus-Hankel systems will be presented that are more stable than classical algorithms.Item Generalized gradient approximations to density functional theory: comparison with exact resultsFilippi, Claudia; Gonze, Xavier; Umrigar, C.J. (Cornell University, 1996-06)In order to assess the accuracy of commonly used approximate exchange-correlation density functionals, we present a comparison of accurate exchange and correlation potentials, exchange energy densities and energy components with the corresponding approximate quantities. Four systems are used as illustrative examples: the model system of two electrons in a harmonic potential and the De, Be and Ne atoms. A new ingredient in the paper is the separation of the exchange-correlation potential into exchange and correlation according to the density functional theory definition.Item Separation of the Exchange-Correlation Potential into Exchange plus Correlation: an Optimized Effective Potential ApproachFilippi, Claudia; Umrigar, C.J.; Gonze, Xavier (Cornell University, 1996-06)Most approximate exchange-correlation functionals used within density functional theory are constructed as the sum of two distinct contributions for exchange and correlation. Separating the exchange component from the entire functional is useful since, for exchange, exact relations exist under uniform density scaling and spin scaling. In the past, accurate exchange-correlation potentials have been generated from essentially exact densities but they have not been correctly decomposed into their separate exchange and correlation components (except for two-electron systems). Using a recently proposed method, equivalent to the solution of an optimized effective potential problem with the corresponding orbitals replaced by the exact Kohn-Sham orbitals, we obtain the separation according to the density functional theory definition. We compare the results for the Ne and Be atoms with those obtained by the previously used appromixate separation scheme.