**3. Application to the boiling two-phase flow analysis in a simulated fuel assembly excited by oscillation acceleration**

Boiling two-phase flow in a simulated fuel assembly excited by oscillation acceleration was performed by the improved ACE-3D in order to investigate how the three-dimensional behavior of boiling two-phase flow in a fuel assembly under oscillation conditions is evaluated by the improved ACE-3D.

### **3.1 Computational condition**

In this analysis, a 7 × 7 fuel assembly in a current BWR core is simulated, as shown in Fig. 9. Fuel rod diameter is 10.8 mm; the narrowest gap between fuel rods is 4.4 mm, and the axial heat length is 3.66 m.

Four subchannels surrounded by nine fuel rods without channel boxes are adopted as the computational domain shown in Fig. 9; this is the smallest domain that can describe the

Development of an Analytical Method on Water-Vapor Boiling Two-Phase

Fig. 11. Time variation in oscillation acceleration

**3.2 Results and discussions** 

directions along the X and Y axes was high.

Fig. 12. Isosurface distribution of void fraction

Flow Characteristics in BWR Fuel Assemblies Under Earthquake Condition 169

In this analysis, in-phase sine wave acceleration was applied in the X and Y directions as shown by the black arrow in Fig. 10. The magnitude and oscillation period of the oscillation acceleration in the X and Y directions were 400 Gal and 0.2 s, respectively, as shown in Fig. 11; these values are based on actual earthquake data measured in the Kashiwazaki-Kariwa nuclear power plant. The computable physical time in this analysis after applying the

After applying the oscillation acceleration, the void fraction distribution fluctuated with the same period as the oscillation acceleration. Figure 12 shows the isosurfaces of the void fraction at *t* = 0.8 s. In the whole area where boiling occurs, the void fraction in the center of the subchannel was relatively low, and the void fraction concentrated in the positive

Figure 13 shows the time variation in the void fraction at Z = 2.3 m in the upstream region of Fig. 12. The oscillation acceleration did not act at *t* = 0.8 s and *t* = 0.9 s and acted in the direction of the black arrow shown in Fig. 13(b). At *t* = 0.8 s, a high void fraction could be seen near the fuel-rod surface in the narrowest region between the fuel rods, as indicated by

oscillation acceleration was 1 s based on the results described in section 2.

three-dimensional behavior of boiling two-phase flow. This computational domain was determined to reflect the basic thermal-hydraulic characteristics in fuel assemblies under earthquake conditions.

In this domain, single-phase water flows in from the bottom of the channel with a mass velocity of 1673 kg/m2s and inlet temperature of 549.15 K. At the exit of the computational domain, pressure was fixed at 7.1 MPa. The mass velocity, inlet temperature, and exit pressure reflect the operating conditions in a current BWR core. The core thermal power is 351.9 W. The axial power distribution of the fuel-rod surfaces is shown in Fig. 9 and it simulates the power distribution in a current BWR core.

Figure 10 shows the boundary conditions and the computational block divisions. Here, the non-slip condition is set for each fuel-rod surface, and the slip condition is set for each symmetric boundary. In this analysis, the computational domain was divided into 9 blocks. The computational grids in each block have 10 and 256 grids in the radial and axial directions, respectively. The number of grids in the peripheral direction is as follows: 30 in block 1, block 3, block 7, and block 9; 60 in block 2, block 4, block 6, and block 8; and 120 in block 5.

In this study, the boiling two-phase flow analysis was performed under steady-state conditions to obtain a steady boiling two-phase flow. Subsequently, oscillation acceleration was applied. The time when the oscillation acceleration was applied is regarded as *t* = 0 s.

Fig. 9. Computational domain and axial power ratio

Fig. 10. Boundary conditions and computational block division

Fig. 11. Time variation in oscillation acceleration

