Direct Monte Carlo Calculations Research Articles

Monte Carlo calculation and similarity evaluation for verification of influential uncertainty parameters to LBLOCA derived by statistical methods

Several studies for identification of influential uncertainty parameters have been conducted, and the most of them used statistical methods utilizing small sample data. On the other hand, the consideration of influential uncertainty variables affects the distribution and the value of the figures of merit (FOMs). Therefore, we can verify influential uncertainty parameter identification results of statistical methods by comparing the distributions and the values of the FOMs. Kang (Kang, 2020) identified important uncertainty variables influential to APR-1400 large break loss of coolant accident (LBLOCA) by utilizing comprehensive correlation and multiple linear regression analysis with hypothesis testing. In this subsequent research, direct Monte Carlo calculations were performed by applying influential uncertainty parameters derived from previous work (Kang, 2020). Then, the similarity of some descriptive statistics, 95-percentile peak cladding temperature (PCT95) with 95% confidence interval, PCT distribution was evaluated. In addition, the PCT95 with 95% tolerance limit by Wilks’ method was also compared to verify previous results obtained by statistical methods.

A HYBRID MONTE-CARLO-DETERMINISTIC METHOD FOR AP1000 EX-CORE DETECTOR RESPONSE SIMULATION

The ex-core detector-response calculation is a typical deep-penetration problem, which is challenging for the Monte Carlo method. The response of the ex-core detector is an important parameter for the safe operation of the nuclear power plants. Meanwhile, evaluation of the ex-core detector response during each step of fuel-loading is used to guide the fuel-loading sequence. The response can also be used to reconstruct core-power distribution for online monitoring of long-term power. The detector used for the ex-core response is the source-range detector which is sensitive to thermal neutrons. For a Monte Carlo shielding calculation of the above detector response, the thermal flux under 0.625eV is needed, which is too small to be tallied by traditional Monte Carlo simulations. In practice, the tally results are close to zero in the detector region under direct Monte Carlo calculation. Even if the number of particles is increased to a significant amount, the statistical variance is still very large. The high variance along with a significant calculation time leads to a small Figure Of Merit (FOM). In order to solve this problem and to improve the tally efficiency of the ex-core detector response, a hybrid Monte-Carlo-deterministic method is employed in this study, and an in-house hybrid Monte-Carlo-deterministic particle transport code, NECP-MCX, is developed in this paper. The method takes the space-energy-dependent adjoint fluxes to generate importance parameters for the mesh-based weight window in the Monte Carlo calculation. Simultaneously, the mesh-based source biasing is performed with the consistent importance parameters to make the starting weight of neutrons matching with the survival weight of the weight windows. As the mesh used in the hybrid Monte-Carlo-deterministic method is superimposed, the mesh of the weight window will not be affected by the complex geometry model. The adjoint flux is obtained by the efficient SN method with the multi-group cross-section data. The whole toolset is convenient to use with single set of the modelling data for both Monte Carlo and deterministic simulations. Compared with the direct Monte Carlo simulation, the hybrid Monte-Carlo-deterministic method has a higher efficiency for a typical deep-penetration problem such as the AP1000 ex-core detector-response simulation.

Determination of backscatter factors based on the quality index for diagnostic kilovoltage x-ray beams.

This study aims to investigate the relationship between backscatter factors and Al-half-value-layers (Al-HVL) by making the quality index (QI) a parameter for diagnostic kilovoltage x-ray beams. Backscatter factors, Bw, for x-ray fluence spectra were calculated from the weighted average of Bw for monoenergetic photons of between 8 and 140keV with field sizes of 10cm×10cm to 40cm×40cm. The value of Bw for monoenergetic photons was calculated from the ratio of the water kerma at the surface of a water phantom and that at the same point free-in-air using the EGSnrc/cavity code. The weighted averaged backscatter factors were validated by comparing them with those of direct Monte Carlo calculations for the x-ray fluence spectra. The Bw for the x-ray fluence spectra were classified by a QI of 0.35, 0.4, 0.5, 0.6, and 0.7 specified by the ratio of the effective energy and maximum energy. The relationship between Bw and Al-HVL was evaluated for the given QI values. The x-ray fluence spectra were generated for tube voltages of 40-140 kVp with Al-HVLs of 0.5-13.2mm using the SpekCalc program. The weighted averaged backscatter factors for x-ray fluence spectra agreed within 0.7% with those of the direct Monte Carlo calculations. The backscatter factors were represented by the fitting curves of R2>0.99 with Al-HVL for the given QI values. It is possible to obtain Bw more accurately by using QI specified by the measured Al-HVL.

