Method and device for calculating radionuclide decay

Abstract

The invention discloses a radionuclide decay calculation method and device, and belongs to the technical field of radionuclides. The method comprises the steps of obtaining a target position index set in a preset coefficient matrix C, wherein the target position index set comprises target nuclides and position indexes of all daughters of the target nuclides; performing repeated pairwise arrangement on the position indexes in the target position index set to obtain a row index and a column index where the target element is located; obtaining the target elements in preset coefficient matrixes C and C-1 according to the row indexes and the column indexes where the target elements are located, and forming coefficient matrixes D and D-1; determining negative decay constants of nuclides corresponding to the position indexes as diagonal elements of a diagonal matrix; and based on the diagonal matrix, the coefficient matrix D and the coefficient matrix D-1, calculating the nuclide concentration field of the target nuclide and all daughters thereof at the decay time. According to the method, the problems of low accuracy of a nuclide concentration calculation result and low calculation efficiency in related technologies can be solved.

Description

Technical Field

[0001] The invention belongs to the technical field of radioactive nuclides, and in particular relates to a method and a device for calculating the decay of radioactive nuclides. Background Art

[0002] After an accident at a nuclear power plant, analyzing, predicting and evaluating the concentration of radionuclides at different distances around the plant site and further evaluating the exposure dose and risk of personnel are important components of the third-level PSA (Probabilistic Safety Assessment) of nuclear power plants. How to quickly and accurately consider the decay of nuclides and calculate the concentration field of radionuclides around the plant site is the premise and key step for reliable calculation of the third-level PSA analysis of nuclear power plants.

[0003] At present, in the three-level PSA analysis of nuclear power plants, the traditional nuclide decay only considers the decay of up to three generations of progeny. In addition, since the calculation of subsequent personnel exposure doses and risks needs to consider the influence of hourly meteorological conditions throughout the year, the calculation results of nuclide concentrations are inaccurate and the calculation efficiency is low. Summary of the invention

[0004] The technical problem to be solved by the present invention is to provide a method and device for calculating the decay of radioactive nuclides in view of the above-mentioned deficiencies in the prior art, which can improve the accuracy of the calculation results of the nuclide concentration and shorten the calculation time.

[0005] In a first aspect, the present invention provides a method for calculating radioactive nuclide decay, comprising: obtaining a target position index set in a preset coefficient matrix C, wherein the target position index set includes position indexes of a target nuclide and all its daughters, the number of position indexes is m, m is a positive integer, and the coefficient matrix C refers to a coefficient matrix associated with a complete decay chain of all nuclides; arranging the m position indexes in the target position index set in pairs with repetition to obtain a row index and a column index where the target element is located; and arranging the target element in the preset coefficient matrices C and C according to the row index and the column index where the target element is located. -1 Get the target element from the matrix and form the coefficient matrices D and D -1 , where C -1 is the inverse matrix of C, D -1 is the inverse matrix of D; the negative decay constants of the nuclides corresponding to the m position indices are determined as the diagonal elements of the diagonal matrix; based on the diagonal matrix, the coefficient matrix D and D -1 , calculate the nuclide concentration field of the target nuclide and all its daughters at the decay time.

[0006] Preferably, before obtaining the target position index set in the preset coefficient matrix C, the calculation method further includes: constructing coefficient matrices C and C related to the complete decay chain of n radionuclides according to the decay constants and branching ratios of the n radionuclides and the position indexes of the daughter nuclei are all after the position index of the mother nucleus. -1 , get the preset coefficient matrices C and C -1 , where n is a positive integer and the n radioactive nuclides include the target nuclide.

[0007] Preferably, the preset coefficient matrix C satisfies formula (1):

[0008]

[0009] Among them, c kj represents the kth row and jth column element of the coefficient matrix C, Λ ik represents the i-th row and k-th column element of the lower triangular matrix Λ,

[0010] The preset C -1 Satisfying formula (2):

[0011]

[0012] Among them, c ik Represents C -1 The element in the i-th row and k-th column of kj Represents C -1 The k-th row and j-th column element of

[0013] The lower triangular matrix Λ satisfies equation (3):

[0014]

[0015] Among them, λ i and λ j are the decay constants of nuclide i and nuclide j, respectively, and b ji Represents the branching ratio of nuclide i and nuclide j.

[0016] Preferably, the method of obtaining a target position index set in a preset coefficient matrix C specifically includes: determining the position index corresponding to the target nuclide in the preset coefficient matrix C; obtaining the position index corresponding to all sub-bodies of the target nuclide in the preset coefficient matrix C; and deduplicating the position index of the target nuclide and the position indexes of all sub-bodies to obtain a target position index set.

[0017] Preferably, the obtaining of the position indexes corresponding to all sub-bodies of the target nuclide in the preset coefficient matrix C specifically includes: for each target nuclide position index, filtering out the target row index from the number of rows greater than the target column index, wherein the element value in the coefficient matrix C corresponding to the target row index is non-zero; marking the target row index as the position index corresponding to the sub-bodies of the target nuclide in the preset coefficient matrix C.

[0018] Preferably, the diagonal matrix, coefficient matrix D and D -1 , calculate the nuclide concentration field of the target nuclide and all its daughters at the decay time, including:

[0019] The nuclide concentration field is calculated according to formula (4):

[0020]

[0021] Where N(0) is the initial concentration of the target nuclide, Satisfying formula (5):

[0022]

[0023] Λ d is a diagonal matrix.

