Poster Session Abstracts

Continuously Varying Material Properties and Tallies for Monte Carlo Calculations
Forrest B. Brown (LANL) & William R. Martin (Univ. Michigan)

Using a high-order Legendre polynomial representation for material density and tallies within each cell, Monte Carlo codes can model continuous variations in material properties and results. We have demonstrated the Monte Carlo techniques for sampling the free-flight distances and performing pathlength flux tallies for this continuous representation. Application to both fixed-source and eigenvalue problems illustrates the benefits of the continuous representation as compared to conventional stepwise approximations. With these new methods, Monte Carlo codes can now be developed which are continuous in energy, angle, space, material properties, and tallied results. [LA-UR-04-0732]

The Integrated TIGER Series version 5.0
Thomas W. Laub
Sandia National Laboratories

The Integrated Tiger Series (ITS) is a powerful and user-friendly software package permitting Monte Carlo solution of linear time-independent coupled electron/photon radiation transport problems, with or without the presence of macroscopic electric and magnetic fields of arbitrary spatial dependence. The package contains programs to perform 1D, 2D, and 3D simulations. Our goal has been to simultaneously maximize operational simplicity and physical accuracy. The ease with which the build system is applied combines with an input scheme based on order-independent descriptive keywords that makes maximum use of defaults and internal error checking to provide experimentalists and theorists alike with a method for the routine but rigorous solution of sophisticated radiation transport problems. Physical rigor is provided by employing accurate cross sections, sampling distributions, and physical models for describing the production and transport of the electron/photon cascade from 1.0 GeV down to 1.0 keV. The availability of source code permits the more sophisticated user to tailor the codes to specific applications and to extend the capabilities of the codes to more complex applications.

Although the last public release of ITS was version 3.0 in 1992, development efforts have continued. Version 5.0, the latest version of ITS, contains:

• Improvements to the ITS 3.0 continuous-energy codes.
• Multigroup codes with adjoint transport capabilities.
• Parallel implementations of all ITS codes.
• More automated subzoning options for combinatorial geometry.
• Additional source distributions, tallies, biasing options, and a torus CG primitive.
• Ability to output subzone energy and charge deposition in a finite-element-like format.
• Ability to use CAD geometry descriptions in the ACIS format.
• Subzoning capabilities for CAD geometries.
• A ray-tracing capability for fast scoping calculations.

Moreover, the general user friendliness of the software has been enhanced through increased internal error checking and improved code portability. In addition, software quality engineering procedures are being implemented.

Ongoing ITS research efforts include:

• Path-length apportioning for photons in subzone structures.
• Improvements in CAD geometry particle tracking efficiency.
• Allowing the use of facet-based geometries.
• Generalized Boltzmann-Fokker-Planck (GBFP) non-analog transport of electrons.
• Conversion of the Fortran77 portions of the codes to Fortran90.
• Extending the transport capability to sub 1-keV energies.

Although ITS V5.0 is not yet generally available, research-only and government-use-only licenses are being pursued. The code is available to NNSA laboratories through the Codes for the Complex.

Mathematical Simulation of the Radiative Transfer in Statistically Inhomogeneous Clouds
Evgueni Kassianov
Pacific Northwest National Laboratory

Cloud inhomogeneity is caused by both cloud field stochastic geometry and inhomogeneous internal structure. The problem of statistical radiative transfer in clouds is aimed at establishing the relationship between the statistical parameters of clouds and radiation. Although the necessity of such statistical treatment was acknowledged long ago and several approaches were suggested, the complex problem as a whole is still far from resolved. The statistical description of radiation transfer in clouds with both types of fluctuations (both geometrical and internal) is very challenging mainly due to complicated mathematical problems and the absence of reliable statistical data about the joint variability of the geometrical and optical properties (e.g., their mutual correlation). Previously, a promising approach to studying the radiative properties of single-layer broken clouds was developed. The term "broken clouds" means that the cloud field has stochastic geometry and deterministic optical parameters inside an individual cloud. This approach was based on analytical or numerical averaging of the stochastic transfer equation over an ensemble of realizations of a cloud field. By using the assumption that broken clouds can be represented as a binary mixture with Markovian statistics, mathematically rigorous methods for radiative calculations have been developed and the influence of the stochastic geometry on solar radiative transfer has been estimated. We generalize this approach for multilayer (statistically inhomogeneous) broken clouds by using the stochastic radiative transfer equation and a new Markovian model of broken clouds. Here, we introduce this generalized approach and show some validation results and applications.