Stochastic polynomial chaos expansion method for random Darcy equation

Abstract A probabilistic collocation based polynomial chaos expansion method is developed for simulation of particle transport in porous medium. The hydraulic conductivity is assumed to be a random field of a given statistical structure. The flow is modeled in a two-dimensional domain with mixed Dirichlet–Neumann boundary conditions. The relevant Karhunen–Loève expansion is constructed by a special randomized singular value decomposition (SVD) of the correlation matrix which makes possible to treat problems of high dimension. The simulation results are compared against a direct Monte Carlo calculation of different Eulerian and Lagrangian statistical characteristics of the solutions. As a byproduct, we suggest an approach to solve an inverse problem of recovering the variance of the log-conductivity.

Monte Carlo algorithm for simulating fermions on Lefschetz thimbles

A possible solution of the notorious sign problem preventing direct Monte Carlo calculations for systems with non-zero chemical potential is to deform the integration region in the complex plane to a Lefschetz thimble. We investigate this approach for a simple fermionic model. We introduce an easy to implement Monte Carlo algorithm to sample the dominant thimble. Our algorithm relies only on the integration of the gradient flow in the numerically stable direction, which gives it a distinct advantage over the other proposed algorithms. We demonstrate the stability and efficiency of the algorithm by applying it to an exactly solvable fermionic model and compare our results with the analytical ones. We report a very good agreement for a certain region in the parameter space where the dominant contribution comes from a single thimble, including a region where standard methods suffer from a severe sign problem. However, we find that there are also regions in the parameter space where the contribution from multiple thimbles is important, even in the continuum limit.

Modeling of inelastic collisions in a multifluid plasma: Excitation and deexcitation

We describe here a model for inelastic collisions for electronic excitation and deexcitation processes in a general, multifluid plasma. The model is derived from kinetic theory, and applicable to any mixture and mass ratio. The principle of detailed balance is strictly enforced, and the model is consistent with all asymptotic limits. The results are verified with direct Monte Carlo calculations, and various numerical tests are conducted for the case of an electron-hydrogen two-fluid system, using a generic, semi-classical model of collision cross sections. We find that in some cases, the contribution of inelastic collisions to the momentum and thermal resistance coefficients is not negligible, in contrast to the assumptions of current multifluid models. This fundamental model is also applied to ionization and recombination processes, the studies on which are currently underway.

Nuclear state densities of odd-mass heavy nuclei in the shell model Monte Carlo approach

While the shell model Monte Carlo approach has been successful in the microscopic calculation of nuclear state densities, it has been difficult to calculate accurately state densities of odd-even heavy nuclei. This is because the projection on an odd number of particles in the shell model Monte Carlo method leads to a sign problem at low temperatures, making it impractical to extract the ground-state energy in direct Monte Carlo calculations. We show that the ground-state energy can be extracted to a good precision by using level counting data at low excitation energies and the neutron resonance data at the neutron threshold energy. This allows us to extend recent applications of the shell-model Monte Carlo method in even-even rare-earth nuclei to the odd-even isotopic chains of $^{149-155}$Sm and $^{143-149}$Nd. We calculate the state densities of the odd-even samarium and neodymium isotopes and find close agreement with the state densities extracted from experimental data.

A fission collision separation method for efficient incident flux response expansion coefficient generation

