Efficient methods for stochastic optimal control

The Hamilton Jacobi Bellman (HJB) equation is central to stochastic optimal control (SOC) theory, yielding the optimal solution to general problems specified by known dynamics and a specified cost functional. Given the assumption of quadratic cost on the control input, it is well known that the HJB reduces to a particular partial differential equation (PDE). While powerful, this reduction is not commonly used as the PDE is of second order, is nonlinear, and examples exist where the problem may not have a solution in a classical sense. Furthermore, each state of the system appears as another dimension of the PDE, giving rise to the curse of dimensionality. Since the number of degrees of freedom required to solve the optimal control problem grows exponentially with dimension, the problem becomes intractable for systems with all but modest dimension.

In the last decade researchers have found that under certain, fairly non-restrictive structural assumptions, the HJB may be transformed into a linear PDE, with an interesting analogue in the discretized domain of Markov Decision Processes (MDP). The work presented in this thesis uses the linearity of this particular form of the HJB PDE to push the computational boundaries of stochastic optimal control.

This is done by crafting together previously disjoint lines of research in computation. The first of these is the use of Sum of Squares (SOS) techniques for synthesis of control policies. A candidate polynomial with variable coefficients is proposed as the solution to the stochastic optimal control problem. An SOS relaxation is then taken to the partial differential constraints, leading to a hierarchy of semidefinite relaxations with improving sub-optimality gap. The resulting approximate solutions are shown to be guaranteed over- and under-approximations for the optimal value function. It is shown that these results extend to arbitrary parabolic and elliptic PDEs, yielding a novel method for Uncertainty Quantification (UQ) of systems governed by partial differential constraints. Domain decomposition techniques are also made available, allowing for such problems to be solved via parallelization and low-order polynomials.

The optimization-based SOS technique is then contrasted with the Separated Representation (SR) approach from the applied mathematics community. The technique allows for systems of equations to be solved through a low-rank decomposition that results in algorithms that scale linearly with dimensionality. Its application in stochastic optimal control allows for previously uncomputable problems to be solved quickly, scaling to such complex systems as the Quadcopter and VTOL aircraft. This technique may be combined with the SOS approach, yielding not only a numerical technique, but also an analytical one that allows for entirely new classes of systems to be studied and for stability properties to be guaranteed.

The analysis of the linear HJB is completed by the study of its implications in application. It is shown that the HJB and a popular technique in robotics, the use of navigation functions, sit on opposite ends of a spectrum of optimization problems, upon which tradeoffs may be made in problem complexity. Analytical solutions to the HJB in these settings are available in simplified domains, yielding guidance towards optimality for approximation schemes. Finally, the use of HJB equations in temporal multi-task planning problems is investigated. It is demonstrated that such problems are reducible to a sequence of SOC problems linked via boundary conditions. The linearity of the PDE allows us to pre-compute control policy primitives and then compose them, at essentially zero cost, to satisfy a complex temporal logic specification.

Discrete Modeling of Granular Media: A NURBS-based Approach

This dissertation is concerned with the development of a new discrete element method (DEM) based on Non-Uniform Rational Basis Splines (NURBS). With NURBS, the new DEM is able to capture sphericity and angularity, the two particle morphological measures used in characterizing real grain geometries. By taking advantage of the parametric nature of NURBS, the Lipschitzian dividing rectangle (DIRECT) global optimization procedure is employed as a solution procedure to the closest-point projection problem, which enables the contact treatment of non-convex particles. A contact dynamics (CD) approach to the NURBS-based discrete method is also formulated. By combining particle shape flexibility, properties of implicit time-integration, and non-penetrating constraints, we target applications in which the classical DEM either performs poorly or simply fails, i.e., in granular systems composed of rigid or highly stiff angular particles and subjected to quasistatic or dynamic flow conditions. The CD implementation is made simple by adopting a variational framework, which enables the resulting discrete problem to be readily solved using off-the-shelf mathematical programming solvers. The capabilities of the NURBS-based DEM are demonstrated through 2D numerical examples that highlight the effects of particle morphology on the macroscopic response of granular assemblies under quasistatic and dynamic flow conditions, and a 3D characterization of material response in the shear band of a real triaxial specimen.

Full and Model-Reduced Structure-Preserving Simulation of Incompressible Fluids

