Publications

We present a fast algorithm for computing the diffracted field from arbitrary binary (sharp-edged) planar apertures and occulters in the scalar Fresnel approximation, for up to moderately high Fresnel numbers ($\lesssim 10^3$). It uses a high-order areal quadrature over the aperture, then exploits a single 2D nonuniform fast Fourier transform (NUFFT) to evaluate rapidly at target points (of order $10^7$ such points per second, independent of aperture complexity). It thus combines the high accuracy of edge integral methods with the high speed of Fourier methods. Its cost is ${\mathcal O}(n^2 \log n)$, where $n$ is the linear resolution required in source and target planes, to be compared with ${\mathcal O}(n^3)$ for edge integral methods. In tests with several aperture shapes, this translates to between 2 and 5 orders of magnitude acceleration. In starshade modeling for exoplanet astronomy, we find that it is roughly $10^4 \times$ faster than the state of the art in accurately computing the set of telescope pupil wavefronts. We provide a documented, tested MATLAB/Octave implementation.
An appendix shows the mathematical equivalence of the boundary diffraction wave, angular integration, and line integral formulae, then analyzes a new non-singular reformulation that eliminates their common difficulties near the geometric shadow edge. This supplies a robust edge integral reference against which to validate the main proposal.

With rheology applications in mind, we present a fast solver for the time-dependent effective viscosity of an infinite lattice containing one or more neutrally buoyant smooth rigid particles per unit cell, in a two-dimensional Stokes fluid with given shear rate. At each time, the mobility problem is reformulated as a 2nd-kind boundary integral equation, then discretized to spectral accuracy by the Nyström method and solved iteratively, giving typically 10 digits of accuracy. Its solution controls the evolution of particle locations and angles in a first-order system of ordinary differential equations. The formulation is placed on a rigorous footing by defining a generalized periodic Green's function for the skew lattice. Numerically, the periodized integral operator is split into a near image sum—applied in linear time via the fast multipole method—plus a correction field solved cheaply via proxy Stokeslets. We use barycentric quadratures to evaluate particle interactions and velocity fields accurately, even at distances much closer than the node spacing. Using first-order time-stepping we simulate, for example, 25 ellipses per unit cell to 3-digit accuracy on a desktop in 1 hour per shear time. Our examples show equilibration at long times, force chains, and two types of blow-ups (jamming) whose power laws match lubrication theory asymptotics.

We introduce a systematically improvable family of variational wave functions for the simulation of strongly correlated fermionic systems. This family consists of Slater determinants in an augmented Hilbert space involving "hidden" additional fermionic degrees of freedom. These determinants are projected onto the physical Hilbert space through a constraint which is optimized, together with the single-particle orbitals, using a neural network parametrization. This construction draws inspiration from the success of hidden particle representations but overcomes the limitations associated with the mean-field treatment of the constraint often used in this context. Our construction provides an extremely expressive family of wave functions, which is proven to be universal. We apply this construction to the ground state properties of the Hubbard model on the square lattice, achieving levels of accuracy which are competitive with and in some cases superior to state-of-the-art computational methods.

The variational quantum Monte Carlo (VQMC) method received significant attention in the recent past because of its ability to overcome the curse of dimensionality inherent in many-body quantum systems. Close parallels exist between VQMC and the emerging hybrid quantum-classical computational paradigm of variational quantum algorithms. VQMC overcomes the curse of dimensionality by performing alternating steps of Monte Carlo sampling from a parametrized quantum state followed by gradient-based optimization. While VQMC has been applied to solve high-dimensional problems, it is known to be difficult to parallelize, primarily owing to the Markov Chain Monte Carlo (MCMC) sampling step. In this work, we explore the scalability of VQMC when autoregressive models, with exact sampling, are used in place of MCMC. This approach can exploit distributed-memory, shared-memory and/or GPU parallelism in the sampling task without any bottlenecks. In particular, we demonstrate the GPU-scalability of VQMC for solving up to ten-thousand dimensional combinatorial optimization problems.