In this analysis, in-phase sine wave acceleration was applied in the X and Y directions as shown by the black arrow in Fig. 10. The magnitude and oscillation period of the oscillation acceleration in the X and Y directions were 400 Gal and 0.2 s, respectively, as shown in Fig. 11; these values are based on actual earthquake data measured in the Kashiwazaki-Kariwa nuclear power plant. The computable physical time in this analysis after applying the oscillation acceleration was 1 s based on the results described in section 2.

#### **3.2 Results and discussions**

168 Nuclear Reactors

three-dimensional behavior of boiling two-phase flow. This computational domain was determined to reflect the basic thermal-hydraulic characteristics in fuel assemblies under

In this domain, single-phase water flows in from the bottom of the channel with a mass velocity of 1673 kg/m2s and inlet temperature of 549.15 K. At the exit of the computational domain, pressure was fixed at 7.1 MPa. The mass velocity, inlet temperature, and exit pressure reflect the operating conditions in a current BWR core. The core thermal power is 351.9 W. The axial power distribution of the fuel-rod surfaces is shown in Fig. 9 and it

Figure 10 shows the boundary conditions and the computational block divisions. Here, the non-slip condition is set for each fuel-rod surface, and the slip condition is set for each symmetric boundary. In this analysis, the computational domain was divided into 9 blocks. The computational grids in each block have 10 and 256 grids in the radial and axial directions, respectively. The number of grids in the peripheral direction is as follows: 30 in block 1, block 3, block 7, and block 9; 60 in block 2, block 4, block 6, and block 8; and 120 in

In this study, the boiling two-phase flow analysis was performed under steady-state conditions to obtain a steady boiling two-phase flow. Subsequently, oscillation acceleration was applied. The time when the oscillation acceleration was applied is regarded as *t* = 0 s.

earthquake conditions.

block 5.

simulates the power distribution in a current BWR core.

Fig. 9. Computational domain and axial power ratio

Fig. 10. Boundary conditions and computational block division

After applying the oscillation acceleration, the void fraction distribution fluctuated with the same period as the oscillation acceleration. Figure 12 shows the isosurfaces of the void fraction at *t* = 0.8 s. In the whole area where boiling occurs, the void fraction in the center of the subchannel was relatively low, and the void fraction concentrated in the positive directions along the X and Y axes was high.

Fig. 12. Isosurface distribution of void fraction

Figure 13 shows the time variation in the void fraction at Z = 2.3 m in the upstream region of Fig. 12. The oscillation acceleration did not act at *t* = 0.8 s and *t* = 0.9 s and acted in the direction of the black arrow shown in Fig. 13(b). At *t* = 0.8 s, a high void fraction could be seen near the fuel-rod surface in the narrowest region between the fuel rods, as indicated by

Development of an Analytical Method on Water-Vapor Boiling Two-Phase

of the void fraction was relatively small.

Flow Characteristics in BWR Fuel Assemblies Under Earthquake Condition 171

void fraction could also be seen away from the fuel rod surface, as shown by the blue circles in Fig. 14(b). The high void fraction in the regions marked by the blue circles in Fig. 14(b) split, and the high void fraction in the blue circles in Fig. 14(c) was formed at *t* = 0.82 s. The high void fraction regions represented by red and blue circles in Fig. 14(c) moved in a direction opposite to the black arrow as shown Fig. 14(d). While the void fraction regions indicated by the red and blue circles in Fig. 14(d) decreased as shown in Fig. 14(e), high void fraction was concentrated in the regions marked by red circles; high void fraction could also be seen in the regions away from the fuel rod surface, such as the regions indicated by the blue circles at *t* = 0.9 s, as shown in Fig. 14(f). Near the fuel rod surface, void fluctuation with a different period to that of the oscillation acceleration was seen while the magnitude

Figure 15 shows the time variation in vapor velocity at Z = 3.4 m. The black arrow shows the direction in which the oscillation acceleration acts at each time. At *t* = 0.78 s, the vapor velocity acted in the direction indicated by red arrows in Fig. 15(a); this direction is opposite, but not parallel to, the black arrow. Between *t* = 0.8 s and *t* = 0.82 s, in spite of the changing direction of the oscillation acceleration, the vapor velocity decreased but still acted in the direction of the red arrows, shown in Fig. 15(b) and Fig. 15(c), because of the effect of inertia. At *t* = 0.8 s, the high void fraction indicated by blue circles in Fig. 14(b) was moved by the vapor velocity. This caused the high void fraction shown in Fig. 15(b) to split, and the