[0024] Preferably, the calculating of the nuclide concentration field specifically includes: calculating the nuclide concentration field in parallel based on a graphics processing unit (GPU).

[0025] In a second aspect, the present invention further provides a device for calculating radioactive nuclide decay, comprising: a first acquisition module, used to obtain a target position index set in a preset coefficient matrix C, wherein the target position index set includes the position index of the target nuclide and all its daughters, the number of position indexes is m, m is a positive integer, and the coefficient matrix C refers to a coefficient matrix associated with the complete decay chain of all nuclides. An arrangement module, connected to the first acquisition module, is used to arrange the m position indexes in the target position index set in pairs with repetitions to obtain the row index and column index of the target element. A second acquisition module, connected to the arrangement module, is used to arrange the row index and column index of the target element in the preset coefficient matrices C and C according to the row index and column index of the target element. -1 Get the target element from the matrix and form the coefficient matrices D and D -1 , where C -1 is the inverse matrix of C, D -1 is the inverse matrix of D. A determination module is connected to the first acquisition module and is used to determine the negative decay constants of the nuclides corresponding to the m position indexes as diagonal elements of the diagonal matrix. A calculation module is connected to the second acquisition module and the determination module and is used to calculate the negative decay constants of the nuclides corresponding to the m position indexes based on the diagonal matrix, the coefficient matrix D and D -1 , calculate the nuclide concentration field of the target nuclide and all its daughters at the decay time.

[0026] Preferably, the preset coefficient matrix C satisfies formula (1):

[0027]

[0028] Among them, c kj represents the kth row and jth column element of the coefficient matrix C, Λ ik represents the i-th row and k-th column element of the lower triangular matrix Λ,

[0029] The preset C -1 Satisfying formula (2):

[0030]

[0031] Among them, c ik Represents C -1 The element in the i-th row and k-th column of kj Represents C -1 The k-th row and j-th column element of

[0032] The lower triangular matrix Λ satisfies equation (3):

[0033]

[0034] Among them, λ i and λ j are the decay constants of nuclide i and nuclide j, respectively, and b ji Represents the branching ratio of nuclide i and nuclide j.

[0035] In a third aspect, the present invention further provides a device for calculating radionuclide decay, comprising: a memory and a processor, wherein the memory and the processor are communicatively connected to each other, the memory stores computer instructions, and the processor executes the method for calculating radionuclide decay described in the first aspect by executing the computer instructions.

[0036] The method and device for calculating radionuclide decay provided by the present invention have the following advantages: -1 Considering the complete decay chain of all nuclides, based on the coefficient matrices C and C -1 The final calculated nuclide concentration field also includes the complete decay chain of the source nuclide, thus improving the accuracy of the calculation results. -1 The coefficient matrices D and D obtained from -1 It is only related to the complete decay chain of the source nuclide, and the dimension of the coefficient matrix D is much smaller than that of the coefficient matrix C, so it can shorten the calculation time of the nuclide concentration and thus improve the calculation efficiency. BRIEF DESCRIPTION OF THE DRAWINGS

[0037] Figure 1 This is a flow chart of a method for calculating radionuclide decay according to Example 1 of the present invention;

[0038] Figure 2 A schematic diagram of a decay chain and branching ratio according to Example 1 of the present invention;

[0039] Figure 3 is a schematic diagram of another decay chain and branching ratio according to Example 1 of the present invention;

[0040] Figure 4 A schematic flow chart of a method for calculating radionuclide decay according to Example 2 of the present invention;

[0041] Figure 5 A schematic diagram of a matrix D according to Embodiment 2 of the present invention;

[0042] Figure 6 A matrix D of Embodiment 2 of the present invention is -1 Schematic diagram;

[0043] Figure 7 This is a schematic diagram of a diagonal matrix in Embodiment 2 of the present invention;

[0044] Figure 8 This is a schematic diagram of the structure of a device for calculating radionuclide decay according to Example 3 of the present invention. DETAILED DESCRIPTION

[0045] In order to enable those skilled in the art to better understand the technical solution of the present invention, the embodiments of the present invention will be further described in detail below with reference to the accompanying drawings.

[0046] It should be understood that the specific embodiments and drawings described herein are only used to explain the present invention rather than to limit the present invention.

[0047] It can be understood that, in the absence of conflict, the various embodiments of the present invention and the various features in the embodiments can be combined with each other.

[0048] It can be understood that, for the convenience of description, the drawings of the present invention only show the parts related to the present invention, while the parts irrelevant to the present invention are not shown in the drawings.

[0049] It can be understood that each unit and module involved in the embodiments of the present invention may correspond to only one physical structure, or may be composed of multiple physical structures, or multiple units and modules may be integrated into one physical structure.

[0050] It can be understood that, in the absence of conflict, the functions and steps marked in the flowcharts and block diagrams of the present invention may occur in an order different from that marked in the drawings.

[0051] It is understood that the flowcharts and block diagrams of the present invention illustrate the possible architectures, functions, and operations of the systems, devices, equipment, and methods according to the various embodiments of the present invention. Each box in the flowchart or block diagram may represent a unit, module, program segment, or code, which contains executable instructions for implementing the specified functions. Moreover, each box or combination of boxes in the block diagram and flowchart may be implemented by a hardware-based system that implements the specified functions, or may be implemented by a combination of hardware and computer instructions.

[0052] It can be understood that the units and modules involved in the embodiments of the present invention can be implemented by software or hardware. For example, the units and modules can be located in a processor.

[0053] Embodiment 1:

