Analysis method for flow-induced vibration of spiral cross-shaped fuel rod

Abstract

The invention provides an analysis method for flow-induced vibration of spiral cross-shaped fuel rods, which comprises the following steps of: establishing a three-dimensional geometric model of a fluid domain by enlarging the distance between rods so as to eliminate the influence of self-supporting contact on grids, and carrying out grid division by adopting an unstructured grid and boundary layer optimization strategy; capturing a transient turbulence pulsating pressure field based on a WALE sub-grid model of a large eddy simulation model; solid domain modeling keeps the original position of a fuel rod and explicitly defines a contact pair, and the contact force is calculated in real time according to a penalty function method; the fluid time-varying pressure load is mapped to the solid surface to execute transient dynamic calculation, and vibration displacement response and contact force time history data are synchronously output; the vibration dominant frequency is extracted through frequency domain power spectrum analysis, and the resonance risk is compared and evaluated by combining the wet mode inherent frequency. According to the method, co-simulation of fluid domain gap processing and solid domain real contact dynamics can be achieved, and an analysis means is provided for prediction of flow-induced vibration and fretting wear of the spiral cross-shaped fuel rod.

Description

Technical Field

[0001] This invention relates to the field of nuclear reactor fuel flow-induced vibration analysis technology, specifically to a method for analyzing the flow-induced vibration of a spiral cross fuel rod based on numerical simulation. Background Technology

[0002] The spiral cross-shaped fuel design represents a major innovation in reactor fuel design. Its cross-shaped cross-section twists axially to form a self-supporting structure, significantly improving core power density and thermal-hydraulic performance. This design abandons the traditional positioning grid, forming periodic contact supports between fuel rods. This increases the fuel surface area, enhances coolant mixing, and reduces flow pressure drop. However, the high-intensity turbulent pulsations generated by the coolant within the complex spiral flow channels can excite unsteady pressure loads on the fuel rod surface. The self-supporting contact points are more likely to amplify the transfer of vibrational energy, ultimately leading to fretting wear and even cladding failure.

[0003] Existing research indicates that over 55% of fuel failures in conventional light water reactors are attributed to flow-induced vibrations between the positioning grid and fuel rods. While helical cross fuels eliminate the grid, their unique geometry makes the fluid-structure interaction mechanism more complex, and their vibration characteristics remain unclear. Current research on flow-induced vibrations in helical cross fuels faces a technological gap. Traditional turbulence models struggle to accurately capture the strongly anisotropic turbulence characteristics within the helical flow channel. While Large Eddy Simulation (LES) possesses high accuracy potential, its application in helical cross fuels lacks systematic usage guidelines. Furthermore, there is no unified standard for simplifying self-supporting structures, and existing gap treatment schemes have fundamental flaws: shortening the straight section of the fuel rod creates gaps but severely reduces the effective length of the fuel rod, distorting the original thermal-hydraulic characteristics. Finally, there are shortcomings in handling fluid-structure interaction vibration analysis of real contact dynamics in the solid domain, failing to accurately simulate the dynamic interactions, energy transfer, and fretting wear between self-supporting structures. The rod spacing extension rule of this invention fully preserves the geometric integrity of the fuel rods, while explicitly modeling the contact in the solid domain and coupling the fluid load for nonlinear transient dynamic solutions, providing a fundamental guarantee for the fidelity of the flow field and the accuracy of vibration response (including contact force) prediction.

[0004] Although a mature framework has been established for the study of flow-induced vibration of conventional fuel rods, its conclusions regarding mechanisms such as lattice excitation force and vortex shedding are not applicable to self-supporting structures without lattice. This lack of technology severely restricts the safety assessment of helical cross fuels for engineering applications, and there is an urgent need to establish an integrated analysis method that combines geometric modeling, fluid simulation, vibration response, and risk assessment. Summary of the Invention

