A fast calculation method of beam transport in degrader section based on transport matrix and Monte Carlo method

Abstract

This invention provides a fast calculation method for beam transport in a de-energizer section based on the transfer matrix and Monte Carlo method, relating to the field of beam physics. The method includes: Step 1, constructing a transfer matrix model of the beam transport section before the de-energizer based on the beam optical parameters and optical element parameters at the starting point of the beam transport line, obtaining the particle phase space distribution parameters at the de-energizer inlet; Step 2, establishing a three-dimensional structural model of the de-energizer based on the particle phase space distribution parameters at the de-energizer inlet and performing Monte Carlo beam transport simulation of the de-energizer section, obtaining the position, momentum, and energy information of the beam particles at the de-energizer outlet; Step 3, performing outlier data removal processing on the position, momentum, and energy information of the beam particles at the de-energizer outlet, obtaining processed particle data. This invention achieves high-precision reconstruction of Monte Carlo particle information into optical parameters, improving the simulation efficiency of beam transport lines containing de-energizer structures.

Description

Technical Field

[0001] This invention relates to the field of beam physics, and in particular to a fast calculation method for beam transmission in a degrader segment based on the transfer matrix and the Monte Carlo method. Background Technology

[0002] In the design of beam transmission lines for hadron therapy devices, the transfer matrix method is often used for preliminary simulation and parameter optimization of the beam transport process due to its high computational efficiency and relatively simple form. This method models linear optical elements such as drift sections and magnets as corresponding transfer matrices, and combined with the optical parameters of the initial beam, it can quickly estimate the beam behavior of the linear transmission section. However, when the beam transmission line contains structures with nonlinear transmission characteristics, such as de-energizers, the traditional transfer matrix method is not applicable to its calculation process. Since the de-energizer involves complex interactions between particles and materials, its beam transport process mostly exhibits nonlinear characteristics, making it difficult to accurately describe the behavior of the beam in the de-energizer using only the linear transfer matrix method.

[0003] To handle such nonlinear transmission processes, the Monte Carlo method is often used in existing technologies to simulate the de-energizer section in order to obtain particle phase space information at the exit. However, there are usually format differences between the particle data output by this method and the optical parameters required by the transfer matrix method, resulting in a less than smooth connection between the two methods. For example, in the beam transmission line design of a certain type of proton therapy system, although the Monte Carlo method can obtain relatively detailed particle motion information when used for system-level simulation, the computational resources consumed are relatively large, which may face efficiency challenges in the rapid iterative optimization at the beginning of the engineering design. Summary of the Invention

[0004] The technical problem to be solved by this invention is to provide a fast calculation method for beam transmission in depressor sections based on the transmission matrix and the Monte Carlo method, so as to realize high-precision reconstruction of Monte Carlo particle information into optical parameters and improve the simulation efficiency of beam transmission lines with depressor structures.

[0005] To solve the above-mentioned technical problems, the technical solution of the present invention is as follows:

[0006] Firstly, a fast calculation method for beam transport in a de-energizer segment based on a transfer matrix and the Monte Carlo method is provided, the method comprising:

[0007] Step 1: Construct a transmission matrix model of the beam transmission section in front of the de-energizer based on the beam optical parameters and optical element parameters at the starting point of the beam transmission line, and obtain the particle phase space distribution parameters at the inlet position of the de-energizer;

[0008] Step 2: Based on the particle phase space distribution parameters at the inlet of the de-energizer, establish a three-dimensional structural model of the de-energizer and perform Monte Carlo beam transport simulation of the de-energizer section to obtain the position, momentum and energy information of the beam particles at the outlet of the de-energizer.

[0009] Step 3: Perform abnormal data removal processing on the position, momentum, and energy information of the beam particles at the outlet of the de-energizer to obtain the processed particle data;

[0010] Step 4: Based on the processed particle data, calculate the mean square position, mean square divergence angle and covariance term of the particles in the beam exiting the de-energizer, and reconstruct the optical energy spectrum parameters of the beam exiting the de-energizer and the parameters of the optical elements behind the de-energizer.

[0011] Step 5: Based on the optical energy spectrum parameters of the beam exiting the de-energizer and the parameters of the optical elements after the de-energizer, construct the transmission matrix model of the beam transmission line after the de-energizer to realize the rapid calculation of the beam transmission process after the de-energizer.

[0012] In a second aspect, a computing device includes:

[0013] One or more processors;

[0014] A storage device for storing one or more programs that, when executed by one or more processors, cause the one or more processors to implement the method.

[0015] Thirdly, a computer-readable storage medium storing a program that, when executed by a processor, implements the method.

[0016] The above-described solution of the present invention has at least the following beneficial effects:

[0017] By eliminating outlier data and reconstructing parameters, the particle position, momentum, and energy information output from Monte Carlo simulations are accurately converted into beam optical parameters and energy spectrum parameters. This reduces the incompatibility issues between the traditional Monte Carlo method and the transfer matrix method in terms of data format. It retains the accuracy of Monte Carlo methods in simulating nonlinear transport processes while leveraging the efficiency of the transfer matrix method, thus enhancing the synergistic capabilities of the two methods. Innovatively, the beam transmission process before and after the depressor is decoupled. The front and rear transmission sections are calculated quickly using the transfer matrix method, while the depressor section is accurately simulated using the Monte Carlo method. This reduces computational resource consumption and time costs, meeting the needs of frequent iterative optimization in beam transmission line design while balancing computational accuracy and efficiency. Attached Figure Description

[0018] Figure 1 This is a flowchart illustrating a fast calculation method for beam transmission in a degrader segment based on a transmission matrix and the Monte Carlo method, provided by an embodiment of the present invention.

[0019] Figure 2 This is a schematic diagram illustrating how an embodiment of the present invention constructs a transmission matrix model of the beam transmission line after the de-energizer based on beam optical parameters and energy spectrum parameters, thereby enabling rapid calculation of the beam transmission process after the de-energizer. Detailed Implementation

[0020] Exemplary embodiments of the present disclosure will now be described in more detail with reference to the accompanying drawings. While exemplary embodiments of the present disclosure are shown in the drawings, it should be understood that the present disclosure may be implemented in various forms and should not be limited to the embodiments set forth herein. Rather, these embodiments are provided so that this disclosure will be thorough and complete, and will fully convey the scope of the disclosure to those skilled in the art.

[0021] like Figure 1 As shown, an embodiment of the present invention proposes a fast calculation method for beam current transmission in a degrader segment based on the transmission matrix and the Monte Carlo method. The method includes the following steps:

[0022] Step 1: Construct a transmission matrix model of the beam transmission section in front of the de-energizer based on the beam optical parameters and optical element parameters at the starting point of the beam transmission line, and obtain the particle phase space distribution parameters at the inlet position of the de-energizer;

[0023] Step 2: Based on the particle phase space distribution parameters at the inlet of the de-energizer, establish a three-dimensional structural model of the de-energizer and perform Monte Carlo beam transport simulation of the de-energizer section to obtain the position, momentum and energy information of the beam particles at the outlet of the de-energizer.

[0024] Step 3: Perform abnormal data removal processing on the position, momentum, and energy information of the beam particles at the outlet of the de-energizer to obtain the processed particle data;

[0025] Step 4: Based on the processed particle data, calculate the mean square position, mean square divergence angle and covariance term of the particles in the beam exiting the de-energizer, and reconstruct the optical energy spectrum parameters of the beam exiting the de-energizer and the parameters of the optical elements behind the de-energizer.

[0026] Step 5: Based on the optical energy spectrum parameters of the beam exiting the de-energizer and the parameters of the optical elements after the de-energizer, construct the transmission matrix model of the beam transmission line after the de-energizer to realize the rapid calculation of the beam transmission process after the de-energizer.

