# Solution of the Static Deflection Mode Shape Function of the Cantilever Beam under Transverse Flow Based on the Boundary Shooting Method

Authors: Xiaowen Wang; Yixian Zhou
Journal: Science and Technology of Nuclear Installations
 (2021)

---

## Abstract

According to the characteristics of the reactor internal structure of nuclear power plants, the vibration of the secondary core support pillar in water can be modeled as the vibration of the cantilever beam structure under the action of transverse flow, and its first beam mode is highly likely to be activated. It is thus necessary to dedicate a separate study on the first-order beam mode. In this work, we study the secondary core support pillar in nuclear reactor AP1000 under the action of transverse flow and focus on the derivation of its static cantilever deflection mode shape function in order to lay a foundation for the calculation of hydrodynamic added mass and frequency for the nuclear reactor internal components and their structural integrity evaluation. First, we proposed a set of nonlinear differential equations for the analysis of the single cantilever beam. Second, to solve the nonlinear differential equations, we used a boundary shooting framework in combination with the Runge–Kutta method. The results of the numerical simulation agree with the analytical solution to a very high degree, which demonstrates the effectiveness of the simulation method. Finally, we solved the static deflection mode shape function of the secondary core support pillar under the normal operating conditions. The nonlinear differential model and simulation method proposed in this paper can be used to solve the static cantilever deflection mode shape function of the equipment support tube.


---

## Body

## 1. Introduction

In the passive advanced nuclear reactor AP1000, the contact between the internal components and the highly energetic flow causes the vibration of equipment, which leads to the fatigue of the equipment immersed in the flowing fluid, such as the secondary support column, the equipment support cylinder, the fuel rod, or the heat exchanger tube. Over time, these effects can result in great damages of the internal components and can even lead to shutdown of the entire power plant [1]. When the scale of the reactor increases or the flow rate of the coolant increases, the internal components are more likely to bend; hence, the flow-induced vibration problem becomes more serious. Since the flow-induced vibration of the internal structure significantly affects the safety of the nuclear reactor, the problem has been studied since 1950. It is difficult to carry out theoretical study because it involves both fluid structure problem and nonlinear science. At present, most of the research studies focus on experimental research and numerical simulation. In the aspect of experimental research, Ferrari et al. [2] studied the vibration of the PWR fuel assembly grid in water and air, and they mainly focused on the nonlinear behavior of the contact between the grid and fuel rod. De Pauw et al. [3] studied the axial flow-induced vibration of a single nuclear fuel rod in a lead-bismuth fast reactor. As for the numerical simulation studies, many works solve the complete equations by using the finite element or finite volume numerical method [1, 4]. We can see that most of the studies focus on different types of fuel assemblies, and boundary conditions are generally fixed at both ends. Only a few works have studied the cantilever structure such as secondary core support pillar mainly by using the experimental and simulation methods, and there still lacks the theoretical research [5]. In the early design stage of engineering problems, the theoretical derivation results are often necessary to provide guidance. According to the engineering characteristics of nuclear power plant internal components, the secondary core support pillar in the reactor can be simplified as the vibration problem of the cantilever beam structure under the action of transverse flow. Due to the large height-diameter ratio, the first-order beam mode is easy to be excited; hence, it is necessary to study the first-order beam mode separately. According to the theory of Au-Yang [6], in order to calculate the hydrodynamic added mass and frequency for the nuclear reactor internal components, it is necessary to obtain the static deflection mode shape function. Generally, many reactor internal components can be simplified as a beam with certain rigidity and elasticity. Therefore, we will focus on a single cantilever elastic beam structure and solve its static deflection mode shape function.The nonlinear deformation of elastic beams has applications in the nuclear industry, aviation, biology, and medicine. There have been many related studies in the field of solid mechanics [7–13]. Since 1970, some numerical and analytical methods have been proposed to solve these problems of elastic beams, such as the finite element method, spectral methods, and DQ method [7]. Later on, a new discrete singular convolution method (DSC), which is simpler and more efficient compared to methods mentioned above, has been widely used in solid mechanics [7, 10]. For example, DSC has been applied to a free vibration analysis of Timoshenko beams by Civalek and Kiracioglu [7]. Clamped, pinned, and sliding boundary conditions are considered in their work. From this work, we can see that the difficulties arise when considering the boundary conditions especially for free edges. Several new methods have been proposed for free boundary conditions, such as local adaptive DQ methods [11], matched interface and boundary method [12], and fictitious domain approach [13]. These methods have great application impact. However, they are often very complicated to implement. The IMB method [13] is an efficient way for implementation. However, they have not obtained efficient results by using the IMB method in their study [7]. Moreover, the fluid-related parameters such as fluid velocity field are not considered in the model and discussed in the paper cited above.The cantilever beam vibration and deflection are also studied in the field of hydraulics [14–17]. As aquatic vegetation is often submerged and highly flexible, it will deflect and interact with the ambient flow environment. It is necessary to study its interaction with the surrounding flow environment in order to determine the flow rate, morphological characteristics, and ecological conditions of the river channel where the vegetation grows [16]. Since the encroachment of vegetation has a destructive effect on the environment, the control of vegetation by hydraulic means is often adopted because it is a nonchemical approach [18]. Duan et al. [14] deduced a set of large deflection nonlinear equations of a single plant and compared the theoretical results with experimental results. The results are in good agreement [14, 15]. Since Duan et al.’s model is applicable for the problem of large-angle deflection of a cantilever beam with certain rigidity and elasticity, the model can be applied to solve the problem of static deflection mode shape function of the internal components of the AP1000 reactor. Based on Duan et al.’s model, this work firstly modeled the deflection problem of a single cantilever beam in the transverse flow. A boundary shooting method in combination with the Runge–Kutta method has been used to solve the nonlinear equations. Comparing with the existing numerical method cited above, the proposed method is efficient and easy to implement for the cantilever problem. Moreover, the method of implementation of the boundary condition is presented and proved in detail. The static cantilever deflection mode shape function of the secondary core support pillar is solved by the proposed method.