high void fraction represented by blue circles in Fig. 14(c) and Fig. 14(d) was formed.

Fig. 15. Time variation in vapor velocity vector at Z = 3.4 m

fluctuation.

(a) *t* = 0.78 s (*fx* = *fy* ≠ 0) (b) *t* = 0.8 s (*fx* = *fy* = 0) c) *t* = 0.82 s (*fx* = *fy* ≠ 0)

Figure 16 shows the time variation in the acceleration vector of the lift force at Z = 3.4 m. The black arrow shows the direction in which the oscillation acceleration acts. The lift force in the red circles in Fig. 16(a) to Fig. 16(c) acted in a direction facing away from the fuel-rod surface. Hence, the vapor velocity was directed along the red arrow, as shown in Fig. 15; this is a direction opposite, but not parallel, to the black arrow. A high void fraction could be seen away from the fuel-rod surface, as shown by blue circles in Fig. 14(b). The lift force in the regions of the blue circles in Fig. 16(b) and Fig. 16(c) also acts in a direction facing away from the fuel rod surface. Hence, a high void fraction could be seen away from the fuel rod surface, as shown by red and blue circles in Fig. 14(c) and Fig. 14(d). Near the fuel rod surface, the magnitude and direction of the lift force were not uniform along the fuel rod and fluctuated with a period different from that of the oscillation acceleration. Consequently, void fraction fluctuation at Z = 3.4 m was significantly dependent on lift force

red circles in Fig. 13(a). At *t* = 0.85 s, a high void fraction moved in a direction opposite to the oscillation acceleration as shown in Fig. 13(b). At *t* = 0.9 s, a high void fraction could be seen near the fuel-rod surface in the regions marked by red circles in Fig. 13(c). This indicates that the magnitude of void fraction fluctuation at Z = 2.3 m is particularly large near the fuel-rod surface. This tendency of void fraction fluctuation is the same between *t* = 0.9 s and 1.0 s, when the oscillation acceleration acted in a direction opposite to that of the black arrow.

Figure 14 shows the time variation in the void fraction at Z = 3.4 m in the downstream region of Fig. 12. The black arrow shows the direction in which the oscillation acceleration acts. At *t* = 0.78 s, the vapor phase moved in a direction opposite to the black arrow and concentrated in the regions marked by red circles, shown in Fig. 14(a). At *t* = 0.8 s, the void fraction in the region marked by the red circles in Fig. 14(b) increased. In addition, a high

Fig. 13. Time variation in void fraction at Z = 2.3m

Fig. 14. Time variation in the void fraction at Z = 3.4 m

red circles in Fig. 13(a). At *t* = 0.85 s, a high void fraction moved in a direction opposite to the oscillation acceleration as shown in Fig. 13(b). At *t* = 0.9 s, a high void fraction could be seen near the fuel-rod surface in the regions marked by red circles in Fig. 13(c). This indicates that the magnitude of void fraction fluctuation at Z = 2.3 m is particularly large near the fuel-rod surface. This tendency of void fraction fluctuation is the same between *t* = 0.9 s and 1.0 s, when the oscillation acceleration acted in a direction opposite to that of the

Figure 14 shows the time variation in the void fraction at Z = 3.4 m in the downstream region of Fig. 12. The black arrow shows the direction in which the oscillation acceleration acts. At *t* = 0.78 s, the vapor phase moved in a direction opposite to the black arrow and concentrated in the regions marked by red circles, shown in Fig. 14(a). At *t* = 0.8 s, the void fraction in the region marked by the red circles in Fig. 14(b) increased. In addition, a high

(a) *t* = 0.8 s (*fx* = *f y* =0) (b) *t* = 0.85 s (*fx* = *fy* ≠ 0) (c) *t* = 0.9 s (*fx* = *fy* = 0)

