T. S. Haut

T. S. Haut, G. Beylkin

The need to compute small con-eigenvalues and the associated con-eigenvectors of positive-definite Cauchy matrices naturally arises when constructing rational approximations with a (near) optimally small $L^{\infty}$ error. Specifically, given a rational function with $n$ poles in the unit disk, a rational approximation with $m\ll n$ poles in the unit disk may be obtained from the $m$th con-eigenvector of an $n\times n$ Cauchy matrix, where the associated con-eigenvalue $λ_{m}>0$ gives the approximation error in the $L^{\infty}$ norm. Unfortunately, standard algorithms do not accurately compute small con-eigenvalues (and the associated con-eigenvectors) and, in particular, yield few or no correct digits for con-eigenvalues smaller than the machine roundoff. We develop a fast and accurate algorithm for computing con-eigenvalues and con-eigenvectors of positive-definite Cauchy matrices, yielding even the tiniest con-eigenvalues with high relative accuracy. The algorithm computes the $m$th con-eigenvalue in $\mathcal{O}(m^{2}n)$ operations and, since the con-eigenvalues of positive-definite Cauchy matrices decay exponentially fast, we obtain (near) optimal rational approximations in $\mathcal{O}(n(\logδ^{-1})^{2})$ operations, where $δ$ is the approximation error in the $L^{\infty}$ norm. We derive error bounds demonstrating high relative accuracy of the computed con-eigenvalues and the high accuracy of the unit con-eigenvectors. We also provide examples of using the algorithm to compute (near) optimal rational approximations of functions with singularities and sharp transitions, where approximation errors close to machine precision are obtained. Finally, we present numerical tests on random (complex-valued) Cauchy matrices to show that the algorithm computes all the con-eigenvalues and con-eigenvectors with nearly full precision.

T. S. Haut, P. G. Maginot, V. Z. Tomov et al.

We propose a graph-based sweep algorithm for solving the steady state, mono-energetic discrete ordinates on meshes of high-order curved mesh elements. Our spatial discretization consists of arbitrarily high-order discontinuous Galerkin finite elements using upwinding at mesh element faces. To determine mesh element sweep ordering, we define a directed, weighted graph whose vertices correspond to mesh elements, and whose edges correspond to mesh element upwind dependencies. This graph is made acyclic by removing select edges in a way that approximately minimizes the sum of removed edge weights. Once the set of removed edges is determined, transport sweeps are performed by lagging the upwind dependency associated with the removed edges. The proposed algorithm is tested on several 2D and 3D meshes composed of high-order curved mesh elements.

T. S. Haut, C. Ahrens, A. Jonko et al.

We present a numerical method for solving the time-independent thermal radiative transfer (TRT) equation or the neutron transport (NT) equation when the opacity or cross-section varies rapidly in energy (frequency). The approach is based on a rigorous homogenization of the TRT/NT equation in the energy (frequency) variable. Discretization of the homogenized TRT/NT equation results in a multigroup-type system, and can therefore be solved by standard methods. We demonstrate the accuracy and efficiency of the approach on three model problems. First we consider the Elsasser band model with constant temperature and a small line spacing. Second, we consider a neutron transport application for fast neutrons incident on iron, where the characteristic resonance spacing necessitates about 16,000 energy discretization parameters if Planck-weighted cross sections are used. Third, we consider an atmospheric TRT problem with an opacity corresponding to water vapor. For all three problems, we demonstrate that we can achieve between 0.1 and 1 percent relative error in the solution, and with several orders of magnitude fewer parameters than a standard multigroup formulation with a comparable accuracy.