Transport Needs in Astrophysics
Chris Fryer
Los Alamos National Laboratory

The effects of transport (from photon radiation to neutrinos) is becoming an increasingly important subject in astrophysics with applications in nearly all facets of astrophysics from star formation, stellar evolution, supernovae and gamma-ray bursts to models of active galactic nuclei and cosmology. The effect of the radiation on these physical problems ranges from the ionization of material where radiation pressure can be neglected to problems where diffusive or optically thin limits apply to full "transport regime" problems where the radiation must be fully coupled to the hydrodynamics. We will review these problems, focusing on the latter case, discussing what techniques astrophysicists currently use. With this clarification, we hope to open up the dialog between the astrophysics and transport communities to take advantage of transport techniques on these astrophysics problems.

Symbolic Implicit Monte Carlo radiation transport in the difference formulation: a piecewise constant discretization
Eugene D. Brooks III, Michael Scott McKinley, Frank Daffin, Abraham Szoke Lawrence Livermore National Laboratory
UCRL-ABS-203979

The equations of radiation transport for thermal photons are notoriously difficult to solve in thick media without resorting to asymptotic approximations such as the diffusion limit. One source of this difficulty is that in thick, absorbing media, thermal emission is almost completely balanced by strong absorption. A new formulation for thermal radiation transport, called the difference formulation, was recently introduced in order to remove the stiff balance between emission and absorption. In the new formulation, thermal emission is replaced by derivative terms that become small in thick media. It was proposed that the difficulties of solving the transport equation in thick media would be ameliorated by the difference formulation, while preserving full rigor and accuracy of the transport solution in the streaming limit. In this paper, the transport equation is solved by the Symbolic Implicit Monte Carlo method and comparisons are made between the standard formulation and the difference formulation. The method was easily adapted to the derivative source terms of the difference formulation, and a remarkable reduction in noise was obtained when the difference formulation is applied to problems involving thick media.

Non-LTE Radiation Transport in High Radiation Plasmas
H. A. Scott
Lawrence Livermore National Laboratory, Livermore, CA 94550, USA

A primary goal of numerical radiation transport is obtaining a self-consistent solution for both the radiation field and plasma properties.  Obtaining such a solution requires consideration of the coupling between the radiation and the plasma.  The different characteristics of this coupling for continuum and line radiation have resulted in two separate sub-disciplines of radiation transport with distinct emphases and computational techniques.  LTE radiation transfer focuses on energy transport and exchange through broad-band radiation, primarily affecting temperature and ionization balance.  Non-LTE line transfer focuses on narrow-band radiation and the response of individual level populations, primarily affecting spectral properties.  Many high-energy density applications, particularly those with high-Z materials, incorporate characteristics of both these regimes.  Applications with large radiation fields including strong line components require a non-LTE broad-band treatment of energy transport and exchange.

We discuss these issues and present a radiation transport treatment which combines features of both types of approaches by explicitly incorporating the dependence of material properties on both temperature and radiation fields.  The additional terms generated by the radiation dependence do not change the character of the system of equations and can easily be added to a numerical transport implementation.  The information needed to characterize the material response to radiation is closely related to that used by the Linear Response Matrix (LRM) approach to near-LTE simulation, and we investigate the use of the LRM for these calculations.  Numerical examples relevant to inertial confinement fusion and Z-pinch applications demonstrate that this method improves both the stability and convergence of the calculations.

This work was performed under the auspices of the U.S. Department of Energy, by the University of California, Lawrence Livermore National Laboratory under contract W-7405-ENG-48.

Finite Element Transport using Wachspress Rational Basis Functions on Quadrilaterals in Diffusive Regions
Gregory Davidson
Oregon State University

Wachspress rational functions are ratios of polynomials having certain properties which make them good candidates for basis functions in finite element methods. It had been theorized that Wachspress rational functions should provide a robust discretization of the neutral particle transport equation in thick diffusive regions, but this had not been confirmed numerically.