[0027] In this embodiment of the invention, by eliminating outlier data and reconstructing parameters, the particle position, momentum, and energy information output by Monte Carlo simulation are accurately converted into beam optical parameters and energy spectrum parameters. This reduces the incompatibility between the data formats of the traditional Monte Carlo method and the transfer matrix method, preserving the accuracy of the Monte Carlo method in simulating nonlinear transport processes while leveraging the efficiency of the transfer matrix method, thus enhancing the synergistic capability of the two methods. Innovatively, the beam transmission process before and after the depressor is decoupled. The before and after transmission sections are rapidly calculated using the transfer matrix method, while the depressor section is accurately simulated using the Monte Carlo method. This reduces computational resource consumption and time costs, meeting the needs of frequent iterative optimization in beam transmission line design while balancing computational accuracy and efficiency.

[0028] In a preferred embodiment of the present invention, step 1 above, which involves constructing a transmission matrix model of the beam transmission section before the de-energizer based on the beam optical parameters at the beam transmission line initiation point and the optical element parameters, to obtain the particle phase space distribution parameters at the de-energizer inlet position, may include:

[0029] In this embodiment of the invention, step 110 involves obtaining the initial Twiss parameters at the starting point of the beam transmission line as input parameters. Specifically, the initial Twiss parameters at the starting point of the beam transmission line are core parameters describing the optical characteristics of the beam at the starting point and are the basic input for constructing the transmission matrix model. These parameters need to be obtained from the initial design parameters of the beam transmission line or the measurement data from previous beam debugging to ensure consistency with the actual initial state of the beam. The initial Twiss parameters include three sets of parameters in the horizontal direction (x-direction) and the vertical direction (y-direction), respectively... , , and , , in For Twiss Parameters that describe the inclination of the beam phase space ellipse; For Twiss Parameters characterize the spatial scale of the beam envelope; For Twiss The parameters are related to the beam divergence angle, and the three sets of parameters satisfy a fixed relationship. This relationship is derived from the geometric properties of the beam phase space ellipse and is used to verify the consistency of the parameters.

[0030] Step 111: Determine the transmission matrix corresponding to each optical element in the beam transmission section before the de-energizer, based on the type and structural parameters of each optical element. Specifically, the optical elements in the beam transmission section before the de-energizer are all linear elements, such as the drift section, diodes, and quadrupoles. Their transmission matrices need to be derived based on the element type and structural parameters to accurately describe the linear transport characteristics of the beam in the elements. The transmission matrix of the drift section is determined. The drift section is a free space without electromagnetic field, where particles move in uniform linear motion, their positions change linearly with the drift length, and the divergence angle remains constant. Let the length of the drift section be L, and its transmission matrix... The form is Substituting the actual drift segment length parameter L into the matrix, the transmission matrix of that drift segment can be determined. The transmission matrix of the dipole magnet is determined because it is used to deflect the beam direction; its transmission matrix needs to distinguish between the deflection directions, horizontal and vertical (there are dipole magnets that deflect vertically). The structural parameters to be obtained include the particle's radius of motion in the magnetic field. Magnetic stiffness related items (Constant related to particle momentum), deflection angle The transfer matrix in the x-direction (deflection plane) for ,in and From deflection angle Calculations show that Curvature that reflects the circular motion of particles The properties of specific particles with specific momentum, along with matrix elements, collectively describe the coupling relationship between the position and divergence angle of the beam during deflection, the deflection direction, and the transmission matrices in the horizontal and vertical directions (in the presence of a vertically deflecting dipole magnet). for ,in For a specific particle with a specific momentum in the y-direction, the matrix describes the transport characteristics of the beam in the vertical direction. Substituting the above structural parameters into the corresponding matrix, the transport matrices of the dipole magnet in both directions can be determined. The transport matrix of the quadrupole magnet is determined; since the quadrupole magnet is used to focus or defocus the beam, the focusing intensity needs to be obtained. Related to the magnetic field gradient and particle momentum, and the magnet length l, the focusing plane, such as the transfer matrix in the x-direction. for ,in and Depend on The product of l and y is used to calculate the focusing characteristics of the beam under the action of a quadrupole magnet, the defocusing plane, such as the transmission moment in the y direction. for ,in and Also by The product of l and is used to calculate the defocusing characteristics of the beam under the action of a quadrupole magnet. By substituting l into the corresponding matrix, the transmission matrix of the quadrupole magnet in different directions can be determined. Through the above process, the type and structural parameters of each linear optical element are transformed into the corresponding transmission matrix, providing a basis for matrix cascading.

[0031] Step 112: Concatenate the transmission matrices of each optical element according to their physical order in the beam transmission line to construct a transmission matrix model for the beam transmission section before the de-energizer. Specifically, the beam's movement in the transmission section before the de-energizer follows a fixed physical path, and the concatenation order of the transmission matrices must be consistent with this path to ensure the transmission matrix model accurately reflects the actual beam transport process. First, the physical order in which the beam passes through each optical element must be determined. For example, starting from the beginning of the transmission line, it sequentially passes through drift section A, quadrupole magnet B, dipole magnet C, drift section D, and finally reaches the de-energizer inlet. This order must be confirmed through the engineering drawings or design documents of the beam transmission line to avoid distortion of the transmission law description due to incorrect order. The concatenation of the transmission matrices follows the matrix multiplication rule, and the multiplication order is the reverse of the physical order in which the beam passes through the elements, i.e., the matrix of the element passed later is multiplied by the matrix of the element passed earlier. Let the beam start from the beginning... Starting from the first component, the position after passing through is... The corresponding transmission matrix is (describe arrive (Transmission); the position after passing through the second element is The corresponding transmission matrix is (describe arrive (transmission); and so on, reaching the energy degrader inlet after passing through the nth element. The corresponding transmission matrix is (describe arrive (Transmission), the cascading process starts from the first element, first calculating the cascading matrix of the first two elements, This matrix describes the beam from After passing through the first two components, it arrives The transmission pattern; then... Matrix with the third element Multiply, we get Description from After passing through the first three components, it arrives The transmission pattern is followed; this process is repeated until the matrix of all n elements is cascaded, ultimately obtaining the transmission section before the de-energizer from the starting point. to the depressor inlet Total transmission matrix The total transmission matrix is ​​the transmission matrix model of the beam transmission section in front of the depressor. Its elements completely include the effects of all linear optical elements on the beam position and divergence angle, and the beam transport process in this section can be quickly calculated through matrix operations.