## 2. Theoretical Analysis

The secondary core support pillar in the reactor is an elastic material that bends or breaks under stress. We consider a cantilever element which has been fixed at the base without support elsewhere, with a load being applied above its base (see Figure1). The cross-section geometry of the cantilever beam affects its rigidity and the bending mechanism [16]. We assume that its cross section does not change along the vertical direction. The maximum bending moment and shear stress are located at the bottom of the cantilever beam according to classical mechanics of material theory.Figure 1
Schematic diagram of cantilever bending.Assuming a single cantilever fully exposed to the flow and no other cantilevers around it to influence the vertical velocity profile, the dominant load acting on the cantilever is the flow-induced drag force and not oscillations due to vortex shedding. The drag force acting on a single cantilever is(1)FD=12ρuz2CDA,where FD is the drag force of the fluid on the immersed part of the cantilever, CD is the drag coefficient, ρ is the density of water, uz is the vertical velocity distribution, and A is the projected area of the immersed part of the object perpendicular to the water flow plane. According to Duan et al.’s model, the differential governing equation of the cantilever beam is as follows:(2)dθds=MEI,(3)dzds=cosθ,(4)dxds=sinθ,(5)dMds=Vcosθ,(6)dVds=wzcosθ,where θ is the deflection angle, s is the distance along the curve, M is the bending moment of a section, I is the moment of inertia, E is the elastic modulus, V is the shear stress acting on the beam section, wz is the distributed load exerted by the water flow on the beam, and x and z, respectively, indicate the horizontal and vertical distances between any section of the beam and its bottom. In the literature, the small deflection hypothesis is usually used to solve the vegetation deflection problem. Since the deflection of the cantilever beam caused by the water flow can be very large, the small deflection hypothesis will cause nonnegligible errors [19]. Therefore, the small deflection hypothesis has not been adopted here. Equations (2)–(6) can be solved when the load distribution is known. Since the solution is the deflection between the location z∗=0 and z∗=1 (where z∗=z/h), in order to facilitate the application of the boundary conditions, we use the water flow depth and other parameters to nondimensionalize equations (2)–(6):(7)dθdz∗=M∗E∗I∗cosθ,(8)ds∗dz∗=1cosθ,(9)dx∗dz∗=tanθ,(10)dM∗dz∗=V∗,(11)dV∗dz∗=w∗z,where z∗=z/h, x∗=x/h, s∗=s/h, V∗=V/0.5ρh2u−2, M∗=M/0.5ρh3u−2, and E∗I∗=V/0.5ρh4u−2. Equations (7)–(11) can be rearranged into 3 dimensionless nonlinear differential equations:(12)dsinθdz∗=M∗E∗I∗,(13)dM∗dz∗=V∗,(14)dV∗dz∗=w∗z.Since the base of the cantilever is fixed, we can obtain(15)sinθz∗=0=0.We consider the situation that the top of the cantilever beam is not immersed in water, so the bending moment and shear force on the top of the cantilever should be zero; thus, the two boundary conditions can be written as(16)M∗z∗=1=0,(17)V∗z∗=1=0.