[0005] To address the complex flow-induced vibration analysis challenges faced by helical cross fuel rods in reactor core environments, this invention provides a method for analyzing the flow-induced vibration of helical cross fuel rods, constructing a complete technical chain encompassing "flow field calculation - load transfer - dynamic response - safety assessment." The method uses a rod spacing expansion method to set the gap, establishes a three-dimensional model, and performs high-precision transient fluid calculations to capture the turbulent pulsating pressure field. The solid domain model maintains the original design position of the fuel rods, explicitly establishing contact pairs between adjacent rods and defining contact mechanical behavior based on penalty functions and Coulomb friction. In transient fluid-structure interaction analysis, not only is the vibration displacement response of the rods caused by fluid loads calculated, but the dynamic normal pressure and tangential friction at the contact interface are also solved in real time. This fully reproduces the nonlinear dynamic response mechanism of the helical cross fuel self-supporting structure under complex flow field excitation, providing mechanical parameters for assessing fretting wear at contact points. Furthermore, time-domain transient dynamic response analysis is conducted, combining frequency-domain power spectral density analysis and wet modal calculations to compare the overlap characteristics of the excitation force's dominant frequency and the structure's natural frequency, quantifying resonance risk.

[0006] To achieve the above objectives, the present invention adopts the following technical solution:

[0007] An analytical method for flow-induced vibration of a helical cross fuel rod includes the following steps:

[0008] Step 1: Geometric Modeling and Gap Setting: Establish a fluid domain model of the spiral cross fuel rod bundle. By increasing the rod spacing and setting a predetermined gap between adjacent fuel rods, the self-supporting contact structure is eliminated, and a three-dimensional geometric model of the spiral cross fuel rod bundle flow channel is obtained.

[0009] Step 2: Fluid Domain Mesh Processing: The fluid domain is then meshed using an unstructured mesh. The specific steps are as follows:

[0010] Step 2-1: Establish a fluid domain model for a single spiral cross fuel rod, mesh the fluid domain model, divide the boundary layer mesh on the wall of the fuel rod in the fluid domain model, adjust the number of boundary layer layers and the thickness of the first boundary layer to divide several different boundary layer mesh generation strategies.

[0011] Step 2-2: The fluid domain calculation uses the same boundary conditions. The transient calculation is performed using the large eddy simulation model until convergence. The wall y+ values ​​are compared under different boundary layer strategies.

[0012] Steps 2-3: Select the group that satisfies the condition that the y+ value is less than 1 and has the fewest grid cells to perform subsequent calculations;

[0013] Steps 2-4: Divide the volume mesh into different sizes, adopt the boundary layer strategy described above, use the Reynolds-Averaged Navier-Stokes turbulence model to carry out steady-state calculations, perform mesh independence analysis, and determine the minimum number of meshes that meet the calculation accuracy.

[0014] Step 3: CFD calculation of the flow channel of the spiral cross fuel rod bundle: Transient fluid dynamics calculations are performed using a large eddy simulation model, specifically in the following steps:

[0015] Step 3-1: The coolant is water. Using the physical property parameters of the corresponding reactor under steady state, set the velocity inlet and pressure outlet. After initialization, use the Reynolds-Averaged Navier-Stokes turbulence model to perform steady-state calculations until convergence, providing an initial flow field environment for transient calculations.

[0016] Step 3-2: The transient calculation adopts the large eddy simulation model, and the subgrid stress model adopts the wall-adapting local eddy-viscosity model. The flow field information of each time step is recorded to obtain the time-varying fluid excitation force load field. The value of the fluid excitation force at each time step is analyzed in the frequency domain by Fourier transform (FFT) to obtain the result of the fluid excitation force power spectral density function.

[0017] Step 4: Solid Domain Modeling and Constraint Setting: Establish a solid domain model of the fuel rod bundle containing the contact structure of adjacent fuel rods. The fuel rods maintain their original design relative positions without increasing the spacing, realistically reflecting the self-supporting contact structure between the rods. Explicitly define contact pairs between the surfaces where adjacent fuel rods may contact. Perform swept mesh generation on each fuel rod, and locally refine the mesh on the surfaces of the areas where adjacent rods are expected to contact to ensure the accuracy of the contact calculation. Perform solid domain mesh independence analysis to determine the appropriate mesh size. Apply fixed constraints at both ends of the fuel rod axis to simulate its actual support boundary conditions in the reactor core.

[0018] Step 5: Unidirectional fluid-structure interaction analysis and contact dynamics calculation: Select a solid domain model containing contact structures, and accurately map the time-varying fluid excitation force load field obtained in Step 3 onto the outer surface of each fuel rod in the selected solid domain model according to the corresponding time step, and perform transient structural dynamics calculation. This calculation solves the vibration displacement response of the fuel rod under fluid pulsating pressure excitation, and at the same time calculates and outputs the dynamically evolving normal contact pressure and tangential friction force on the contact interface in real time.