[0032] Step 113: Input the initial Twiss parameters into the transmission matrix model of the beam transmission section before the de-energizer for calculation to obtain the particle phase space distribution parameters at the de-energizer inlet position; specifically, the initial Twiss parameters need to be converted into the beam at the start of the transmission line. The phase space vector, which includes the horizontal position x and the horizontal divergence angle. Vertical position y, vertical divergence angle ,Right now ,in and The unit is meters. and This is a dimensionless angular quantity (radians), and the vector is related to the initial Twiss parameters through the beam phase space elliptic equation, for example, the horizontal direction satisfies... ( (For horizontal emittance), ensuring that the phase space vector accurately reflects the beam state described by the initial Twiss parameters, and using the phase space vector... Substitute the total transmission matrix obtained in step 112 The energy depressor inlet is calculated using matrix multiplication. phase space vector ,Right now = × The calculation needs to be performed separately for the horizontal and vertical directions. The horizontal submatrix is... Let the top left 2×2 part be its elements. The horizontal position of the entrance Horizontal divergence angle at the entrance = The vertical submatrix is Let the bottom right 2×2 part have the following elements. The vertical position of the entrance Vertical divergence angle at the entrance Calculated =[ , , , ]ᵀ is the particle phase space distribution parameter at the depressor inlet position. This parameter fully reflects the position and divergence angle distribution of the beam when it arrives at the depressor inlet, providing accurate initial beam conditions for Monte Carlo simulation of the depressor section.

[0033] By standardizing the entire process of initial Twiss parameter acquisition, optical element transfer matrix determination, matrix cascading, and phase space parameter calculation, the efficiency advantage of the transfer matrix method in processing linear elements is fully utilized, avoiding inefficiency or error accumulation caused by non-standard calculation processes.

[0034] In a preferred embodiment of the present invention, step 2 above, which involves establishing a three-dimensional structural model of the de-energizer based on the particle phase space distribution parameters at the de-energizer inlet and performing Monte Carlo beam transport simulation of the de-energizer section to obtain the position, momentum, and energy information of the beam particles at the de-energizer outlet, may include:

[0035] In this embodiment of the invention, step 220 involves using the particle phase space distribution parameters at the inlet of the de-energizer as the initial beam condition. Specifically, this includes: the particle phase space distribution parameters at the inlet of the de-energizer are the core initial input for Monte Carlo beam transport simulation, and the complete composition of these parameters needs to be determined first, including the horizontal position. Horizontal divergence angle Vertical position Vertical divergence angle These parameters directly reflect the spatial distribution and motion trend of the beam when it reaches the de-energizer inlet. First, these parameters need to be organized into a data format recognizable by the Monte Carlo simulation program, and then the inlet state of each particle should be sequentially associated with its particle number, i.e., each particle corresponds to a set of (…). , , , The parameters are then used to derive the initial momentum components of the particles through the relationship between beam momentum and divergence angle. Let the initial momentum of the beam along the transmission line axis (z-direction) be... Determined by the initial beam energy and particle mass, satisfying = ,in The total energy of the inlet beam. The rest mass of the particle. Given the speed of light in a vacuum, the horizontal momentum component... = × Vertical momentum component = × momentum component in the z-direction = (Since the beam mainly propagates along the z-direction and the divergence angle is small, the momentum in the z-direction is approximately constant.) Through this derivation process, the phase space distribution parameters are transformed into the initial position and momentum information of the particles required for Monte Carlo simulation, forming complete initial beam conditions.

[0036] Step 221: Based on the initial beam conditions, establish a complete three-dimensional geometric model of the energy degrader. Specifically, the initial beam conditions determine the key dimensions and structural layout of the energy degrader's three-dimensional geometric model. It is necessary to ensure that the model completely covers the beam's transmission path within the energy degrader to avoid simulated particle leakage due to missing geometric boundaries. First, determine the constituent components of the energy degrader's three-dimensional geometric model, including the wedge-shaped energy degrading device, vacuum cavity, and mechanical support structure. The structural parameters of each component need to be determined in conjunction with beam transmission requirements and engineering design specifications. The core parameter of the wedge-shaped scattering device is the wedge angle. The material thickness d needs to be designed according to the beam energy modulation requirements to ensure that the beam can reach the target energy range after scattering; the size of the vacuum cavity needs to be expanded according to the outline of the wedge scattering device to reserve sufficient beam transmission space, and the spacing between the inner wall and the scattering device needs to be set to avoid particle collisions with the cavity wall from interfering with the simulation results; the mechanical support structure needs to be designed around the vacuum cavity and the scattering device, and its position needs to avoid the beam transmission path to ensure that it does not affect the beam transport.

[0037] When constructing the geometric structure of each component, to ensure the accuracy and closure of the structural boundaries, a two-dimensional graphic of each component is constructed, and then stretched along the z-direction to form a three-dimensional structure. The specific implementation process is as follows: First, obtain the two-dimensional contour vertex coordinate set of a single component, such as the wedge-shaped scattering device. The coordinates are based on the center of the energy degrader inlet as the origin, with the x-axis representing the horizontal direction and the y-axis representing the vertical direction. The coordinates of each vertex are represented as (…). , () (where the vertex index is used), sort the vertex coordinate set in ascending order of x-coordinate, and if the x-coordinates are the same, sort them in ascending order of y-coordinate to obtain the sorted vertex sequence. Then construct the lower convex hull. Starting from the starting point of the sorted vertex sequence, traverse the vertices one by one. If the cross product of the vectors formed by the current three vertices is greater than or equal to zero (i.e., the vector...), then... with vector cross product Then remove the middle vertex. The process continues until the lower convex hull vertex sequence is obtained after traversal. Then, the upper convex hull is constructed. Starting from the end point of the sorted vertex sequence, the vertices are traversed in reverse, and the same cross product judgment logic as for the lower convex hull is executed to obtain the upper convex hull vertex sequence. The lower and upper convex hull vertex sequences are merged, and duplicate vertices are removed to obtain the two-dimensional convex hull contour of the component. Finally, the two-dimensional convex hull contour is stretched along the z-direction to the designed thickness to form a closed three-dimensional geometric structure. The stretching thickness is determined according to the functional requirements of each component. After completing the three-dimensional modeling of the wedge-shaped scattering device, vacuum cavity, and mechanical support structure in sequence according to the above method, the components are combined according to the actual assembly relationship to ensure that the positions of each component are accurately aligned. The vacuum cavity completely encloses the wedge-shaped scattering device, and the mechanical support structure is symmetrically distributed on the outside of the vacuum cavity without obstructing the beam transmission channel, thus forming a complete three-dimensional geometric model of the energy degrader.

[0038] Step 222 involves converting the 3D geometric model of the energy degrader into an input file for constructing a solid geometry format, obtaining a geometric description suitable for Monte Carlo simulation. Specifically, this includes: first, exporting the 3D geometric model of the energy degrader constructed in Step 221 as a STEP format file; analyzing the topology of the 3D geometric model to extract the boundary representation information of each component, including the geometric features and relationships of faces, edges, and vertices; then, decomposing the complex solid structure of each component into basic CSG fundamental surfaces recognizable by the Monte Carlo simulation program. Commonly used fundamental surfaces include planes, cylindrical surfaces, and spheres. The decomposition process must follow the principle of breaking down complex structures into simpler sub-bodies, for example, decomposing the cuboid of the vacuum cavity... The structure is decomposed into a closed region enclosed by six planes. The cylindrical support structure is further decomposed into a solid composed of one cylindrical surface and two circular planes. For the decomposed basic sub-body, Boolean logic operations (intersection, union, complement) are used to define their combination relationships to restore the geometry of the original components. For example, when constructing a wedge-shaped scattering device, the multiple planes constituting the wedge are first combined using the union operation to form the basic outline, and then the area overlapping with the vacuum cavity is removed using the complement operation. When constructing the geometry of the overall energy degrader, the intersection operation is used to ensure the accurate positional relationship of each component and avoid geometric interference. During the conversion process, a unique identifier must be assigned to each basic surface and the combined body, and its geometric parameters, such as the equation parameters of the planes, must be recorded. ,in , , For plane normal vector components, The constant term, the equation of the axis of the cylindrical surface, and the radius are given. The axis is a straight line in space, and the radius is the cross-sectional radius of the cylinder. These parameters will be directly written into the CSG format input file. The final generated CSG format file must contain the geometric description of all components, Boolean operation logic, and identifier information to ensure that the Monte Carlo simulation program can accurately identify the three-dimensional geometry of the energy reducer.

[0039] Step 223: Using the initial beam conditions and the input file containing the constructed solid geometry as input parameters, configure the running parameters of the Monte Carlo simulation program. Specifically, this includes: first, importing the initial beam condition data compiled in Step 220 into the input card of the Monte Carlo simulation program, determining the type of particle source and the corresponding beam particle types, such as protons, heavy ions, and the spatial and momentum distribution of particle emission, ensuring that the particle source parameters completely match the initial beam conditions; then, importing the CSG format input file generated in Step 222 as the geometric input parameter into the program, specifying the reading path and parsing method of the geometric file, ensuring that the program can correctly load the three-dimensional geometric model of the de-energizer; based on this, configure the core parameters for simulation operation, one of which is the number of particles. The simulation is determined based on several factors: First, the statistical accuracy requirements are set to ensure the simulation results are statistically significant and to avoid excessive data fluctuations due to insufficient particle count. Second, the physical process mechanism is selected, including the elastic scattering mechanism, inelastic scattering model, and energy loss mechanism of the particles, to ensure accurate simulation of particle scattering and energy loss processes within the depressor. Third, simulation accuracy control parameters are set, including an upper limit for the particle tracking step size to avoid missed collisions and energy cutoff thresholds caused by excessively large step sizes; particles below the threshold are stopped from tracking, reducing computational resource consumption. Fourth, output parameters are set, specifying the types of particle information to be recorded, including the position, momentum, and energy of particles at the depressor exit, and determining the data output format and storage path.

[0040] Step 224: Execute the Monte Carlo simulation program to simulate the complete transport process of the beam in the three-dimensional structure of the de-energizer, obtaining the position, momentum, and energy information of each beam particle at the de-energizer exit. Specifically, this includes: starting the configured Monte Carlo simulation program; the program first reads the initial beam conditions from the input card, generates a particle source that meets the parameter settings, and each particle is emitted from the de-energizer inlet according to the preset position, angular distribution, and momentum distribution. Then, the program loads the CSG format geometry mechanism, tracks the trajectory of each particle within the de-energizer based on the geometric description, and calculates the interaction between the particle and the material's atomic nuclei and electrons, including changes in scattering direction, based on the configured physical process mechanism when the particle moves to the de-energizer material region. The program calculates the change in momentum direction (from the scattering angle) and energy loss (from the energy loss model to the kinetic energy decay). When a particle moves to the geometric boundary, the program determines whether the particle has passed through the de-energizer (reached the exit). If it has not passed through, the program continues to track its movement in the subsequent geometric region until the particle passes through the de-energizer or the energy is lower than the cutoff threshold. For particles that pass through the de-energizer exit, the program records the three-dimensional coordinates of the exit position, the three-dimensional components of momentum, and the total energy in real time. After the simulation program finishes running, the recorded information of all exit particles is stored in a file in a specified path according to a preset format. Each particle corresponds to one line of data, which includes position, momentum, and energy parameters in sequence, forming a complete dataset of de-energizer exit particle information.

[0041] By standardizing the model construction, format conversion, and parameter configuration process, the nonlinear transport environment of the beam within the de-energizer is realistically reproduced, avoiding simulation distortion caused by geometric errors or parameter mismatches. This achieves the standardization and completeness of Monte Carlo simulation input parameters. Through the transformation of initial beam conditions and the adaptation of geometric formats, a reliable input basis is provided for the simulation program, ensuring that the acquired exit particle information accurately reflects the state of the beam after passing through the de-energizer.

[0042] In a preferred embodiment of the present invention, step 3 above, which involves processing the position, momentum, and energy information of the beam particles at the outlet of the de-energizer to remove abnormal data and obtain processed particle data, may include:

[0043] In this embodiment of the invention, step 330 involves statistically analyzing the position, momentum, and energy information of each beam particle at the de-energizer outlet, calculating the average and standard deviation of energy, horizontal position, horizontal divergence angle, vertical position, and vertical divergence angle, respectively. Specifically, this includes: first, determining the statistical object as the raw data of all beam particles at the de-energizer outlet output by the Monte Carlo simulation. This data includes the horizontal position, horizontal divergence angle, vertical position, vertical divergence angle, and energy of each particle. Since the Monte Carlo simulation outputs a large amount of particle data to ensure statistical representativeness, the total number of particles participating in the statistics needs to be determined first, denoted as N; next, calculating the average and standard deviation of the five physical quantities, the average... The calculation method is to determine the horizontal positions of all N particles. , Indicates the first After summing the horizontal positions of each particle, divide by the total number of particles N, the formula is: Standard deviation The calculation method is to first calculate the horizontal position of each particle. Compared with the average The square of the difference, summed from all squares, divided by the total number of particles N, and then the arithmetic square root of the result is calculated using the following formula: ,in It reflects the degree of dispersion of all particles in the horizontal direction.

[0044] average value The calculation method is to calculate the horizontal divergence angle of all N particles. ( =1, 2, Indicates the first After summing the horizontal divergence angles of each particle and dividing by the total number of particles N, the formula is: Standard deviation The calculation method is to first calculate the horizontal divergence angle of each particle. Compared with the average The square of the difference, summed from all squares, divided by the total number of particles N, and then the arithmetic square root of the result is calculated using the following formula: ,in It reflects the degree of dispersion of all particles in the horizontal direction of motion.

[0045] average value The calculation method is to calculate the vertical positions of all N particles. Indicates the first After summing the vertical positions of each particle, divide by the total number of particles N, the formula is: Standard deviation The calculation method is to first calculate the vertical position of each particle. Compared with the average The square of the difference, summed from all squares, divided by the total number of particles N, and then the arithmetic square root of the result is calculated using the following formula: ,in It reflects the degree of dispersion of all particles in the vertical direction.

[0046] average value The calculation method is to calculate the vertical divergence angle of all N particles. Indicates the first After summing the vertical divergence angles of each particle and dividing by the total number of particles N, the formula is: Standard deviation The calculation method is to first calculate the vertical divergence angle of each particle. Compared with the average The square of the difference, summed from all squares, divided by the total number of particles N, and then the arithmetic square root of the result is calculated using the following formula: ,in It reflects the degree of dispersion of all particles in the vertical direction of motion.

[0047] average value The calculation method is to use the energy of all N particles. Indicates the first The formula is: (The sum of the energies of each particle) divided by the total number of particles N, The formula is: Standard deviation The calculation method is to first calculate the energy of each particle. Compared with the average The square of the difference, summed from all squares, divided by the total number of particles N, and then the arithmetic square root of the result is calculated using the following formula: ,in It reflects the degree of energy dispersion of all particles.

[0048] Step 331: Based on the mean and standard deviation, set the normal data range for each component; identify and remove abnormal particle data that exceed the corresponding normal data range in each component, and retain valid particle data within the normal data range; specifically, this includes: firstly, based on the mean and standard deviation of each physical quantity calculated in step 330, set the normal data range for each physical quantity; combining the conventional standards for statistical data processing in beam physics and the physical characteristics of particle transport within the de-energizer, uniformly set the normal data range as [μ−3σ, μ+3σ], where μ represents the mean of the corresponding physical quantity, used to determine the distribution center of the physical quantity, and σ represents the standard deviation of the corresponding physical quantity, used to measure the distribution width of the physical quantity; 3σ is chosen as the boundary because, under the normal distribution statistical law, the data within this range can cover approximately 99.73% of normal particles, maximizing the retention of statistically representative particle data, while excluding extreme outlier data caused by abnormal scattering or simulation errors.

[0049] Next, anomaly identification and removal are performed on each of the five physical quantity components of each particle. Anomaly identification and removal are performed on the horizontal position component. For the first... A particle, if its horizontal position satisfy or If the particle is found to be abnormal in its horizontal position component, the entire particle should be marked as an anomalous particle; if In If the particle is within the normal range, then its horizontal position component conforms to the normal range; the identification and elimination of anomalies in the horizontal divergence angle component, for the first... A particle, if its horizontal divergence angle is... satisfy or If the particle exhibits an anomaly in its horizontal divergence angular component, the entire particle must be marked as an anomalous particle; if In If the particle is within the normal range, then its horizontal divergence angle component conforms to the normal range; the identification and elimination of anomalies in the vertical position component, for the first... A particle, if its vertical position satisfy or If the particle exhibits an anomaly in its vertical position component, the entire particle must be marked as an anomalous particle; if In If the particle is within the normal range, then its vertical position component conforms to the normal range; the identification and elimination of anomalies in the vertical divergence angle component, for the first... A particle, if its vertical divergence angle is... satisfy or If the particle exhibits an anomaly in its vertical divergence angular component, the entire particle must be marked as an anomalous particle; if In Within the range, the particle conforms to the normal range in the vertical divergence angular component; the identification and elimination of abnormal energy components, for the first... A particle, if its energy satisfy or If the particle exhibits an anomaly in its energy component, then the entire particle must be marked as an anomalous particle; if In If the particle is within the specified range, then its energy component falls within the normal range.

[0050] After identifying component anomalies for all particles, all particle data marked as anomalous are removed from the original dataset, retaining only particle data where all components fall within the corresponding normal data range. These retained particle data are considered valid particle data. The key to this step is to comprehensively eliminate non-representative particles caused by extreme scattering within the de-energizer, such as abnormal deflection due to particle impacts on the mechanical support structure of the de-energizer, or statistical errors in Monte Carlo simulations, through multi-component collaborative screening. This avoids interference from such particles in the calculation of beam optical parameters and ensures that the reconstructed Twiss parameters accurately reflect the overall optical characteristics of the exit beam of the de-energizer.

[0051] Step 332 involves integrating the effective particle data to form processed particle data; specifically, this includes grouping the five physical quantity components of all effective particles according to their categories to form five independent effective data sequences, namely, effective horizontal position sequences. For the effective number of particles, , Indicates the first Sequence of horizontal position and effective horizontal divergence angle of each effective particle , Indicates the first The horizontal divergence angle and effective vertical position sequence of each effective particle , Indicates the first The sequence of vertical positions and effective vertical divergence angles of each effective particle. , Indicates the first Vertical divergence angle and effective energy sequence of each effective particle , Indicates the first The energy of an effective particle.

[0052] During the classification and aggregation process, it is necessary to strictly maintain the correlation between the five physical quantity components of each effective particle, that is, for the same effective particle (denoted as the first...). (one), its , , , , Must correspond to the same particle number To avoid deviations in subsequent covariance calculations due to data misalignment, since covariance calculations require paired operations based on the position and divergence angle data of the same particle, the five collected valid data sequences are standardized according to the conventional format of beam physics calculations. This includes standardizing the units of physical quantities and data precision to ensure consistency in units and accuracy during the calculation process. This prevents unit conversion errors or inconsistencies in precision from affecting the accuracy of the calculation results. The standardized five valid data sequences and the number of valid particles are then... They are collectively packaged into a processed particle dataset.

[0053] By first calculating the average value and standard deviation of each physical quantity, then removing outlier data within the range of [μ−3σ, μ+3σ], and finally integrating the valid data, non-representative particles can be eliminated to the greatest extent possible. This ensures that the processed particle data conforms to the normal physical laws of beam transport, and that the reconstructed Twiss parameters and energy spectrum parameters can truly reflect the actual optical characteristics of the beam exiting the de-energizer, thus avoiding calculation errors caused by interference from outlier data.

[0054] In a preferred embodiment of the present invention, step 4 above, which involves calculating the mean square position, mean square divergence angle, and covariance term of the particles in the de-energizer exit beam based on the processed particle data, and reconstructing the optical energy spectrum parameters of the de-energizer exit beam and the parameters of the optical elements behind the de-energizer, may include:

[0055] In this embodiment of the invention, step 440 involves calculating the mean square position, mean square divergence angle, and covariance term of position and divergence angle in the horizontal direction based on the processed particle data, and simultaneously calculating the mean square position, mean square divergence angle, and covariance term of position and divergence angle in the vertical direction. Specifically, this includes: first, determining the composition of the processed particle data, which is the effective particle dataset obtained in step 332, containing the number of effective particles. Effective particle horizontal position sequence ( Indicates the first (Horizontal position of each effective particle), horizontal divergence angle sequence ( Indicates the first The horizontal divergence angle of an effective particle, i.e., the tangent of the angle between the particle's direction of motion and the transmission line axis (dimensionless), and the vertical position sequence. ( Indicates the first Vertical position of each effective particle), vertical divergence angle sequence ( Indicates the first (Vertical divergence angle of an effective particle, dimensionless).

[0056] Calculate the three core statistics in the horizontal direction, and the mean square position in the horizontal direction. The average level of the squares of the positions of all effective particles in the horizontal direction is the basis for calculating optical parameters. It is calculated by taking the horizontal position of each effective particle... Perform a squaring operation, sum all the squared results, and divide by the effective particle number M. The formula is: (from =1 to M) ;in It is the first The square of the horizontal position of each effective particle. (from =1 to M) This represents the sum of the squares of the horizontal positions of all effective particles, which, when divided by M, yields the square of the average horizontal position.

[0057] Horizontal mean square divergence angle It reflects the average squared divergence angle of all effective particles in the horizontal direction, and is calculated by taking the horizontal divergence angle of each effective particle as the average value. Perform a squaring operation, sum all the squared results, and divide by the effective particle number M. The formula is: (from =1 to M) ;in It is the first The square of the horizontal divergence angle of each effective particle. (From j=1 to M) This represents the sum of the squares of the horizontal divergence angles of all effective particles, which, when divided by M, yields the square of the average horizontal divergence angle.

[0058] Covariance term of horizontal position and divergence angle This reflects the correlation between the horizontal particle position and the divergence angle, playing a crucial role in calculating geometric emittance. The calculation method involves setting the horizontal position of each effective particle... With horizontal divergence angle Multiply, sum all the products, and divide by the effective particle number M, using the following formula: (from =1 to M) ( ) / M; where It is the first The product of the horizontal position of each effective particle and its divergence angle. (from =1 to M) ( ) represents the sum of the products of all effective particles, which, when divided by M, yields the covariance.

[0059] Simultaneously calculate three core statistics in the vertical direction, following the same calculation logic as in the horizontal direction, with the vertical mean square position... The calculation method involves determining the vertical position of each effective particle. Perform a squaring operation, sum all the squared results, and divide by the effective particle number M. The formula is: (from =1 to M) ; It is the first The vertical position of each effective particle, this parameter reflects the average squared value of the particle position in the vertical direction; the mean square divergence angle in the vertical direction. The calculation method is to calculate the vertical divergence angle of each effective particle. Perform a squaring operation, sum all the squared results, and divide by the effective particle number M. The formula is: (from =1 to M) ; It is the first The vertical divergence angle of each effective particle, which reflects the average squared value of the vertical divergence angle of the particles, and the covariance term of the vertical position and divergence angle. The calculation method involves determining the vertical position of each effective particle. With vertical divergence angle Multiply, sum all the products, and divide by the effective particle number M, using the following formula: (from =1 to M) ( ) / M; This parameter reflects the correlation between the vertical particle position and the divergence angle.

[0060] Step 441: Calculate the geometrical reactivity in the horizontal direction using the calculated mean square position, mean square divergence angle, and covariance term. Specifically, this includes inputting the mean square position in the horizontal direction obtained in step 440 as the input parameter. Horizontal mean square divergence angle Horizontal covariance term Geometric reactivity in the horizontal direction It is a core parameter characterizing the volume of the beam in the horizontal phase space (position-divergence angle space), directly reflecting the beam quality. Its calculation follows the standard logic of beam physics. Combining the formula given in the document, the calculation formula is: = The specific calculation process consists of two steps. The first step is to calculate the radicand term within the square root. First, the mean square position in the horizontal direction is... Mean square divergence angle with horizontal direction Multiply, we get × Then calculate the horizontal covariance term. Squared, we get Then subtract the latter from the former, that is × - The second step is to perform an arithmetic square root operation on the difference obtained in the first step; the result is the geometric radii in the horizontal direction. The unit is mm·mrad.

[0061] Step 442: Calculate the geometrical reactivity in the vertical direction using the calculated vertical mean square position, mean square divergence angle, and covariance term; specifically, this includes: the input parameter being the vertical mean square position obtained in step 440. Vertical mean square divergence angle Vertical covariance term Vertical geometrical reactivity The volume of the beam in the vertical phase space, together with the emittance in the horizontal direction, constitutes a complete description of the beam's phase space quality. The calculation formula is as follows: = The calculation logic is consistent with the horizontal direction; the specific calculation process is as follows: First, calculate the radicand term within the square root, first squaring the vertical mean square position. Mean square divergence angle with vertical direction Multiply, we get × Then calculate the vertical covariance term. Squared, we get Then subtract the latter from the former, that is × - The second step is to perform an arithmetic square root operation on the difference, and the result is the geometric radiance in the vertical direction. The unit is also mm·mrad.

[0062] Step 443: Based on the geometrical emittance, mean square position, mean square divergence angle, and covariance term in the horizontal direction, the Twiss parameters in the horizontal direction are reconstructed; specifically, this includes: the input parameter being the mean square position in the horizontal direction from step 440. Horizontal mean square divergence angle Horizontal covariance term and the horizontal geometrical emittance of step 441 The Twiss parameter is the core input of the transfer matrix method, including the horizontal direction. function( ), function( ), Function ( ), the three jointly describe the focusing characteristics and transmission law of the beam in the horizontal direction. The reconstruction logic meets the core requirement of converting particle position and momentum information into Twiss parameters. For the specific calculation processes of the three parameters, in the horizontal direction Function ( ), which reflects the focusing degree of the beam in the horizontal direction, and the calculation formula is ; Divide the obtained in step 440 by the obtained in step 441, and the quotient is , in the horizontal direction Function ( ), which reflects the coupling degree between the position and divergence angle of the beam in the horizontal direction, and the calculation formula is ; Take the negative value of the obtained in step 440, and then divide it by , and the quotient is , is a dimensionless parameter; in the horizontal direction Function ( ), which reflects the distribution characteristics of the divergence angle of the beam in the horizontal direction, and the calculation formula is ; Divide the obtained in step 440 by , and the quotient is .

[0063] Step 444, based on the geometric emittance, mean square position, mean square divergence angle and covariance term in the vertical direction, reconstruct the Twiss parameters in the vertical direction; specifically including: the input parameters are the mean square position in the vertical direction of step 440, the mean square divergence angle in the vertical direction, the covariance term in the vertical direction, and the geometric emittance in the vertical direction of step 442. The Twiss parameters in the vertical direction include (vertical function), (vertical function), (vertical function). Their physical meanings correspond to those in the horizontal direction, only describing the focusing characteristics and transmission law of the beam in the vertical direction. The reconstruction logic is the same as that in the horizontal direction.

[0064] For the specific calculation processes of the three parameters, the vertical function , which reflects the focusing degree of the beam in the vertical direction, and the calculation formula is = / ; Divide the Divided by the result obtained in step 442 , the quotient is , in the vertical direction function , which reflects the coupling degree between the vertical position and the divergence angle of the beam, and the calculation formula is =- / ; Take the negative value of the obtained in step 440, and then divide it by , the quotient is , which is a dimensionless parameter; in the vertical direction function , which reflects the distribution characteristics of the vertical divergence angle of the beam, and the calculation formula is = / ; Divide the obtained in step 440 by , the quotient is .

[0065] Step 445, perform statistical analysis on the energy information in the processed particle data, and obtain the average energy and energy dispersion parameter of the beam through normal distribution fitting, which are used as energy spectrum parameters; specifically including: the energy information in the processed particle data is the effective particle energy sequence ( represents the th effective particle energy), the energy spectrum parameter describes the energy distribution characteristics of the beam, and the core function of the energy degrader is to modulate the beam energy. Therefore, the energy spectrum parameter is crucial for the calculation of the transmission line after the energy degrader. For example, the magnetic rigidity of the magnet is related to the energy, and the energy distribution affects the transmission trajectory. The specific processing process is divided into two steps. The first step is to perform basic statistics on the energy sequence and calculate the arithmetic mean energy (from = 1 to M) , and the energy standard deviation<000|0534]]= , and the two are used as the initial parameters for subsequent fitting; the second step is to fit the energy sequence with a normal distribution. The probability density function of the normal distribution is , where E is the energy variable, is the fitting average energy, [[ID=|2]] is the fitting energy dispersion. When fitting, use and as the initial values, and adjust and through the least squares method or maximum likelihood estimation to make the deviation between the fitting curve and the energy histogram the smallest. The finally determined and are the energy spectrum parameters.

[0066] By transforming the raw particle data into the inputs required by the transfer matrix method, such as Twiss parameters and energy spectrum parameters, the nonlinear simulation accuracy of Monte Carlo simulation is combined with the linear computational efficiency of the transfer matrix method, ensuring the accuracy of beam parameter calculation and providing reliable support for transmission simulation after the de-energizer.

[0067] In a preferred embodiment of the present invention, step 5 above, which involves constructing a transmission matrix model of the beam transmission line after the de-energizer based on the optical energy spectrum parameters of the exit beam and the parameters of the optical elements after the de-energizer, to achieve rapid calculation of the beam transmission process after the de-energizer, may include:

[0068] In this embodiment of the invention, step 550 involves using the reconstructed horizontal Twiss parameters and vertical Twiss parameters as the beam optical input parameters for calculating the starting point of the beam transmission line after the energy degrader; specifically, this includes: first determining the composition of the reconstructed Twiss parameters, wherein the horizontal Twiss parameters include horizontal... function (denoted as) ),level function (denoted as) ),level function (denoted as) Vertical Twiss parameters include vertical... function (denoted as) ), vertical function (denoted as) ), vertical function (denoted as) These parameters are all reconstructed through steps 443 and 44 and are the core indicators describing the optical characteristics of the beam. The starting point for calculating the beam transmission line after the de-energizer is the exit position of the de-energizer. The transfer matrix method requires the beam optical parameters at this starting point as input in order to derive the beam information at subsequent positions through matrix operations. The core reason for choosing the Twiss parameter as input is that the Twiss parameter can completely characterize the distribution characteristics of the beam in the phase space (position-divergence angle space). The transmission line after the de-energizer is composed of linear optical elements (such as drift segments and magnets). The calculation of the linear segment by the transfer matrix method needs to be based on the distribution characteristics of the phase space.