[0054] like Figure 1 As shown, this embodiment provides a method for calculating radionuclide decay, which can be applied to scenarios such as three-level PSA analysis, accident consequence assessment, safety analysis, and environmental impact assessment in nuclear power plants. The calculation method includes:

[0055] Step 101, obtain a target position index set in a preset coefficient matrix C, wherein the target position index set includes position indexes of the target nuclide and all its daughters, the number of position indexes is m, m is a positive integer, and the coefficient matrix C refers to a coefficient matrix associated with the complete decay chains of all nuclides.

[0056] In this embodiment, the target nuclide refers to the source nuclide, such as the parent nuclide involved in the environmental release source of a nuclear accident. The position index refers to an index that is consistent with the row index and the column index, and is identified by the row index or the column index. It can be understood that the elements located on the diagonal of the matrix have a position index, and the positions of the elements located on the non-diagonal of the matrix need to be identified by a combination of the row index and the column index. Since the coefficient matrix C is a coefficient matrix related to the complete decay chain of all nuclides, the nuclide concentration field calculated based on the coefficient matrix C also takes into account the complete decay chain, thereby improving the accuracy of the calculation results.

[0057] Among them, the preset coefficient matrix C satisfies formula (1):

[0058]

[0059] c kj represents the kth row and jth column element of the coefficient matrix C, Λ ik Represents the i-th row and k-th column element of the lower triangular matrix Λ.

[0060] Default C -1 Satisfying formula (2):

[0061]

[0062] Among them, c ik Represents C -1 The element in the i-th row and k-th column of kj Represents C -1 The k-th row and j-th column element of .

[0063] The lower triangular matrix Λ satisfies equation (3):

[0064]

[0065] Among them, λ i and λ j are the decay constants of radionuclide i and nuclide j, respectively, and b ji Represents the branching ratio of radionuclide i and nuclide j.

[0066] The following describes the derivation process of the preset coefficient matrix C:

[0067] Consider a system with n radionuclides, where all daughter nuclei (i.e., daughter bodies) are sorted after the parent nucleus (i.e., target nuclide), where the vector N j represents the number of atoms of nuclide j at a certain time t, and the nuclide decay satisfies the following differential equation:

[0068]

[0069] Among them, b ji represents the branching ratio of radionuclide i decaying to j, λ j is the decay constant of radionuclide j. The differential equation of the above radionuclide decay chain can be expressed in matrix notation as:

[0070]

[0071] where Λ is a lower triangular matrix whose elements are:

[0072]

[0073] In the above formula (3), λ i and λ j are the decay constants of radionuclide i and nuclide j, respectively, and b ji represents the branching ratio of radionuclide i and nuclide j. The decay constants of n radionuclides are different, that is, the diagonal elements of Λ are different, and Λ is a lower triangular matrix, so Λ has n different eigenvalues, and the n eigenvectors corresponding to the n eigenvalues ​​are linearly independent, so Λ is a diagonalizable matrix.

[0074] By diagonalizing Λ, equation (7) can be written as:

[0075]

[0076] Here d is a diagonal matrix whose main diagonal elements are negative decay constants Λ dii =-λ j , are n different eigenvalues ​​of Λ, the n column vectors of the coefficient matrix C are n linearly independent eigenvectors of Λ, and the eigenvalues ​​corresponding to these n eigenvectors are the diagonal matrix Λ d The elements on the main diagonal of . Among them, C -1 is the inverse matrix of matrix C.

[0077] The coefficient matrices C and C are given below -1 The solution process of . Let us use formula (9) to represent the lower triangular matrix Λ.

[0078]

[0079] Since the decay constant λ of the nuclide j ,λ i The branch ratio b of nuclide i decaying to j ji are all known, so all elements in the lower triangular matrix Λ are known.

[0080] From the above analysis, we can see that solving the coefficient matrix C is to solve the n eigenvectors of the matrix Λ. jj The corresponding eigenvector c j It can be solved by equation (10), where I is the unit matrix.

[0081] (Λ-Λ jj I) c j =0 (10)

[0082] For ease of explanation, equation (11) is used instead of equation (10).

[0083]

[0084] According to formula (11), when i≤j, (Λ ii -Λ jj )c ij =0, and when i≠j, Λ jj ≠Λ ii Therefore, when i<j, c ij = 0. When i = j, c ij It can take any value, let c ij =1.

[0085] When i=j+1, according to formula (11), Λ j+1,j c jj +(Λj+1,j+1 -Λ jj )c j+1,j = 0. Therefore,

[0086] When i=j+2, according to formula (11), we have Λ j+2,j c jj +Λ j+2,j+1 c j+1,j +(Λ j+2,j+2 -Λ jj )c j+2,j = 0, so we have

[0087] In general, when i>j, according to formula (11), Λ i,j c jj +Λ i,j+1 c j+1,j +…+Λ i,i-1 c i-1,j +(Λ ii -Λ jj )c ij = 0. Therefore,

[0088] In summary, the coefficient matrix C is shown in formula (1).

[0089]

[0090] Might as well use d ij Represents C -1 The matrix element in row i and column j is C, because C is a lower triangular matrix, so C -1 It is also a lower triangular matrix. So when i<j, d ij =0.

[0091] Because CC -1 =I, so we have The expansion can be written as formula (12).

[0092]

[0093] From the above conditions, we can see that when k<i, d ki =0, when i<k, there is c ik =0, so c ii d ii = 1. And because c ii =1, so d ii =1, for i=1,2,…,n.