(a) *t* = 0.78 s (*fx* = *fy* ≠ 0) (b) *t* = 0.8 s (*fx* = *fy* = 0) (c) *t* = 0.82 s (*fx* = *fy* ≠ 0)

(d) *t* = 0.84 s (*fx* = *fy* ≠ 0) (e) *t* = 0.86 s (*fx* = *fy* ≠ 0) (f) *t* = 0.9 s (*fx* = *fy* = 0)

Fig. 13. Time variation in void fraction at Z = 2.3m

Fig. 14. Time variation in the void fraction at Z = 3.4 m

black arrow.

void fraction could also be seen away from the fuel rod surface, as shown by the blue circles in Fig. 14(b). The high void fraction in the regions marked by the blue circles in Fig. 14(b) split, and the high void fraction in the blue circles in Fig. 14(c) was formed at *t* = 0.82 s. The high void fraction regions represented by red and blue circles in Fig. 14(c) moved in a direction opposite to the black arrow as shown Fig. 14(d). While the void fraction regions indicated by the red and blue circles in Fig. 14(d) decreased as shown in Fig. 14(e), high void fraction was concentrated in the regions marked by red circles; high void fraction could also be seen in the regions away from the fuel rod surface, such as the regions indicated by the blue circles at *t* = 0.9 s, as shown in Fig. 14(f). Near the fuel rod surface, void fluctuation with a different period to that of the oscillation acceleration was seen while the magnitude of the void fraction was relatively small.

Figure 15 shows the time variation in vapor velocity at Z = 3.4 m. The black arrow shows the direction in which the oscillation acceleration acts at each time. At *t* = 0.78 s, the vapor velocity acted in the direction indicated by red arrows in Fig. 15(a); this direction is opposite, but not parallel to, the black arrow. Between *t* = 0.8 s and *t* = 0.82 s, in spite of the changing direction of the oscillation acceleration, the vapor velocity decreased but still acted in the direction of the red arrows, shown in Fig. 15(b) and Fig. 15(c), because of the effect of inertia. At *t* = 0.8 s, the high void fraction indicated by blue circles in Fig. 14(b) was moved by the vapor velocity. This caused the high void fraction shown in Fig. 15(b) to split, and the high void fraction represented by blue circles in Fig. 14(c) and Fig. 14(d) was formed.

Fig. 15. Time variation in vapor velocity vector at Z = 3.4 m

Figure 16 shows the time variation in the acceleration vector of the lift force at Z = 3.4 m. The black arrow shows the direction in which the oscillation acceleration acts. The lift force in the red circles in Fig. 16(a) to Fig. 16(c) acted in a direction facing away from the fuel-rod surface. Hence, the vapor velocity was directed along the red arrow, as shown in Fig. 15; this is a direction opposite, but not parallel, to the black arrow. A high void fraction could be seen away from the fuel-rod surface, as shown by blue circles in Fig. 14(b). The lift force in the regions of the blue circles in Fig. 16(b) and Fig. 16(c) also acts in a direction facing away from the fuel rod surface. Hence, a high void fraction could be seen away from the fuel rod surface, as shown by red and blue circles in Fig. 14(c) and Fig. 14(d). Near the fuel rod surface, the magnitude and direction of the lift force were not uniform along the fuel rod and fluctuated with a period different from that of the oscillation acceleration. Consequently, void fraction fluctuation at Z = 3.4 m was significantly dependent on lift force fluctuation.

Development of an Analytical Method on Water-Vapor Boiling Two-Phase

Fig. 18. Time variation in the bubble diameter at Z = 3.4 m

analysis code ACE-3D under earthquake conditions.

fraction fluctuation.

oscillation conditions.

distribution.

**4. Conclusion** 

oscillation acceleration.

time during the earthquake.

Flow Characteristics in BWR Fuel Assemblies Under Earthquake Condition 173