[0019] Step 6: Vibration characteristics and contact force extraction: Obtain time-domain data of vibration displacement at different axial positions of the fuel rod, perform Fourier transform (FFT) frequency domain analysis on it to obtain the power spectral density distribution of the vibration displacement, and obtain the dominant vibration frequency of the fuel rod; extract time-varying normal contact pressure and tangential friction force data at key positions in the contact area, perform time-domain statistics and frequency domain analysis on these contact force signals to obtain the dominant frequency characteristics and energy distribution of contact force pulsation;

[0020] Step 7: Resonance Risk Assessment: Construct a solid domain model and an outer fluid domain model of a single helical cross fuel rod. Use a common-node mesh coupling at the solid-fluid interface to conduct wet modal analysis and obtain the first six underwater natural frequencies of the fuel rod. Compare the vibration dominant frequency with the natural frequency. If the relative deviation between any natural frequency and the vibration dominant frequency is less than or equal to 5%, a resonance risk is determined to exist; otherwise, the risk is considered acceptable.

[0021] The method of increasing the spacing between the fuel rods in step 1 is to keep the original geometric features of the fuel rods unchanged, and to form a predetermined gap by increasing the overall size of the rod bundle array. This fully preserves the curved surface profile and torsional features of the spiral cross fuel, and while eliminating the self-supporting contact structure, ensures that the flow channel morphology change rate is less than 5%, thus maintaining the original flow field characteristics to the maximum extent.

[0022] The selection of the turbulence model in step 3 is determined by comparing the ability of different models to capture wall pressure fluctuations.

[0023] In step 4, a solid domain model of the rod bundle containing the contact structure of adjacent fuel rods is established, in which the fuel rods maintain their original relative positions and contact pairs are explicitly defined in the contact area of ​​adjacent rods; the contact algorithm adopts surface-to-surface contact based on the penalty function method, defines the contact stiffness and friction coefficient, and allows the contact state to change dynamically during the calculation process.

[0024] In step 5:

[0025] Normal contact pressure F n The calculation formula is:

[0026]

[0027] Tangential friction force f t The calculation formula is:

[0028] f t =μF n (2)

[0029] In the formula, g is the real-time gap between the contact surfaces, and a negative value indicates the penetration depth |g|; k n is the normal contact stiffness coefficient; μ is the friction coefficient.

[0030] The convergence condition in step 2-2 is that the fluid transient calculation residual reaches 0.001. This convergence value ensures the statistical stability of transient turbulent pulsating pressure by balancing computational efficiency and accuracy requirements. Under the premise of meeting the analytical accuracy of the wall y+ value, the residual threshold of 0.001 can effectively capture the high-frequency pulsating component dominated by eddy shedding, while avoiding excessive consumption of computational resources.

[0031] The convergence condition in step 3-1 is that the steady-state calculation residual reaches 0.001. This convergence value provides a physically reasonable initial flow field for transient large eddy simulation, and by controlling the residual to suppress numerical oscillations, it ensures that the evolution of the turbulent kinetic energy spectrum in the subsequent LES calculation conforms to the actual physical process.

[0032] Compared with the prior art, the present invention has the following advantages:

[0033] 1. This invention provides an analysis method for flow-induced vibration of helical cross fuel rods, which can provide a simulation benchmark for experimental verification.

[0034] 2. This invention uses the rod spacing expansion method to set the gap and adjust the spatial arrangement of the rod bundle array, overcoming the limitations brought about by shortening the length of the straight section of the fuel rod, and achieving a unity of geometric fidelity and mesh computability.

[0035] 3. The fluid domain calculation of this invention is based on large eddy simulation and fine mesh, which accurately analyzes the transient fluid excitation force on the wall.

[0036] 4. Because the solid domain model completely preserves the geometric features of the fuel rod and the fluid pressure load is mapped to cover the entire outer surface, this invention can realize time-frequency domain analysis of fuel rod vibration at any axial position, and combine modal calculation to quantify resonance risk.

[0037] 5. This invention explicitly establishes and solves the dynamic contact problem between fuel rods in the solid domain. By coupling fluid load and contact dynamics, it can directly output the time history data of normal pressure and tangential friction in the contact area. Attached Figure Description

[0038] Figure 1 This is a flowchart of the method of the present invention.

[0039] Figure 2 This is a schematic diagram of a certain spiral cross fuel rod bundle structure.

[0040] Figure 3 This is a schematic diagram of the boundary layer division on the wall of the fuel rod in the fluid domain.

[0041] Figure 4 This is a schematic diagram of the spiral cross fuel contact point.