A new fission collision (generation) separation method has been developed for generating incident flux response expansion coefficients and implemented into the hybrid stochastic/deterministic coarse mesh transport (COMET) method for whole core criticality (eigenvalue) calculations. During the stochastic pre-computation of the response library, all coefficients associated with each fission neutron generation are tallied separately. This enables COMET to directly calculate the eigenvalue (k) dependent response coefficients on-the-fly as a superposition of contributions from the nth fission neutron generation scaled by the factor 1kn. The new method is as accurate as or more accurate and significantly more efficient than the existing response methods, regardless of whether the fission source is treated explicitly or implicitly as in COMET. The method also removes the eigenvalue range restriction on the response library.The 3D C5G7 benchmark problem was used to test the accuracy and efficiency of the new response coefficient generation method at both lattice (local) and core (global) levels. For lattice calculations, the response coefficients computed by the new method were compared with those computed by the direct Monte Carlo method. It was found that the truncation errors of the fission neutron separation method at expansion order 3 are much lower than the stochastic uncertainties of the response coefficients and consequently can be ignored. The average truncation errors of the new method ranged from 0% to 0.004% and 0% to 0.08% for the surface-to-surface and fission density response coefficients, respectively. For core level calculations, the COMET whole-core solutions based on the response coefficient libraries generated by the new method and the interpolation method were compared with direct (reference) Monte Carlo calculations using the MCNP code. It was found that the error in the eigenvalues predicted by COMET based on the new method is within 11–26pcm for three core configurations with and without control rods. These eigenvalue results are about 15pcm more accurate than those based on the interpolation method, although still within the overall uncertainty. The corresponding assembly and pin fission densities are almost identical for the two methods. These comparisons indicate that the new method is as accurate as or slightly more accurate than the interpolation method, depending on the comparison parameter. The new method was found to be about 5 times faster than the interpolation method for the generation of the response function library.

Three dosimetry models of lipoma arborescens treated by 90Y synovectomy.

Lipoma arborescens (LA) is a benign intra-articular lipomatous proliferation of the synovial membrane. This extremely rare condition has previously been treated by intra-articular (90)Y radiosynoviorthesis but dosimetry literature on this form of radionuclide therapy is nonexistent. The authors detail methodology for successful treatment of LA and provide for the first time estimates of radiation dosimetry. The authors also analyze the biodistribution of the radiopharmaceutical over the course of the patient's treatment through sequential imaging. A patient with bilateral LA underwent intracavity injection of (90)Y citrate colloid to the right and left knee joint spaces (181 and 198 MBq, respectively). SPECT/CT datasets were acquired over 9 days to quantify the biodistribution and kinetics of the radiopharmaceutical. Radiation dosimetry was performed using the MIRD schema (through OLINDA software), a custom voxel-based method, and a direct Monte Carlo calculation (OEDIPE). Follow-up MRI showed marked reduction in LA size in both knees. Mean absorbed doses to the LA were 21.2 ± 0.8 and 42.9 ± 2.3 Gy using OLINDA, 8.1 ± 0.3 and 16.7 ± 0.5 Gy using voxel based methodology, and 8.2 ± 0.3 and 15.7 ± 0.5 Gy for OEDIPE in the right and left LA, respectively. Distribution of the radiopharmaceutical within the joint space alters over the imaging period, with less than 1% of the remaining activity having moved posteriorly in the knee cavity. No uptake was detected outside of the joint space after assessment with whole-body scintigraphy. An activity of approximately 185 MBq successfully relieved clinical symptoms of LA. There was good correlation between direct Monte Carlo and voxel based techniques, but OLINDA was shown to overestimate the absorbed dose to the tumor. Accurate dosimetry may help select an activity more tailored to the specific size and location of the LA.

Direct MC conversion of absorbed dose to graphite to absorbed dose to water for 6060Co radiation