[0094] Because CC -1 =I, so when i>j, The expansion can be written as formula (13).

[0095]

[0096] When i<k, c ik =0, and c ii =1, so when i>j,

[0097]

[0098] As analyzed above, when k<j, there is d kj =0, so when i>j, formula (14) can be simplified to:

[0099] In summary, the matrix C -1 As shown in formula (2):

[0100]

[0101] Optionally, before obtaining the target position index set in the preset coefficient matrix C, the calculation of the radionuclide decay further includes: constructing coefficient matrices C and C related to the complete decay chain of n radionuclides according to the decay constants and branching ratios of the n radionuclides and the position indexes of the daughter nuclei are all after the position index of the parent nucleus. -1 , get the preset coefficient matrices C and C -1 , where n is a positive integer and the n radioactive nuclides include the target nuclide.

[0102] In this embodiment, it is assumed that the nuclide library has 6 nuclides (Te-131, I-131, Te-131m, Kr-85, Kr-85m, Rb-85), that is, n=6, and the coefficient matrix C refers to the coefficient matrix associated with the complete decay chains of the 6 nuclides. The position indexes and decay constants corresponding to the 6 nuclides are shown in Table 1. The decay chains and branches corresponding to the nuclides are as follows: Figure 2 and Figure 3 shown.

[0103] Table 1 Examples of position indexes and decay constants of six nuclides

[0104] Nuclide Position Index Decay constant (1 / s) Kr-85m 1 4.3E-5 Kr-85 2 2.04E-9 Rb-85 3 0 Te-131m 4 6.42E-6 Te-131 5 4.62E-4 I-131 6 1.00E-6

[0105] According to the decay constants and branching ratios of the six nuclides, and formula (3), the lower triangular matrix is ​​obtained:

[0106]

[0107] According to the position index of the child core is after the position index of the parent core, the lower triangular matrix and formula (1) are obtained:

[0108]

[0109] Then according to formula (2), we get:

[0110]

[0111] Optionally, step 101: obtaining a target position index set in a preset coefficient matrix C includes steps 1011 to 1013:

[0112] Step 1011, determining the position index corresponding to the target nuclide in the preset coefficient matrix C.

[0113] In this embodiment, it is assumed that the source nuclides of the third-level PSA analysis of the nuclear power plant are Kr-85m and Kr-85, that is, two target nuclides (s=2), and the position index corresponding to the target nuclide Kr-85m in the preset coefficient matrix C is 1, and the position index corresponding to the target nuclide Kr-85 in the preset coefficient matrix C is 2.

[0114] Step 1012: Obtain the corresponding position indexes of all daughter bodies of the target nuclide in the preset coefficient matrix C.

[0115] In this embodiment, all daughters of the target nuclide are not limited to the decay of the third-generation daughters, but refer to all daughters produced by the complete decay chain.

[0116] Specifically, obtaining the position indexes corresponding to all sub-bodies of the target nuclide in the preset coefficient matrix C includes: for each target nuclide position index, filtering out the target row index from the number of rows greater than the target column index, wherein the element value in the coefficient matrix C corresponding to the target row index is non-zero; marking the target row index as the position index corresponding to the sub-bodies of the target nuclide in the preset coefficient matrix C.

[0117] In this embodiment, for the target column index 1 where the position index 1 of the target nuclide Kr-85m is located, all the cases where the number of rows in the matrix C is greater than 1 are queried, all the elements whose corresponding position data are not 0 are found, and the row index corresponding to the element is marked as the position index corresponding to all the daughter nuclei of the target nuclide Kr-85m. For example, the position index corresponding to the target nuclide Kr-85m is 1. For all the cases where the number of columns in C is 1 and the number of rows is greater than 1, all the elements whose position data are not 0, i.e. -0.214 and -0.786, are found, and the corresponding target row indexes are 2 and 3, and 2 and 3 are marked as the corresponding position indexes of the daughters of the target nuclide Kr-85m in C. For the target nuclide Kr-85, the corresponding position index is 2. For all the cases where the number of columns in C is 2 and the number of rows is greater than 2, all the elements whose position data are not 0, i.e. -1, are found, and the corresponding target row index is 3, and 3 is marked as the corresponding position index of the daughters of the target nuclide Kr-85 in C.

[0118] Step 1013, deduplicate the position index of the target nuclide and the position indexes of all sub-bodies to obtain a target position index set.

[0119] In this embodiment, the position index (1, 2) of the target nuclide and the position index (2, 3, 3) of all daughters are deduplicated to obtain the target position index set (1, 2, 3), that is, the nuclide corresponding to the position index (1, 2, 3) in C is the position index of all nuclides after decay (including the parent and all daughters). Among them, the number of position indexes is 3, that is, m=3.

[0120] Step 102: Arrange the m position indexes in the target position index set in pairs with repetition to obtain the row index and column index where the target element is located.

[0121] In this embodiment, the three position indexes (1, 2, 3) in the target position index set are arranged in pairs with repetition to obtain 9 arrangement results, that is, the row index and column index of the 9 target elements are: (1, 1), (1, 2), (1, 3), (2, 1), (2, 2), (2, 3), (3, 1), (3, 2), (3, 3). If m = 4, 16 arrangement results are obtained, that is, there are 16 target elements.

[0122] Step 103: add the target element to the preset coefficient matrices C and C according to the row index and column index of the target element. -1 Get the target element from the matrix and form the coefficient matrices D and D -1 , where C -1 is the inverse matrix of C, D -1 is the inverse matrix of D.