[0042] Figure 5 The curve shows the vibration displacement response at the axial half-point of the fuel rod at the center position. Detailed Implementation

[0043] The method of the present invention will be further described in detail below with reference to the accompanying drawings and specific embodiments:

[0044] like Figure 1 As shown, the present invention provides an analytical method for flow-induced vibration of a spiral cross fuel rod, comprising the following steps:

[0045] Step 1: Geometric Modeling and Gap Setting: Establish a fluid domain model for the helical cross fuel rod bundle. By increasing the rod spacing, obtain a three-dimensional geometric model of the helical cross fuel rod bundle flow channel. The method of increasing the rod spacing is to maintain the original geometric features of the fuel rods and form a predetermined gap by increasing the overall size of the rod bundle array. This fully preserves the curved surface profile and torsional features of the helical cross fuel, ensuring that the flow channel morphology change rate is less than 5% while eliminating self-supporting contact structures, thus maximizing the maintenance of the original flow field characteristics. Figure 2 The image shows a schematic diagram of the fluid domain model with 7 spiral cross fuel rods arranged in a hexagonal pattern.

[0046] Step 2: Fluid Domain Mesh Processing: The fluid domain is then meshed using an unstructured mesh. The specific steps are as follows:

[0047] Step 2-1: Establish a fluid domain model for a single helical cross fuel rod. Mesh the fluid domain model, creating boundary layer meshes on the walls of the fuel rod. Adjust the number of boundary layer layers and the thickness of the first boundary layer to generate several different boundary layer meshes; for example... Figure 3 As shown, this is a set of generation strategies for the boundary layer mesh on the fuel rod wall in the fluid domain model. From this enlarged schematic diagram, it can be clearly seen that the prism layer mesh grows gradually along the normal of the complex helical surface. This strategy ensures that the y+ value in the near-wall region is less than 1 throughout the entire domain, which meets the analytical requirements of the large eddy simulation LES model and lays a geometrically faithful foundation for accurately capturing the turbulent fluctuations in the helical channel in step 3-2.

[0048] Step 2-2: The fluid domain calculation uses the same boundary conditions. The transient calculation is performed using the Large Eddy Simulation (LES) model until the fluid transient calculation residual reaches the convergence condition of 0.001. The wall y+ values ​​are compared under different boundary layer strategies.

[0049] Steps 2-3: Select the group that satisfies the condition that the y+ value is less than 1 and has the fewest grid cells to perform subsequent calculations;

[0050] Steps 2-4: Divide the volume mesh into different sizes, adopt the boundary layer strategy described above, use the RANS (Reynolds-Averaged Navier-Stokes) turbulence model to carry out steady-state calculations, perform mesh independence analysis, and determine the minimum number of meshes that meet the calculation accuracy.

[0051] Step 3: CFD calculation of the flow channel of the spiral cross fuel rod bundle: Transient fluid dynamics calculations are performed using the Large Eddy Simulation (LES) model, specifically in the following steps:

[0052] Step 3-1: The coolant is water. Using the physical property parameters of the corresponding reactor under steady state, set the velocity inlet and pressure outlet. First, perform initialization, and then use the RANS turbulence model to perform steady-state calculations until the steady-state calculation residual reaches the convergence condition of 0.001, so as to provide an initial flow field environment for transient calculations.

[0053] Step 3-2: Transient calculations use the Large Eddy Simulation (LES) model, and the subgrid stress model is the Wall-Adapting Local Eddy-viscosity (WALE) model. The flow field information at each time step is recorded to obtain the time-varying fluid excitation force load field. The value of the fluid excitation force at each time step is analyzed in the frequency domain by Fourier Transform (FFT) to obtain the result of the fluid excitation force power spectral density function.

[0054] Step 4: Solid Domain Modeling and Constraint Setting: Establish a solid domain model of the fuel rod bundle containing the contact structure of adjacent fuel rods. The fuel rods maintain their original design relative positions without widening the spacing, realistically reflecting the self-supporting contact structure between the rods; for example... Figure 4 As shown, the contact points of two helical cross fuel rods under a certain pitch are displayed, with a total of 5 contact points. Contact pairs are explicitly defined between the surfaces of adjacent fuel rods that may come into contact. The contact algorithm adopts a surface-to-surface contact method based on the penalty function method, defining the contact stiffness and friction coefficient, and allowing the contact state to change dynamically during the calculation. A swept mesh is generated for each fuel rod, and the mesh is locally refined on the surface of the area where adjacent rods are expected to come into contact to ensure the accuracy of the contact calculation. Solid domain mesh independence analysis is performed to determine a suitable mesh size, and fixed constraints are applied to both ends of the fuel rods along the axial direction to simulate their actual support boundary conditions in the reactor core.