The ARPANSA calibration service for (60)Co gamma rays is based on a primary standard graphite calorimeter that measures absorbed dose to graphite. Measurements with the calorimeter are converted to the absorbed dose to water using the calculation of the ratio of the absorbed dose in the calorimeter to the absorbed dose in a water phantom. ARPANSA has recently changed the basis of this calculation from a photon fluence scaling method to a direct Monte Carlo (MC) calculation. The MC conversion uses an EGSnrc model of the cobalt source that has been validated against water tank and graphite phantom measurements, a step that is required to quantify uncertainties in the underlying interaction coefficients in the MC code. A comparison with the Bureau International des Poids et Mesures (BIPM) as part of the key comparison BIPM.RI(I)-K4 showed an agreement of 0.9973 (53).

Adjoint-Weighted Tallies for k-Eigenvalue Calculations with Continuous-Energy Monte Carlo

A Monte Carlo method is developed that performs adjoint-weighted tallies in continuous-energy k-eigenvalue calculations. Each contribution to a tally score is weighted by an estimate of the relative magnitude of the fundamental adjoint mode, by way of the iterated fission probability, at the phase-space location of the contribution. The method is designed around the power iteration method such that no additional random walks are necessary, resulting in a minimal increase in computational time. The method is implemented in the Monte Carlo N-Particle (MCNP) code. These adjoint-weighted tallies are used to calculate adjoint-weighted fluxes, point reactor kinetics parameters, and reactivity changes from first-order perturbation theory. The results are benchmarked against discrete ordinates calculations, experimental measurements, and direct Monte Carlo calculations.

Accuracy of out-of-field dose calculations by a commercial treatment planning system

The dosimetric accuracy of treatment planning systems (TPSs) decreases for locations outside the treatment field borders. However, the true accuracy of specific TPSs for locations beyond the treatment field borders is not well documented. Our objective was to quantify the accuracy of out-of-field dose predicted by the commercially available Eclipse version 8.6 TPS (Varian Medical Systems, Palo Alto, CA) for a clinical treatment delivered on a Varian Clinac 2100. We calculated (in the TPS) and determined (with thermoluminescent dosimeters) doses at a total of 238 points of measurement (with distance from the field edge ranging from 3.75 to 11.25 cm). Our comparisons determined that the Eclipse TPS underestimated out-of-field doses by an average of 40% over the range of distances examined. As the distance from the treatment field increased, the TPS underestimated the dose with increasing magnitude-–up to 55% at 11.25 cm from the treatment field border. These data confirm that accuracy beyond the treatment border is inadequate, and out-of-field data from TPSs should be used only with a clear understanding of this limitation. Studies that require accurate out-of-field dose should use other dose reconstruction methods, such as direct measurements or Monte Carlo calculations.

Dose reduction in CT using bismuth shielding: measurements and Monte Carlo simulations.

In this research, using direct measurements and Monte Carlo calculations, the potential dose reduction achieved by bismuth shielding in computed tomography was evaluated. The patient dose was measured using an ionisation chamber in a polymethylmethacrylate (PMMA) phantom that had five measurement points at the centre and periphery. Simulations were performed using the MCNPX code. For both the bare and the bismuth-shielded phantom, the differences of dose values between experiment and simulation were within 9%. The dose reductions due to the bismuth shielding were 1.2-55% depending on the measurement points, X-ray tube voltage and the type of shielding. The amount of dose reduction was significant for the positions covered by the bismuth shielding (34 - 46% for head and 41 - 55% for body phantom on average) and negligible for other peripheral positions. The artefact on the reconstructed images were minimal when the distance between the shielding and the organs was >1 cm, and hence the shielding should be selectively located to protect critical organs such as the eye lens, thyroid and breast. The simulation results using the PMMA phantom was compared with those using a realistically voxelised phantom (KTMAN-2). For eye and breast, the simulation results using the PMMA and KTMAN-2 phantoms were similar with each other, while for thyroid the simulation results were different due to the discrepancy of locations and the sizes of the phantoms. The dose reductions achieved by bismuth and lead shielding were compared with each other and the results showed that the difference of the dose reductions achieved by the two materials was less than 2-3%.