## 3. Methods

In order to solve the differential equations using the boundary shooting method and Runge–Kutta method, we have introduced three variablesx1,x2,andx3:(18)x1=sinθ,x2=dsinθdz∗,x3=d2sinθdz∗2.Then, equations (12)–(14) can be written as(19)dx1dz∗=x2,dx2dz∗=x3,dx3dz∗=1/2ρh3uz2CDDEI.Three boundary conditions can be written as(20)x1z∗=0=0,x2z∗=1=0,x3z∗=1=0.The boundary shooting method uses the initial conditions to replace the known boundary conditions to solve the problem. By constantly changing the initial value, it gradually approaches the given boundary condition. Therefore, the selection of the formula for adjusting the boundary condition is very important. We set two adjustment equations:(21)x2i0=x2i−10−αx2i−11−b2,(22)x3i0=x3i−10−αx3i−11−b3,where b2 and b3 represent the boundary conditions at the top of the cantilever z∗=1, α represents the shooting coefficient, whose value is smaller than 1 and depends on the load exerted on the free end of the beam, and i represents the iteration step number. Next, we prove that equations (21) and (22) are necessary conditions to obtain a stable solution.Firstly, according to the fourth-order Runge–Kutta method, we havex10=0. Considering arbitrary initial condition x20, we can obtain the first iteration:(23)x211=x210+l6K1x110,0+2K2x110,0+2K3x110,0+K4x110,0,where l indicates the step length and Kii=1,2,3,4 represents the fourth-order Runge–Kutta coefficient. Therefore, the iterative process at the top of the cantilever x21N is(24)x21N=x210+∑k=0N−1fx11k,k.Secondly, by introducingF1=∑k=0N−1fx11k,k, according to equation (21), the second iteration can be written as(25)x220=x210+α⋅x21N=1−α⋅x210−α⋅F1+αb2.From equations (24) and (25), the second iteration value at z∗=1 can be written as(26)x22N=x220+F2=1−α⋅x210−α⋅F1+F2+αb2,and the third iteration at z∗=1 can be written as(27)x23N=x230+F3=1−α2⋅x210−α⋅1−α⋅F1+α⋅F2+F3−b21−α2+b2.Then, we obtain the final iteration:(28)x2jN=x2j0+Fj=1−αj−1⋅x210−α⋅Fj−1−α1−α⋅Fj−2+α21−α⋅Fj−3−α21−α2⋅Fj−4−⋯+Fj−b21−αj−1+b2,where Fj=∑k=0N−1fx11k,k is the end-to-end increment calculated by the Runge–Kutta method in each iteration. When the increment is small and tends to be 1, equation (28) can be written as(29)x2jN≈−α⋅Fj−1+Fj+b2.Since the values ofFj and Fj−1 are similar after multiple iterations, it can be seen from equation (29) that x2jN will tend to be b2. x3jN will tend to be b3, thus corresponding to the correct boundary conditions. Because the magnitude of Fj is affected by the load acting on the beam, for beams with smaller loads, the value of α should be close to 1. The accurate solution can be obtained with a small number of iterations. However, if the load on the beam is large, it is necessary to reduce the value of α to avoid oscillations of the solution.