This thesis outlines the construction of several types of structured integrators for incompressible fluids. We first present a vorticity integrator, which is the Hamiltonian counterpart of the existing Lagrangian-based fluid integrator. We next present a model-reduced variational Eulerian integrator for incompressible fluids, which combines the efficiency gains of dimension reduction, the qualitative robustness to coarse spatial and temporal resolutions of geometric integrators, and the simplicity of homogenized boundary conditions on regular grids to deal with arbitrarily-shaped domains with sub-grid accuracy.

Both these numerical methods involve approximating the Lie group of volume-preserving diffeomorphisms by a finite-dimensional Lie-group and then restricting the resulting variational principle by means of a non-holonomic constraint. Advantages and limitations of this discretization method will be outlined. It will be seen that these derivation techniques are unable to yield symplectic integrators, but that energy conservation is easily obtained, as is a discretized version of Kelvin's circulation theorem.

Finally, we outline the basis of a spectral discrete exterior calculus, which may be a useful element in producing structured numerical methods for fluids in the future.

Representing measures on the Royden boundary for solutions of Δu=Pu on a Riemannian manifold

Consider the Royden compactification R* of a Riemannian n-manifold R, Γ = R*\R its Royden boundary, Δ its harmonic boundary and the elliptic differential equation Δu = Pu, P ≥ 0 on R. A regular Borel measure m^P can be constructed on Γ with support equal to the closure of Δ^P = {q ϵ Δ : q has a neighborhood U in R* with _U^ʃ_ᴖR^{P ˂ ∞} }. Every enegy-finite solution to u (i.e. E(u) = D(u) + ^ʃ_Ru²P ˂ ∞, where D(u) is the Dirichlet integral of u) can be represented by u(z) = ^ʃ_Γu(q)K(z,q)dm^P(q) where K(z,q) is a continuous function on ^Rx Γ . A _P^~_E-function is a nonnegative solution which is the infimum of a downward directed family of energy-finite solutions. A nonzero _P^~_E-function is called _P^~_E-minimal if it is a constant multiple of every nonzero _P^~_E-function dominated by it. THEOREM. There exists a _P^~_E-minimal function if and only if there exists a point in q ϵ Γ such that m^P(q) > 0. THEOREM. For q ϵ Δ^P , m^P(q) > 0 if and only if m⁰(q) > 0 .

Eigenvalue problems associated with Korn's inequalities in the theory of elasticity

Interest in the possible applications of a priori inequalities in linear elasticity theory motivated the present investigation. Korn's inequality under various side conditions is considered, with emphasis on the Korn's constant. In the "second case" of Korn's inequality, a variational approach leads to an eigenvalue problem; it is shown that, for simply-connected two-dimensional regions, the problem of determining the spectrum of this eigenvalue problem is equivalent to finding the values of Poisson's ratio for which the displacement boundary-value problem of linear homogeneous isotropic elastostatics has a non-unique solution.

Previous work on the uniqueness and non-uniqueness issue for the latter problem is examined and the results applied to the spectrum of the Korn eigenvalue problem. In this way, further information on the Korn constant for general regions is obtained.

A generalization of the "main case" of Korn's inequality is introduced and the associated eigenvalue problem is a gain related to the displacement boundary-value problem of linear elastostatics in two dimensions.

Fast numerical methods for mixed, singular Helmholtz boundary value problems and Laplace eigenvalue problems - with applications to antenna design, sloshing, electromagnetic scattering and spectral geometry

This thesis presents a novel class of algorithms for the solution of scattering and eigenvalue problems on general two-dimensional domains under a variety of boundary conditions, including non-smooth domains and certain "Zaremba" boundary conditions - for which Dirichlet and Neumann conditions are specified on various portions of the domain boundary. The theoretical basis of the methods for the Zaremba problems on smooth domains concern detailed information, which is put forth for the first time in this thesis, about the singularity structure of solutions of the Laplace operator under boundary conditions of Zaremba type. The new methods, which are based on use of Green functions and integral equations, incorporate a number of algorithmic innovations, including a fast and robust eigenvalue-search algorithm, use of the Fourier Continuation method for regularization of all smooth-domain Zaremba singularities, and newly derived quadrature rules which give rise to high-order convergence even around singular points for the Zaremba problem. The resulting algorithms enjoy high-order convergence, and they can tackle a variety of elliptic problems under general boundary conditions, including, for example, eigenvalue problems, scattering problems, and, in particular, eigenfunction expansion for time-domain problems in non-separable physical domains with mixed boundary conditions.