[0069] Step 551: Using the beam optical input parameters and the fitted energy spectrum parameters as input conditions, construct the corresponding transmission matrix according to the type and parameters of each optical element in the beam transmission line after the de-energizer. Specifically, this includes: First, determining the complete composition of the input conditions. In addition to the beam optical input parameters determined in step 550, it is also necessary to combine the energy spectrum parameters fitted in step 445, because the beam energy directly affects the transmission characteristics of the optical elements in the transmission line. The optical elements in the beam transmission line after the de-energizer mainly include the drift section, the diode magnet, and the quadrupole magnet. According to the type and specific structural parameters of each element, a corresponding transmission matrix needs to be constructed. The specific process is as follows: Construct the transmission matrix of the drift section. The drift section is a free space without electromagnetic field, which only realizes the straight-line transmission of the beam. Its core parameter is the length L. The transmission matrix of the drift section with length L is... The matrix is ​​structured as follows: the first row and first column contain 1s, the first row and second column contain L, the second row and first column contain 0s, and the second row and second column contain 1s. This matrix reflects the physical law that the beam's position increases linearly with the transmission distance and the divergence angle remains constant during the drift segment. When constructing the matrix, simply substitute the actual drift segment length L to determine the specific values. A transmission matrix is ​​then constructed using a dipole magnet, which is used to deflect the beam. Key parameters include the particle trajectories. Magnetic stiffness related items Deflection angle The transmission matrix of the dipole magnet in the x-direction needs to be combined with , , Calculation, matrix elements contain , - / , Items, etc., should be constructed first according to... Sure Substitute and The actual values ​​are used to obtain the transmission matrix in the x-direction; the beam is undeflected in the y-direction, and its transmission matrix is ​​similar in form to that of the drift section, only requiring adjustment of length-related elements in the matrix according to the equivalent length of the magnet in the y-direction; the transmission matrix of the quadrupole magnet is constructed, which is used to achieve beam focusing or defocusing, and the core parameters include focusing intensity. The length is denoted as l, and the quadrupole magnet is in the focal plane ( The transfer matrix (with positive values) is in the form of a matrix where the first row and first column elements are... ( l), The element in the first row and second column is ( l) / The element in the second row and first column is - ( l), The element in the second row and second column is ( l); in the defocus plane ( If it is negative, record it as negative. =- , The transfer matrix (where the integers are positive) is in the form of a matrix with the first row and first column containing the following elements: h( l), the element in the first row and second column is sinh( l) / The element in the second row and first column is sinh( l), The element in the second row and second column is cosh ( l), during construction, first determine based on the function of the magnet. The positive and negative signs, substituted into reality The transmission matrix of the corresponding plane is calculated using the value of l.