## 4. Results and Discussion

### 4.1. Validation of the Numerical Method

In order to validate the numerical method, we use the value of parameters provided by Li and Xie [16] and compare the numerical results with the analytical solution. The parameters of the cantilever and field of flow are as follows: stiffness parameter EI=2.5×10−7N⋅m2, diameter of the cantilever D=2.4×10−5m, water depth h=0.05m, water density ρ=1×103kg/m3, average velocity of flow u¯=0.5m/s, and drag coefficient CD=1.2. The vertical velocity profile is assumed to satisfy the power law [14]:(30)uzu¯=m+1mzh1/m,where m = 3.156. By solving equations (12)–(14), we can obtain the distribution of deflection angle, bending moment, and shear force along the cantilever beam. The analytical solution is expressed as follows:(31)sinθz∗=0.20z∗3.63−0.96z∗2+1.19z∗,M∗z∗=3.2×10−40.73z∗2.63−1.91z∗+1.19,V∗z∗=3.2×10−41.91−1.91z∗1.63.The numerical and analytical solutions of the deflection curve, the sine of the deflection angle, the bending moment, and the shear are shown in Figure2. We can see that the numerical solution and the analytical solution are in good agreement. The maximum relative errors in Figures 2(a)–2(d) are 3.02 × 10−6, 3.79 × 10−6, 8.84 × 10−4, and 1.95 × 10−4, respectively. Therefore, the method proposed in this study is accurate to a very high degree.Figure 2
Deflection curve, sine of the deflection angle sinθ, bending moment M∗, and shear V∗ of the cantilever beam versus z∗.
(a)(b)(c)(d)


### 4.2. Influence of the Iteration Number

The number of iterations significantly affects the calculation time and accuracy. The results of angle of the cantilever beam with various iteration numbers are shown in Figure3, while the parameters of calculation are kept the same as in Figure 2. The results show that the numerical solution with 20 iterations and 100 iterations and analytical solutions are in good agreement. Among them, the maximum relative error after 20 iterations can be less than 0.5%. We can see that an accurate numerical solution can be obtained after 20 iterations. For the calculation of this paper, the number of iterations is set to Ni = 100.Figure 3
Evolution of sinθ as a function of z∗ for different iterations.


### 4.3. Influence ofα

According to equation (29), the value of α will affect the number of iterations. Figure 4 shows the variation of Ni (number of iterations needed for the numerical solution to converge) with respect to α under different average velocities of the incoming flow. The convergence criterion is that the error between the numerical and the analytical solution of the shear force is less than 1%. We can see that, for different flow velocities, when α is greater than 0.075, the required number of iterations Ni is less than 100. When α is greater than 0.2, the required number of iterations Ni is less than 20. This result is consistent with the theoretical prediction. Convergence gets faster while α approaches more closely to 1. However, when the load is large, numerical oscillations could occur. And when 0.2≤α≤1, the convergence is relatively fast. The shooting coefficient is set to α = 0.5 in this work.Figure 4
The variation ofNi (number of iterations needed for the numerical solution to converge) with respect to α under different average velocities of the incoming flow.


### 4.4. Influence of the Profile of Flow Field

The flow field profile affects the deflection. We compared the sine of the deflection angle of the cantilever beam with the power law flow field profile and the following three profiles:(1)
Uniformuz/u¯=1(2)
Linearuz/u¯=2z/h(3)
Open channels [20], uz/u¯=1−0.522z/h2/1+1+501−z/hz/h2Table1 shows the numerical solution, the analytical results, and the relative error of the deflection angle when the dimensionless height is z∗=0.3. We can see that the velocity distribution has a great impact on the deflection. The power law distribution and the velocity distribution in open channels are closer to the actual situation. If a uniform distribution is used, a larger deflection angle will be obtained; if a linear distribution is used, a smaller deflection angle will be obtained. We can conclude that our method can effectively solve the governing equations of cantilever beams subjected to complex loads.Table 1
The numerical solution, the analytical results, and the relative error of the deflection angle for various flow fields.
Flow profileNumerical solutionAnalytical solutionRelative errorUniform0.29740.29930.0063Linear0.08100.08200.0121Open channel0.1729——Power law0.21360.21250.0052