We have derived a discontinuous finite element discretization on quadrilaterals, and we have performed an asymptotic analysis of the transport discretization in the interior of thick diffusive regions. Using this, we derived an asymptotic-P1 diffusion synthetic acceleration preconditioner. These derivations were implemented to provide numerical results demonstrating the behavior of Wachspress rational functions in the thick diffusive limit.

Analytical and computational results indicate that Wachspress rational functions provide a robust spatial discretization in the thick diffusive limit.

The modern level of the knowledge of spatial structure of clouds does not allow building optical models of the corresponding scattering media without introducing some approximation and this undoubtedly leads to the accuracy of the prediction being limited by such lack of information. That is one of the basic reasons why different approximations are successfully used to study the radiation effects of 3 dimensional structure of real clouds. However, even such a simplification assumes that a differential equation with nonlinear coefficients has to be solved. Further simplification can be achieved if it is possible to assume the horizontal inhomogeneity being weak with respect to the corresponding average characteristics of the medium. This assumption leads to a much simpler solution and recently a perturbation approach based on an introduction of a small parameter into the radiative transfer equation has been developed to study the effect of the extinction coefficient variation.

A different type of perturbation approach may be formulated on the basis of the joint consideration of both the direct and corresponding adjoint problem. The most obvious advantage of such an approach is that it allows one to study the effect of different variations of the medium characteristics (including the phase function and single scattering albedo) on the radiance characteristics. The general equations relevant to such a perturbation calculation are formulated. This approach was applied to study the effect of the horizontal inhomogeneity on the radiation propagation through the realistic cloud and the results obtained are demonstrated and discussed.

We will discuss the Implicit Monte Carlo Diffusion method (JCP 172, 534). Implicit Monte Carlo (IMC) is often employed to numerically simulate radiative transfer. In problems with regions that are characterized by a small mean free path, IMC can take a prohibitive amount of time. because many steps must be simulated to advance the particle through the time step. Problems containing regions with a small mean free path can frequently be accurately simulated much more quickly by employing the diffusion equation as an approximation. However, the diffusion approximation is not accurate in regions where  the mean free path is large.

We present a method for accelerating time-dependent Monte Carlo radiative transfer calculations by using a descretization fo the diffusion equation to calculate probabilities that are used to advance particles in regions with small mean free paths. The method is demonstrated on problems with one and two dimensional orthogonal grids. It results in decreases in run time or more than an order of magnitude on these problems, while producing answers with accuracy comparable to pure IMC simulations.

An Evaluation of the Difference Formulation for Photon Transport in a Two Level Atomic System in Slab Geometry
F. C. Daffin, M. S. McKinley, E. D. Brooks III and A. Szöke
Lawrence Livermore National Laboratory

The difference formulation of photon transport, introduced for thermal radiation by some of the authors, has the unique property that the transport equation is written in terms that become small for optically thick systems. In this presentation we extend the difference formulation for radiation transport to the case of a single atomic line in a strongly absorbing and emitting medium of twolevel atoms. We examine the accuracy, performance and stability of the difference formulation within the framework of the Symbolic Implicit Monte Carlo method. We find that the difference formulation offers a significant computational advantage over the standard formulation for a thick system. The correct treatment of the line profile, however, requires that the difference formulation in the core of the line, where the photon field is strong, be mixed with the standard formulation in the wings, where the photon field is weak. This may limit the computational advantage of the method. We bypass this problem by using the gray approximation, and develop three Monte Carlo solution methods based on different degrees of implicitness for the treatment of the source terms. We find only conditional stability unless the source terms are treated fully implicitly.