[0123] In this embodiment, the coefficient matrices D and D -1 is a matrix of size mxm. Since the coefficient matrices D and D -1 is the coefficient matrix C and C -1 The matrix is ​​composed of some target elements, thus reducing the matrix dimension and simplifying the subsequent calculation complexity, thus improving the subsequent calculation speed and efficiency. For example, according to the above 9 arrangement results, the coefficient matrices C and C are respectively -1 Get the target element from the matrix and get the coefficient matrices D and D -1 :

[0124]

[0125] Step 104: determine the negative decay constants of the nuclides corresponding to the m position indices as diagonal elements of a diagonal matrix.

[0126] In this embodiment, the diagonal matrix is ​​also a matrix of size mxm, and its diagonal elements are respectively the negative decay constants of the nuclides corresponding to the m position indices. d for:

[0127]

[0128] Step 105, based on the diagonal matrix, coefficient matrix D and D -1 , calculate the nuclide concentration field of the target nuclide and all its daughters at the decay time.

[0129] Specifically, the nuclide concentration field is calculated according to formula (4):

[0130]

[0131] Where N(0) is the initial concentration of the target nuclide, Satisfying formula (5):

[0132]

[0133] Λ d is a diagonal matrix.

[0134] In this embodiment, if the initial concentration field of the nuclide is N(0), based on the parallel calculation on the graphics processing unit GPU based on the cupy library, assuming that the decay time is t, the nuclide concentration field after decay is Wherein, N(0) can be a one-dimensional, two-dimensional or three-dimensional concentration field. In this embodiment, a one-hour decay is considered, that is, t=3600 seconds, and the calculation results of N(t) are shown in Table 2.

[0135] Table 2 Concentration field calculation results

[0136] Nuclide N(0) N(t) Kr-85m 1.00E+10 8.5665E+09 Kr-85 0 1.4576E+04 Rb-85 0 0

[0137] Since GPU can process thousands of threads at the same time, it is very suitable for the task of nuclide concentration field, which involves a lot of calculations. GPU can greatly improve the calculation speed, and through parallel calculation, it can more effectively utilize computing resources, reduce calculation time, and improve research efficiency. Therefore, this embodiment is based on the parallel calculation of nuclide concentration field by graphics processing unit GPU, which can further improve the calculation speed and efficiency. This embodiment is based on Python's cupy library to perform parallel calculations on GPU, and the total calculation time is 1.047 seconds, while the calculation time based on the prior art is as long as more than ten seconds, so the calculation efficiency of this embodiment is high and the calculation speed is fast.

[0138] The calculation method of radioactive nuclide decay provided by this embodiment obtains the coefficient matrix D from the coefficient matrix C constructed by considering the complete decay chains of all nuclides. The coefficient matrix D only involves the complete decay chains of the target nuclide and all its daughters, and has a lower matrix dimension. Then, the nuclide concentration field is calculated according to the coefficient matrix D, which can improve the accuracy of the nuclide concentration calculation result and shorten the calculation time. Specifically, since the coefficient matrix C is a coefficient matrix related to the complete decay chains of all nuclides, the nuclide concentration field calculated based on the coefficient matrix C also considers the complete decay chains, rather than only considering three generations of decay, thereby improving the accuracy of the calculation result. Since the coefficient matrix D and D -1 are selected from the coefficient matrix C and C -1 by screening out some target elements, the matrix dimension is reduced, and the subsequent calculation complexity is simplified, so the subsequent calculation speed and calculation efficiency can be improved. In addition, parallel computing of the nuclide concentration field based on the graphics processing unit GPU can further improve the calculation speed and calculation efficiency.

[0139] Embodiment 2:

[0140] As Figure 4 shown, this embodiment provides a calculation method of radioactive nuclide decay, including:

[0141] S1. For all nuclides in the ICRS Report No. 107 dataset (assuming there are n nuclides, if n = 1498), combining the decay constants and branching ratios of each nuclide, and sorting them according to the position index of all daughter nuclei being behind the parent nucleus, the lower triangular matrix Λ is obtained after sorting:

[0142]

[0143] Among them, λ i and λ j are the decay constants of nuclide i and nuclide j respectively, and b ji represents the branching ratio of nuclide i and nuclide j.

[0144] S2. According to Equation (1) and Equation (2), calculate its coefficient matrix C and its inverse C -1 , including S21 - S22.

[0145] S21. Construct the coefficient matrix C, and the specific steps are as follows:

[0146]

[0147] Among them, c kj represents the element in the k-th row and j-th column of the matrix C, and Λ ik represents the element in the i-th row and k-th column of the triangular matrix Λ.

[0148] For all cases where i < j in formula (1), cij = 0;

[0149] For all cases where i = j in formula (1), c ij = 1;

[0150] For all cases where i > j in formula (1), solve it using the following iterative method:

[0151] For example, for the element in the second row and first column of the coefficient matrix C,

[0152] For the element in the third row and first column of the coefficient matrix C,

[0153] And so on, the element in the nth row and first column of the coefficient matrix C can be obtained,

[0154] Similarly, for the element in the third row and second column of the coefficient matrix C,

[0155] And so on, the element in the nth row and second column of the coefficient matrix C can be obtained,

[0156] ……

[0157] And so on, the element in the nth row and qth column of the coefficient matrix C can be obtained,

[0158] The element in the nth row and (n - 1)th column of the coefficient matrix C,