[0055] Step 5: Unidirectional fluid-structure interaction analysis and contact dynamics calculation: Select a solid domain model containing contact structures, and accurately map the time-varying fluid excitation force load field obtained from the fluid domain calculation in Step 3 onto the outer surface of each fuel rod in the selected solid domain model according to the corresponding time step, and perform transient structural dynamics calculation. This calculation solves the vibration displacement response of the fuel rod under fluid pulsating pressure excitation, and at the same time calculates and outputs the dynamically evolving normal contact pressure and tangential friction force on the contact interface in real time.

[0056] Normal contact pressure F n The calculation formula is:

[0057]

[0058] Tangential friction force f t The calculation formula is:

[0059] f t =μF n (2)

[0060] In the formula, g is the real-time gap between the contact surfaces, and a negative value indicates the penetration depth |g|; k n is the normal contact stiffness coefficient; μ is the friction coefficient.

[0061] Step 6: Vibration Characteristics and Contact Force Extraction: Obtain time-domain data of vibration displacement at different axial positions of the fuel rod, perform Fourier Transform (FFT) frequency-domain analysis on it to obtain the power spectral density distribution of the vibration displacement, and obtain the dominant vibration frequency of the fuel rod; for example... Figure 5 As shown, the vibration displacement response curve at the axial half of the fuel rod at the center position is obtained within 0.5 s. This curve fully records the continuous evolution of the displacement at this point over time under the scouring of the coolant. The time-domain curve shows that the maximum amplitude caused by the fluid excitation force is 2.3 micrometers. Based on this vibration time history data, FFT analysis can be performed to accurately extract the vibration dominant frequency distribution characteristics. The time-varying normal contact pressure and tangential friction force data at key locations in the contact area are extracted. Time-domain statistics and frequency-domain analysis are performed on these contact force signals to obtain the dominant frequency characteristics and energy distribution of the contact force pulsation.

[0062] Step 7: Resonance Risk Assessment: Construct a solid domain model and an outer fluid domain model of a single helical cross fuel rod. Use a common-node mesh coupling at the solid-fluid interface to conduct wet modal analysis and obtain the first six underwater natural frequencies of the fuel rod. Compare the vibration dominant frequency with the natural frequency. If the relative deviation between any natural frequency and the vibration dominant frequency is less than or equal to 5%, a resonance risk is determined to exist; otherwise, the risk is considered acceptable.

Claims