We are exploring the use of the SIMC method for the simulation of an ICF capsule implosion. First, we have had to introduce variance reduction techniques to perform a satisfying calculation. Among them, the photons direction biasing is important to ensure uniform radiation at the ablation front. We present numerical examples compared with another Monte Carlo method (Fleck's method).

Having implemented deterministic transport within RAGE (an AMR Eulerian code with cell-centered data structures), once (at some cost of storage) with a ``staggered-mesh'' $P_1$ scheme, and once with an explicit cell-centered $P_1$ method, we now describe a new cell-centered, Riemann-solver based Implicit (gray) $P_1$ package. This package couples both nearest neighbor energy densities and flux densities in approximate Riemann solver expressions for face-centered values of flux and pressure; these in turn are summed around the zonal control volume to update energy and flux via coupled matrix equations. Because this would require us to invert a $(2.mesh)^2$ matrix in 1d, and a $(4 .mesh)^2$ matrix in 3d, we propose to implement an ADI method instead. We propose to wrap this inside a Jacobi-free Newton method in the future in order to reduce residuals and allow larger than explicit timesteps without loss of symmetry. This poster will present early 1d results.

Nonlinear Yvon-Mertens Method for the Transport Equation
Dmitriy Anistratov & Loren Roberts (NC State University)

The Yvon-Mertens (YM) method is a nonlinear projective-iteration method. It is defined by a system of equations consisting of two parts: the transport (high-order) and moments (low-order) equations. It possesses a combination of certain valuable features of the quasidiffusion and flux methods. There exist different ways to formulate the low-order problem of the YM method in multidimensional geometries, depending on projection operators and definitions of linear- fractional functionals. We will present the results of our studies of its convergence properies and various discretizations of the low-order equations of the YM method.

Positive Discrete Ordinates Techniques
Kirk Mathews, Ph.D.
Professor of Nuclear Engineering Air Force Institute of Technology

Many practitioners avoid using discrete ordinates methods because they can produce unphysical, negative values for particle fluxes. This can arise through the combination of Gauss-like angular quadrature sets with truncated spherical-harmonic expansions of the group-to-group cross sections. This cause is particularly difficult to avoid (with current production software) when high energy resolution is demanded. Negative fluxes can also result from the use of higher-than-first-order linear spatial differencing (spatial quadrature) methods. To eliminate these artifacts, my students and I have developed
   (a) nonnegative spatial quadratures,
   (b) nonnegative representations of anisotropic, group-to-group scattering cross sections, and
   (c) angular quadratures that use these cross sections.
      The combination of these three techniques yields a discrete ordinate method that cannot produce negative fluxes, given nonnegative sources and incident fluxes at the boundaries.

SERRANO: A Linear Discontinuous Finite Element Sn Transport Code
Kent Budge
Los Alamos National Laboratory

SERRANO is a library of Sn transport software components that provides a parallel sweep capability for continuous and discontinuous finite element discretizations of the uncollided Boltzmann equation, as well as physics-based source iteration preconditioners (such as DSA or TSA) for a Krylov solution of the collided Boltzmann equation. We compare the two discretizations for problems that include both optically thin and optically thick regions to investigate numerical behavior at material discontinuities, which may introduce unresolved internal boundary layers. We also present performance results for various sweeper settings and choices of preconditioner and Krylov method. In particular, we report on the effectiveness of using a continuous form of the DSA preconditioner for the acceleration of the discontinuous discretization.

Backward Error Compensation for Transportation Equations with Application to the Level Set Method
Yingjie Liu, Georgia Inst of Tech
Collaborator: Todd Dupont, The Univ of Chicago

Level set method uses a level set function, usually an approximate signed distance function, Phi, to represent the interface as the zero set of Phi. When Phi is advanced to the next time level by a transportation equation, its new zero level set will represent the new interface position. We update the level set function Phi forward in time and then backward to get another copy of the level set function, say Phi_1. Phi_1 and Phi should have been equal if there were no numerical error. Therefore Phi-Phi_1 provides us the information of error and this information can be used to compensate Phi before updating Phi forward again in time. One nice property is that it has the convenience of possibly improving the temporal and spatial order of an odd order scheme simultaneously. We found that when applying this idea to semi-Lagrangian schemes, e.g., CIR scheme(which has no CFL restriction, a nice feature for local refinement), the property is still valid (MacCormack scheme has the same benefit but may not be easily applied here). Numerical results for interface movements with level set equation computed by the new methods are promising. Also we have found some interesting theoretical results for applying this idea to a general linear scheme.

Central Schemes on Overlapping Cells
Yingjie Liu, School of Math, Georgia Inst. of Tech.

Nessyahu and Tadmor's central scheme (J. Comput. Phys, 87(1990)) has the benefit of not using Riemann solvers or characteristic decomposition for solving hyperbolic conservation laws and related convection diffusion equations. But the staggered averaging causes large dissipation when the time step size is small comparing to the mesh size. The recent work of Kurganov and Tadmor (J. Comput. Phys, 160(2000)) overcomes the problem by use of a variable control volume and obtains a semi-discrete non-staggered central scheme. Motivated by this work, we introduce overlapping cell averages of the solution at the same discrete time level, and develop a simple alternative technique to control the $O(1/\Delta t)$ dependence of the dissipation. Semi-discrete form of the central scheme can also be obtained to which
the TVD Runge-Kutta time discretization of Osher and Shu (J. Comput. Phys, 77(1988)) can be applied. This technique is essentially independent of the reconstruction and the shape of the mesh, thus could also be useful for unstructured mesh. The overlapping cell representation of the solution also opens new possibilities for reconstructions. Generally more compact reconstruction can be achieved. We demonstrate through numerical examples that combining two classes of the overlapping cells in the reconstruction can achieve higher resolution.

A Piecewise Linear Finite Element Discretization of the Diffusion Equation
Teresa Bailey

A piecewise linear (PWL) finite element spatial discretization has been developed for the multi-dimensional diffusion equation.  It uses piecewise linear weight and basis functions in the finite element approximation.  This method can solve the diffusion equation on arbitrary polygonal (2D) or polyhedral (3D) grids, which allows for the solution of problems with complex shapes.  My presentation will describe the implementation of the PWL method into an existing diffusion code that is part of the KULL project at Lawrence Livermore National Laboratory.  The new method, which generates a symmetric positive definite coefficient matrix, will be compared against the existing method, which is a vertex-centered discretization with an asymmetric coefficient matrix. 

"Three-Dimensional Radiative Transfer in Cloudy Atmospheres"
- A new Volume for Springer-Verlag
Alexander Marshak
NASA/GSFC, Code 913, Greenbelt, MD 20771
and
Anthony B. Davis,
Los Alamos National Laboratory, MS B-244, Los Alamos, NM 87545

Three-dimensional cloud radiative community has matured enough to prepare a volume on 3D radiative transfer in cloudy atmosphere that will be published by Springer-Verlag this year. Many leading 3D radiative transfer scientists are amongst the coauthors of the book. The book starts with the basic 3D radiative transfer problem, describes its computational solutions and phenomenological models, discusses the effects of cloud inhomogeneity for remote sensing, addresses climate problems in realistic atmospheres, and studies cloud-vegetation interactions. In our poster presentation, we will discuss the outline of the book, give examples from different chapters, mostly focusing on broken clouds and cloud-vegetation interactions. A hard-copy of the book in galley-proof format will be on hand for perusal.

A Moment-Preserving Nonanalog Method for Charged Particle Transport
Anil K. Prinja
The University of New Mexico

Extremely short collision mean free paths and near-singular elastic and inelastic differential cross sections (DCS) make analog Monte Carlo simulation an impractical tool for charged particle transport. The widely used alternative, the condensed history method, while efficient, also suffers from several limitations  arising from the use of precomputed smooth distributions for sampling. There is much interest in developing computationally efficient algorithms that implement the correct transport mechanics. Here we present a nonanalog transport-based method that incorporates the correct transport mechanics and is computationally efficient for implementation in single event Monte Carlo codes. Our method systematically preserves important physics and is mathematically rigorous. It builds on higher order Fokker-Planck and Boltzmann Fokker-Planck representations of the scattering and energy-loss process, and we accordingly refer to it as a Generalized Boltzmann Fokker-Planck (GBFP) approach.

We postulate the existence of nonanalog single collision scattering and energy-loss distributions (differential cross sections) and impose the constraint that the first few momentum transfer and energy-loss moments be identical to corresponding analog values. This is effected through a decomposition or hybridizing scheme wherein the singular forward peaked, small energy-transfer collisions are isolated and de-singularized using different moment-preserving strategies, while the large angle, large energy-transfer collisions are described by the exact (analog) DCS or approximated to a high degree of accuracy. The inclusion of the latter component allows the higher angle and energy-loss moments to be accurately captured. This procedure yields a regularized transport model characterized by longer mean free paths and smoother scattering and energy transfer kernels than analog. In practice, acceptable accuracy is achieved with two rigorously preserved moments, but accuracy can be systematically increased to analog level by preserving successively higher moments with almost no change to the algorithm. Details of specific moment-preserving strategies will be described and results presented for dose in heterogeneous media due to a pencil beam and a line source of monoenergetic electrons. Error and runtimes of our nonanalog formulations will be contrasted against condensed history implementations.

Fully Implicit Solution of Large-Scale Non-Equilibrium Radiation Diffusion with High Order Time Integration
Carol S. Woodward
Center for Applied Scientific Computing
Lawrence Livermore National Laboratory
P.O. Box 808, L-561
Livermore, CA  94551

Modeling radiation diffusion processes has traditionally been accomplished through inaccurate and nonscalable solution methods based on decoupling linearized equations.  We present a solution method for fully implicit radiation diffusion problems with material energy coupling and highly nonlinear fusion source terms discretized on meshes having millions of spatial zones. This solution method makes use of high order in time integration techniques, inexact Newton--Krylov nonlinear solvers, and multigrid preconditioners.  Our method fully converges the nonlinear solution.  We explore the advantages and disadvantages of high order time integration methods for the fully implicit formulation on both one- and three-dimensional problems with tabulated opacities and compare results to a typical semi-implicit method.

This work was done in collaboration with Peter N. Brown and Dan E. Shumaker, Lawrence Livermore National Laboratory.

This work was performed under the auspices of the U.S. Department of Energy by University of California Lawrence Livermore National Laboratory under contract No. W-7405-Eng-48.

Implicit Riemann Solvers for the Pn equations
Ryan McClarren, Univ. of Michigan / Sandia National Labs
Thomas Brunner, Sandia National Labs
James Paul Holloway, Univ. of Michigan
Thomas Mehlhorn, Sandia National Labs

The spherical harmonics (Pn) approximation to the transport equation for time dependent problems
has previously been treated using Riemann solvers and explicit time integration. Here we present an
implicit time integration method for the Pn equations using Riemann solvers. Both first-order and highresolution
spatial discretization schemes are detailed. One facet of the high-resolution scheme is that
a system of nonlinear equations must be solved at each time step. This nonlinearity is the result of
slope reconstruction techniques necessary to avoid the introduction of artifical extrema in the numerical
solution. Results are presented that show auspicious agreement with analytical solutions using time steps
well beyond the CFL limit.

AGILE-BOLTZTRAN: A neutrino radiation hydrodynamics code for core-collapse supernovae simulations
Bronson Messer
University of Chicago

AGILE-BOLTZTRAN is a spherically symmetric neutrino radiation hydrodynamics code developed specifically for core-collapse supernova simulations. The heart of the code consists of two primary parts: a fully implicit Boltzmann solver for neutrino transport via a discrete ordinates method and a Lagrangian hydrodynamics module with adaptive mesh redistribution. I will describe the general structure of the code, with particular attention paid to some unique problems presented by calculating fermion transport in an astrophysical setting. Comparisons to other neutrino transport algorithms currently employed in core-collapse SNe simulations, with their attendant advantages and disadvantages, will also be discussed.

Parallel Implementation of Block-Structured Adaptive Mesh Refinement for the Discrete Ordinates Method
Louis H. Howell
Lawrence Livermore National Laboratory

The BoxLib framework for structured adaptive mesh refinement (AMR) decomposes each refinement level of a mesh into rectangular grids, and distributes one or more of these grids to each active processor. In time-dependent calculations coarse levels are advanced at larger timesteps than fine levels. Some of the issues that must be addressed for efficient implementation of the discrete ordinates method in this framework are: parallel sweeping strategies that take into account the dependencies between grids, conservative coupling between refinement levels, and convergence acceleration for problems with short mean free paths. I will present these issues in the context of radiative transport coupled in a fully-implicit fashion to fluid energy. Other topics for discussion include the additional dependencies between ordinate directions introduced by 2D axisymmetric discretization, the closed dependency loops that become possible with some 3D AMR grid layouts, and parallel scaling results for both 2D and 3D adaptive meshes.