KR20100066572A - Reactor dosimetry applications using a parallel 3-d radiation transport code - Google Patents

Abstract

The invention relates generally to a method for the calculation of radiation field distributions employing a new parallel 3-D radiation transport code and, a multi-processor computer architecture. The code solves algorithms using a domain decomposition approach. For example, angular and spatial domains can be partitioned into subsets and, the subsets can be independently allocated and processed.

Description

FIELD CIRCUIT Calculation Method and Computer Program

Cross Reference to Related Applications

This provisional patent application claims priority to US Provisional Patent Application No. 60 / 975,525, filed September 27, 2007.

The present invention relates generally to the calculation of the radiation field distribution and is particularly useful for predicting neutron-dose response and internal components for the reactor cavity.

Various methods can be used to obtain numerical solutions to the Linearized Boltzmann Equation (LBE) for neutron and gamma radiation transmission applications. Discrete ordinates method (S _N ) is one such method, especially in the field of atomic engineering. The numerical solution of the S _N equation is achieved through simultaneous discretization of the phase space, ie the angle, space and energy domains. Simultaneous discretization of the phase space leads to many unknowns in the S _N equation, so extensive computational resources are required to solve this problem.

For large scale 3-D neutron and gamma transmission applications, the main memory required to generate the numerical solution of the LBE using the S _N equation may exceed the current computing capacity of a typical single processor workstation. For example, in a conventional two-loop pressurized water reactor (PWR) characterized by approximately 1.5 million spatial meshes, an S ₈ quadrature set, a P ₃ extension of the diffusion kernel, and a 47-neutron energy group. The solution to the full 3-D neutron transfer problem for the above can lead to a main memory requirement of approximately 45 GByte. The critical computational resources required can preclude the use of single processor workstations to solve this problem.

It is required to overcome these difficulties by developing new solution algorithms for the S _N equations in order to take advantage of the multiprocessor computational architecture, ie distributed memory architecture. For example, it would be desirable to configure a number of physically dependent workstations that are linked to each other via a network backbone to establish what is generally referred to as a cluster computing environment. This kind of computing platform has become widespread in recent years, especially in scientific computing and large scale numerical simulation. However, it is necessary to develop specialized algorithms to exploit the characteristics of the cluster environment.

Thus, there is room for improvement in the solution algorithm set of the S _N equation to take advantage of the multi-processor computational architecture. There is also room for improvement in methods for obtaining numerical solutions of LBEs for the calculation of radiation field distributions such as neutron and gamma radiation field distributions. There is also room for improvement in how to predict the dosimetric response in an accurate and efficient manner for the application of the reactor.

A method of calculating a radiation field distribution, comprising applying a 3-D radiation transmission computer code, the code comprising a domain decomposition algorithm comprising a domain selected from the group consisting of angle and spatial domains, This domain is assigned and processed independently on a multiprocessor computer architecture.

A computer program for calculating a radiation field distribution, comprising code segments that, when executed, divide the angle and spatial domains into subsets, assign these subsets independently, and process the subsets on a multi-processor architecture.

1A shows the geometric and physical distribution of a 3-D transmission model for a two-loop PWR. 1B shows the 2-D section of the model on the xy plane at z = 0.0 for the 2-loop PWR.
FIG. 2 shows the measured to calculated (M / C) ratio of dose measurement data calculated using the directional theta weighted adaptive differential scheme.
3A, 3B and 3C show M / for a modified reactor pressure vessel with a modified thickness as compared to the unmodified thickness for core upper position 3a, core midplane position 3b and core lower position 3c. The C ratio is shown.
4 shows the acceleration obtained for a given (up to 20) processors using different domain decomposition strategies.
5 shows a flowchart of an embodiment of the invention, wherein a 3-D radiation transmission computer program is implemented to generate a radiation field distribution.

The present invention relates to a method of calculating a radiation field distribution for a system. The radiation field may include neutron and gamma radiation. In one aspect of the invention, the ex-vessel neutron dosimetry response of the cavity of the reactor can be calculated. The type of reactor is not limited and may include various commercial designs known in the art. Suitable reactors may include, but are not limited to, pressurized water reactors (PWRs) and boiling water reactors (BWRs). For simplicity of explanation, this aspect of the invention will be described with reference to a two-loop commercial PWR. Computational modeling of the PWR to generate a dosimetric response is used in the design and operation of the PWR.

