Publications

With the funding provided by this award, we developed numerical codes for the study of magnetically confined plasmas for fusion applications. Accordingly, our work can be divided into two separate categories: 1) the design and analysis of novel numerical methods providing high accuracy and high efficiency; 2) the study of the equilibrium and stability of magnetically confined plasmas with some of these numerical codes, as well as the study of the nature of the turbulent behavior which may arise in the presence of instabilities. We first developed new numerical schemes based on integral equation methods for the computation of steady-state magnetic configurations in fusion experiments, providing high accuracy for the magnetic field and its derivatives, which are required for stability and turbulence calculations. We employed different integral formulations depending on the application of interest: axisymmetric or non-axisymmetric equilibria, force-free or magnetohydrodynamic equilibria, fixed-boundary equilibria or free-boundary equilibria. While efficient, these methods do not yet apply to plasma boundaries which are not smooth, a situation which is fairly common in magnetic confinement experiments. To address this temporary weakness, we also constructed a new steady-state solver based on the Hybridizable Discontinuous Galerkin (HDG) method, which provides full geometric flexibility. In addition to these numerical tools focused on steady-states, we also contributed to the improvement of the speed and accuracy of codes simulating the plasma dynamics of fusion plasmas, by developing a novel velocity space representation for the efficient solution of kinetic equations, which most accurately describe the time evolution of hot plasmas in fusion experiments. Using the tools discussed above, we studied several questions pertaining to the equilibrium and stability of magnetically confined plasmas. In particular, we derived a new simple model for axisymmetric devices called tokamaks, to predict how elongated a fusion plasma can be before it becomes unstable and collapses. We also looked at the effect of the shape of the outer plasma surface on key properties of the steady-state magnetic configurations, and how these properties impact turbulence in fusion plasmas, and the corresponding transport of momentum. Likewise, we studied the role of large localized flows on the steady-state magnetic configurations, and how they may influence plasma stability and turbulence. Non-axisymmetric steady-state magnetic configurations are inherently more complex than axisymmetric steady-state configurations, and the subject of ongoing controversies regarding the regularity of the equations determining such steady-states, and their solutions. Implementing an existing NYU code in a new geometry, we studied the nature of the singularity of the solutions observed in the code, and methods to eliminate them. Our main conclusion is that by appropriately tailoring the plasma boundary, it is possible to eliminate the singularities otherwise appearing in our simulations, and to obtain steady-states which appear to be smooth. To gain further insights on incompletely understood turbulence phenomena, we proposed a new reduced model capturing most of these phenomena, which is simple enough to not require expensive numerical simulations on massive supercomputers to investigate them. We demonstrated the strong similarity between our simulations and published results obtained from computationally expensive simulations, and plan to rely on our reduced model to identify the key mechanisms determining the evolution and strength turbulent driven transport in fusion plasmas. Finally, we proposed a new framework for tokamak reactor design studies, enabling us to consider the relative merits of steady-state versus pulsed fusion reactors. We found that pulsed fusion reactors may benefit most from recent advances in magnet technology, and the availability of very high field magnets. As such, they may become more desirable than steady-state tokamak reactors for cost efficient electricity generation.

A new scheme is proposed to construct a $\mathcal{C}^n\$ function extension for smooth functions defined on a smooth domain $D\in \mathbb{R}^d$. Unlike the PUX scheme, which requires the extrapolation of the volume grid via an expensive ill-conditioned least squares fitting, the scheme relies on an explicit formula consisting of a linear combination of function values in $D,$ which only extends the function along the normal direction. To be more precise, the $\mathcal{C}^n$ extension requires only $n+1$ function values along the normal directions in the original domain and ensures $\mathcal{C}^n$ smoothness by construction. When combined with a shrinking function and a smooth window function, the scheme can be made stable and robust for a broad class of domains with complex smooth boundary.

In this paper, the random iterative method is introduced to massive multiple-input multiple-output (MIMO) systems for the efficient downlink linear precoding. By adopting the random sampling into the traditional iterative methods, the matrix inversion within the linear precoding schemes can be approximated statistically, which not only achieves a faster exponential convergence with low complexity but also experiences a global convergence without suffering from the various convergence requirements. Specifically, based on the random iterative method, the randomized iterative precoding algorithm (RIPA) is firstly proposed and we show its approximation error decays exponentially and globally along with the number of iterations. Then, with respect to the derived convergence rate, the concept of conditional sampling is introduced, so that further optimization and enhancement are carried out to improve both the convergence and the efficiency of the randomized iterations. After that, based on the equivalent iteration transformation, the modified randomized iterative precoding algorithm (MRIPA) is presented, which achieves a better precoding performance with low-complexity for various scenarios of massive MIMO. Finally, simulation results based on downlink precoding in massive MIMO systems are given to show the system gains of RIPA and MRIPA in terms of performance and complexity.

We present an efficient basis for imaginary time Green's functions based on a low rank decomposition of the spectral Lehmann representation. The basis functions are simply a set of well-chosen exponentials, so the corresponding expansion may be thought of as a discrete form of the Lehmann representation using an effective spectral density which is a sum of δ functions. The basis is determined only by an upper bound on the product βωmax, with β the inverse temperature and ωmax an energy cutoff, and a user-defined error tolerance ϵ. The number r of basis functions scales as (log(βωmax)log(1/ϵ)). The discrete Lehmann representation of a particular imaginary time Green's function can be recovered by interpolation at a set of r imaginary time nodes. Both the basis functions and the interpolation nodes can be obtained rapidly using standard numerical linear algebra routines. Due to the simple form of the basis, the discrete Lehmann representation of a Green's function can be explicitly transformed to the Matsubara frequency domain, or obtained directly by interpolation on a Matsubara frequency grid. We benchmark the efficiency of the representation on simple cases, and with a high precision solution of the Sachdev-Ye-Kitaev equation at low temperature. We compare our approach with the related intermediate representation method, and introduce an improved algorithm to build the intermediate representation basis and a corresponding sampling grid.