[0159] In summary, all elements in all rows and all columns of the coefficient matrix C can be obtained.

[0160] S22. Construct matrix C -1 , and the specific steps are as follows:

[0161]

[0162] where c ik represents the element in the i-th row and k-th column of matrix C -1 , and d kj represents the element in the k-th row and j-th column of matrix C -1 .

[0163] For all cases where i < j in formula (2), d ij = 0;

[0164] For all cases where i = j in formula (2), d ij = 1;

[0165] For all cases where i>j in formula (2), the following iterative method is used to solve:

[0166] For matrix C -1 The element in row 2 and column 1, d 21 =-c 21 d 11 ,

[0167] For the matrix C -1 The 3rd row and 1st column element in d 31 =-(c 31 d 11 +c 32 d 21 ),

[0168] By analogy, we can find the matrix C -1 The nth row and 1st column element in

[0169]

[0170] Similarly, for the matrix C -1 The element in row 3 and column 2, d 32 =-c 32 d 22 ,

[0171] By analogy, we can find the matrix C -1 The element at the nth row and the second column in

[0172]

[0173] …

[0174] Matrix C -1 The element in the nth row and n-1th column, d n,n-1 =-c n,n-1 d n-1,n-1 ;

[0175] In summary, the matrix C can be obtained -1 The elements of all rows and all columns in .

[0176] S3, according to the environmental release source term of the nuclear accident (assuming the number of target nuclides is s), obtain the position indexes corresponding to all parent nuclei and all their daughters, and the specific steps include S31-S35.

[0177] S31, obtain the position index corresponding to each target nuclide in the coefficient matrix C, such as index = id1, id2, ... id s ;

[0178] S32, for id1, query the coefficient matrix C for all cases where the number of columns is id1 and the number of rows is greater than id1, find all elements whose corresponding position data is not 0, mark the target row index corresponding to the element as the position index corresponding to all sub-cores of id1, that is, find all c in the coefficient matrix C u,id1 ≠0, (where u=id1+1,…,n) record the above target row indexes in sequence, assuming that the position indexes corresponding to all sub-cores of id1 are: u1, u2, u5…;

[0179] S33, for id2, query the coefficient matrix C for all cases where the number of columns is id2 and the number of rows is greater than id2, find all elements whose corresponding position data is not 0, mark the target row index corresponding to the element as the position index corresponding to all sub-cores of id2, that is, find all c in the coefficient matrix C u,id2 ≠0, (where u=id2+1,…,n), record the above target row indexes in sequence, assuming that the position indexes corresponding to all sub-cores of id2 are: u1, u3, u4…;

[0180] S34, and so on, for id s , query all columns in the coefficient matrix C as id s And the number of rows is greater than id s In the case of , find all elements whose corresponding position data is not 0, and mark the target row index corresponding to the element as id s The position indexes corresponding to all sub-cores, that is, find all c in the coefficient matrix C u,id ≠0, (where u = ids+1, ..., n), record the above target row indexes in sequence, assuming that id s The position indexes corresponding to all sub-cores are: u1, u4, u7…;

[0181] S35, for all the above-mentioned mother core position indexes index = id1, id2, ... id s And the position indexes u1,u2,u5…,u1,u3,u4…,u1,u4,u7… of all its sub-cores are removed after duplicate operation to obtain m position indexes. For the above assumption, the repeated position indexes u1,u4 are merged, that is, index_new = u1,u2,u3,u4,u5,u7…,id1,id2,…id s , which is the position index of all nuclides after decay (including the parent and all daughter nuclides, (assuming the number of index_new is m)).

[0182] In this embodiment, Am-241 is taken as an example, that is, it is assumed that the source nuclide (i.e., target nuclide) input into the three-level PSA analysis atmospheric diffusion model is Am-241 (i.e., s=1), and its initial activity is 1E+10Bq. A total of 14 position indexes corresponding to the parent nucleus Am-241 and all its daughters are obtained, i.e., m=14, as shown in Table 3.

[0183] Table 3 Position index corresponding to Am-241 and all its descendants

[0184] Nuclide Position Index Am-241 1194 Np-237 1200 Pa-233 1207 U-233 1211 Th-229 1215 Ra-225 1218 Ac-225 1222 Fr-221 1229 At-217 1233 Bi-213 1236 Tl-209 1240 Po-213 1247 Pb-209 1251 Bi-209 1497

[0185] S4, based on the position index obtained in S3, adds the coefficient matrices C and C -1 In the above example, we can obtain all the elements whose rows and columns are index_new in sequence to form new coefficient matrices D and D -1 , D and D -1 A matrix of size mxm.

[0186] In this embodiment, the coefficient matrix D of size 14x14 is as follows Figure 5 As shown, D -1 like Figure 6 As shown, both are lower triangular matrices.

[0187] S5, determine Λ d The value of Λ d is a diagonal matrix of size mxm, whose diagonal elements are the negative decay constants of the nuclides corresponding to the position index index_new, respectively.

[0188]

[0189] In this embodiment, a 14x14 diagonal matrix is ​​as follows: Figure 7 shown.

[0190] S6, if the initial concentration field of the target nuclide is N(0), based on the cupy library in parallel calculation on the GPU, assuming that the decay time is t, the nuclide concentration field after decay is Here, N(0) can be a one-dimensional, two-dimensional or three-dimensional concentration field.

[0191] In this embodiment, a 1-hour decay is considered, that is, t=3600s, and the calculation results of N(t) are shown in Table 4.

[0192] Table 4 Nuclide concentration calculation results