The method of the present invention involves the application of a 3-D parallel radiation transmission code referred to as Rapid Parallel Transport Of Radiation-Multiple 3D Geometries (RAPTOR-M3G). The transmission code provides a set of parallel algorithms for solving the S _N equation. This method is based on a domain decomposition algorithm, where the space, angle, and / or energy domains are divided into subsets that can be independently allocated and processed in a multi-processor architecture. Examples of suitable 3-D parallel deterministic transmission codes known include PENTRAN ^™ and PARTISN. (Sjoden GE and Haghighat A., "PENTRAN-Parallel Enviroment Neutralparticle TRANport in 3-D Cartesian Geometry", Proceedings of the Joint International Conference on Mathematicla Methods and Supercomputing for Nuclear Application, Vol. 1, pp. 232-234, Saratoga Springs, New York (1997). Compared to conventional uniprocessor applications, the method of the present invention reduces the memory requirements per processor in addition to the computational load, resulting in an efficient solution to large scale 3-D problems.

RAPTOR-M3G computer code is developed in Fortran 90 using the Message Transfer Interface (MPI) parallel library. (Gropp W., Lusk E., and Skellum A., Using MPI Portable Paralle Programming with the Message Passing Interface, The MIT Press, Cambridge, Massachusetts (1999)). Several features of the RAPTOR-M3G include:

Multi-group S _N equations for 3-D Cartesian coordinates (RAPTOR-XYZ) and cylindrical geometry (RAPTOR-RTZ) for non-uniform orthogonal structured mesh (MA Hunter, G. Longoni, and SL Anderson, "Extension of RAPTOR-M3G to r-θ-z geometry for use in reactor dosimetry applications ", Proceedings of the 13th International Symposium on Reactor Dosimetry, The Metherlands (2008));

Spatial, angular and combined spatial / angle domain decomposition algorithms;

Positive finite weighting differential schemes: Zero / Theta Weighted, and Directional Theta Weighted;

Automatic generation of level-symmetric orthogonality set to approximately 20 (Longoni G. et al., "Investigation of New Quadrature Sets for Discrete Ordinates Method with Application to Non-Conventional Problems," Transaction of the American Nuclear Society, Vol. 84, pages 224-226 (2001);

Parallel memory: allows local allocation of spatial and angular sub-domains, thus reducing memory requirements per processor;

Parallel tasking: Simultaneous solution of S _N equations for multiple processors to reduce computation time when compared to single processor technology;

Parallel I / O: Each processor accesses its storage locally to reduce I / O time.

ㆍ Compatibility and Integration with BOT3P (R. Orsi, "Potential Enhanced Performances in Radiation Transport Analysis on Structured Mesh Grids Made Available by BOT3P", Nuclear Science and Engineering, Vol.157, 110-116 (2007), Automated Mesh Generator , And GIP, a multi-group cross-section preprocessor.

Figure 1A illustrates the geometry of a 3-D transmission model, such as a two-loop PWR, for an embodiment of the present invention. The PWR may include a 12-foot atomic core, a thermal protection design and a 3-inch reactor cavity air gap. Model geometries include water-core mixtures, core shrouds, core barrels, heat shields, reactor pressure vessels (RPVs) including stainless steel-steel liners, and reflective separators. It is usually common to have an RPV cylindrical shape of a PWR, for example closed at both ends by a lower head and a controllable upper head. The upper and lower inner regions above and below the reactor are modeled using steel-water mixtures. The bottom interior of the RPV includes a core barrel (ie core support structure). The core barrel is surrounded by a heat shield between the core barrel and the inner wall of the RPV. In some embodiments, neutron pads are used instead of heat shields. The core cover sets the inside of the core barrel. Each downcomer surrounds the reactor core barrel. Cooling solution, typically water, is circulated to the downcomer.

FIG. 1B shows the 2-D cross section of the model on the x-y plane at z = 0.0 for the 2-loop PWR. Also shown is the mass distribution of the PWRs. Model geometry and mesh discretization are generated using BOT3P code, version 5.2. This model extends from 0.0 cm to 245.0 cm along the x and y axes and from -200.0 cm to 200.0 cm along the z axis. A uniform mesh is applied to this model conductor, with a mesh size of 2.0 x 2.0 x 4.0 cm determined along the x, y and z axes, respectively, yielding a total of 1,464,100 meshes.

Cross sections for material mixtures in the transfer model can be found in the RSGC Data Library Collection BUGLE-96, "Coupled 47 Neutron, 20 Gamma-Ray Group Cross Section Library Derived from ENDF / B-VI for LWR Shielding and Pressure Vessel Dosimetry Application, "Oak Ridge National Laboratory, Oak Ridge, TN (1999) and GIP computer code as part of the DOORS package (RSICC Computer Code Collection DOORS 3.2a," One-, Two- and Three-Dimensional Discrete Ordinates Nuetron / Phton Transport Code System, "Oak Ridge National Laboratory, Oak Ridge, TN (2003). S ₈ level symmetric orthogonal set and spherical harmonic extension of P ₃ spreading kernel are used for the transmission calculation. Passive neutron The detector's system can be installed in the reactor cavity air gap between the reflective separation and the pressure vessel, and the dosimetry system accurately measures the fast neutron exposure on the reactor vessel's beltline. A pure metal film can be placed in the reactor cavity surrounded by aluminum cells, which minimizes the distortion of the fast neutron spectrum, resulting in efficient free-field measurements. The neutron dosimeter being defined is not explicitly defined in the transmission model.

One aspect of the present invention includes a domain decomposition algorithm for discretization of the S _N equation, which divides the spatial and / or angular energy domains into subsets that can be independently allocated and processed on a multi-processor architecture. use.

The spatial and angular discretization of the S _N equation and the angular domain decomposition algorithm described herein are specified for the 3D Cartesian XYZ version of the code. The formation of the S _N equation developed for RAPTOR-RTZ differs from RAPTOR-XYZ due to the presence of the diffusion redistribution period.

The phase space of the S _N equation is discretized into angle, space and energy. Thus, the result set of linear algebraic equations is suitable for solutions on digital computers. The energy domain is discretized into multiple discrete intervals, g = 1 ... G, using a multigroup approach, starting with the highest energy particle (g = 1) and down to the lowest energy particle (g = G). It ends. The transfer equation (ie LBE) of the multigroup approximation is formulated in equation (1).

Angular domains are discretized by considering a limited set of directions and applying a suitable orthogonal integration scheme. Each discrete direction can be visualized as a point on the surface of a unit sphere having an associated surface area mathematically corresponding to the weight in the orthogonal direction. The combination of discrete directions and corresponding weights is referred to as an orthogonal set. In general, an orthogonal set must satisfy a number of conditions in order to be determined accurately and mathematically, for example, a level-symmetric orthogonal set (LQn) and a Legendre polygonal orthogonal set (Longoni G and Haghighat A., "DEvelopment of New Quadrature Sets with the Ordinate Splitting Technique, "Proceedings of the ANS International Meeting on Mathematical Methods for Nuclear Application (M & C 2001), Salt Lake City, UT, 9-13 September 2001, American Nuclear Society, Inc., Several approaches can be used, such as: Grange Park, IL (2001) The orthogonal set developed and implemented in RAPTOR-M3G is based on the LQn method.

Spatial variables can be discretized by various techniques, such as finite differential and finite element methods. The formula developed in RAPTOR-M3G is based on a finite differential scheme that divides the spatial domain into computational cells, such as finite meshes, where the cross section is an assumed constant within each cell. In the 3D Cartesian geometry, the angular flux at the cell-center position is obtained from equation (2).

In equation (2), the angle and energy dependence are represented by the indices m and g, respectively. The term q _{i , j, k} represents the sum of diffusion, division and external sources at the cell center. The indices i, j, k represent cell-centric values and the weights a _{i , j, k, m, g} , b _{i , j, k, m, g} and c _{i , k, m, g} are between 0.5 and 1.0 RAPTOR-M3G calculates weights during transmission sweeps using theta-weighted (TW), zero-weighted (ZW) or adaptive directional theta weighted (DTW) differential methods (B. Petrovic and Haghighat, "New Directional Theta-Weighted S _N Differencing Scheme and Its Application to Pressure Vessel Fluence Calculations," Proceedings of the 1996 Radiation Protection and Shielding Topical Meeting, Falmouth, MA, pp. 1-10 (1996)) .

The S _N equation is solved by starting from the boundary of the problem domain and traveling through each direction. This solution process is also referred to as transmission sweep. The angular flux defined at the center cell location is obtained from the boundary condition or boundary angle flux previously calculated in the adjacent cell. The cell center angle flux is calculated using equation (2). The angular flux from the calculating cell is calculated using an additional relationship called "differential solution".

Transfer sweeps are performed within an iterative process called source iterations, also called fixed point iterations or Richardson iterations. This process continues until a suitable hunting criterion is met, ie, the relative error for any average Stala flux between two iterations is below a predetermined cutoff value (Adams ML and Larsen EW, "Fast Iterative Methods for Discrete-Ordinates Particle Transoirt Caculations, "Progress in Nuclear Energy, Vol. 40, n.1 (2002). For radiation protection calculations, this cutoff value is set to 1.0e ^-3 or 1.0e ^-4 . It is common.

The parallel algorithm developed in RAPTOR-M3G is based on the decomposition of the angle and / or spatial domains on the processor network. The RAPTOR-M3G creates _a number of processors based on virtual topology that are assigned to angle and spatial domains, respectively, designated P _a and P _s . The total number of processors required for any decomposition is P _n = P _a · P _s . Based on this information, the processor network is mapped onto the spatial and angular domains that create a virtual topology that associates each processor to its local sub-domain.

The angular domain is divided on an octant basis, with processors assigned on the angular domain assigned sequentially to the local octets. The local number of eight elements allocated per processor is given by equation (3).

Transmission sweeps are performed locally on the N _{loct octet} on the P _a processor, and the MPI forwarders for the angular domain are used to synchronize the angular flux among the processors and identify the reflective boundary conditions.

The spatial domain is partitioned along the z axis by sequentially assigning P _s processors to multiple xy planes. The total number of fine meshes along the z axis, ie km, is split on the P _s processor, and a mapping array, ie kmloc, is used to assign the xy plane to the P _s processor. The number of xy planes allocated to the P _s processor is arbitrary, but the condition of equation (4) needs to be satisfied to define a spatial decomposition that is topologically consistent with the problem geometry.

The flexibility of mapping these processors on the spatial domain to any number of x-y planes may depend on the fact that the number of z planes may not be exactly partitionable by the number of processors on the spatial domain. Uneven partitioning of the x-y plane on the Ps processor can lead to processor load imbalance with the resulting loss of performance. In the present invention, hybrid angle / spatial decomposition schemes can be applied to overcome this difficulty. Hybrid decomposition involves a combination of angular and spatial domains that include simultaneous division of these domains. Hybrid domain degradation is further described in Example 1 below.

5 shows a flowchart of an embodiment of the present invention, in which a 3-D radiation computer program for generating a radiation field distribution is implemented. This embodiment includes obtaining geometry and material information for the system to be modeled. This information may be obtained from various sources such as, but not limited to, reactor drawings. The S _N order required for proper calculation can be selected. In addition, the Pn extension order for the angular domain (ie, diffusion kernel and angular flux) may be selected. Suitable differential measures (ie, TW, ZW or DTW) may be selected. A 3-D model of the system can be generated and discretized using a suitable mesh generator for Cartesian XYZ or RTZ geometry. Cross section tables for each material are generated by mixing the materials using a suitable data set such as BUGLE-96. The number of processors and the corresponding decomposition scheme for the problem solved are selected. Thereafter, the calculation is performed. In the next calculation, post-processing and analysis of the results produced can be performed.

The RAPTOR-M3G computer code used in the present invention provides an accurate and efficient (eg, reduced computation time) solution to the radiation transmission problem. In one aspect of the invention, the vessel outer neutron dose calculation response is calculated for the cavity of the reactor vessel. Compared with the actual measurements, the fast neutron response of the reactor cavity air gap of the 2-loop PWR calculated by RAPTOR-M3G was on average 96% accurate. In addition, a solution to the transmission problem was obtained in approximately 106 minutes with clock time on a 20-processor computer cluster using a hybrid angle / space domain decomposition scheme.

While particular embodiments of the present invention have been described, those skilled in the art will recognize that various modifications and substitutions may be made in detail in light of the overall disclosure herein. For example, one aspect of the invention disclosed herein relates to the nuclear industry, in particular nuclear reactors. However, the present invention may be used in a wide variety of other fields, such as the medical field. For example, the present invention can be used to determine the dose of radiation delivered to a patient for cancer treatment and / or management. Accordingly, the specific devices disclosed are exemplary and are not intended to limit the scope of the invention, which is defined by the full scope of the appended claims and all equivalents thereof.

Yes

Example 1- RAPTOR - M3G  Parallel performance analysis

The transfer calculation described in the example was performed using a 20 processor computer cluster, RAPTOR-M3G running on EAGLE-1. The cluster consists of five nodes with two dual-core dual processor AMD Opteron 64-bit architectures. The total cluster memory available, or RAM, was 40 GBytes, and the network interconnect was characterized by 1 GBit / s bandwidth. Using this hardware configuration, the RAPTOR-M3G completed the entire 3D transfer calculation for a two-loop PWR in approximately 106 minutes on 20 processors. No significant difference in performance was observed using the DTW, TW or ZW differential schemes.

In addition, a simple test problem was set up to analyze the parallel performance of the code. A test problem consisting of a 50 x 50 x 50 cm box using a uniformly distributed fixed source was discretized into a 1 cm uniform mesh. An S8 orthogonal set and P0 isotropic diffusion were used, with one set of energy group cross sections. Wall clock periods, acceleration and parallel efficiency were used to evaluate the parallel performance of the RAPTOR-M3G. Acceleration and parallel efficiencies are defined in equations 5 and 6, respectively.

Where T _s and T _p are wall clock times required by uniprocessor and multiprocessor calculations, respectively. N _p is the number of processors used to achieve wall clock time. 4 shows a comparison of the acceleration obtained for 20 processors using different decomposition schemes.

Acceleration obtained using spatial decomposition gradually decreased as the number of processors increased. This behavior is thought to be due to a finer computational granularity per processor, as the number of operations per processor decreases as the spatial domain breaks down into smaller sub-domains, increasing inter-processor communication time. Performance is reduced. Network data transfer between nodes has generally been a limiting factor in distributed memory architecture. The larger number of iterations required to converge on the problem was due to further reducing the performance of the spatial decomposition scheme. However, hybrid decomposition with simultaneous division of angular and spatial domains yielded better results. This behavior is thought to be due to the sparse computational granularity induced by this decomposition, even for hybrid decomposition, the number of iterations required to hunt the problem does not increase by spatial decomposition.

Example 2-Measured Dose Response RAPTOR - M3G  Comparison of responses

A comparison was made between the measured linear response and the corresponding prediction obtained using RAPTOR-M3G. IRDF-2002 Dosimetry Library (I. Kodeli and A. Trkov, "Validation of the IRDF-2002 Dosimetry Library," Nuclear Instruments and Methods in Physics Research Section A: Accelerators, Spectrometers, Detectors and Associated Equipment, Volume 57, Publication 3, Pages 664-681 (2007) were used to generate the wrapped dosimetric responses to the neutron responses listed in Table 1.

The measured dose response for the responses listed in Table 1 was compared with the response calculated by RAPTOR-M3G.

Since the reactions listed in Table 1 were measured using a Cadmium protective metal film, the thermal component of the neutron spectrum was suppressed.

The measured response was obtained at four azimuth angles, 0 °, 15 °, 30 ° and 45 ° of the core midplane of the reactor cavity air gap. Because the two-loop PWR reactor is generally characterized by peak fast neutron fluctuations at the 0 ° position, additional measurements were obtained at this location due to the proximity of the neutron fuel with RPV. Specifically, at 0 °, measurements were obtained axially symmetric on the top and bottom of the active core. The initially calculated dosimetric response was found to continue to over-predict the measured data.

Further investigation revealed that the RPV thickness used in the transmission model was smaller than that measured during the in-service inspection of the reactor pressure vessel (RPV) and confirmed the initial findings. The new RPV thickness was derived into the transmission model, and the accuracy of the calculated dosimetry data improved by an average of ~ 8%. The measured to calculated value (M / C) ratio of the dosimetric data calculated using the DTW adaptive differential scheme is shown in FIG. 2. As shown in FIG. 2, the M / C ratio continued to be located at each position and within the 10% range for all dosimetric materials. Excessive predictions at 30 ° and 45 ° azimuth locations can be reduced by using non-uniform mesh refinement at these locations, where the curve of the system becomes more relevant. The average M / C ratio for all dosimetry sites was 0.96.

3A-3C provide M / C values obtained using ISI modified RPV thickness when compared to values obtained without thickness modification. This comparison was performed and provided for all dosimetry samples positioned at 0 ° azimuth positions. The modified RPV thickness using ISI measurements improved the accuracy of the calculated response at all dosimetry locations. Similar results were obtained at 15 °, 30 ° and 45 ° azimuth positions.

The M / C ratios using the RPC thicknesses modified using the measured and calculated response rate and ISI measurements at each dosimetry site are shown in Table 2. The average M / C ratio for the reactions listed in Table 2 across all dosimetry sites was 0.96.

Claims (20)

A method of calculating the radiation field distribution,
Applying a 3-D radiation transmission computer code,
The code includes a domain decomposition algorithm including a domain selected from the group consisting of an angular and spatial domain,
The domain is allocated and processed on a multi-processor computer architecture.
Method of calculating the radiation field.
The method of claim 1,
The radiation field distribution comprises a radiation field selected from the group consisting of neutron radiation and gamma radiation.
Method of calculating the radiation field.
The method of claim 1,
The radiation field distribution is for a cavity of the reactor
Method of calculating the radiation field.
The method of claim 3, wherein
The reactor is a 2-loop pressurized water reactor
Method of calculating the radiation field.
The method of claim 1,
The algorithm is structured to numerically solve the linear Boltzmann equation
Method of calculating the radiation field.
The method of claim 1,
The multi-processor computer architecture includes a plurality of physically independent workstations linked to each other by a network connection.
Method of calculating the radiation field.
The method of claim 1,
Dividing the domain into subsets that are independently allocated and processed on the multi-processor architecture;
Method of calculating the radiation field.
The method of claim 5, wherein
Partitioning the domain includes using parallel memory on a multiprocessor.
Method of calculating the radiation field.
The method of claim 1,
The radiation field distribution includes an ex-vessel neutron dosimetry response outside the vessel in the reactor cavity.
Method of calculating the radiation field.
The method of claim 7, wherein
The dosimetric response is obtained on 3-D Cartesian geometry.
Method of calculating the radiation field.
The method of claim 1,
The algorithm is structured to solve the S _N equations simultaneously using a parallel algorithm.
Method of calculating the radiation field.
The method of claim 1,
Each processor in the multi-processor computer architecture has local access to its own device selected from the group consisting of storage and local memory.
Method of calculating the radiation field.
The method of claim 1,
Further discretizing the energy domain into a plurality of discrete intervals.
Method of calculating the radiation field.
The method of claim 1,
Dividing the spatial domain into computed cells;
Method of calculating the radiation field.
The method of claim 14,
The assumed cross section within each cell is constant
Method of calculating the radiation field.
The method of claim 1,
The domain decomposition is a hybrid decomposition in which the angle and spatial domains are divided simultaneously.
Method of calculating the radiation field.
17. The method of claim 16,
The number of processors is up to 20
Method of calculating the radiation field.
The method of claim 17,
The problem is solved in less than 2 hours
Method of calculating the radiation field.
The method of claim 17,
Solutions derived using this method have an accuracy of at least 90% when compared with the calculated values.
Method of calculating the radiation field.
A computer program for calculating a radiation field distribution,
A code segment that, when executed, divides the angle and spatial domains into subsets, assigns the subsets independently, and processes the subsets on a multi-processor architecture.
Computer programs.

Images