[0070] Step 552: Combine the transmission matrices according to the actual transmission order of the beam in the transmission line to form a complete transmission matrix model of the beam transmission line after the de-energizer. Specifically, this includes: First, it is necessary to determine the actual transmission order of the beam in the transmission line after the de-energizer. This order is determined by the physical design of the transmission line. For example, a typical transmission order is drift section, quadrupole magnet (focusing), drift section, dipole magnet (deflection), quadrupole magnet (defocusing), drift section. It is necessary to strictly follow the principle that the elements that the beam passes through first are in front, and the elements that pass through last are in the back. Since the multiplication of the transmission matrix does not satisfy the commutative law, an incorrect matrix combination order will cause the final total matrix to deviate completely from the actual transmission law. Therefore, the accuracy of the order is a key prerequisite for model construction.

[0071] The core logic of matrix combination is based on the transmission matrix principle of linear combination systems. If the beam originates from position... Transmit to location The transmission matrix is ​​M( | ), from position Transmit to location The transmission matrix is ​​M( | ), then the beam from position Transmit to location The total transmission matrix is ​​M( | )=M( | )×M( | That is, the subsequent transmission matrix is ​​multiplied by the preceding transmission matrix, and so on. The transmission line after the de-energizer starts from the calculation point. (Energy reducer outlet) to transmission line terminal The total transmission matrix is ​​obtained by multiplying the transmission matrices of all components in the actual transmission order. This total matrix is ​​the complete transmission matrix model of the beam transmission line after the de-energizer. The specific combination process is as follows: according to the actual beam transmission order, all optical components in the transmission line after the de-energizer are numbered sequentially from 1 to n, where component 1 is the first component the beam passes through, and component n is the last component the beam passes through. The transmission matrix corresponding to each component is denoted as follows: (Transmission matrix of component 1) (Transmission matrix of element 2)... (Transmission matrix of element n), calculate beam current from To the outlet position of component 1 The transmission matrix, M( | )= This matrix contains only the transmission characteristics of the first element, calculating the beam current from... To the outlet position of component 2 The transmission matrix: M( | )= ×M( | )= × This matrix contains the combined transmission characteristics of the first two elements, and so on, progressively calculating the beam current from... The transmission matrix is ​​then used to calculate the beam current from the exit positions of each subsequent component. To the terminal Total transmission matrix Its expression is = ×……× × This matrix covers the transmission characteristics of all linear optical elements in the transmission line after the de-energizer, and can completely describe the linear transport process of the beam from the de-energizer outlet to the transmission line terminal.

[0072] Step 553: Calculate the beam optical input parameters using the complete beam transmission line transmission matrix model after the de-energizer to obtain the beam optical parameter distribution at any position in the beam transmission line after the de-energizer. Specifically, this includes: first, determining the mathematical expression of the beam optical input parameters. In beam physics, the beam position and divergence angle are usually represented by phase space vectors, where the horizontal phase space vector is denoted as... x represents the horizontal position, and x' represents the horizontal divergence angle; the vertical phase space vector is denoted as... Where y is the vertical position and y' is the vertical divergence angle, the Twiss parameters determined in step 550 can be transformed into a calculation starting point using the fundamental relationships of beam optics. The statistical properties of the phase space vector provide an initial statistical basis for transmission calculations. The core logic of transmission calculations is that, for any position in the transmission line after the de-energizer... First determine the beam from to The transmission matrix M ( | ), then The phase space vector at point M ( | Perform matrix multiplication to obtain the position. The phase space vector at that location is used to derive the beam optical parameters at that position. The specific calculation process, taking the horizontal direction as an example, is detailed below: Determining the position... The corresponding transmission matrix M ( | Assume the matrix has the following specific form: the element in the first row and first column is 'a', the element in the first row and second column is 'b', the element in the second row and first column is 'c', and the element in the second row and second column is 'd', i.e., M( | )=[[a, b], [c, d]], to get the starting point of the calculation. horizontal phase space vector at the location ( )=( , )ᵀ, among which for The horizontal position at that location for Horizontal divergence angle at the location, calculate position horizontal phase space vector at the location ( According to the rules of matrix multiplication, ( )=M( | )× ( ), after unfolding, the horizontal position is obtained. ( )=a× +b× Horizontal divergence angle ( )=c× +d× The position can be directly obtained through these two expressions. The horizontal position and divergence angle of the beam are used to derive the position. The horizontal Twiss parameter at that location is based on ( )and ( Based on the statistical properties of the data, combined with the calculation logic of geometric emissivity in step 441 and the reconstruction logic of Twiss parameters in step 443, the position can be further calculated. level function ( ), function ( ), function ( The vertical transmission calculation process is exactly the same as the horizontal one; only the horizontal phase space vector needs to be changed. Replace with vertical phase space vector The horizontal transmission matrix M ( | By replacing the matrix with the vertical transfer matrix, the position can be obtained. The vertical position y( ), vertical divergence angle y' ( and vertical Twiss parameters ( ), ( ), ( ).

[0073] Step 554: Based on the beam optical parameter distribution and energy spectrum parameters, a rapid simulation calculation of the entire beam transport process from the de-energizer outlet to the beam transmission line terminal is achieved. Specifically, this includes: firstly, determining that the core inputs for the full-process simulation calculation are two parts: one is the beam optical parameter distribution obtained in step 553, which covers all positions from the de-energizer outlet to the transmission line terminal, including the beam position, divergence angle, and Twiss parameters at each position, describing the spatial and focusing characteristics of the beam; the other is the energy spectrum parameters obtained in step 445, describing the energy characteristics of the beam. Both need to work together to fully reflect the physical behavior of the beam during transmission, because the spatial transmission characteristics of the beam are affected by energy. The specific simulation calculation process is as follows: the correlation and integration of optical parameters and energy spectrum parameters is performed, integrating the optical parameters at each position in the beam optical parameter distribution... The optical parameters are correlated with the energy spectrum parameters; for example, for energies higher than... The particle whose radius of motion in the diode magnet is slightly larger than that of the average energy particle will cause the particle to move at a different position. The horizontal deflection position is slightly greater than the average energy of the particle x ( At this point, it is necessary to base it on energy dispersion. , regarding position The transmission matrix is ​​fine-tuned to ensure that the optical parameters at that location cover the transmission characteristics of particles with different energies, thus ensuring that the simulation results closely resemble reality. When the amplitude is small, the fine-tuning amplitude can be ignored to balance simulation accuracy and computational efficiency. The rationality verification of the transmission process is based on the integrated optical parameter distribution. The focus is on verifying whether the beam parameters at key locations in the transmission line meet the design requirements. For example, check whether the horizontal position of the beam at the terminal is at the center of the target area and whether the horizontal divergence angle is less than the maximum threshold allowed by the design. If the parameters do not meet the requirements, the structural parameters of the optical components in the transmission line need to be adjusted in reverse. Then, the optical parameter distribution is recalculated through steps 551 to 553 until the parameters at key locations meet the design requirements. This achieves iterative optimization of the beam transmission line. The beam parameters at all locations from the degrader outlet to the transmission line terminal are organized in the transmission sequence to form a full-process beam transport simulation report. The report should clearly indicate the parameter values, design target values, and deviation ranges at each key location to provide an intuitive basis for the verification and adjustment of the transmission line design scheme.

[0074] By using the reconstructed Twiss parameters as the starting point for the transfer matrix method, and combining the energy spectrum parameters to construct the element transfer matrix and complete model, the Monte Carlo method's accurate description of the nonlinear process within the depressor is combined with the transfer matrix method's efficient calculation of the linear segment, forming a seamless full-process simulation link and improving the simulation efficiency of transmission lines containing depressors.

[0075] Embodiments of the present invention also provide a computing device, including: a processor and a memory storing a computer program, wherein the computer program, when executed by the processor, performs the method described above. All implementations in the above method embodiments are applicable to this embodiment and can achieve the same technical effects.

[0076] Embodiments of the present invention also provide a computer-readable storage medium storing instructions that, when executed on a computer, cause the computer to perform the method described above. All implementations in the above method embodiments are applicable to this embodiment and can achieve the same technical effects.

[0077] The above description represents the preferred embodiments of the present invention. It should be noted that those skilled in the art can make various improvements and modifications without departing from the principles of the present invention, and these improvements and modifications should also be considered within the scope of protection of the present invention.

Claims

1. A fast calculation method for beam transmission in a degrader segment based on the transfer matrix and Monte Carlo method, characterized in that, The method includes: Step 1: Construct a transmission matrix model of the beam transmission section in front of the de-energizer based on the beam optical parameters and optical element parameters at the starting point of the beam transmission line, and obtain the particle phase space distribution parameters at the inlet position of the de-energizer; Step 2: Based on the particle phase space distribution parameters at the inlet of the de-energizer, establish a three-dimensional structural model of the de-energizer and perform Monte Carlo beam transport simulation of the de-energizer section to obtain the position, momentum and energy information of the beam particles at the outlet of the de-energizer. Step 3: Perform abnormal data removal processing on the position, momentum, and energy information of the beam particles at the outlet of the de-energizer to obtain the processed particle data; Step 4: Based on the processed particle data, calculate the mean square position, mean square divergence angle and covariance term of the particles in the beam exiting the de-energizer, and reconstruct the optical energy spectrum parameters of the beam exiting the de-energizer and the parameters of the optical elements behind the de-energizer. Step 5: Based on the optical energy spectrum parameters of the beam exiting the de-energizer and the parameters of the optical elements after the de-energizer, construct the transmission matrix model of the beam transmission line after the de-energizer to realize the rapid calculation of the beam transmission process after the de-energizer.

2. The method of claim 1, wherein the method is a fast calculation method of beam transport in a degrader section based on transport matrix and Monte Carlo method. Based on the beam optical parameters at the beam transmission line initiation point and the optical element parameters, a transmission matrix model of the beam transmission section before the depressor is constructed, yielding the particle phase space distribution parameters at the depressor inlet position, including: Obtain the initial Twiss parameters of the beam transmission line start point as input parameters; Based on the type and structural parameters of each optical element in the beam current transmission section before the depletor, determine the transmission matrix corresponding to each optical element; The transmission matrices of each optical element are cascaded and combined according to their physical order in the beam transmission line to construct the transmission matrix model of the beam transmission section in front of the depressor. The initial Twiss parameters are input into the transmission matrix model of the beam transmission section in front of the depressor for calculation, and the particle phase space distribution parameters at the depressor inlet position are obtained.

3. The fast calculation method for beam transmission in the degrader segment based on the transfer matrix and Monte Carlo method according to claim 2, characterized in that, Based on the particle phase space distribution parameters at the depressor inlet, a three-dimensional structural model of the depressor is established, and Monte Carlo beam transport simulation of the depressor section is performed to obtain the position, momentum, and energy information of the beam particles at the depressor outlet, including: The particle phase space distribution parameters at the inlet of the depressor are used as the initial beam conditions; Based on the initial beam conditions, a complete three-dimensional geometric model of the de-energizer is established; The three-dimensional geometric model of the degrader is converted into an input file in the format of a constructed solid geometry to obtain a geometric description suitable for Monte Carlo simulation. The Monte Carlo simulation program's running parameters are configured using the input files of the initial beam conditions and the constructed solid geometry as input parameters. The Monte Carlo simulation program was executed to simulate the complete transport process of the beam in the three-dimensional structure of the de-energizer, and the position, momentum and energy information of each beam particle at the de-energizer outlet were obtained.

4. The fast calculation method for beam transmission in the degrader segment based on the transfer matrix and Monte Carlo method according to claim 3, characterized in that, The three-dimensional geometric model of the energy degrader includes a wedge-shaped scattering device, a vacuum cavity, and a mechanical support structure.

5. The fast calculation method for beam transmission in the degrader segment based on the transfer matrix and Monte Carlo method according to claim 4, characterized in that, The position, momentum, and energy information of the beam particles at the depressor outlet are processed to remove outliers, resulting in processed particle data, including: Statistical analysis was performed on the position, momentum, and energy information of each beam particle at the outlet of the de-energizer, and the average value and standard deviation of energy, horizontal position, horizontal divergence angle, vertical position, and vertical divergence angle were calculated respectively. Based on the mean and standard deviation, the normal data range corresponding to each component is set; abnormal particle data that exceeds the corresponding normal data range in each component is identified and removed, while valid particle data within the normal data range is retained. The effective particle data is integrated to form processed particle data.

6. The method of claim 5, wherein the method is characterized by, The normal data range is set to [μ-3σ, μ+3σ], where μ represents the average value, used to determine the center value of each physical quantity in the beam, and σ represents the standard deviation, used to measure the dispersion of the physical quantities.

7. The fast calculation method for beam transmission in the degrader segment based on the transfer matrix and Monte Carlo method according to claim 6, characterized in that, Based on the processed particle data, the mean square position, mean square divergence angle, and covariance term of the particles in the de-energizer exit beam are calculated. The optical energy spectrum parameters of the de-energizer exit beam and the parameters of the optical elements behind the de-energizer are then reconstructed, including: Based on the processed particle data, the mean square position, mean square divergence angle, and covariance term of position and divergence angle in the horizontal direction are calculated respectively. At the same time, the mean square position, mean square divergence angle, and covariance term of position and divergence angle in the vertical direction are also calculated. Using the calculated mean square position, mean square divergence angle, and covariance term in the horizontal direction, the geometrical reactivity in the horizontal direction is calculated. The geometrical reactivity in the vertical direction is calculated using the calculated vertical mean square position, mean square divergence angle, and covariance term. Based on the geometric emittance, mean square position, mean square divergence angle and covariance term in the horizontal direction, the Twiss parameters in the horizontal direction are reconstructed. Based on the geometric emittance, mean square position, mean square divergence angle and covariance term in the vertical direction, the Twiss parameters in the vertical direction are reconstructed. Statistical analysis was performed on the energy information in the processed particle data, and the average energy and energy dispersion parameters of the beam were obtained by fitting a normal distribution, which were used as energy spectrum parameters.

8. The method of claim 7, wherein the method is characterized by, Based on the optical energy spectrum parameters of the exit beam of the de-energizer and the parameters of the optical elements after the de-energizer, a transmission matrix model of the beam transmission line after the de-energizer is constructed to achieve rapid calculation of the beam transmission process after the de-energizer, including: The horizontal and vertical Twiss parameters obtained from the reconstruction are used as the beam optical input parameters for calculating the beam transmission line after the de-energizer. Using the beam optics input parameters and the fitted energy spectrum parameters as input conditions, a corresponding transmission matrix is ​​constructed based on the type and parameters of each optical element in the beam transmission line after the de-energizer. The transmission matrices are combined according to the actual transmission order of the beam in the transmission line to form a complete transmission matrix model of the beam transmission line after the de-energizer. The beam optical input parameters are transmitted using a complete beam transmission line transmission matrix model after the de-energizer, and the beam optical parameter distribution at any position in the beam transmission line after the de-energizer is obtained. Based on the beam optical parameter distribution and energy spectrum parameters, a rapid simulation calculation of the entire beam transport process from the de-energizer outlet to the beam transmission line terminal is achieved.

9. A computing device, comprising: include: One or more processors; A storage device for storing one or more programs, which, when executed by one or more processors, cause the one or more processors to implement the method as described in any one of claims 1 to 8.

10. A computer-readable storage medium, characterized in that, The computer-readable storage medium stores a program that, when executed by a processor, implements the method as described in any one of claims 1 to 8.

Images