Uncovering the Lagrangian from observations of trajectories

We approach the problem of automatically modeling a mechanical system from data about its dynamics, using a method motivated by variational integrators. We write the discrete Lagrangian as a quadratic polynomial with varying coefficients, and then use the discrete Euler-Lagrange equations to numerically solve for the values of these coefficients near the data points. This method correctly modeled the Lagrangian of a simple harmonic oscillator and a simple pendulum, even with significant measurement noise added to the trajectories.

I. Electronic processes in α-sulfur. II. Polaron motion in a D.C. electric field

Part I: The mobilities of photo-generated electrons and holes in orthorhombic sulfur are determined by drift mobility techniques. At room temperature electron mobilities between 0.4 cm²/V-sec and 4.8 cm²/V-sec and hole mobilities of about 5.0 cm²/V-sec are reported. The temperature dependence of the electron mobility is attributed to a level of traps whose effective depth is about 0.12 eV. This value is further supported by both the voltage dependence of the space-charge-limited, D.C. photocurrents and the photocurrent versus photon energy measurements.

As the field is increased from 10 kV/cm to 30 kV/cm a second mechanism for electron transport becomes appreciable and eventually dominates. Evidence that this is due to impurity band conduction at an appreciably lower mobility (4.10^-4 cm²/V-sec) is presented. No low mobility hole current could be detected. When fields exceeding 30 kV/cm for electron transport and 35 kV/cm for hole transport are applied, avalanche phenomena are observed. The results obtained are consistent with recent energy gap studies in sulfur.

The theory of the transport of photo-generated carriers is modified to include the case of appreciable thermos-regeneration from the traps in one transit time.

Part II: An explicit formula for the electric field E necessary to accelerate an electron to a steady-state velocity v in a polarizable crystal at arbitrary temperature is determined via two methods utilizing Feynman Path Integrals. No approximation is made regarding the magnitude of the velocity or the strength of the field. However, the actual electron-lattice Coulombic interaction is approximated by a distribution of harmonic oscillator potentials. One may be able to find the “best possible” distribution of oscillators using a variational principle, but we have not been able to find the expected criterion. However, our result is relatively insensitive to the actual distribution of oscillators used, and our E-v relationship exhibits the physical behavior expected for the polaron. Threshold fields for ejecting the electron for the polaron state are calculated for several substances using numerical results for a simple oscillator distribution.

Energy inequalities and error estimates for axisymmetric torsion of thin elastic shells of revolution

The problem motivating this investigation is that of pure axisymmetric torsion of an elastic shell of revolution. The analysis is carried out within the framework of the three-dimensional linear theory of elastic equilibrium for homogeneous, isotropic solids. The objective is the rigorous estimation of errors involved in the use of approximations based on thin shell theory.

The underlying boundary value problem is one of Neumann type for a second order elliptic operator. A systematic procedure for constructing pointwise estimates for the solution and its first derivatives is given for a general class of second-order elliptic boundary-value problems which includes the torsion problem as a special case.

The method used here rests on the construction of “energy inequalities” and on the subsequent deduction of pointwise estimates from the energy inequalities. This method removes certain drawbacks characteristic of pointwise estimates derived in some investigations of related areas.

Special interest is directed towards thin shells of constant thickness. The method enables us to estimate the error involved in a stress analysis in which the exact solution is replaced by an approximate one, and thus provides us with a means of assessing the quality of approximate solutions for axisymmetric torsion of thin shells.

Finally, the results of the present study are applied to the stress analysis of a circular cylindrical shell, and the quality of stress estimates derived here and those from a previous related publication are discussed.

I. Quantum mechanics of one-dimensional two-particle models. II. A quantum mechanical treatment of inelastic collisions

Part I

