The power iteration method is the standard Monte Carlo approach for obtaining the eigenfunctions of a nuclear system, but the power method sometimes converges very slowly. Most discussions give a mathematical reason for the slow convergence of the Monte Carlo power method using the same concepts and terminology as when the power method is applied to a deterministic problem. This note first looks at why the convergence is slow from an intuitive Monte Carlo neutron perspective. Second, this note proposes building an eigenfunction intuitively in a cumulative (and noniterative) neutron by neutron manner that tends to better direct neutrons to where the neutrons need to be. Third, a very similar method for building the second eigenfunction is speculatively proposed.

Improvements incorporated into MCNP variance reduction methodology and code releases since 2000 are discussed. Some of the improvements are modifications or generalizations of older techniques, and some are entirely new. In particular, pulse-height-tally variance reduction is now possible in MCNP, and the dxtran technique has been generalized to allow an arbitrary nesting of dxtran spheres. A new precollision, next-event estimator is discussed along with flux-at-a-point image tallies. Additionally, the event log analyzer is a tool designed to help the user understand what causes the variance in the user’s particular MCNP calculation.

Monte Carlo particle transport is often introduced primarily as a method to solve linear integral equations such as the Boltzmann transport equation. This paper discusses some common misconceptions about Monte Carlo methods that are often associated with an equation-based focus. Many of the misconceptions apply directly to standard Monte Carlo codes such as MCNP and some are worth noting so that one does not unnecessarily restrict future methods.

Particlefluxes on surfaces are difficult to calculate with Monte Carlo codes because the score requires a division by the surface-crossing angle cosine, and grazing angles lead to inaccuracies. We revisit the standard practice of dividing by half of a cosine "cutoff" for particles whose surface-crossing cosines are below the cutoff We concentrate on the flux crossing an external boundary, deriving the standard approach in a manner that explicitly points out three assumptions: (a) that the external boundary surface flux is isotropic or mostly isotropic, (b) that the cosine cutoff is small, and (c) that the minimum possible surface-crossing cosine is 0. We find that the requirement for accuracy of the standard surface flux estimate is more restrictive for external boundaries (a very isotropic surface flux) than for internal surfaces (an isotropic or linearly anisotropic surface flux). Numerical demonstrations involve analytic and semianalytic solutions for monoenergetic point sources irradiating surfaces with no scattering. We conclude with a discussion of potentially more robust approaches.

The usual Monte Carlo approach to the thermal radiative transfer problem is to view Monte Carlo as a solution technique for the nonlinear thermal radiative transfer equations. The equations contain time derivatives which are approximated by introducing small time steps. An alternative approach avoids time steps by using Monte Carlo to directly sample the time at which the next event occurs. That is, the time is advanced on a natural event-by-event basis rather than by introducing an artificial time step.

This technical note shows that it is possible and effective to use Monte Carlo variance-reduction methods for the probability of initiation problem. The benefits are threefold. First, the proper use of variance reduction obviates using an arbitrary definition of a "divergent chain." Second, because chains of all lengths are allowed, there is no bias introduced by ignoring some long chains because they meet the divergent chain definition. Third, variance-reduction methods might drastically increase the efficiency of some of these calculations.

Recently, we proposed a modified power iteration method that simultaneously determines the dominant and subdominant eigenvalues and eigenfunctions of a matrix or a continuous operator. One advantage of this method is the convergence rate to the dominant eigenfunction being vertical bar k(3)vertical bar/k(1) instead of vertical bar k(2)vertical bar/k(1), a potentially significant acceleration. One challenge for a Monte Carlo implementation of this method is that the second eigenfunction is represented by particles of both positive and negative weights that somehow must sum (cancel) to estimate the second eigenfunction faithfully. Our previous Monte Carlo work has demonstrated the improved convergence rate by using a point flux estimator method and a binning method to effect this cancellation. This paper presents an exact method that cancels over a region instead of at points or in small bins and has the potential of being significantly more efficient than the other two.

Monte Carlo simulations may involve skewed, heavy-tailed distributions. When variances of those distributions exist, statistically valid confidence intervals can be obtained using the central limit theorem, providing that the simulation is run ''long enough.'' If variances do not exist, however, valid confidence intervals are difficult or impossible to obtain. The main result in this paper establishes that upon replacing ordinary Monte Carlo sampling of such heavy-tailed distributions with ex post facto sampling, estimates having finite moments of all orders are ensured for the most common class of infinite variance distributions. We conjecture that this phenomenon applies to all distributions (having finite means) when the ex post facto process is iterated.

A method to provide an unbiased Monte Carlo estimate of the reciprocal of an integral is described. In Monte Carlo transport calculations, one often uses a single sample as an estimate of an integral. This paper shows that a similar situation exists with respect to a single sample for an unbiased estimate of the reciprocal of an integral. If an appropriate approximation to the integrand is known, then obtaining a single unbiased estimate of the reciprocal of an integral will not be much more time consuming than obtaining a single unbiased estimate of the integral itself.