### 4.5. Solution of the Beam Mode Function

The model and method proposed in this paper are used to solve the mode shape function of the secondary core support pillar in the AP1000 reactor. It can be considered as a cantilever beam; thus, we solve equation (19). According to the literature [21], it can be assumed that the incoming flow acting on the support pillar is uniform. Other parameters are set as follows: the average flow velocity of the coolant in the reactor u¯=16m/s, the support pillar diameter D=90mm, the water support column height h=500mm, and the drag force coefficient CD=1.2. The bending stiffness EI=6.4×105N⋅m2 (since the material is SA-479 stainless steel bar). Figure 5 shows the numerical solution under the uniform incoming flow condition. According to the relationship x∗z∗=∫0z∗tanθdz∗, the following static deflection mode function relationship is obtained (see the red dashed line in Figure 5):(32)x∗z∗=1.1×10−4z∗4−4.5×10−4z∗3+6.7×10−4z∗2+7×10−7z∗.Figure 5
The deflection of the secondary core support pillar under the operating condition of the reactor, where the red dotted line represents equation (32).


## 4.1. Validation of the Numerical Method

In order to validate the numerical method, we use the value of parameters provided by Li and Xie [16] and compare the numerical results with the analytical solution. The parameters of the cantilever and field of flow are as follows: stiffness parameter EI=2.5×10−7N⋅m2, diameter of the cantilever D=2.4×10−5m, water depth h=0.05m, water density ρ=1×103kg/m3, average velocity of flow u¯=0.5m/s, and drag coefficient CD=1.2. The vertical velocity profile is assumed to satisfy the power law [14]:(30)uzu¯=m+1mzh1/m,where m = 3.156. By solving equations (12)–(14), we can obtain the distribution of deflection angle, bending moment, and shear force along the cantilever beam. The analytical solution is expressed as follows:(31)sinθz∗=0.20z∗3.63−0.96z∗2+1.19z∗,M∗z∗=3.2×10−40.73z∗2.63−1.91z∗+1.19,V∗z∗=3.2×10−41.91−1.91z∗1.63.The numerical and analytical solutions of the deflection curve, the sine of the deflection angle, the bending moment, and the shear are shown in Figure2. We can see that the numerical solution and the analytical solution are in good agreement. The maximum relative errors in Figures 2(a)–2(d) are 3.02 × 10−6, 3.79 × 10−6, 8.84 × 10−4, and 1.95 × 10−4, respectively. Therefore, the method proposed in this study is accurate to a very high degree.Figure 2
Deflection curve, sine of the deflection angle sinθ, bending moment M∗, and shear V∗ of the cantilever beam versus z∗.
(a)(b)(c)(d)


## 4.2. Influence of the Iteration Number

The number of iterations significantly affects the calculation time and accuracy. The results of angle of the cantilever beam with various iteration numbers are shown in Figure3, while the parameters of calculation are kept the same as in Figure 2. The results show that the numerical solution with 20 iterations and 100 iterations and analytical solutions are in good agreement. Among them, the maximum relative error after 20 iterations can be less than 0.5%. We can see that an accurate numerical solution can be obtained after 20 iterations. For the calculation of this paper, the number of iterations is set to Ni = 100.Figure 3
Evolution of sinθ as a function of z∗ for different iterations.


## 4.3. Influence ofα

According to equation (29), the value of α will affect the number of iterations. Figure 4 shows the variation of Ni (number of iterations needed for the numerical solution to converge) with respect to α under different average velocities of the incoming flow. The convergence criterion is that the error between the numerical and the analytical solution of the shear force is less than 1%. We can see that, for different flow velocities, when α is greater than 0.075, the required number of iterations Ni is less than 100. When α is greater than 0.2, the required number of iterations Ni is less than 20. This result is consistent with the theoretical prediction. Convergence gets faster while α approaches more closely to 1. However, when the load is large, numerical oscillations could occur. And when 0.2≤α≤1, the convergence is relatively fast. The shooting coefficient is set to α = 0.5 in this work.Figure 4
The variation ofNi (number of iterations needed for the numerical solution to converge) with respect to α under different average velocities of the incoming flow.