We initiate the study of neural-network quantum state algorithms for analyzing continuous-variable lattice quantum systems in first quantization. A simple family of continuous-variable trial wavefunctons is introduced which naturally generalizes the restricted Boltzmann machine (RBM) wavefunction introduced for analyzing quantum spin systems. By virtue of its simplicity, the same variational Monte Carlo training algorithms that have been developed for ground state determination and time evolution of spin systems have natural analogues in the continuum. We offer a proof of principle demonstration in the context of ground state determination of a stoquastic quantum rotor Hamiltonian. Results are compared against those obtained from partial differential equation (PDE) based scalable eigensolvers. This study serves as a benchmark against which future investigation of continuous-variable neural quantum states can be compared, and points to the need to consider deep network architectures and more sophisticated training algorithms.

Understanding the physics of large cosmological surveys down to small (nonlinear) scales will significantly improve our knowledge of the Universe. Large N-body simulations have been built to obtain predictions in the non-linear regime. However, N-body simulations are computationally expensive and generate large amount of data, putting burdens on storage. These data are snapshots of the simulated Universe at different times, and fine sampling is necessary to accurately save its whole history. We employ a deep neural network model to predict the nonlinear N-body simulation at an intermediate time step given two widely separated snapshots. Our results outperform the cubic Hermite interpolation benchmark method in interpolating N-body simulations. This work can greatly reduce the storage requirement and allow us to reconstruct the cosmic history from far fewer snapshots of the universe.

Cosmologists aim to model the evolution of initially low amplitude Gaussian density fluctuations into the highly non-linear "cosmic web" of galaxies and clusters. They aim to compare simulations of this structure formation process with observations of large-scale structure traced by galaxies and infer the properties of the dark energy and dark matter that make up 95% of the universe. These ensembles of simulations of billions of galaxies are computationally demanding, so that more efficient approaches to tracing the non-linear growth of structure are needed. We build a V-Net based model that transforms fast linear predictions into fully nonlinear predictions from numerical simulations. Our NN model learns to emulate the simulations down to small scales and is both faster and more accurate than the current state-of-the-art approximate methods. It also achieves comparable accuracy when tested on universes of significantly different cosmological parameters from the one used in training. This suggests that our model generalizes well beyond our training set.

An identification is found between meta-learning and the problem of determining the ground state of a randomly generated Hamiltonian drawn from a known ensemble. A model-agnostic meta-learning approach is proposed to solve the associated learning problem and a preliminary experimental study of random Max-Cut problems indicates that the resulting Meta Variational Monte Carlo accelerates training and improves convergence.

Multivariate regression techniques are commonly applied to explore the associations between large numbers of outcomes and predictors. In real-world applications, the outcomes are often of mixed types, including continuous measurements, binary indicators, and counts, and the observations may also be incomplete. Building upon the recent advances in mixed-outcome modeling and sparse matrix factorization, generalized co-sparse factor regression (GOFAR) is proposed, which utilizes the flexible vector generalized linear model framework and encodes the outcome dependency through a sparse singular value decomposition (SSVD) of the integrated natural parameter matrix. To avoid the estimation of the notoriously difficult joint SSVD, GOFAR proposes both sequential and parallel unit-rank estimation procedures. By combining the ideas of alternating convex search and majorization–minimization, an efficient algorithm is developed to solve the sparse unit-rank problem and implemented in the R package gofar. Extensive simulation studies and two real-world applications demonstrate the effectiveness of the proposed approach.

The Bayesian approach to data analysis provides a powerful way to handle uncertainty in all observations, model parameters, and model structure using probability theory. Probabilistic programming languages make it easier to specify and fit Bayesian models, but this still leaves us with many options regarding constructing, evaluating, and using these models, along with many remaining challenges in computation. Using Bayesian inference to solve real-world problems requires not only statistical skills, subject matter knowledge, and programming, but also awareness of the decisions made in the process of data analysis. All of these aspects can be understood as part of a tangled workflow of applied Bayesian statistics. Beyond inference, the workflow also includes iterative model building, model checking, validation and troubleshooting of computational problems, model understanding, and model comparison. We review all these aspects of workflow in the context of several examples, keeping in mind that in practice we will be fitting many models for any given problem, even if only a subset of them will ultimately be relevant for our conclusions.