A method for simultaneously obtaining the two largest eigenvalues and their associated eigenfunctions is demonstrated mathematically and empirically. The method uses estimates of the eigenvalue in two different regions rather than the single estimate traditionally used. The method can be generalized to obtain the several largest eigenfunctions, if those are desired as well. Additionally, it is shown that using multiple estimates of the eigenvalues accelerates the convergence of the eigenfunctions.

a sequence of fission generations (cycles), during each of which M particles per cycle are processed. After a series of "inactive" cycles during which the fission source "converges," a series of "active" cycles are performed. For each active cycle, the eigenvalue and eigenfunction are estimated; after N >> 1 active cycles are performed, the results are averaged to obtain estimates of the eigenvalue and eigenfunction and their standard deviations. This method has several disadvantages: (i) the estimate of k depends on the number M of particles per cycle, (iii) for optically thick systems, the eigenfunction estimate may not converge due to undersampling of the fission source, and (iii) since the fission source in any cycle depends on the estimated fission source from the previous cycle (the fission sources in different cycles are correlated), the estimated variance in k is smaller than the real variance. For an acceptably large number M of particles per cycle, the estimate of k is nearly independent of M; this essentially takes care of item (i). Item (ii) can be addressed by taking M sufficiently large, but for optically thick systems a sufficiently large M can easily be unrealistic. Item (iii) cannot be accounted for by taking M or N sufficiently large; it is an inherent deficiency due to the correlations between fission source estimates. In the new FMC method, the eigenvalue problem (expressed in terms of the Boltzmann equation) is integrated over the energy and direction variables. Then these equations are multiplied by J special "tent" functions in space and integrated over the spatial variable. This yields J equations that are exactly satisfied by the eigenvalue k and J space-angle-energy moments of the eigenfunction. Multiplying and dividing by suitable integrals of the eigenfunction, one obtains J algebraic equations for k and the space-angle-energy moments of the eigenfunction, which contain nonlinear functionals that depend weakly on the eigenfunction. In the FMC method, information from the standard Monte Carlo solution for each active cycle is used to estimate the functionals, and at the end of each cycle the J equations for k and the space-angle-energy moments of the eigenfunction are solved. Finally, these results are averaged over N active cycles to obtain estimated means and standard deviations for k and the space-angle-energy moments of the eigenfunction. Our limited testing shows that for large single fissile systems such as a commercial reactor core, (i) the FMC estimate of the eigenvalue is at least one order of magnitude more accurate than estimates obtained from the standard Monte Carlo approach, (ii) the FMC estimate of the eigenfunction converges and is several orders of magnitude more accurate than the standard estimate, and (iii) the FMC estimate of the standard deviation in k is at least one order of magnitude closer to the correct standard deviation than the standard estimate. These advances occur because: (i) the Monte Carlo estimates of the nonlinear functionals are much more accurate than the direct Monte Carlo estimates of the eigenfunction, (ii) the system of discrete equations that determines the FMC estimates of k is robust, and (iii) the functionals are only very weakly correlated between different fission generations. The FMC method was developed only late in the project and has to date received limited testing. Current work, which is taking place after the conclusion of this project, involves further development and testing for more complex and realistic problems. We expect that the FMC method will become a practical tool for Monte Carlo criticality simulations.

AVATAR (Automatic Variance And Time of Analysis Reduction), accessed through the graphical user interface application Justine, is a superset of MCNP that automatically invokes THREEDANT for a three-dimensional deterministic adjoint calculation on a mesh independent of the Monte Carlo geometry, calculates weight windows, and runs MCNP. Computational e ciency increases by a factor of 2 to 5 for a three-detector oil well logging tool model. User e ciency increases dramatically, since AVATAR eliminates the need for deep intuition and hours of tedious data manipulation.

This paper applies a recently developed variance reduction technique to the first-principles calculations of photon detector responses. This technique makes possible the direct comparison of pulse height calculations with measurements without the need for unfolding techniques. Comparisons are made between several experiments and the calculations to demonstrate the utility of the supertrack Monte Carlo technique for reproducing and interpreting experimental count rate spectra.

It is shown that a zero-variance solution exists for any linear Monte Carlo problem and that this solution can be obtained by sampling the random numbers proportional to the expected score subsequently produced by using these random numbers in the random walk.

A continuous random walk procedure for solving some elliptic partial differential equations at at a single point is generalized to estimate the solution everywhere. The Monte Carlo method described here is exact (except at the boundary) in the sense that the only error is the statistical sampling error that tends to zero as the sample size increases. A method to estimate the error introduced at the boundary is provided so that the boundary error can always be made less than the statistical error.