## 4.4. Influence of the Profile of Flow Field

The flow field profile affects the deflection. We compared the sine of the deflection angle of the cantilever beam with the power law flow field profile and the following three profiles:(1)
Uniformuz/u¯=1(2)
Linearuz/u¯=2z/h(3)
Open channels [20], uz/u¯=1−0.522z/h2/1+1+501−z/hz/h2Table1 shows the numerical solution, the analytical results, and the relative error of the deflection angle when the dimensionless height is z∗=0.3. We can see that the velocity distribution has a great impact on the deflection. The power law distribution and the velocity distribution in open channels are closer to the actual situation. If a uniform distribution is used, a larger deflection angle will be obtained; if a linear distribution is used, a smaller deflection angle will be obtained. We can conclude that our method can effectively solve the governing equations of cantilever beams subjected to complex loads.Table 1
The numerical solution, the analytical results, and the relative error of the deflection angle for various flow fields.
Flow profileNumerical solutionAnalytical solutionRelative errorUniform0.29740.29930.0063Linear0.08100.08200.0121Open channel0.1729——Power law0.21360.21250.0052


## 4.5. Solution of the Beam Mode Function

The model and method proposed in this paper are used to solve the mode shape function of the secondary core support pillar in the AP1000 reactor. It can be considered as a cantilever beam; thus, we solve equation (19). According to the literature [21], it can be assumed that the incoming flow acting on the support pillar is uniform. Other parameters are set as follows: the average flow velocity of the coolant in the reactor u¯=16m/s, the support pillar diameter D=90mm, the water support column height h=500mm, and the drag force coefficient CD=1.2. The bending stiffness EI=6.4×105N⋅m2 (since the material is SA-479 stainless steel bar). Figure 5 shows the numerical solution under the uniform incoming flow condition. According to the relationship x∗z∗=∫0z∗tanθdz∗, the following static deflection mode function relationship is obtained (see the red dashed line in Figure 5):(32)x∗z∗=1.1×10−4z∗4−4.5×10−4z∗3+6.7×10−4z∗2+7×10−7z∗.Figure 5
The deflection of the secondary core support pillar under the operating condition of the reactor, where the red dotted line represents equation (32).


## 5. Conclusions

This work studies the mode shape function of the secondary core support pillar in the AP1000 reactor under the action of transverse flow by using the theoretical analysis and numerical calculation. The following conclusions are drawn:(1)
In light of the large deflection nonlinear differential equations proposed by Duan et al. [14], a set of nonlinear differential equations for the analysis of the single cantilever beam was proposed, and the applicability of the model to the cantilever beam deflection problem was demonstrated.(2)
We use a boundary shooting framework in combination with the Runge–Kutta method to solve the aforementioned nonlinear differential control equations and obtain the numerical solutions of the deflection curve, deflection angle, bending moment, and shear force of the cantilever beam. The results of the numerical simulation agree with the analytical solution to a very high degree, which demonstrates the effectiveness of the simulation method. Meanwhile, the influence of the number of iterations and the shooting coefficient is studied, which provides a basis for the selection of parameters.(3)
We have quantitatively studied the influence of velocity profile on the deflection angle of the cantilever beam. Comparing with the actual situation, if the uniform profile is adopted, a larger deflection angle will be obtained; if the linear profile is adopted, a smaller deflection angle will be obtained. Our method can solve the governing equations with cantilever subjected to complex loads. In this case, an analytical solution cannot be obtained.(4)
The static shape function of the secondary core support pillar in the AP1000 reactor is solved. This work lays a foundation for the calculation of hydrodynamic added mass and frequency for the nuclear reactor internal components and their structural integrity evaluation.



---

*Source: 4861312-2021-11-30.xml*