This case study is a reimplementation of the algorithm described in Barmherzig and Sun (2022) [1] as a Stan model. This requires the new features available in Stan 2.30

In this paper, we introduce a randomized iterative method for signal detection in uplink large-scale multiple-input multiple-output (MIMO) systems, which not only achieves a low computational complexity but also enjoys a global and exponentially fast convergence. First of all, by adopting the random sampling into the iterations, the randomized iterative detection algorithm (RIDA) is proposed for large-scale MIMO systems. We show that RIDA converges exponentially fast in terms of mean squared error (MSE). Furthermore, this global convergence always holds, and does not depend on the standard requirements such as N≫K , where N and K denote the numbers of antennas at the sides of base station and users. This broadly extends the applications of low-complexity detection in uplink large-scale MIMO systems. Then, based on a new conditional sampling, optimization and enhancements are given to further improve both the convergence and efficiency of RIDA, resulting in the modified randomized iterative detection algorithm (MRIDA). Meanwhile, with respect to MRIDA, further complexity reduction by exploiting the matrix structure is given while its implementation by deep neural networks (DNN) is also presented for a better detection performance.

Marked power spectra are two-point statistics of a marked field obtained by weighting each location with a function that depends on the local density around that point. We consider marked power spectra of the galaxy field in redshift space that up-weight low density regions, and perform a Fisher matrix analysis to assess the information content of this type of statistics using the Molino mock catalogs built upon the Quijote simulations. We identify four different ways to up-weight the galaxy field, and compare the Fisher information contained in their marked power spectra to the one of the standard galaxy power spectrum, when considering monopole and quadrupole of each statistic. Our results show that each of the four marked power spectra can tighten the standard power spectrum constraints on the cosmological parameters Om, Ob, h, ns, Mν by 15−25 and on s8 by a factor of 2. The same analysis performed by combining the standard and four marked power spectra shows a substantial improvement compared to the power spectrum constraints that is equal to a factor of 6 for σ8 and 2.5−3 for the other parameters. Our constraints may be conservative, since the galaxy number density in the Molino catalogs is much lower than the ones in future galaxy surveys, which will allow them to probe lower density regions of the large-scale structure.

We introduce a family of neural quantum states for the simulation of strongly interacting systems in the presence of spatial periodicity. Our variational state is parametrized in terms of a permutationally invariant part described by the Deep Sets neural-network architecture. The input coordinates to the Deep Sets are periodically transformed such that they are suitable to directly describe periodic bosonic systems. We show example applications to both one- and two-dimensional interacting quantum gases with Gaussian interactions, as well as to $4He$ confined in a one-dimensional geometry. For the one-dimensional systems we find very precise estimations of the ground-state energies and the radial distribution functions of the particles. In two dimensions we obtain good estimations of the ground-state energies, comparable to results obtained from more conventional methods.

Quaternions are important for a wide variety of rotation-related problems in computer graphics, machine vision, and robotics. We study the nontrivial geometry of the relationship between quaternions and rotation matrices by exploiting the adjugate matrix of the characteristic equation of a related eigenvalue problem to obtain the manifold of the space of a quaternion eigenvector. We argue that quaternions parameterized by their corresponding rotation matrices cannot be expressed, for example, in machine learning tasks, as single-valued functions: the quaternion solution must instead be treated as a manifold, with different algebraic solutions for each of several single-valued sectors represented by the adjugate matrix. We conclude with novel constructions exploiting the quaternion adjugate variables to revisit several classic pose estimation applications: 2D point-cloud matching, 2D point-cloud-to-projection matching, 3D point-cloud matching, 3D orthographic point-cloud-to-projection matching, and 3D perspective point-cloud-to-projection matching. We find an exact solution to the 3D orthographic least squares pose extraction problem, and apply it successfully also to the perspective pose extraction problem with results that improve on existing methods.

We propose Pathfinder, a variational method for approximately sampling from differentiable log densities. Starting from a random initialization, Pathfinder locates normal approximations to the target density along a quasi-Newton optimization path, with local covariance estimated using the inverse Hessian estimates produced by the optimizer. Pathfinder returns draws from the approximation with the lowest estimated Kullback-Leibler (KL) divergence to the true posterior. We evaluate Pathfinder on a wide range of posterior distributions, demonstrating that its approximate draws are better than those from automatic differentiation variational inference (ADVI) and comparable to those produced by short chains of dynamic Hamiltonian Monte Carlo (HMC), as measured by 1-Wasserstein distance. Compared to ADVI and short dynamic HMC runs, Pathfinder requires one to two orders of magnitude fewer log density and gradient evaluations, with greater reductions for more challenging posteriors. Importance resampling over multiple runs of Pathfinder improves the diversity of approximate draws, reducing 1-Wasserstein distance further and providing a measure of robustness to optimization failures on plateaus, saddle points, or in minor modes. The Monte Carlo KL divergence estimates are embarrassingly parallelizable in the core Pathfinder algorithm, as are multiple runs in the resampling version, further increasing Pathfinder's speed advantage with multiple cores.