[0193] Nuclide N(0) N(t) Am-241 1.00E+10 1.00E+10 Np-237 0 6.15E-02 Pa-233 0 5.49E-06 U-233 0 3.38E-10 Th-229 0 0.00E+00 Ra-225 0 7.42E-10 Ac-225 0 6.61E-11 Fr-221 0 4.54E-10 At-217 0 7.56E-11 Bi-213 0 2.12E-10 Tl-209 0 8.05E-12 Po-213 0 9.44E-11 Pb-209 0 0.00E+00 Bi-209 0 0.00E+00

[0194] Based on the Python cupy library, parallel computing is performed on the GPU, with a total of 10000x10000 calculations and a total computing time of 28.388 seconds. The computing time of this embodiment is much shorter than that of the prior art.

[0195] The calculation method of radionuclide decay provided in this embodiment proposes a set of differential equations for calculating radionuclide decay chains, calculates radionuclide decay based on matrix analysis, and is coupled with the Python-based cupy calculation library, which can assign calculation tasks to GPU for acceleration, greatly improving the calculation speed of nuclide decay and saving a lot of time for subsequent three-level PSA analysis of nuclear power plants. The calculation method of this embodiment is not affected and restricted by different calculation modes such as atmospheric diffusion model and dose calculation model in the three-level PSA analysis of nuclear power plants. It has universality and can be directly embedded in the three-level PSA analysis platform of nuclear power plants for nuclear power plant accident consequence analysis and public radiation risk assessment.

[0196] Embodiment 3:

[0197] like Figure 8 As shown, this embodiment provides a device for calculating radionuclide decay, including:

[0198] The first acquisition module 31 is used to obtain a target position index set in a preset coefficient matrix C, wherein the target position index set includes the position index of the target nuclide and all its daughters, the number of position indexes is m, m is a positive integer, and the coefficient matrix C refers to the coefficient matrix associated with the complete decay chain of all nuclides.

[0199] The arrangement module 32 is connected to the first acquisition module 31 and is used to arrange the m position indexes in the target position index set in pairs with repetition to obtain the row index and column index where the target element is located.

[0200] The second acquisition module 33 is connected to the arrangement module 32 and is used to respectively select the row index and column index of the target element in the preset coefficient matrix C and C -1 Get the target element from the matrix and form the coefficient matrices D and D -1 , where C -1 is the inverse matrix of C, D -1 is the inverse matrix of D.

[0201] The determination module 34 is connected to the first acquisition module 31 and is used to determine the negative decay constants of the nuclides corresponding to the m position indexes as diagonal elements of the diagonal matrix.

[0202] The calculation module 35 is connected to the second acquisition module 33 and the determination module 34, and is used to calculate the diagonal matrix, the coefficient matrix D and D -1 , calculate the nuclide concentration field of the target nuclide and all its daughters at the decay time.

[0203] Optionally, the preset coefficient matrix C satisfies formula (1):

[0204]

[0205] Among them, c kj represents the kth row and jth column element of the coefficient matrix C, Λ ik Represents the i-th row and k-th column element of the lower triangular matrix Λ.

[0206] Default C -1 Satisfying formula (2):

[0207]

[0208] Among them, c ik Represents C -1 The element in the i-th row and k-th column of kj Represents C -1 The k-th row and j-th column element of .

[0209] The lower triangular matrix Λ satisfies equation (3):

[0210]

[0211] Among them, λ i and λ j are the decay constants of nuclide i and nuclide j, respectively, and b ji Represents the branching ratio of nuclide i and nuclide j.

[0212] Optionally, the computing device further includes a construction module. The construction module is used to construct coefficient matrices C and C related to the complete decay chain of n radionuclides according to the decay constants and branching ratios of the n radionuclides and the position indexes of the daughter nuclei are all after the position index of the parent nucleus. -1 , get the preset coefficient matrices C and C -1 , where n is a positive integer and the n radioactive nuclides include the target nuclide.

[0213] Optionally, the first acquisition module is used to determine the position index corresponding to the target nuclide in the preset coefficient matrix C, and obtain the position index corresponding to all sub-bodies of the target nuclide in the preset coefficient matrix C, and is also used to deduplicate the position index of the target nuclide and the position index of all sub-bodies to obtain a target position index set.

[0214] Optionally, the first acquisition module is used to filter out the target row index from the number of rows greater than the target column index in the target column index where the position index of each target nuclide is located, wherein the element value in the coefficient matrix C corresponding to the target row index is non-zero, and is also used to mark the target row index as the corresponding position index of the sub-body of the target nuclide in the preset coefficient matrix C.

[0215] Optionally, the calculation module is used to calculate the nuclide concentration field according to formula (4):

[0216]

[0217] Where N(0) is the initial concentration of the target nuclide, Satisfying formula (5):

[0218]

[0219] Λ d is a diagonal matrix.

[0220] Optionally, the calculation module is used to calculate the nuclide concentration field in parallel based on a graphics processing unit GPU.

[0221] Embodiment 4:

[0222] This embodiment provides a device for calculating radionuclide decay, including: a memory and a processor, the memory and the processor are communicatively connected to each other, the memory stores computer instructions, and the processor executes the method for calculating radionuclide decay described in Example 1 or Example 2 by executing the computer instructions.

[0223] It is to be understood that the above embodiments are merely exemplary embodiments used to illustrate the principles of the present invention, but the present invention is not limited thereto. For those of ordinary skill in the art, various modifications and improvements can be made without departing from the spirit and essence of the present invention, and these modifications and improvements are also considered to be within the scope of protection of the present invention.