Solutions of Schrödinger’s equation for system of two particles bound in various stationary one-dimensional potential wells and repelling each other with a Coulomb force are obtained by the method of finite differences. The general properties of such systems are worked out in detail for the case of two electrons in an infinite square well. For small well widths (1-10 a.u.) the energy levels lie above those of the noninteresting particle model by as much as a factor of 4, although excitation energies are only half again as great. The analytical form of the solutions is obtained and it is shown that every eigenstate is doubly degenerate due to the “pathological” nature of the one-dimensional Coulomb potential. This degeneracy is verified numerically by the finite-difference method. The properties of the square-well system are compared with those of the free-electron and hard-sphere models; perturbation and variational treatments are also carried out using the hard-sphere Hamiltonian as a zeroth-order approximation. The lowest several finite-difference eigenvalues converge from below with decreasing mesh size to energies below those of the “best” linear variational function consisting of hard-sphere eigenfunctions. The finite-difference solutions in general yield expectation values and matrix elements as accurate as those obtained using the “best” variational function.

The system of two electrons in a parabolic well is also treated by finite differences. In this system it is possible to separate the center-of-mass motion and hence to effect a considerable numerical simplification. It is shown that the pathological one-dimensional Coulomb potential gives rise to doubly degenerate eigenstates for the parabolic well in exactly the same manner as for the infinite square well.

Part II

A general method of treating inelastic collisions quantum mechanically is developed and applied to several one-dimensional models. The formalism is first developed for nonreactive “vibrational” excitations of a bound system by an incident free particle. It is then extended to treat simple exchange reactions of the form A + BC →AB + C. The method consists essentially of finding a set of linearly independent solutions of the Schrödinger equation such that each solution of the set satisfies a distinct, yet arbitrary boundary condition specified in the asymptotic region. These linearly independent solutions are then combined to form a total scattering wavefunction having the correct asymptotic form. The method of finite differences is used to determine the linearly independent functions.

The theory is applied to the impulsive collision of a free particle with a particle bound in (1) an infinite square well and (2) a parabolic well. Calculated transition probabilities agree well with previously obtained values.

Several models for the exchange reaction involving three identical particles are also treated: (1) infinite-square-well potential surface, in which all three particles interact as hard spheres and each two-particle subsystem (i.e. BC and AB) is bound by an attractive infinite-square-well potential; (2) truncated parabolic potential surface, in which the two-particle subsystems are bound by a harmonic oscillator potential which becomes infinite for interparticle separations greater than a certain value; (3) parabolic (untruncated) surface. Although there are no published values with which to compare our reaction probabilities, several independent checks on internal consistency indicate that the results are reliable.

A Continuum Model for Slip-Twinning Interactions in Magnesium and Magnesium Alloys

Due to their high specific strength and low density, magnesium and magnesium-based alloys have gained great technological importance in recent years. However, their underlying hexagonal crystal structure furnishes Mg and its alloys with a complex mechanical behavior because of their comparably smaller number of energetically favorable slip systems. Besides the commonly studied slip mechanism, another way to accomplish general deformation is through the additional mechanism of deformation-induced twinning. The main aim of this thesis research is to develop an efficient continuum model to understand and ultimately predict the material response resulting from the interaction between these two mechanisms.

The constitutive model we present is based on variational constitutive updates of plastic slips and twin volume fractions and accounts for the related lattice reorientation mechanisms. The model is applied to single- and polycrystalline pure magnesium. We outline the finite-deformation plasticity model combining basal, pyramidal, and prismatic dislocation activity as well as a convexification based approach for deformation twinning. A comparison with experimental data from single-crystal tension-compression experiments validates the model and serves for parameter identification. The extension to polycrystals via both Taylor-type modeling and finite element simulations shows a characteristic stress-strain response that agrees well with experimental observations for polycrystalline magnesium. The presented continuum model does not aim to represent the full details of individual twin-dislocation interactions, yet it is sufficiently efficient to allow for finite element simulations while qualitatively capturing the underlying microstructural deformation mechanisms.

I. Non-linear effects in a self-consistent calculation of SU_3 symmetry breaking in strong interaction coupling constants. II. Direct-channel reggeization of strong interaction scattering amplitudes

The thesis is divided into two parts. Part I generalizes a self-consistent calculation of residue shifts from SU₃ symmetry, originally performed by Dashen, Dothan, Frautschi, and Sharp, to include the effects of non-linear terms. Residue factorizability is used to transform an overdetermined set of equations into a variational problem, which is designed to take advantage of the redundancy of the mathematical system. The solution of this problem automatically satisfies the requirement of factorizability and comes close to satisfying all the original equations.