1. A method for analyzing flow-induced vibration of a helical cross fuel rod, characterized in that: Includes the following steps: Step 1: Geometric Modeling and Gap Setting: Establish a fluid domain model of the spiral cross fuel rod bundle. By increasing the rod spacing and setting a predetermined gap between adjacent fuel rods, the self-supporting contact structure is eliminated, and a three-dimensional geometric model of the spiral cross fuel rod bundle flow channel is obtained. Step 2: Fluid Domain Mesh Processing: The fluid domain is then meshed using an unstructured mesh. The specific steps are as follows: Step 2-1: Establish a fluid domain model for a single spiral cross fuel rod, mesh the fluid domain model, divide the boundary layer mesh on the wall of the fuel rod in the fluid domain model, adjust the number of boundary layer layers and the thickness of the first boundary layer to divide several different boundary layer mesh generation strategies. Step 2-2: The fluid domain calculation uses the same boundary conditions. The transient calculation is performed using the large eddy simulation model until convergence. The wall y+ values ​​are compared under different boundary layer strategies. Steps 2-3: Select the group that satisfies the condition that the y+ value is less than 1 and has the fewest grid cells to perform subsequent calculations; Steps 2-4: Divide the volume mesh into different sizes, adopt the boundary layer strategy described above, use the Reynolds-Averaged Navier-Stokes turbulence model to carry out steady-state calculations, perform mesh independence analysis, and determine the minimum number of meshes that meet the calculation accuracy. Step 3: CFD calculation of the flow channel of the spiral cross fuel rod bundle: Transient fluid dynamics calculations are performed using a large eddy simulation model, specifically in the following steps: Step 3-1: The coolant is water. Using the physical property parameters of the corresponding reactor under steady state, set the velocity inlet and pressure outlet. After initialization, use the Reynolds-Averaged Navier-Stokes turbulence model to perform steady-state calculations until convergence, providing an initial flow field environment for transient calculations. Step 3-2: The transient calculation adopts the large eddy simulation model, and the subgrid stress model adopts the wall-adapting local eddy-viscosity model. The flow field information of each time step is recorded to obtain the time-varying fluid excitation force load field. The value of the fluid excitation force at each time step is analyzed in the frequency domain by Fourier transform (FFT) to obtain the result of the fluid excitation force power spectral density function. Step 4: Solid Domain Modeling and Constraint Setting: Establish a solid domain model of the fuel rod bundle containing the contact structure of adjacent fuel rods. The fuel rods maintain their original design relative positions without increasing the spacing, realistically reflecting the self-supporting contact structure between the rods. Explicitly define contact pairs between the surfaces where adjacent fuel rods may contact. Perform swept mesh generation on each fuel rod, and locally refine the mesh on the surfaces of the areas where adjacent rods are expected to contact to ensure the accuracy of the contact calculation. Perform solid domain mesh independence analysis to determine the appropriate mesh size. Apply fixed constraints at both ends of the fuel rod axis to simulate its actual support boundary conditions in the reactor core. Step 5: Unidirectional fluid-structure interaction analysis and contact dynamics calculation: Select a solid domain model containing contact structures, and accurately map the time-varying fluid excitation force load field obtained in Step 3 onto the outer surface of each fuel rod in the selected solid domain model according to the corresponding time step, and perform transient structural dynamics calculation. This calculation solves the vibration displacement response of the fuel rod under fluid pulsating pressure excitation, and at the same time calculates and outputs the dynamically evolving normal contact pressure and tangential friction force on the contact interface in real time. Step 6: Vibration characteristics and contact force extraction: Obtain time-domain data of vibration displacement at different axial positions of the fuel rod, perform Fourier transform (FFT) frequency domain analysis on it to obtain the power spectral density distribution of the vibration displacement, and obtain the dominant vibration frequency of the fuel rod; extract time-varying normal contact pressure and tangential friction force data at key positions in the contact area, perform time-domain statistics and frequency domain analysis on these contact force signals to obtain the dominant frequency characteristics and energy distribution of contact force pulsation; Step 7: Resonance Risk Assessment: Construct a solid domain model and an outer fluid domain model of a single helical cross fuel rod. Use a common-node mesh coupling at the solid-fluid interface to conduct wet modal analysis and obtain the first six underwater natural frequencies of the fuel rod. Compare the vibration dominant frequency with the natural frequency. If the relative deviation between any natural frequency and the vibration dominant frequency is less than or equal to 5%, a resonance risk is determined to exist; otherwise, the risk is considered acceptable.

2. The method according to claim 1, characterized in that: The method of increasing the spacing between the fuel rods in step 1 is to keep the original geometric features of the fuel rods unchanged, and to form a predetermined gap by increasing the overall size of the rod bundle array. This fully preserves the curved surface profile and torsional features of the spiral cross fuel, and while eliminating the self-supporting contact structure, ensures that the flow channel morphology change rate is less than 5%, thus maintaining the original flow field characteristics to the maximum extent.

3. The method according to claim 1, characterized in that: The selection of the turbulence model in step 3 is determined by comparing the ability of different models to capture wall pressure fluctuations.

4. The method according to claim 1, characterized in that: In step 4, a solid domain model of the rod bundle containing the contact structure of adjacent fuel rods is established, wherein the fuel rods maintain their original relative positions and contact pairs are explicitly defined in the contact region of adjacent rods; The contact algorithm uses a surface-to-surface contact method based on the penalty function method, defines the contact stiffness and friction coefficient, and allows the contact state to change dynamically during the calculation process.

5. The method according to claim 1, characterized in that: In step 5: Normal contact pressure F n The calculation formula is: Tangential friction force f t The calculation formula is: f t =μF n (2) In the formula, g is the real-time gap between the contact surfaces, and a negative value indicates the penetration depth |g|; k n is the normal contact stiffness coefficient; μ is the friction coefficient.

6. The method according to claim 1, characterized in that: The convergence condition in step 2-2 is that the fluid transient calculation residual reaches 0.

001.

7. The method according to claim 1, characterized in that: The convergence condition in step 3-1 is that the steady-state calculation residual reaches 0.001.

Images