Dosimetric and radiobiological evaluation of dose distribution perturbation due to head heterogeneities for Linac and Gamma Knife stereotactic radiotherapy

Introduction. In SRT/SRS, dedicated treatment planning systems are used for the calculation of the dose distribution. The majority of these systems utilize the standard TMR/OAR formalism for dose calculation as well as they usually neglect any perturbation due to head heterogeneities. The aim of this study is to examine the errors due to head heterogeneities for both absolute and relative dose distributions in stereotactic radiotherapy. Materials and methods. Dosimetric measurements in phantoms have been made for linac stereotactic irradiation. CT-based phantoms have been used for Monte Carlo simulations for both linac-based stereotactic system and Gamma Knife unit. Absolute and relative dose distributions have been compared between homogeneous and heterogeneous media. DVH and TCP results are presented for all cases. Results. The maximum absolute dose difference at the isocenter was 2.2% and 6.9% for the linac and Gamma Knife respectively. The impact of heterogeneity in the target DVH was minor for the linac technique whereas considerable difference was observed for the Gamma Knife treatment. This was reflected also to the radiobiological evaluation, where the maximum TCP difference for the linac system was 2.7% and for the Gamma Knife was 4%. Discussion and conclusions. The errors rising from the existence of head heterogeneities are not negligible especially for the Gamma Knife which uses lower energy beams. The errors of the absolute dose calculation could be easily eliminated by implementing a simple heterogeneity correction algorithm at the TPS. Nevertheless, the errors for not taking into account the lateral electron transport would require a more sophisticated approach and even direct Monte Carlo calculation.

Monte Carlo efficiency transfer method for full energy peak efficiency calibration of three type HPGe detectors: A coaxial N-type, a coaxial P-type and four BEGe detectors

Monte Carlo efficiency transfer method was used to determine full energy peak efficiency of three type high-purity germanium (HPGe) detectors: a coaxial N-type, a coaxial P-type and four broad energy germanium detectors (BEGe). The comparisons between calculations and point source experiments showed that the relative deviations were mostly within ±5% for BEGe detectors in the energy range 13.9–1332 keV, within ±3% for the coaxial N-type detector in the energy range 17.5–662 keV, and within ±3% for the coaxial P-type detector in the energy range 59.5–1332 keV, respectively. An evaluation of Monte Carlo efficiency transfer method at different source-to-detector geometry configurations indicates that peak efficiency may be underestimated when point sources are placed nearer than the reference position to detector on the detector axis, and may be overestimated when point sources are placed far from the detector axis, especially for low-energy photons. Some discussions about direct Monte Carlo calculation method and Monte Carlo efficiency transfer method for full-energy peak efficiency calculation are also presented in this paper.

Trajectory reversal approach for electron backscattering from solid surfaces

The backscattering of medium energy electrons from solid surfaces is investigated by analysis of a linearized Boltzmann-type kinetic equation. A closed expression is derived for the Green's function in an infinite medium valid for a spherically symmetric potential describing the interaction with the ionic subsystem. The solution is expressed in terms of fluctuations of the energy loss and scattering angles and the collision statistics associated with them. Since the fluctuation part is independent of the boundary conditions of the considered problem, solution of the backscattering problem requires an appropriate treatment of the collision statistics. In this context, the exact solution for the Oswald--Kasper--Gaukler model is derived and its limitations are analyzed. An exact approach is presented and implemented in an efficient Monte Carlo scheme based on the trajectory reversal technique. The resulting procedure is faster than the conventional Monte Carlo algorithm by several orders of magnitude. Results for the angular distribution are compared with conventional Monte Carlo calculations and perfectly agree with the latter within their statistical uncertainty. A second approximate algorithm is also given. The approximation involved in this second procedure turns out to be very reasonable: deviations from direct Monte Carlo calculations remain below $\ensuremath{\sim}5%$ for energies exceeding $200\phantom{\rule{0.3em}{0ex}}\mathrm{eV}$. The integral elastic-backscattering coefficient for normal incidence for a large number of materials in the energy range $50\phantom{\rule{0.3em}{0ex}}\mathrm{eV}--10\phantom{\rule{0.3em}{0ex}}\mathrm{keV}$ is found to approximately exhibit a universal dependence on the ratio of the inelastic and the transport mean free paths, the so-called scattering parameter.