Part II investigates some consequences of direct channel Regge poles and treats the problem of relating Reggeized partial wave expansions made in different reaction channels. An analytic method is introduced which can be used to determine the crossed-channel discontinuity for a large class of direct-channel Regge representations, and this method is applied to some specific representations.

It is demonstrated that the multi-sheeted analytic structure of the Regge trajectory function can be used to resolve apparent difficulties arising from infinitely rising Regge trajectories. Also discussed are the implications of large collections of "daughter trajectories."

Two things are of particular interest: first, the threshold behavior in direct and crossed channels; second, the potentialities of Reggeized representations for us in self-consistent calculations. A new representation is introduced which surpasses previous formulations in these two areas, automatically satisfying direct-channel threshold constraints while being capable of reproducing a reasonable crossed channel discontinuity. A scalar model is investigated for low energies, and a relation is obtained between the mass of the lowest bound state and the slope of the Regge trajectory.

Minimum drag profiles in infinite cavity flows

The problem considered is that of minimizing the drag of a symmetric plate in infinite cavity flow under the constraints of fixed arclength and fixed chord. The flow is assumed to be steady, irrotational, and incompressible. The effects of gravity and viscosity are ignored.

Using complex variables, expressions for the drag, arclength, and chord, are derived in terms of two hodograph variables, Γ (the logarithm of the speed) and β (the flow angle), and two real parameters, a magnification factor and a parameter which determines how much of the plate is a free-streamline.

Two methods are employed for optimization:

(1) The parameter method. Γ and β are expanded in finite orthogonal series of N terms. Optimization is performed with respect to the N coefficients in these series and the magnification and free-streamline parameters. This method is carried out for the case N = 1 and minimum drag profiles and drag coefficients are found for all values of the ratio of arclength to chord.

(2) The variational method. A variational calculus method for minimizing integral functionals of a function and its finite Hilbert transform is introduced, This method is applied to functionals of quadratic form and a necessary condition for the existence of a minimum solution is derived. The variational method is applied to the minimum drag problem and a nonlinear integral equation is derived but not solved.

Large plane deformations of thin elastic sheets of neo-hookean material

Large plane deformations of thin elastic sheets of neo-Hookean material are considered and a method of successive substitutions is developed to solve problems within the two-dimensional theory of finite plane stress. The first approximation is determined by linear boundary value problems on two harmonic functions, and it is approached asymptotically at very large extensions in the plane of the sheet. The second and higher approximations are obtained by solving Poisson equations. The method requires modification when the membrane has a traction-free edge.

Several problems are treated involving infinite sheets under uniform biaxial stretching at infinity. First approximations are obtained when a circular or elliptic inclusion is present and when the sheet has a circular or elliptic hole, including the limiting cases of a line inclusion and a straight crack or slit. Good agreement with exact solutions is found for circularly symmetric deformations. Other examples discuss the stretching of a short wide strip, the deformation near a boundary corner which is traction-free, and the application of a concentrated load to a boundary point.

Relativistic velocity - potential hydrodynamics and stellar stability

The equations of relativistic, perfect-fluid hydrodynamics are cast in Eulerian form using six scalar "velocity-potential" fields, each of which has an equation of evolution. These equations determine the motion of the fluid through the equation

U_ʋ=µ^-1 (ø,_ʋ + αβ,_ʋ + ƟS,_ʋ).

Einstein's equations and the velocity-potential hydrodynamical equations follow from a variational principle whose action is

I = (R + 16π p) (-g)^1/2 d⁴x,

where R is the scalar curvature of spacetime and p is the pressure of the fluid. These equations are also cast into Hamiltonian form, with Hamiltonian density –T₀⁰ (-g^oo)^-1/2.

The second variation of the action is used as the Lagrangian governing the evolution of small perturbations of differentially rotating stellar models. In Newtonian gravity this leads to linear dynamical stability criteria already known. In general relativity it leads to a new sufficient condition for the stability of such models against arbitrary perturbations.

By introducing three scalar fields defined by

ρ ᵴ = ∇λ + ∇x(xi + ∇xɣi)

(where ᵴ is the vector displacement of the perturbed fluid element, ρ is the mass-density, and i, is an arbitrary vector), the Newtonian stability criteria are greatly simplified for the purpose of practical applications. The relativistic stability criterion is not yet in a form that permits practical calculations, but ways to place it in such a form are discussed.