dependent on the void fraction, and a local high void fraction results in a local large bubble diameter. Thus, a strongly inhomogeneous bubble diameter distribution results from void

It is necessary to adequately evaluate the influence of the void fraction upon bubble diameter in order to avoid a strongly inhomogeneous bubble diameter distribution under

According to our results, void fraction fluctuation in the downstream region is significantly dependent on the lift force caused by a strongly inhomogeneous bubble diameter

(a) *t* = 0.78 s (*fx* = *fy* ≠ 0) (b) *t* = 0.8 s (*fx* = *fy* = 0) (c) *t* = 0.82 s (*fx* = *fy* ≠ 0)

A new external force term, which can simulate the oscillation acceleration, was added to the momentum conservation equations in order to apply the three-dimensional two-fluid model

A boiling two-phase flow excited by applying vertical and horizontal oscillation acceleration was simulated in order to confirm that the simulation can be performed under oscillation conditions. It was confirmed that the void fraction fluctuation with the same period as that of the oscillation acceleration could be calculated in the case of both horizontal and vertical

The influence of the oscillation period of the oscillation acceleration on the boiling twophase flow behavior in a fuel assembly was investigated in order to evaluate the highest frequency necessary for the improved method to be consistent with the time-series data of oscillation acceleration and the shortest period of oscillation acceleration for which the boiling two-phase flow shows quasi-steady time variation. It was confirmed that a boiling two-phase flow analysis consistent with the time-series data of oscillation acceleration and with a time interval greater than 0.01 s, can be performed. It was also shown that an effective analysis can be performed by extracting an earthquake motion of about 1 s at any

The three-dimensional behavior of boiling two-phase flow in a fuel assembly under oscillation conditions was evaluated using a simulated fuel assembly excited by oscillation acceleration. On the basis of this evaluation, it was confirmed that void fraction fluctuation

Figure 17 shows the time variation in the Eotvos number at Z = 3.4 m and also shows a range of Eotvos number from 4 to 10 for which the effect of bubble deformation upon the lift force is dependent upon Eotvos number, as shown in Eq. (7). The black arrow shows the direction in which the oscillation acceleration acts. The red and blue circles in Fig. 17 correspond to regions where the magnitude of the lift force was large; the lift force acted in a direction facing away from the fuel rod surface, as shown in Fig. 16. In these regions, the effect of bubble deformation on the lift force was dominant because the Eotvos number exhibited high values. Near the fuel rod surface, the Eotvos numbers less than 4 and greater than 10 were mixed, indicating that the magnitude and direction of the lift force were not uniform near the fuel rod surface.

Fig. 16. Time variation in lift force vector at Z = 3.4 m

Fig. 17. Time variation in the Eotvos number at Z = 3.4 m

Figure 18 shows the variation in bubble diameter with time at Z = 3.4 m. The black arrow shows the direction in which the oscillation acceleration acts. Bubble diameters greater than 7 mm are distributed in the region where the Eotvos number is greater than 10, as shown in Fig. 17. The bubble diameter distribution shown in Fig. 18 is strongly inhomogeneous and physically invalid because large bubble diameters are mainly observed in small regions in the subchannel, while small bubble diameters of less than 3 mm are observed in the center of the subchannel. This strongly inhomogeneous bubble diameter distribution resulted in locally high Eotvos numbers and fluctuation in the direction of the lift force vectors.

The region where large bubble diameters are seen corresponds to the region of high void fraction, as shown in Fig. 14. According to Eq. (8), the bubble diameter is significantly dependent on the void fraction, and a local high void fraction results in a local large bubble diameter. Thus, a strongly inhomogeneous bubble diameter distribution results from void fraction fluctuation.

It is necessary to adequately evaluate the influence of the void fraction upon bubble diameter in order to avoid a strongly inhomogeneous bubble diameter distribution under oscillation conditions.

According to our results, void fraction fluctuation in the downstream region is significantly dependent on the lift force caused by a strongly inhomogeneous bubble diameter distribution.

Fig. 18. Time variation in the bubble diameter at Z = 3.4 m