Claims

1. A method for calculating radionuclide decay, characterized in that: include: Obtain a target position index set in a preset coefficient matrix C, wherein the target position index set includes position indexes of the target nuclide and all its daughters, the number of position indexes is m, m is a positive integer, and the coefficient matrix C refers to a coefficient matrix associated with the complete decay chains of all nuclides; Arrange the m position indexes in the target position index set in pairs with repetition to obtain the row index and column index of the target element; According to the row index and column index of the target element, the preset coefficient matrices C and C -1 Get the target element from the matrix and form the coefficient matrices D and D -1 , where C -1 is the inverse matrix of C, D -1 is the inverse matrix of D; The negative decay constants of the nuclides corresponding to the m position indices are determined as the diagonal elements of the diagonal matrix; Based on the diagonal matrix, coefficient matrix D and D -1 , calculate the nuclide concentration field of the target nuclide and all its daughters at the decay time.

2. The method according to claim 1, characterized in that Before obtaining the target position index set in the preset coefficient matrix C, the method further includes: According to the decay constants and branching ratios of n radionuclides, and the position index of the daughter nucleus is after the position index of the parent nucleus, the coefficient matrices C and C related to the complete decay chain of n radionuclides are constructed. -1 , get the preset coefficient matrices C and C -1 , where n is a positive integer and the n radioactive nuclides include the target nuclide.

3. The method according to claim 1, characterized in that The preset coefficient matrix C satisfies formula (1): Among them, c kj represents the kth row and jth column element of the coefficient matrix C, Λ ik represents the i-th row and k-th column element of the lower triangular matrix Λ, The preset C -1 Satisfying formula (2): Among them, c ik Represents C -1 The element in the i-th row and k-th column of kj Represents C -1 The k-th row and j-th column element of The lower triangular matrix Λ satisfies equation (3): Among them, λ i and λ j are the decay constants of nuclide i and nuclide j, respectively, and b ji Represents the branching ratio of nuclide i and nuclide j.

4. The method according to claim 1, characterized in that: The step of obtaining a target position index set in a preset coefficient matrix C specifically includes: Determine the position index corresponding to the target nuclide in the preset coefficient matrix C; Obtain the corresponding position indexes of all daughter bodies of the target nuclide in the preset coefficient matrix C; The position index of the target nuclide and the position indexes of all its daughter bodies are deduplicated to obtain a target position index set.

5. The method according to claim 4, characterized in that The step of obtaining the position indexes corresponding to all daughter bodies of the target nuclide in the preset coefficient matrix C specifically includes: For each target nuclide position index in the target column index, select the target row index from the number of rows greater than the target column index, where the element value in the coefficient matrix C corresponding to the target row index is non-zero; The target row index is marked as the corresponding position index of the sub-body of the target nuclide in the preset coefficient matrix C.

6. The method according to claim 1, characterized in that The diagonal matrix, coefficient matrix D and D -1 , calculate the nuclide concentration field of the target nuclide and all its daughters at the decay time, including: The nuclide concentration field is calculated according to formula (4): Where N(0) is the initial concentration of the target nuclide, Satisfying formula (5): Λ d is a diagonal matrix.

7. The method according to claim 1, characterized in that The calculation of the nuclide concentration field specifically includes: The nuclide concentration field is calculated in parallel based on the graphics processing unit GPU.

8. A device for calculating radionuclide decay, characterized in that: include: The first acquisition module is used to obtain a target position index set in a preset coefficient matrix C, wherein the target position index set includes position indexes of the target nuclide and all its daughters, the number of position indexes is m, m is a positive integer, and the coefficient matrix C refers to a coefficient matrix associated with the complete decay chain of all nuclides, An arrangement module, connected to the first acquisition module, is used to arrange the m position indexes in the target position index set in pairs with repetition to obtain the row index and column index of the target element. The second acquisition module is connected to the arrangement module and is used to select the target element in the preset coefficient matrix C and C according to the row index and column index of the target element. -1 Get the target element from the matrix and form the coefficient matrices D and D -1 , where C -1 is the inverse matrix of C, D -1 is the inverse matrix of D, A determination module, connected to the first acquisition module, is used to determine the negative decay constants of the nuclides corresponding to the m position indexes as diagonal elements of the diagonal matrix, A calculation module is connected to the second acquisition module and the determination module, and is used to calculate the diagonal matrix, the coefficient matrix D and D -1 , calculate the nuclide concentration field of the target nuclide and all its daughters at the decay time.

9. The device according to claim 8, characterized in that The preset coefficient matrix C satisfies formula (1): Among them, c kj represents the kth row and jth column element of the coefficient matrix C, Λ ik represents the i-th row and k-th column element of the lower triangular matrix Λ, The preset C -1 Satisfying formula (2): Among them, c ik Represents C -1 The element in the i-th row and k-th column of kj Represents C -1 The k-th row and j-th column element of The lower triangular matrix Λ satisfies equation (3): Among them, λ i and λ j are the decay constants of nuclide i and nuclide j, respectively, and b ji Represents the branching ratio of nuclide i and nuclide j.

10. A device for calculating radionuclide decay, characterized in that: include: A memory and a processor, wherein the memory and the processor are communicatively connected to each other, the memory stores computer instructions, and the processor executes the method for calculating the decay of a radionuclide according to any one of claims 1 to 7 by executing the computer instructions.

Images