Multi-Assembly Analysis with an Advanced Correlated Sampling Method

The Advanced Correlated Sampling (ACS) method, which can greatly enhance the calculation efficiency of Monte Carlo (MC) methods using a perturbation technique and a pseudo-scattering technique, is improved so as to be able to handle multi-assembly analysis. Two boundary connection methods—the neutron path connection and neutron current connection methods—are developed, and their accuracies for a 2×2 multi-assembly system are confirmed. The speedup factor against the corresponding direct MC calculation is a maximum of about 100. These boundary connection methods are also applied in pin-wise 2-D quarter core calculations as an example of a large-scale system. A newly developed perturbed neutron source iteration technique allows for very smooth source convergence. Furthermore, by incorporating the Direct Response Matrix (DRM) method as the acceleration means for the source convergence, the number of source iteration cycles is significantly shortened. The results confirm that this improved ACS method has the ability to perform pin-wise whole core calculations very accurately and efficiently.

Multi-Assembly Analysis with an Advanced Correlated Sampling Method

The Advanced Correlated Sampling (ACS) method, which can greatly enhance the calculation efficiency of Monte Carlo (MC) methods using a perturbation technique and a pseudo-scattering technique, is improved so as to be able to handle multi-assembly analysis. Two boundary connection methods—the neutron path connection and neutron current connection methods—are developed, and their accuracies for a 2×2 multi-assembly system are confirmed. The speedup factor against the corresponding direct MC calculation is a maximum of about 100. These boundary connection methods are also applied in pin-wise 2-D quarter core calculations as an example of a large-scale system. A newly developed perturbed neutron source iteration technique allows for very smooth source convergence. Furthermore, by incorporating the Direct Response Matrix (DRM) method as the acceleration means for the source convergence, the number of source iteration cycles is significantly shortened. The results confirm that this improved ACS method has the ability to perform pin-wise whole core calculations very accurately and efficiently.

Determination of the Ion Sputtering-induced In-depth Distribution by Means of Elastic Peak Electron Spectroscopy

Elastic peak depth profiling was carried out on an Mo/Si multilayer system using a rotating specimen, grazing angle of incidence (86° with respect to the surface normal) and 0.5 keV Ar ion energy. The depth profiling was simulated by dynamic TRIM (T-DYN) code. The T-DYN code provided the in-depth distribution of the elements in any moment of the depth profiling process. Direct Monte Carlo calculation was developed to calculate the intensity of the elastic peak emitted from the distribution determined by the T-DYN code. Good agreement was found between the calculated and measured elastic peak depth profiles, proving that the Monte Carlo routine developed is readily applicable to inhomogeneous layers and the mixing of the Mo/Si layer is well described by the T-DYN code.

Determination of the Ion Sputteringinduced Indepth Distribution by Means of Elastic Peak Electron Spectroscopy

Elastic peak depth profiling was carried out on an Mo/Si multilayer system using a rotating specimen, grazing angle of incidence (86° with respect to the surface normal) and 0.5 keV Ar ion energy. The depth profiling was simulated by dynamic TRIM (T-DYN) code. The T-DYN code provided the in-depth distribution of the elements in any moment of the depth profiling process. Direct Monte Carlo calculation was developed to calculate the intensity of the elastic peak emitted from the distribution determined by the T-DYN code. Good agreement was found between the calculated and measured elastic peak depth profiles, proving that the Monte Carlo routine developed is readily applicable to inhomogeneous layers and the mixing of the Mo/Si layer is well described by the T-DYN code. © 1997 by John Wiley & Sons, Ltd.