Tag: flow3d

Schematic diagram of overlapping powder particles.

Because the particle size of the powder used for HP-LPBF is much smaller than 100 μm, the effect of van der Waals forces must be considered. Therefore, the cohesive force �jkr of the Hertz–Mindlin model was used instead of van der Waals forces,⁷² with the following expression:�jkr=−4��0�*�1.5+4�*3�*�3,

(4)1�*=(1−��2)��+(1−��2)��,

(5)1�*=1��+1��,

(6)

where �* is the equivalent Young’s modulus in GPa; �* is the equivalent particle radius in m; �0 is the surface energy of the powder particles in J/m²; α is the contact radius in m; �� and �� are the Young’s modulus of the unit particles � and �⁠, respectively, in GPa; and �� and �� are the Poisson’s ratio of the unit particles � and �⁠, respectively.

2. Model building

Figure 2 shows a 3D powder bed model generated using DEM with a coarse powder geometry of 1000 × 400 × 30 μm³. The powder layer thickness is 30 μm, and the powder bed porosity is 40%. The average particle size of this spherical powder is 31.7 μm and is normally distributed in the range of 15–53 μm. The geometry of the fine powder was 1000 × 400 × 20 μm³, with a layer thickness of 20 μm, and the powder bed porosity of 40%. The average particle size of this spherical powder is 11.5 μm and is normally distributed in the range of 5–25 μm. After the 3D powder bed model is generated, it needs to be imported into the CFD simulation software for calculation, and the imported geometric model is shown in Fig. 3. This geometric model is mainly composed of three parts: protective gas, powder bed, and substrate. Under the premise of ensuring the accuracy of the calculation, the mesh size is set to 3 μm, and the total number of coarse powder meshes is 1 704 940. The total number of fine powder meshes is 3 982 250.

FIG. 2.

Three-dimensional powder bed model: (a) coarse powder, (b) fine powder.

FIG. 3.

Geometric modeling of the powder bed computational domain: (a) coarse powder, (b) fine powder.

B. Modeling of fluid mechanics simulation

In order to solve the flow, melting, and solidification problems involved in HP-LPBF molten pool, the study must follow the three governing equations of conservation of mass, conservation of energy, and conservation of momentum.⁷³ The VOF method, which is the most widely used in fluid dynamics, is used to solve the molten pool dynamics model.

1. VOF

VOF is a method for tracking the free interface between the gas and liquid phases on the molten pool surface. The core idea of the method is to define a volume fraction function F within each grid, indicating the proportion of the grid space occupied by the material, 0 ≤ F ≤ 1 in Fig. 4. Specifically, when F = 0, the grid is empty and belongs to the gas-phase region; when F = 1, the grid is completely filled with material and belongs to the liquid-phase region; and when 0 < F < 1, the grid contains free surfaces and belongs to the mixed region. The direction normal to the free surface is the direction of the fastest change in the volume fraction F (the direction of the gradient of the volume fraction), and the direction of the gradient of the volume fraction can be calculated from the values of the volume fractions in the neighboring grids.⁷⁴ The equations controlling the VOF are expressed as follows:𝛻��+�⋅(��→)=0,

(7)

where t is the time in s and �→ is the liquid velocity in m/s.

FIG. 4.

Schematic diagram of VOF.

The material parameters of the mixing zone are altered due to the inclusion of both the gas and liquid phases. Therefore, in order to represent the density of the mixing zone, the average density �¯ is used, which is expressed as follows:⁷²�¯=(1−�1)�gas+�1�metal,

(8)

where �1 is the proportion of liquid phase, �gas is the density of protective gas in kg/m³, and �metal is the density of metal in kg/m³.

2. Control equations and boundary conditions

Figure 5 is a schematic diagram of the HP-LPBF melting process. First, the laser light strikes a localized area of the material and rapidly heats up the area. Next, the energy absorbed in the region is diffused through a variety of pathways (heat conduction, heat convection, and surface radiation), and this process triggers complex phase transition phenomena (melting, evaporation, and solidification). In metals undergoing melting, the driving forces include surface tension and the Marangoni effect, recoil due to evaporation, and buoyancy due to gravity and uneven density. The above physical phenomena interact with each other and do not occur independently.

FIG. 5.

Schematic diagram of HP-LPBF melting process.

Laser heat sourceThe Gaussian surface heat source model is used as the laser heat source model with the following expression:�=2�0��2exp(−2�12��2),(9)where � is the heat flow density in W/m², �0 is the absorption rate of SS316L, �� is the radius of the laser focal spot in m, and �1 is the radial distance from the center of the laser focal spot in m. The laser focal spot can be used for a wide range of applications.
Energy absorptionThe formula for calculating the laser absorption �0 of SS316L is as follows:�0=0.365(�0[1+�0(�−20)]/�)0.5,(10)where �0 is the direct current resistivity of SS316L at 20 °C in Ω m, �0 is the resistance temperature coefficient in ppm/°C, � is the temperature in °C, and � is the laser wavelength in m.
Heat transferThe basic principle of heat transfer is conservation of energy, which is expressed as follows:𝛻𝛻𝛻�(��)��+�·(��→�)=�·(�0��)+��,(11)where � is the density of liquid phase SS316L in kg/m³, �� is the specific heat capacity of SS316L in J/(kg K), 𝛻� is the gradient operator, t is the time in s, T is the temperature in K, 𝛻�� is the temperature gradient, �→ is the velocity vector, �0 is the coefficient of thermal conduction of SS316L in W/(m K), and �� is the thermal energy dissipation term in the molten pool.
Molten pool flowThe following three conditions need to be satisfied for the molten pool to flow:
- Conservation of mass with the following expression:𝛻�·(��→)=0.(12)
- Conservation of momentum (Navier–Stokes equation) with the following expression:𝛻𝛻𝛻𝛻��→��+�(�→·�)�→=�·[−pI+�(��→+(��→)�)]+�,(13)where � is the pressure in Pa exerted on the liquid phase SS316L microelement, � is the unit matrix, � is the fluid viscosity in N s/m², and � is the volumetric force (gravity, atmospheric pressure, surface tension, vapor recoil, and the Marangoni effect).
- Conservation of energy, see Eq. (11)
Surface tension and the Marangoni effectThe effect of temperature on the surface tension coefficient is considered and set as a linear relationship with the following expression:�=�0−��dT(�−��),(14)where � is the surface tension of the molten pool at temperature T in N/m, �� is the melting temperature of SS316L in K, �0 is the surface tension of the molten pool at temperature �� in Pa, and σdσ/ dT is the surface tension temperature coefficient in N/(m K).In general, surface tension decreases with increasing temperature. A temperature gradient causes a gradient in surface tension that drives the liquid to flow, known as the Marangoni effect.
Metal vapor recoilAt higher input energy densities, the maximum temperature of the molten pool surface reaches the evaporation temperature of the material, and a gasification recoil pressure occurs vertically downward toward the molten pool surface, which will be the dominant driving force for the molten pool flow.⁷⁵ The expression is as follows:��=0.54�� exp ��−��0��,(15)where �� is the gasification recoil pressure in Pa, �� is the ambient pressure in kPa, �� is the latent heat of evaporation in J/kg, �0 is the gas constant in J/(mol K), T is the surface temperature of the molten pool in K, and T_e is the evaporation temperature in K.
Solid–liquid–gas phase transitionWhen the laser hits the powder layer, the powder goes through three stages: heating, melting, and solidification. During the solidification phase, mutual transformations between solid, liquid, and gaseous states occur. At this point, the latent heat of phase transition absorbed or released during the phase transition needs to be considered.⁶⁸ The phase transition is represented based on the relationship between energy and temperature with the following expression:�=��,(�<��),�(��)+�−��−��,(��<�<��)�(��)+(�−��)��,(��<�),,(16)where �� and �� are solid and liquid phase density, respectively, of SS316L in kg/m³. �� and �� unit volume of solid and liquid phase-specific heat capacity, respectively, of SS316L in J/(kg K). �� and ��⁠, respectively, are the solidification temperature and melting temperature of SS316L in K. �� is the latent heat of the phase transition of SS316L melting in J/kg.

3. Assumptions

The CFD model was computed using the commercial software package FLOW-3D.⁷⁶ In order to simplify the calculation and solution process while ensuring the accuracy of the results, the model makes the following assumptions:

It is assumed that the effects of thermal stress and material solid-phase thermal expansion on the calculation results are negligible.
The molten pool flow is assumed to be a Newtonian incompressible laminar flow, while the effects of liquid thermal expansion and density on the results are neglected.
It is assumed that the surface tension can be simplified to an equivalent pressure acting on the free surface of the molten pool, and the effect of chemical composition on the results is negligible.
Neglecting the effect of the gas flow field on the molten pool.
The mass loss due to evaporation of the liquid metal is not considered.
The influence of the plasma effect of the molten metal on the calculation results is neglected.

It is worth noting that the formulation of assumptions requires a trade-off between accuracy and computational efficiency. In the above models, some physical phenomena that have a small effect or high difficulty on the calculation results are simplified or ignored. Such simplifications make numerical simulations more efficient and computationally tractable, while still yielding accurate results.

4. Initial conditions

The preheating temperature of the substrate was set to 393 K, at which time all materials were in the solid state and the flow rate was zero.

5. Material parameters

The material used is SS316L and the relevant parameters required for numerical simulations are shown in Table I.^46,77,78

TABLE I.

SS316L-related parameters.

Property	Symbol	Value
Density of solid metal (kg/m³)	�metal	7980
Solid phase line temperature (K)	��	1658
Liquid phase line temperature (K)	��	1723
Vaporization temperature (K)	��	3090
Latent heat of melting (⁠ J/kg⁠)	��	2.60×105
Latent heat of evaporation (⁠ J/kg⁠)	��	7.45×106
Surface tension of liquid phase (N /m⁠)	�	1.60
Liquid metal viscosity (kg/m s)	��	6×10−3
Gaseous metal viscosity (kg/m s)	�gas	1.85×10−5
Temperature coefficient of surface tension (N/m K)	��/�T	0.80×10−3
Molar mass (⁠ kg/mol⁠)	M	0.05 593
Emissivity	�	0.26
Laser absorption	�0	0.35
Ambient pressure (kPa)	��	101 325
Ambient temperature (K)	�0	300
Stefan–Boltzmann constant (W/m² K⁴)	�	5.67×10−8
Thermal conductivity of metals (⁠ W/m K⁠)	�	24.55
Density of protective gas (kg/m³)	�gas	1.25
Coefficient of thermal expansion (/K)	��	16×10−6
Generalized gas constant (⁠ J/mol K⁠)	R	8.314

III. RESULTS AND DISCUSSION

With the objective of studying in depth the evolutionary patterns of single-track and double-track molten pool development, detailed observations were made for certain specific locations in the model, as shown in Fig. 6. In this figure, P1 and P2 represent the longitudinal tangents to the centers of the two melt tracks in the XZ plane, while L1 is the transverse profile in the YZ plane. The scanning direction is positive and negative along the X axis. Points A and B are the locations of the centers of the molten pool of the first and second melt tracks, respectively (x = 1.995 × 10⁻⁴, y = 5 × 10⁻⁷, and z = −4.85 × 10⁻⁵).

FIG. 6.

Schematic diagram of observation position.

A. Single-track simulation

A series of single-track molten pool simulation experiments were carried out in order to investigate the influence law of laser power as well as scanning speed on the HP-LPBF process. Figure 7 demonstrates the evolution of the 3D morphology and temperature field of the single-track molten pool in the time period of 50–500 μs under a laser power of 100 W and a scanning speed of 800 mm/s. The powder bed is in the natural cooling state. When t = 50 μs, the powder is heated by the laser heat and rapidly melts and settles to form the initial molten pool. This process is accompanied by partial melting of the substrate and solidification together with the melted powder. The molten pool rapidly expands with increasing width, depth, length, and temperature, as shown in Fig. 7(a). When t = 150 μs, the molten pool expands more obviously, and the temperature starts to transfer to the surrounding area, forming a heat-affected zone. At this point, the width of the molten pool tends to stabilize, and the temperature in the center of the molten pool has reached its peak and remains largely stable. However, the phenomenon of molten pool spatter was also observed in this process, as shown in Fig. 7(b). As time advances, when t = 300 μs, solidification begins to occur at the tail of the molten pool, and tiny ripples are produced on the solidified surface. This is due to the fact that the melt flows toward the region with large temperature gradient under the influence of Marangoni convection and solidifies together with the melt at the end of the bath. At this point, the temperature gradient at the front of the bath is significantly larger than at the end. While the width of the molten pool was gradually reduced, the shape of the molten pool was gradually changed to a “comet” shape. In addition, a slight depression was observed at the top of the bath because the peak temperature at the surface of the bath reached the evaporation temperature, which resulted in a recoil pressure perpendicular to the surface of the bath downward, creating a depressed region. As the laser focal spot moves and is paired with the Marangoni convection of the melt, these recessed areas will be filled in as shown in Fig. 7(c). It has been shown that the depressed regions are the result of the coupled effect of Marangoni convection, recoil pressure, and surface tension.⁷⁹ By t = 500 μs, the width and height of the molten pool stabilize and show a “comet” shape in Fig. 7(d).

FIG. 7.

Single-track molten pool process: (a) t = 50 ��⁠, (b) t = 150 ��⁠, (c) t = 300 ��⁠, (d) t = 500 ��⁠.

Figure 8 depicts the velocity vector diagram of the P1 profile in a single-track molten pool, the length of the arrows represents the magnitude of the velocity, and the maximum velocity is about 2.36 m/s. When t = 50 μs, the molten pool takes shape, and the velocities at the two ends of the pool are the largest. The variation of the velocities at the front end is especially more significant in Fig. 8(a). As the time advances to t = 150 μs, the molten pool expands rapidly, in which the velocity at the tail increases and changes more significantly, while the velocity at the front is relatively small. At this stage, the melt moves backward from the center of the molten pool, which in turn expands the molten pool area. The melt at the back end of the molten pool center flows backward along the edge of the molten pool surface and then converges along the edge of the molten pool to the bottom center, rising to form a closed loop. Similarly, a similar closed loop is formed at the front end of the center of the bath, but with a shorter path. However, a large portion of the melt in the center of the closed loop formed at the front end of the bath is in a nearly stationary state. The main cause of this melt flow phenomenon is the effect of temperature gradient and surface tension (the Marangoni effect), as shown in Figs. 8(b) and 8(e). This dynamic behavior of the melt tends to form an “elliptical” pool. At t = 300 μs, the tendency of the above two melt flows to close the loop is more prominent and faster in Fig. 8(c). When t = 500 μs, the velocity vector of the molten pool shows a stable trend, and the closed loop of melt flow also remains stable. With the gradual laser focal spot movement, the melt is gradually solidified at its tail, and finally, a continuous and stable single track is formed in Fig. 8(d).

FIG. 8.

Vector plot of single-track molten pool velocity in XZ longitudinal section: (a) t = 50 ��⁠, (b) t = 150 ��⁠, (c) t = 300 ��⁠, (d) t = 500 ��⁠, (e) molten pool flow.

In order to explore in depth the transient evolution of the molten pool, the evolution of the single-track temperature field and the melt flow was monitored in the YZ cross section. Figure 9(a) shows the state of the powder bed at the initial moment. When t = 250 μs, the laser focal spot acts on the powder bed and the powder starts to melt and gradually collects in the molten pool. At this time, the substrate will also start to melt, and the melt flow mainly moves in the downward and outward directions and the velocity is maximum at the edges in Fig. 9(b). When t = 300 μs, the width and depth of the molten pool increase due to the recoil pressure. At this time, the melt flows more slowly at the center, but the direction of motion is still downward in Fig. 9(c). When t = 350 μs, the width and depth of the molten pool further increase, at which time the intensity of the melt flow reaches its peak and the direction of motion remains the same in Fig. 9(d). When t = 400 μs, the melt starts to move upward, and the surrounding powder or molten material gradually fills up, causing the surface of the molten pool to begin to flatten. At this time, the maximum velocity of the melt is at the center of the bath, while the velocity at the edge is close to zero, and the edge of the melt starts to solidify in Fig. 9(e). When t = 450 μs, the melt continues to move upward, forming a convex surface of the melt track. However, the melt movement slows down, as shown in Fig. 9(f). When t = 500 μs, the melt further moves upward and its speed gradually becomes smaller. At the same time, the melt solidifies further, as shown in Fig. 9(g). When t = 550 μs, the melt track is basically formed into a single track with a similar “mountain” shape. At this stage, the velocity is close to zero only at the center of the molten pool, and the flow behavior of the melt is poor in Fig. 9(h). At t = 600 μs, the melt stops moving and solidification is rapidly completed. Up to this point, a single track is formed in Fig. 9(i). During the laser action on the powder bed, the substrate melts and combines with the molten state powder. The powder-to-powder fusion is like the convergence of water droplets, which are rapidly fused by surface tension. However, the fusion between the molten state powder and the substrate occurs driven by surface tension, and the molten powder around the molten pool is pulled toward the substrate (a wetting effect occurs), which ultimately results in the formation of a monolithic whole.^38,80,81

FIG. 9.

Evolution of single-track molten pool temperature and melt flow in the YZ cross section: (a) t = 0 ��⁠, (b) t = 250 ��⁠, (c) t = 300 ��⁠, (d) t = 350 ��⁠, (e) t = 400 ��⁠, (f) t = 450 ��⁠, (g) t = 500 ��⁠, (h) t = 550 ��⁠, (i) t = 600 ��⁠.

The wetting ability between the liquid metal and the solid substrate in the molten pool directly affects the degree of balling of the melt,^82,83 and the wetting ability can be measured by the contact angle of a single track in Fig. 10. A smaller value of contact angle represents better wettability. The contact angle α can be calculated by�=�1−�22,

(17)

where �1 and �2 are the contact angles of the left and right regions, respectively.

FIG. 10.

Schematic of contact angle.

Relevant studies have confirmed that the wettability is better at a contact angle α around or below 40°.⁸⁴ After measurement, a single-track contact angle α of about 33° was obtained under this process parameter, which further confirms the good wettability.

B. Double-track simulation

In order to deeply investigate the influence of hatch spacing on the characteristics of the HP-LPBF process, a series of double-track molten pool simulation experiments were systematically carried out. Figure 11 shows in detail the dynamic changes of the 3D morphology and temperature field of the double-track molten pool in the time period of 2050–2500 μs under the conditions of laser power of 100 W, scanning speed of 800 mm/s, and hatch spacing of 0.06 mm. By comparing the study with Fig. 7, it is observed that the basic characteristics of the 3D morphology and temperature field of the second track are similar to those of the first track. However, there are subtle differences between them. The first track exhibits a basically symmetric shape, but the second track morphology shows a slight deviation influenced by the difference in thermal diffusion rate between the solidified metal and the powder. Otherwise, the other characteristic information is almost the same as that of the first track. Figure 12 shows the velocity vector plot of the P2 profile in the double-track molten pool, with a maximum velocity of about 2.63 m/s. The melt dynamics at both ends of the pool are more stable at t = 2050 μs, where the maximum rate of the second track is only 1/3 of that of the first one. Other than that, the rest of the information is almost no significant difference from the characteristic information of the first track. Figure 13 demonstrates a detailed observation of the double-track temperature field and melts flow in the YZ cross section, and a comparative study with Fig. 9 reveals that the width of the second track is slightly wider. In addition, after the melt direction shifts from bottom to top, the first track undergoes four time periods (50 μs) to reach full solidification, while the second track takes five time periods. This is due to the presence of significant heat buildup in the powder bed after the forming of the first track, resulting in a longer dynamic time of the melt and an increased molten pool lifetime. In conclusion, the level of specimen forming can be significantly optimized by adjusting the laser power and hatch spacing.

FIG. 11.

Double-track molten pool process: (a) t = 2050 ��⁠, (b) t = 2150 ��⁠, (c) t = 2300 ��⁠, (d) t = 2500 ��⁠.

FIG. 12.

Vector plot of double-track molten pool velocity in XZ longitudinal section: (a) t = 2050 ��⁠, (b) t = 2150 ��⁠, (c) t = 2300 ��⁠, (d) t = 2500 ��⁠.

FIG. 13.

Evolution of double-track molten pool temperature and melt flow in the YZ cross section: (a) t = 2250 ��⁠, (b) t = 2300 ��⁠, (c) t = 2350 ��⁠, (d) t = 2400 ��⁠, (e) t = 2450 ��⁠, (f) t = 2500 ��⁠, (g) t = 2550 ��⁠, (h) t = 2600 ��⁠, (i) t = 2650 ��⁠.

In order to quantitatively detect the molten pool dimensions as well as the remolten region dimensions, the molten pool characterization information in Fig. 14 is constructed by drawing the boundary on the YZ cross section based on the isothermal surface of the liquid phase line. It can be observed that the heights of the first track and second track are basically the same, but the depth of the second track increases relative to the first track. The molten pool width is mainly positively correlated with the laser power as well as the scanning speed (the laser line energy density �⁠). However, the remelted zone width is negatively correlated with the hatch spacing (the overlapping ratio). Overall, the forming quality of the specimens can be directly influenced by adjusting the laser power, scanning speed, and hatch spacing.

FIG. 14.

Double-track molten pool characterization information on YZ cross section.

In order to study the variation rule of the temperature in the center of the molten pool with time, Fig. 15 demonstrates the temperature variation curves with time for two reference points, A and B. Among them, the red dotted line indicates the liquid phase line temperature of SS316L. From the figure, it can be seen that the maximum temperature at the center of the molten pool in the first track is lower than that in the second track, which is mainly due to the heat accumulation generated after passing through the first track. The maximum temperature gradient was calculated to be 1.69 × 10⁸ K/s. When the laser scanned the first track, the temperature in the center of the molten pool of the second track increased slightly. Similarly, when the laser scanned the second track, a similar situation existed in the first track. Since the temperature gradient in the second track is larger than that in the first track, the residence time of the liquid phase in the molten pool of the first track is longer than that of the second track.

FIG. 15.

Temperature profiles as a function of time for two reference points A and B.

C. Simulation analysis of molten pool under different process parameters

In order to deeply investigate the effects of various process parameters on the mesoscopic-scale temperature field, molten pool characteristic information and defects of HP-LPBF, numerical simulation experiments on mesoscopic-scale laser power, scanning speed, and hatch spacing of double-track molten pools were carried out.

1. Laser power

Figure 16 shows the effects of different laser power on the morphology and temperature field of the double-track molten pool at a scanning speed of 800 mm/s and a hatch spacing of 0.06 mm. When P = 50 W, a smaller molten pool is formed due to the lower heat generated by the Gaussian light source per unit time. This leads to a smaller track width, which results in adjacent track not lapping properly and the presence of a large number of unmelted powder particles, resulting in an increase in the number of defects, such as pores in the specimen. The surface of the track is relatively flat, and the depth is small. In addition, the temperature gradient before and after the molten pool was large, and the depression location appeared at the biased front end in Fig. 16(a). When P = 100 W, the surface of the track is flat and smooth with excellent lap. Due to the Marangoni effect, the velocity field of the molten pool is in the form of “vortex,” and the melt has good fluidity, and the maximum velocity reaches 2.15 m/s in Fig. 16(b). When P = 200 W, the heat generated by the Gaussian light source per unit time is too large, resulting in the melt rapidly reaching the evaporation temperature, generating a huge recoil pressure, forming a large molten pool, and the surface of the track is obviously raised. The melt movement is intense, especially the closed loop at the center end of the molten pool. At this time, the depth and width of the molten pool are large, leading to the expansion of the remolten region and the increased chance of the appearance of porosity defects in Fig. 16(c). The results show that at low laser power, the surface tension in the molten pool is dominant. At high laser power, recoil pressure is its main role.

FIG. 16.

Simulation results of double-track molten pool under different laser powers: (a) P = 50 W, (b) P = 100 W, (c) P = 200 W.

Table II shows the effect of different laser powers on the characteristic information of the double-track molten pool at a scanning speed of 800 mm/s and a hatch spacing of 0.06 mm. The negative overlapping ratio in the table indicates that the melt tracks are not lapped, and 26/29 indicates the melt depth of the first track/second track. It can be seen that with the increase in laser power, the melt depth, melt width, melt height, and remelted zone show a gradual increase. At the same time, the overlapping ratio also increases. Especially in the process of laser power from 50 to 200 W, the melting depth and melting width increased the most, which increased nearly 2 and 1.5 times, respectively. Meanwhile, the overlapping ratio also increases with the increase in laser power, which indicates that the melting and fusion of materials are better at high laser power. On the other hand, the dimensions of the molten pool did not change uniformly with the change of laser power. Specifically, the depth-to-width ratio of the molten pool increased from about 0.30 to 0.39 during the increase from 50 to 120 W, which further indicates that the effective heat transfer in the vertical direction is greater than that in the horizontal direction with the increase in laser power. This dimensional response to laser power is mainly affected by the recoil pressure and also by the difference in the densification degree between the powder layer and the metal substrate. In addition, according to the experimental results, the contact angle shows a tendency to increase and then decrease during the process of laser power increase, and always stays within the range of less than 33°. Therefore, in practical applications, it is necessary to select the appropriate laser power according to the specific needs in order to achieve the best processing results.

TABLE II.

Double-track molten pool characterization information at different laser powers.

Laser power (W)	Depth (μm)	Width (μm)	Height (μm)	Remolten region (μm)	Overlapping ratio (%)	Contact angle (°)
50	16	54	11	/	−10	23
100	26/29	74	14	18	23.33	33
200	37/45	116	21	52	93.33	28

2. Scanning speed

Figure 17 demonstrates the effect of different scanning speeds on the morphology and temperature field of the double-track molten pool at a laser power of 100 W and a hatch spacing of 0.06 mm. With the gradual increase in scanning speed, the surface morphology of the molten pool evolves from circular to elliptical. When � = 200 mm/s, the slow scanning speed causes the material to absorb too much heat, which is very easy to trigger the overburning phenomenon. At this point, the molten pool is larger and the surface morphology is uneven. This situation is consistent with the previously discussed scenario with high laser power in Fig. 17(a). However, when � = 1600 mm/s, the scanning speed is too fast, resulting in the material not being able to absorb sufficient heat, which triggers the powder particles that fail to melt completely to have a direct effect on the bonding of the melt to the substrate. At this time, the molten pool volume is relatively small and the neighboring melt track cannot lap properly. This result is consistent with the previously discussed case of low laser power in Fig. 17(b). Overall, the ratio of the laser power to the scanning speed (the line energy density �⁠) has a direct effect on the temperature field and surface morphology of the molten pool.

FIG. 17.

Simulation results of double-track molten pool under different scanning speed: (a) � = 200 mm/s, (b) � = 1600 mm/s.

Table III shows the effects of different scanning speed on the characteristic information of the double-track molten pool under the condition of laser power of 100 W and hatch spacing of 0.06 mm. It can be seen that the scanning speed has a significant effect on the melt depth, melt width, melt height, remolten region, and overlapping ratio. With the increase in scanning speed, the melt depth, melt width, melt height, remelted zone, and overlapping ratio show a gradual decreasing trend. Among them, the melt depth and melt width decreased faster, while the melt height and remolten region decreased relatively slowly. In addition, when the scanning speed was increased from 200 to 800 mm/s, the decreasing speeds of melt depth and melt width were significantly accelerated, while the decreasing speeds of overlapping ratio were relatively slow. When the scanning speed was further increased to 1600 mm/s, the decreasing speeds of melt depth and melt width were further accelerated, and the un-lapped condition of the melt channel also appeared. In addition, the contact angle increases and then decreases with the scanning speed, and both are lower than 33°. Therefore, when selecting the scanning speed, it is necessary to make reasonable trade-offs according to the specific situation, and take into account the factors of melt depth, melt width, melt height, remolten region, and overlapping ratio, in order to achieve the best processing results.

TABLE III.

Double-track molten pool characterization information at different scanning speeds.

Scanning speed (mm/s)	Depth (μm)	Width (μm)	Height (μm)	Remolten region (μm)	Overlapping ratio (%)	Contact angle (°)
200	55/68	182	19/32	124	203.33	22
1600	13	50	11	/	−16.67	31

3. Hatch spacing

Figure 18 shows the effect of different hatch spacing on the morphology and temperature field of the double-track molten pool under the condition of laser power of 100 W and scanning speed of 800 mm/s. The surface morphology and temperature field of the first track and second track are basically the same, but slightly different. The first track shows a basically symmetric morphology along the scanning direction, while the second track shows a slight offset due to the difference in the heat transfer rate between the solidified material and the powder particles. When the hatch spacing is too small, the overlapping ratio increases and the probability of defects caused by remelting phenomenon grows. When the hatch spacing is too large, the neighboring melt track cannot overlap properly, and the powder particles are not completely melted, leading to an increase in the number of holes. In conclusion, the ratio of the line energy density � to the hatch spacing (the volume energy density E) has a significant effect on the temperature field and surface morphology of the molten pool.

FIG. 18.

Simulation results of double-track molten pool under different hatch spacings: (a) H = 0.03 mm, (b) H = 0.12 mm.

Table IV shows the effects of different hatch spacing on the characteristic information of the double-track molten pool under the condition of laser power of 100 W and scanning speed of 800 mm/s. It can be seen that the hatch spacing has little effect on the melt depth, melt width, and melt height, but has some effect on the remolten region. With the gradual expansion of hatch spacing, the remolten region shows a gradual decrease. At the same time, the overlapping ratio also decreased with the increase in hatch spacing. In addition, it is observed that the contact angle shows a tendency to increase and then remain stable when the hatch spacing increases, which has a more limited effect on it. Therefore, trade-offs and decisions need to be made on a case-by-case basis when selecting the hatch spacing.

TABLE IV.

Double-track molten pool characterization information at different hatch spacings.

Hatch spacing (mm)	Depth (μm)	Width (μm)	Height (μm)	Remolten region (μm)	Overlapping ratio (%)	Contact angle (°)
0.03	25/27	82	14	59	173.33	30
0.12	26	78	14	/	−35	33

In summary, the laser power, scanning speed, and hatch spacing have a significant effect on the formation of the molten pool, and the correct selection of these three process parameters is crucial to ensure the forming quality. In addition, the melt depth of the second track is slightly larger than that of the first track at higher line energy density � and volume energy density E. This is mainly due to the fact that a large amount of heat accumulation is generated after the first track, forming a larger molten pool volume, which leads to an increase in the melt depth.

D. Simulation analysis of molten pool with powder particle size and laser focal spot diameter

Figure 19 demonstrates the effect of different powder particle sizes and laser focal spot diameters on the morphology and temperature field of the double-track molten pool under a laser power of 100 W, a scanning speed of 800 mm/s, and a hatch spacing of 0.06 mm. In the process of melting coarse powder with small laser focal spot diameter, the laser energy cannot completely melt the larger powder particles, resulting in their partial melting and further generating excessive pore defects. The larger powder particles tend to generate zigzag molten pool edges, which cause an increase in the roughness of the melt track surface. In addition, the molten pool is also prone to generate the present spatter phenomenon, which can directly affect the quality of forming. The volume of the formed molten pool is relatively small, while the melt depth, melt width, and melt height are all smaller relative to the fine powder in Fig. 19(a). In the process of melting fine powders with a large laser focal spot diameter, the laser energy is able to melt the fine powder particles sufficiently, even to the point of overmelting. This results in a large number of fine spatters being generated at the edge of the molten pool, which causes porosity defects in the melt track in Fig. 19(b). In addition, the maximum velocity of the molten pool is larger for large powder particle sizes compared to small powder particle sizes, which indicates that the temperature gradient in the molten pool is larger for large powder particle sizes and the melt motion is more intense. However, the size of the laser focal spot diameter has a relatively small effect on the melt motion. However, a larger focal spot diameter induces a larger melt volume with greater depth, width, and height. In conclusion, a small powder size helps to reduce the surface roughness of the specimen, and a small laser spot diameter reduces the minimum forming size of a single track.

FIG. 19.

Arabian Journal for Science and Engineering Aims and scope Submit manuscript

Simulation results of double-track molten pool with different powder particle size and laser focal spot diameter: (a) focal spot = 25 μm, coarse powder, (b) focal spot = 80 μm, fine powder.

Table V shows the maximum temperature gradient at the reference point for different powder sizes and laser focal spot diameters. As can be seen from the table, the maximum temperature gradient is lower than that of HP-LPBF for both coarse powders with a small laser spot diameter and fine powders with a large spot diameter, a phenomenon that leads to an increase in the heat transfer rate of HP-LPBF, which in turn leads to a corresponding increase in the cooling rate and, ultimately, to the formation of finer microstructures.

TABLE V.

Maximum temperature gradient at the reference point for different powder particle sizes and laser focal spot diameters.

Laser power (W)	Scanning speed (mm/s)	Hatch spacing (mm)	Average powder size (μm)	Laser focal spot diameter (μm)	Maximum temperature gradient (×10⁷ K/s)
100	800	0.06	31.7	25	7.89
100	800	0.06	11.5	80	7.11

IV. CONCLUSIONS

In this study, the geometrical characteristics of 3D coarse and fine powder particles were first calculated using DEM and then numerical simulations of single track and double track in the process of forming SS316L from monolayer HP-LPBF at mesoscopic scale were developed using CFD method. The effects of Marangoni convection, surface tension, recoil pressure, gravity, thermal convection, thermal radiation, and evaporative heat dissipation on the heat and mass transfer in the molten pool were considered in this model. The effects of laser power, scanning speed, and hatch spacing on the dynamics of the single-track and double-track molten pools, as well as on other characteristic information, were investigated. The effects of the powder particle size on the molten pool were investigated comparatively with the laser focal spot diameter. The main conclusions are as follows:

The results show that the temperature gradient at the front of the molten pool is significantly larger than that at the tail, and the molten pool exhibits a “comet” morphology. At the top of the molten pool, there is a slightly concave region, which is the result of the coupling of Marangoni convection, recoil pressure, and surface tension. The melt flow forms two closed loops, which are mainly influenced by temperature gradients and surface tension. This special dynamic behavior of the melt tends to form an “elliptical” molten pool and an almost “mountain” shape in single-track forming.
The basic characteristics of the three-dimensional morphology and temperature field of the second track are similar to those of the first track, but there are subtle differences. The first track exhibits a basically symmetrical shape; however, due to the difference in thermal diffusion rates between the solidified metal and the powder, a slight asymmetry in the molten pool morphology of the second track occurs. After forming through the first track, there is a significant heat buildup in the powder bed, resulting in a longer dynamic time of the melt, which increases the life of the molten pool. The heights of the first track and second track remained essentially the same, but the depth of the second track was greater relative to the first track. In addition, the maximum temperature gradient was 1.69 × 10⁸ K/s during HP-LPBF forming.
At low laser power, the surface tension in the molten pool plays a dominant role. At high laser power, recoil pressure becomes the main influencing factor. With the increase of laser power, the effective heat transfer in the vertical direction is superior to that in the horizontal direction. With the gradual increase of scanning speed, the surface morphology of the molten pool evolves from circular to elliptical. In addition, the scanning speed has a significant effect on the melt depth, melt width, melt height, remolten region, and overlapping ratio. Too large or too small hatch spacing will lead to remelting or non-lap phenomenon, which in turn causes the formation of defects.
When using a small laser focal spot diameter, it is difficult to completely melt large powder particle sizes, resulting in partial melting and excessive porosity generation. At the same time, large powder particles produce curved edges of the molten pool, resulting in increased surface roughness of the melt track. In addition, spatter occurs, which directly affects the forming quality. At small focal spot diameters, the molten pool volume is relatively small, and the melt depth, the melt width, and the melt height are correspondingly small. Taken together, the small powder particle size helps to reduce surface roughness, while the small spot diameter reduces the forming size.

REFERENCES

S. L. Sing and W. Y. Yeong , “ Laser powder bed fusion for metal additive manufacturing: Perspectives on recent developments,” Virtual Phys. Prototyping. 15, 359–370 (2020).https://doi.org/10.1080/17452759.2020.1779999
Google ScholarCrossref
A. M. Khorasani , I. G. Jithin , J. K. Veetil , and A. H. Ghasemi , “ A review of technological improvements in laser-based powder bed fusion of metal printers,” Int. J. Adv. Manuf. Technol. 108, 191–209 (2020).https://doi.org/10.1007/s00170-020-05361-3
Google ScholarCrossref
Y. Qin , A. Brockett , Y. Ma , A. Razali , J. Zhao , C. Harrison , W. Pan , X. Dai , and D. Loziak , “ Micro-manufacturing: Research, technology outcomes and development issues,” Int. J. Adv. Manuf. Technol. 47, 821–837 (2010).https://doi.org/10.1007/s00170-009-2411-2
Google ScholarCrossref
B. Nagarajan , Z. Hu , X. Song , W. Zhai , and J. Wei , “ Development of micro selective laser melting: The state of the art and future perspectives,” Engineering. 5, 702–720 (2019).https://doi.org/10.1016/j.eng.2019.07.002
Google ScholarCrossref
Y. Wei , G. Chen , W. Li , Y. Zhou , Z. Nie , J. Xu , and W. Zhou , “ Micro selective laser melting of SS316L: Single tracks, defects, microstructures and thermal/mechanical properties,” Opt. Laser Technol. 145, 107469 (2022).https://doi.org/10.1016/j.optlastec.2021.107469
Google ScholarCrossref
Y. Wei , G. Chen , W. Li , M. Li , Y. Zhou , Z. Nie , and J. Xu , “ Process optimization of micro selective laser melting and comparison of different laser diameter for forming different powder,” Opt. Laser Technol. 150, 107953 (2022).https://doi.org/10.1016/j.optlastec.2022.107953
Google ScholarCrossref
H. Zhiheng , B. Nagarajan , X. Song , R. Huang , W. Zhai , and J. Wei , “ Formation of SS316L single tracks in micro selective laser melting: Surface, geometry, and defects,” Adv. Mater. Sci. Eng. 2019, 9451406.https://doi.org/10.1155/2019/9451406
Crossref
B. Nagarajan , Z. Hu , S. Gao , X. Song , R. Huang , M. Seita , and J. Wei , “ Effect of in-situ laser remelting on the microstructure of SS316L fabricated by micro selective laser melting,” in Advanced Surface Enhancement, edited by Sho Itoh and Shashwat Shukla , Lecture Notes in Mechanical Engineering ( Springer Singapore, Singapore, 2020), pp. 330–336.
Google ScholarCrossref
H. Zhiheng , B. Nagarajan , X. Song , R. Huang , W. Zhai , and J. Wei , “ Tailoring surface roughness of micro selective laser melted SS316L by in-situ laser remelting,” in Advanced Surface Enhancement, edited by Sho Itoh and Shashwat Shukla , Lecture Notes in Mechanical Engineering ( Springer Singapore, Singapore, 2020), pp. 337–343.
Google Scholar
J. Fu , Z. Hu , X. Song , W. Zhai , Y. Long , H. Li , and M. Fu , “ Micro selective laser melting of NiTi shape memory alloy: Defects, microstructures and thermal/mechanical properties,” Opt. Laser Technol. 131, 106374 (2020).https://doi.org/10.1016/j.optlastec.2020.106374
Google ScholarCrossref
E. Abele and M. Kniepkamp , “ Analysis and optimisation of vertical surface roughness in micro selective laser melting,” Surf. Topogr.: Metrol. Prop. 3, 034007 (2015).https://doi.org/10.1088/2051-672X/3/3/034007
Google ScholarCrossref
S. Qu , J. Ding , J. Fu , M. Fu , B. Zhang , and X. Song , “ High-precision laser powder bed fusion processing of pure copper,” Addit. Manuf. 48, 102417 (2021).https://doi.org/10.1016/j.addma.2021.102417
Google ScholarCrossref
Y. Wei , G. Chen , M. Li , W. Li , Y. Zhou , J. Xu , and Z. wei , “ High-precision laser powder bed fusion of 18Ni300 maraging steel and its SiC reinforcement composite materials,” J. Manuf. Process. 84, 750–763 (2022).https://doi.org/10.1016/j.jmapro.2022.10.049
Google ScholarCrossref
B. Liu , R. Wildman , T. Christopher , I. Ashcroft , and H. Richard , “ Investigation the effect of particle size distribution on processing parameters optimisation in selective laser melting process,” in 2011 International Solid Freeform Fabrication Symposium ( University of Texas at Austin, 2011).
Google Scholar
T. D. McLouth , G. E. Bean , D. B. Witkin , S. D. Sitzman , P. M. Adams , D. N. Patel , W. Park , J.-M. Yang , and R. J. Zaldivar , “ The effect of laser focus shift on microstructural variation of Inconel 718 produced by selective laser melting,” Mater. Des. 149, 205–213 (2018).https://doi.org/10.1016/j.matdes.2018.04.019
Google ScholarCrossref
Y. Qian , Y. Wentao , and L. Feng , “ Mesoscopic simulations of powder bed fusion: Research progresses and conditions,” Electromachining Mould 06, 46–52 (2017).https://doi.org/10.3969/j.issn.1009-279X.2017.06.012
Google Scholar
J. Fu , S. Qu , J. Ding , X. Song , and M. W. Fu , “ Comparison of the microstructure, mechanical properties and distortion of stainless Steel 316L fabricated by micro and conventional laser powder bed fusion,” Addit. Manuf. 44, 102067 (2021).https://doi.org/10.1016/j.addma.2021.102067
Google ScholarCrossref
N. T. Aboulkhair , I. Maskery , C. Tuck , I. Ashcroft , and N. M. Everitt , “ The microstructure and mechanical properties of selectively laser Melted AlSi10Mg: The effect of a conventional T6-like heat treatment,” Mater. Sci. Eng. A 667, 139–146 (2016).https://doi.org/10.1016/j.msea.2016.04.092
Google ScholarCrossref
S. Y. Chen , J. C. Huang , C. T. Pan , C. H. Lin , T. L. Yang , Y. S. Huang , C. H. Ou , L. Y. Chen , D. Y. Lin , H. K. Lin , T. H. Li , J. S. C. Jang , and C. C. Yang , “ Microstructure and mechanical properties of open-cell porous Ti-6Al-4V fabricated by selective laser melting,” J. Alloys Compd. 713, 248–254 (2017).https://doi.org/10.1016/j.jallcom.2017.04.190
Google ScholarCrossref
Y. Bai , Y. Yang , D. Wang , and M. Zhang , “ Influence mechanism of parameters process and mechanical properties evolution mechanism of Maraging steel 300 by selective laser melting,” Mater. Sci. Eng. A 703, 116–123 (2017).https://doi.org/10.1016/j.msea.2017.06.033
Google ScholarCrossref
Y. Bai , Y. Yang , Z. Xiao , M. Zhang , and D. Wang , “ Process optimization and mechanical property evolution of AlSiMg0.75 by selective laser melting,” Mater. Des. 140, 257–266 (2018).https://doi.org/10.1016/j.matdes.2017.11.045
Google ScholarCrossref
Y. Liu , M. Zhang , W. Shi , Y. Ma , and J. Yang , “ Study on performance optimization of 316L stainless steel parts by high-efficiency selective laser melting,” Opt. Laser Technol. 138, 106872 (2021).https://doi.org/10.1016/j.optlastec.2020.106872
Google ScholarCrossref
D. Gu , Y.-C. Hagedorn , W. Meiners , G. Meng , R. J. S. Batista , K. Wissenbach , and R. Poprawe , “ Densification behavior, microstructure evolution, and wear performance of selective laser melting processed commercially pure titanium,” Acta Mater. 60, 3849–3860 (2012).https://doi.org/10.1016/j.actamat.2012.04.006
Google ScholarCrossref
N. Read , W. Wang , K. Essa , and M. M. Attallah , “ Selective laser melting of AlSi10Mg alloy: Process optimisation and mechanical properties development,” Mater. Des. 65, 417–424 (2015).https://doi.org/10.1016/j.matdes.2014.09.044
Google ScholarCrossref
I. A. Roberts , C. J. Wang , R. Esterlein , M. Stanford , and D. J. Mynors , “ A three-dimensional finite element analysis of the temperature field during laser melting of metal powders in additive layer manufacturing,” Int. J. Mach. Tools Manuf. 49(12–13), 916–923 (2009).https://doi.org/10.1016/j.ijmachtools.2009.07.004
Google ScholarCrossref
K. Dai and L. Shaw , “ Finite element analysis of the effect of volume shrinkage during laser densification,” Acta Mater. 53(18), 4743–4754 (2005).https://doi.org/10.1016/j.actamat.2005.06.014
Google ScholarCrossref
K. Carolin , E. Attar , and P. Heinl , “ Mesoscopic simulation of selective beam melting processes,” J. Mater. Process. Technol. 211(6), 978–987 (2011).https://doi.org/10.1016/j.jmatprotec.2010.12.016
Google ScholarCrossref
F.-J. Gürtler , M. Karg , K.-H. Leitz , and M. Schmidt , “ Simulation of laser beam melting of steel powders using the three-dimensional volume of fluid method,” Phys. Procedia 41, 881–886 (2013).https://doi.org/10.1016/j.phpro.2013.03.162
Google ScholarCrossref
P. Meakin and R. Jullien , “ Restructuring effects in the rain model for random deposition,” J. Phys. France 48(10), 1651–1662 (1987).https://doi.org/10.1051/jphys:0198700480100165100
Google ScholarCrossref
J-m Wang , G-h Liu , Y-l Fang , and W-k Li , “ Marangoni effect in nonequilibrium multiphase system of material processing,” Rev. Chem. Eng. 32(5), 551–585 (2016).https://doi.org/10.1515/revce-2015-0067
Google ScholarCrossref
W. Ye , S. Zhang , L. L. Mendez , M. Farias , J. Li , B. Xu , P. Li , and Y. Zhang , “ Numerical simulation of the melting and alloying processes of elemental titanium and boron powders using selective laser alloying,” J. Manuf. Process. 64, 1235–1247 (2021).https://doi.org/10.1016/j.jmapro.2021.02.044
Google ScholarCrossref
U. S. Bertoli , A. J. Wolfer , M. J. Matthews , J.-P. R. Delplanque , and J. M. Schoenung , “ On the limitations of volumetric energy density as a design parameter for selective laser melting,” Mater. Des. 113, 331–340 (2017).https://doi.org/10.1016/j.matdes.2016.10.037
Google ScholarCrossref
W. E. King , H. D. Barth , V. M. Castillo , G. F. Gallegos , J. W. Gibbs , D. E. Hahn , C. Kamath , and A. M. Rubenchik , “ Observation of keyhole-mode laser melting in laser powder-bed fusion additive manufacturing,” J. Mater. Process. Technol. 214(12), 2915–2925 (2014).https://doi.org/10.1016/j.jmatprotec.2014.06.005
Google ScholarCrossref
L. Cao , “ Numerical simulation of the impact of laying powder on selective laser melting single-pass formation,” Int. J. Heat Mass Transfer 141, 1036–1048 (2019).https://doi.org/10.1016/j.ijheatmasstransfer.2019.07.053
Google ScholarCrossref
L. Huang , X. Hua , D. Wu , and F. Li , “ Numerical study of keyhole instability and porosity formation mechanism in laser welding of aluminum alloy and steel,” J. Mater. Process. Technol. 252, 421–431 (2018).https://doi.org/10.1016/j.jmatprotec.2017.10.011
Google ScholarCrossref
K. Q. Le , C. Tang , and C. H. Wong , “ On the study of keyhole-mode melting in selective laser melting process,” Int. J. Therm. Sci. 145, 105992 (2019).https://doi.org/10.1016/j.ijthermalsci.2019.105992
Google ScholarCrossref
J.-H. Cho and S.-J. Na , “ Theoretical analysis of keyhole dynamics in polarized laser drilling,” J. Phys. D: Appl. Phys. 40(24), 7638 (2007).https://doi.org/10.1088/0022-3727/40/24/007
Google ScholarCrossref
W. Ye , “ Mechanism analysis of selective laser melting and metallurgy process based on base element powder of titanium and boron,” Ph.D. dissertation ( Nanchang University, 2021).
Google Scholar
R. Ammer , M. Markl , U. Ljungblad , C. Körner , and U. Rüde , “ Simulating fast electron beam melting with a parallel thermal free surface lattice Boltzmann method,” Comput. Math. Appl. 67(2), 318–330 (2014).https://doi.org/10.1016/j.camwa.2013.10.001
Google ScholarCrossref
H. Chen , Q. Wei , S. Wen , Z. Li , and Y. Shi , “ Flow behavior of powder particles in layering process of selective laser melting: Numerical modeling and experimental verification based on discrete element method,” Int. J. Mach. Tools Manuf. 123, 146–159 (2017).https://doi.org/10.1016/j.ijmachtools.2017.08.004
Google ScholarCrossref
F. Verhaeghe , T. Craeghs , J. Heulens , and L. Pandelaers , “ A pragmatic model for selective laser melting with evaporation,” Acta Mater. 57(20), 6006–6012 (2009).https://doi.org/10.1016/j.actamat.2009.08.027
Google ScholarCrossref
C. H. Fu and Y. B. Guo , “ Three-dimensional temperature gradient mechanism in selective laser melting of Ti-6Al-4V,” J. Manuf. Sci. Eng. 136(6), 061004 (2014).https://doi.org/10.1115/1.4028539
Google ScholarCrossref
Y. Xiang , Z. Shuzhe , L. Junfeng , W. Zhengying , Y. Lixiang , and J. Lihao , “ Numerical simulation and experimental verification for selective laser single track melting forming of Ti6Al4V,” J. Zhejiang Univ. (Eng. Sci.) 53(11), 2102–2109 + 2117 (2019).https://doi.org/10.3785/j.issn.1008-973X.2019.11.007
Google Scholar
Q. He , H. Xia , J. Liu , X. Ao , and S. Lin , “ Modeling and numerical studies of selective laser melting: Multiphase flow, solidification and heat transfer,” Mater. Des. 196, 109115 (2020).https://doi.org/10.1016/j.matdes.2020.109115
Google ScholarCrossref
L. Cao , “ Mesoscopic-scale numerical simulation including the influence of process parameters on SLM single-layer multi-pass formation,” Metall. Mater. Trans. A 51, 4130–4145 (2020).https://doi.org/10.1007/s11661-020-05831-z
Google ScholarCrossref
L. Cao , “ Mesoscopic-scale numerical investigation including the influence of process parameters on LPBF multi-layer multi-path formation,” Comput. Model. Eng. Sci. 126(1), 5–23 (2021).https://doi.org/10.32604/cmes.2021.014693
Google ScholarCrossref
H. Yin and S. D. Felicelli , “ Dendrite growth simulation during solidification in the LENS process,” Acta Mater. 58(4), 1455–1465 (2010).https://doi.org/10.1016/j.actamat.2009.10.053
Google ScholarCrossref
P. Nie , O. A. Ojo , and Z. Li , “ Numerical modeling of microstructure evolution during laser additive manufacturing of a nickel-based superalloy,” Acta Mater. 77, 85–95 (2014).https://doi.org/10.1016/j.actamat.2014.05.039
Google ScholarCrossref
Z. Liu and H. Qi , “ Effects of substrate crystallographic orientations on crystal growth and microstructure formation in laser powder deposition of nickel-based superalloy,” Acta Mater. 87, 248–258 (2015).https://doi.org/10.1016/j.actamat.2014.12.046
Google ScholarCrossref
L. Wei , L. Xin , W. Meng , and H. Weidong , “ Cellular automaton simulation of the molten pool of laser solid forming process,” Acta Phys. Sin. 64(01), 018103–018363 (2015).https://doi.org/10.7498/aps.64.018103
Google ScholarCrossref
R. Acharya , J. A. Sharon , and A. Staroselsky , “ Prediction of microstructure in laser powder bed fusion process,” Acta Mater. 124, 360–371 (2017).https://doi.org/10.1016/j.actamat.2016.11.018
Google ScholarCrossref
M. R. Rolchigo and R. LeSar , “ Modeling of binary alloy solidification under conditions representative of additive manufacturing,” Comput. Mater. Sci. 150, 535–545 (2018).https://doi.org/10.1016/j.commatsci.2018.04.004
Google ScholarCrossref
S. Geng , P. Jiang , L. Guo , X. Gao , and G. Mi , “ Multi-scale simulation of grain/sub-grain structure evolution during solidification in laser welding of aluminum alloys,” Int. J. Heat Mass Transfer 149, 119252 (2020).https://doi.org/10.1016/j.ijheatmasstransfer.2019.119252
Google ScholarCrossref
W. L. Wang , W. Q. Liu , X. Yang , R. R. Xu , and Q. Y. Dai , “ Multi-scale simulation of columnar-to-equiaxed transition during laser selective melting of rare earth magnesium alloy,” J. Mater. Sci. Technol. 119, 11–24 (2022).https://doi.org/10.1016/j.jmst.2021.12.029
Google ScholarCrossref
Q. Xia , J. Yang , and Y. Li , “ On the conservative phase-field method with the N-component incompressible flows,” Phys. Fluids 35, 012120 (2023).https://doi.org/10.1063/5.0135490
Google ScholarCrossref
Q. Xia , G. Sun , J. Kim , and Y. Li , “ Multi-scale modeling and simulation of additive manufacturing based on fused deposition technique,” Phys. Fluids 35, 034116 (2023).https://doi.org/10.1063/5.0141316
Google ScholarCrossref
A. Hussein , L. Hao , C. Yan , and R. Everson , “ Finite element simulation of the temperature and stress fields in single layers built without-support in selective laser melting,” Mater. Des. 52, 638–647 (2013).https://doi.org/10.1016/j.matdes.2013.05.070
Google ScholarCrossref
J. Ding , P. Colegrove , J. Mehnen , S. Ganguly , P. M. Sequeira Almeida , F. Wang , and S. Williams , “ Thermo-mechanical analysis of wire and arc additive layer manufacturing process on large multi-layer parts,” Comput. Mater. Sci. 50(12), 3315–3322 (2011).https://doi.org/10.1016/j.commatsci.2011.06.023
Google ScholarCrossref
Y. Du , X. You , F. Qiao , L. Guo , and Z. Liu , “ A model for predicting the temperature field during selective laser melting,” Results Phys. 12, 52–60 (2019).https://doi.org/10.1016/j.rinp.2018.11.031
Google ScholarCrossref
X. Luo , M. Liu , L. Zhenhua , H. Li , and J. Shen , “ Effect of different heat-source models on calculated temperature field of selective laser melted 18Ni300,” Chin. J. Lasers 48(14), 1402005–1402062 (2021).https://doi.org/10.3788/CJL202148.1402005
Google ScholarCrossref
J. F. Li , L. Li , and F. H. Stott , “ Thermal stresses and their implication on cracking during laser melting of ceramic materials,” Acta Mater. 52(14), 4385–4398 (2004).https://doi.org/10.1016/j.actamat.2004.06.005
Google ScholarCrossref
P. Aggarangsi and J. L. Beuth , “ Localized preheating approaches for reducing residual stress in additive manufacturing,” paper presented at the 2006 International Solid Freeform Fabrication Symposium, The University of Texas in Austin on August 14–16, 2006.
K. Dai and L. Shaw , “ Thermal and mechanical finite element modeling of laser forming from metal and ceramic powders,” Acta Mater. 52(1), 69–80 (2004).https://doi.org/10.1016/j.actamat.2003.08.028
Google ScholarCrossref
A. H. Nickel , D. M. Barnett , and F. B. Prinz , “ Thermal stresses and deposition patterns in layered manufacturing,” Mater. Sci. Eng. A 317(1–2), 59–64 (2001).https://doi.org/10.1016/S0921-5093(01)01179-0
Google ScholarCrossref
M. F. Zaeh and G. Branner , “ Investigations on residual stresses and deformations in selective laser melting,” Prod. Eng. 4(1), 35–45 (2010).https://doi.org/10.1007/s11740-009-0192-y
Google ScholarCrossref
P. Bian , J. Shi , Y. Liu , and Y. Xie , “ Influence of laser power and scanning strategy on residual stress distribution in additively manufactured 316L steel,” Opt. Laser Technol. 132, 106477 (2020).https://doi.org/10.1016/j.optlastec.2020.106477
Google ScholarCrossref
B. M. Marques , C. M. Andrade , D. M. Neto , M. C. Oliveira , J. L. Alves , and L. F. Menezes , “ Numerical analysis of residual stresses in parts produced by selective laser melting process,” Procedia Manuf. 47, 1170–1177 (2020).https://doi.org/10.1016/j.promfg.2020.04.167
Google ScholarCrossref
W. Mu , “ Numerical simulation of SLM forming process and research and prediction of forming properties,” MA thesis ( Anhui Jianzhu University, 2022).
Google Scholar
Y. Zhang , “ Multi-scale multi-physics modeling of laser powder bed fusion process of metallic materials with experiment validation,” Ph.D. dissertation ( Purdue University, 2018).
Google Scholar
Y. Qian , “ Mesoscopic simulation studies of key processing issues for powder bed fusion technology,” Ph.D. dissertation ( Tsinghua University, 2019).
Google Scholar
N. V. Brilliantov , S. Frank , J.-M. Hertzsch , and T. Pöschel , “ Model for collisions in granular gases,” Phys. Rev. E 53(5), 5382–5392 (1996).https://doi.org/10.1103/PhysRevE.53.5382
Google ScholarCrossref
Z. Xiao , “ Research on microscale selective laser melting process of high strength pure copper specimens,” MA thesis ( Hunan University, 2022).
Google Scholar
Z. Li , K. Mukai , M. Zeze , and K. C. Mills , “ Determination of the surface tension of liquid stainless steel,” J. Mater. Sci. 40(9–10), 2191–2195 (2005).https://doi.org/10.1007/s10853-005-1931-x
Google ScholarCrossref
R. Scardovelli and S. Zaleski , “ Analytical relations connecting linear interfaces and volume fractions in rectangular grids,” J. Comput. Phys. 164(1), 228–237 (2000).https://doi.org/10.1006/jcph.2000.6567
Google ScholarCrossref
D.-W. Cho , W.-I. Cho , and S.-J. Na , “ Modeling and simulation of arc: Laser and hybrid welding process,” J. Manuf. Process. 16(1), 26–55 (2014).https://doi.org/10.1016/j.jmapro.2013.06.012
Google ScholarCrossref
76.Flow3D. Version 11.1.0: User Manual ( FlowScience, Santa Fe, NM, USA, 2015).
Y. Tian , L. Yang , D. Zhao , Y. Huang , and J. Pan , “ Numerical analysis of powder bed generation and single track forming for selective laser melting of ss316l stainless steel,” J. Manuf. Process. 58, 964–974 (2020).https://doi.org/10.1016/j.jmapro.2020.09.002
Google ScholarCrossref
C. Tang , K. Q. Le , and C. H. Wong , “ Physics of humping formation in laser powder bed fusion,” Int. J. Heat Mass Transfer 149, 119172 (2020).https://doi.org/10.1016/j.ijheatmasstransfer.2019.119172
Google ScholarCrossref
L. Cao , “ Mesoscopic-scale simulation of pore evolution during laser powder bed fusion process,” Comput. Mater. Sci. 179, 109686 (2020).https://doi.org/10.1016/j.commatsci.2020.109686
Google ScholarCrossref
R. Li , J. Liu , Y. Shi , W. Li , and W. Jiang , “ Balling behavior of stainless steel and nickel powder during selective laser melting process,” Int. J. Adv. Manuf. Technol. 59(9–12), 1025–1035 (2012).https://doi.org/10.1007/s00170-011-3566-1
Google ScholarCrossref
S. A. Khairallah and A. Anderson , “ Mesoscopic simulation model of selective laser melting of stainless steel powder,” J. Mater. Process. Technol. 214(11), 2627–2636 (2014).https://doi.org/10.1016/j.jmatprotec.2014.06.001
Google ScholarCrossref
J. Liu , D. Gu , H. Chen , D. Dai , and H. Zhang , “ Influence of substrate surface morphology on wetting behavior of tracks during selective laser melting of aluminum-based alloys,” J. Zhejiang Univ. Sci. A 19(2), 111–121 (2018).https://doi.org/10.1631/jzus.A1700599
Google ScholarCrossref
L. Li , J. Li , and T. Fan , “ Phase-field modeling of wetting and balling dynamics in powder bed fusion process,” Phys. Fluids 33, 042116 (2021).https://doi.org/10.1063/5.0046771
Google ScholarCrossref
X. Nie , Z. Hu , H. Zhu , Z. Hu , L. Ke , and X. Zeng , “ Analysis of processing parameters and characteristics of selective laser melted high strength Al-Cu-Mg alloys: from single tracks to cubic samples,” J. Mater. Process. Technol. 256, 69–77 (2018).https://doi.org/10.1016/j.jmatprotec.2018.01.030
Google ScholarCrossref

원본 바로가기

에스티아이씨앤디가 대한민국 최대 생산제조기술 전회(SIMTOS 2024)에 참여합니다.

주식회사 에스티아이씨앤디에서는 대한민국 최대 생산제조기술 전회(SIMTOS 2024)에 참가하여 FLOW-3D를 주조품 연구 개발에 활용하는 다양한 사례와 제품 정보 소개를 드릴 예정입니다.
2024년 4월 1일(월)~4월 5일(금)까지 KINTEX 제1, 2전시장에서 만나뵐 수 있습니다.
전시장에서는 FLOW-3D를 이용한 다양한 연구 개발 주제에 대해 전시, 홍보 예정이오니 여러분들의 많은 참여를 부탁드립니다.

* 전시회의 부스명은 ‘에스티아이씨앤디’입니다. 부스를 찾으실 때 ‘에스티아이씨앤디’으로 찾아주세요!
주식회사 에스티아이씨앤디 Booth No. : 08E310

무료 관람을 원하시는 분께서는 저희에게 “회사명 / 성함 / 휴대번호”를 알려주시면, 초청장을 보내드리겠습니다.

연락처 : 02-2026-0455
이메일 : flow3d@stikorea.co.kr

감사합니다.

FLOW-3D 교육사례

고객 맞춤형 교육

FLOW-3D를 이용한 연구개발을 검토하시거나 구매하신 고객에게, 당사에서는 현업에서 해결해야 하는 문제들에 대한 FLOW-3D 사용법과 맞춤형 해석 교육을 지원합니다.

고객 맞춤 교육과정을 통해 수치해석에 대한 부담 없이 업무 지식을 기반으로 FLOW-3D 를 쉽게 이해하고 연구 개발에 활용할 수 있습니다.

Education Customer List

한국생산기술연구원 UST FLOW-3D AM 교육
한국해양과학기술원 해양분야 교육(세굴,GMO)
현대자동차 FLOW-3D CAST 업그레이드 버전 교육
SK 이노베이션 FLOW-3D 교육
세메스 FLOW-3D 교육
LG 이노텍 FLOW-3D 교육
디알엑시온 FLOW-3D CAST 교육
지오시스템리서치 FLOW-3D 교육
수엔지니어링 FLOW-3D 교육
서울대학교 FLOW-3D WELD 교육
한양대학교 FLOW-3D WELD 교육
한국수자원공사 FLOW-3D 교육
율촌화학 FLOW-3D 교육
세종대학교 FLOW-3D 교육
동명대학교 FLOW-3D WELD 교육
한경대학교 FLOW-3D 해양 교육
한국기계연구원 FLOW-3D 교육
한국생산기술연구원 FLOW-3D 교육
한국원자력연구원 WELD 교육
LS 전선 FLOW-3D 해석 교육
LS메탈 FLOW-3D 해석 교육
아람코퍼레이션 수자원분야 FLOW-3D 해석 교육
성균관대학교 FLOW-3D 교육
삼성 SDI WELD 교육
농어촌공사 안전진단본부 측수로 불완전원류 해석 교육
동아대학교 신소재공학과 WELD 교육
한경대학교 환경수리학 연구실 해석 교육
LG에너지 솔루션 WELD 교육
현대 트랜시스 CAST 교육

이 외에도 수 많은 고객들이 교육을 통해 FLOW-3D 수치해석을 경험하고 연구개발에 활용하고 있습니다.
연구개발 업무에 어떻게 활용할 수 있는지에 대해 문의 사항이 있으시면 언제든지 아래 연락처로 연락 주십시오.
감사합니다.

고객센터 및 교육 담당자

전화 : 02)2026-0450, 02)2026-0455
이메일 : flow3d@stikorea.co.kr

교육 장소 안내

지하철 1호선/가산디지털단지역 (8번출구), 지하철 7호선/가산디지털단지역 (5번출구)
우림라이온스밸리 B동 302호 또는 교육장
당사 건물에 주차할 경우 무료 주차 1시간만 지원되오니, 가능하면 대중교통을 이용해 주시기 바랍니다.

오시는 길

Figura 1 – Mapa de localização da PCH Salto Paraopeba

하천 저수지 물리적 모델의 침적 과정에 대한 전산 유체 역학 모델링(CFD) 기준

Natália Melo da Silva1¹; Jorge Luis Zegarra Tarqui²,Edna Maria de Faria Viana ³

Abstract

저수지 침전은 수력 발전의 지속 가능한 발전을 위한 주요 문제 중 하나이며 브라질에 매우 중요합니다. 브라질의 주요 에너지원은 수력발전소에서 나옵니다. 소규모 수력 발전소(SHP)는 재생 에너지의 보완적 발전을 위한 중요한 대안입니다.

이들의 설계, 건설, 운영 및 재동력을 최적화하기 위해 저수지 내 퇴적물의 유체 역학 및 이동을 연구하는 것이 매우 중요합니다.

3차원 전산유체역학 – CFD 3D 모델링은 복잡한 흐름 문제에 가장 적합한 방법입니다. 제안된 방법은 MG Jeceaba 자치구에 위치한 PCH Salto Paraopeba의 유체 역학 및 퇴적물 이동 현상을 재현하고 평가하는 것을 목표로 하며, 취수구의 완전한 침전으로 인해 작동이 중단되었습니다.

모델의 검증은 미나스제라이스 연방대학교의 수력학 연구 센터(CPH)에 구축된 축소된 물리적 모델의 실험 데이터를 사용하여 수행됩니다.

Abstract: The reservoir silting is one of the main problems for sustainable development in the
generation of hydroelectric energy and it is of great significance for Brazil. The main source of energy
in Brazil comes from hydroelectric power plant. The Small Hydroelectric Power Plant (SHP) are an
important alternative for complementary generation of renewable energy.
Seeking to optimize the design, construction, operation, and repowering of these, it is extremely
important to study the hydrodynamics and transport of sediments in their reservoirs. Threedimensional Computational Fluid Dynamics – CFD 3D modeling is the most appropriate method for
complex flow problems. The proposed method aims to reproduce and evaluate the hydrodynamic and
sediment transport phenomena of the PCH Salto Paraopeba, located in the municipality of Jeceaba,
MG, which stopped working due to the complete silting up of its water intake. The validation of the
model will be done using experimental data from the reduced physical model, built at the Hydraulic
Research Center (CPH) at the Federal University of Minas Gerais.

Keywords

퇴적물 수송, 물리적 모델, 소규모 수력 발전소, Sediment transport, physical model, Small Hydroelectric Power Plant.

Figura 2 – PCH Salto Paraopeba e modelo reduzido.

REFERÊNCIAS

ALBARELLO, L. “Guia para implantação de pequenas centrais hidrelétricas- PCHs”. Dissertação de Mestrado. Trabalho de Conclusão de Curso de Especialista. Programa de Pós-Graduação em Eficiência Energética Aplicada aos Processos Produtivos. UFSM, Panambi /RS, 2014.

ANEEL, SIGA – Sistema de Informações de Geração da ANEEL. Disponível em: . Acesso em 10 de maio de 2023.

CAMPELLO, B.S.C. “Estudo Da Velocidade de Queda e do Início do Movimento das Partículas de Borracha e Areia”. Dissertação de Mestrado. Programa de Pós-Graduação em Saneamento, Meio Ambiente e Recursos Hídricos. UFMG, Belo Horizonte /MG, 2017.

CAMPOS, A.S. “A Importância do Coeficiente de Rugosidade de Manning na Avaliação Numérica
do Assoreamento de Reservatórios A Fio D’água”. Dissertação de Mestrado. Programa de PósGraduação em Saneamento, Meio Ambiente e Recursos Hídricos. UFMG, Belo Horizonte /MG,
2018.

CARVALHO, N. O. et al. Guia de avaliação de assoreamento de reservatórios, ANEEL, Brasília,
2000.

CARVALHO, N. O. Hidrossedimentologia prática, CPRM, Rio de Janeiro, 1994.

CARVALHO, N. O. Hidrossedimentologia prática, CPRM, Rio de Janeiro, 2008.

EMIG, PCH Salto do Paraopeba. Disponível em: < https://www.cemig.com.br/usina/pch-salto-doparaopeba/>. Acesso em 12 de maio de 2023.

ELETROBRÁS; Instituto de Pesquisas Hidráulicas – IPH. Diagnóstico das Condições sedimentológicas dos Principais Rios Brasileiros. UFRGS, Rio de Janeiro, 1992.

FLOW-3D®. FLOW-3D® 2022R2 – User Manual. Disponível em: < https://users.flow3d.com/flow3d/> FORTUNA, A.O. (2000). Técnicas Computacionais para Dinâmica dos Fluidos – Conceitos Básicos e Aplicações. São Paulo – SP.

GONÇALVES, M.O. “Análise Comparativa Entre Modelo Reduzido e Modelos Computacionais Uni
e Bidimensionais”. Dissertação de Mestrado. Programa de Pós-Graduação em Engenharia de
Recursos Hídricos e Ambiental. UFRP, Curitiba/PR, 2017.

HILLEBRAND, G.; KLASSEN, I.; OLSEN, N. R. B. (2017). “3D CFD modelling of velocities and
sediment transport in the Iffezheim hydropower reservoir”. Hydrology Research 48 (1), pp. 147–159.
JULIEN, P. Y. (2010). Erosion and sedimentation, Cambridge University Press, Cambridge, UK, 2nd
edn.

MAHMOOD, K. (1987). Reservoir sedimentation – impact, extent and mitigation. World Bank Tech.
Paper No. 71. Washington, DC.

MIRANDA, R.B. “A influência do assoreamento na geração de energia hidrelétrica: estudo de caso
na Usina de Três Irmãos – SP”. Dissertação de Mestrado. Programa de Pós-Graduação em Ciências
da Engenharia Ambiental. USP, São Carlos/SP, 2011.

MOHAMMAD, M.E.; AL-ANSARI, N.; KNUTSSON, S.; LAUE, J. (2020). “A Computational
Fluid Dynamics Simulation Model of Sediment Deposition in a Storage Reservoir Subject to Water
Withdrawal”. Water, 12, 959.

OLIVEIRA, M. A.“Repotenciação de pequenas centrais hidrelétricas: Avaliação técnica e
econômica”. Dissertação de Mestrado. Programa de Pós-Graduação em Engenharia de Energia.
UNIFEI, Itajubá/MG, 2012.

SALIBA, A.P.M. Notas de aula, Modelos fundo móvel, Disciplina Introdução a Modelagem Física
em Engenharia, Universidade Federal de Minas Gerais. Belo Horizonte, 2020.

SOARES, W.S. “Taxa de Assoreamento no Reservatório da Usina Hidrelétrica do Funil, MG”.
Dissertação de Mestrado. Programa de Pós-Graduação em Tecnologias e Inovações Ambientais.
UFLA, Lavras/MG, 2015.

원문보기

비선형 파력의 영향에 따른 잔해 언덕 방파제 형상의 효과에 대한 수치 분석

Numerical Analysis of the Effects of Rubble Mound Breakwater Geometry Under the Effect of Nonlinear Wave Force

Research Article-Civil Engineering
Published: 07 December 2023
(2023)

22 Accesses
Explore all metrics

Cite this article

Abstract

Assessing the interaction of waves and porous offshore structures such as rubble mound breakwaters plays a critical role in designing such structures optimally. This study focused on the effect of the geometric parameters of a sloped rubble mound breakwater, including the shape of the armour, method of its arrangement, and the breakwater slope. Thus, three main design criteria, including the wave reflection coefficient (K_r), transmission coefficient (K_t), and depreciation wave energy coefficient (K_d), are discussed. Based on the results, a decrease in wavelength reduced the K_r and increased the K_t and K_d. The rubble mound breakwater with the Coreloc armour layer could exhibit the lowest K_r compared to other armour geometries. In addition, a decrease in the breakwater slope reduced the K_r and K_d by 3.4 and 1.25%, respectively. In addition, a decrease in the breakwater slope from 33 to 25° increased the wave breaking height by 6.1% on average. Further, a decrease in the breakwater slope reduced the intensity of turbulence depreciation. Finally, the armour geometry and arrangement of armour layers on the breakwater with its different slopes affect the wave behaviour and interaction between the wave and breakwater. Thus, layering on the breakwater and the correct use of the geometric shapes of the armour should be considered when designing such structures.

파도와 잔해 더미 방파제와 같은 다공성 해양 구조물의 상호 작용을 평가하는 것은 이러한 구조물을 최적으로 설계하는 데 중요한 역할을 합니다. 본 연구는 경사진 잔해 둔덕 방파제의 기하학적 매개변수의 효과에 초점을 맞추었는데, 여기에는 갑옷의 형태, 배치 방법, 방파제 경사 등이 포함된다. 따라서 파동 반사 계수(K_r), 투과 계수(K_t) 및 감가상각파 에너지 계수(K_d)에 대해 논의합니다. 결과에 따르면 파장이 감소하면 K가 감소합니다._r그리고 K를 증가시켰습니다_t 및 K_d. Coreloc 장갑 층이 있는 잔해 언덕 방파제는 가장 낮은 K를 나타낼 수 있습니다._r 다른 갑옷 형상과 비교했습니다. 또한 방파제 경사가 감소하여 K가 감소했습니다._r 및 K_d 각각 3.4%, 1.25% 증가했다. 또한 방파제 경사가 33°에서 25°로 감소하여 파도 파쇄 높이가 평균 6.1% 증가했습니다. 또한, 방파제 경사의 감소는 난류 감가상각의 강도를 감소시켰다. 마지막으로, 경사가 다른 방파제의 장갑 형상과 장갑 층의 배열은 파도 거동과 파도와 방파제 사이의 상호 작용에 영향을 미칩니다. 따라서 이러한 구조를 설계 할 때 방파제에 층을 쌓고 갑옷의 기하학적 모양을 올바르게 사용하는 것을 고려해야합니다.

Keywords

Rubble mound breakwater
Computational fluid dynamics
Armour layer
Wave reflection coefficient
Wave transmission coefficient
Wave energy dissipation coefficient

References

Sollitt, C.K.; Cross, R.H.: Wave transmission through permeable breakwaters. In Coastal Engineering. pp. 1827–1846. (1973)
Sulisz, W.: Wave reflection and transmission at permeable breakwaters of arbitrary cross-section. Coast. Eng. 9(4), 371–386 (1985)Article Google Scholar
Kobayashi, N.; Wurjanto, A.: Numerical model for waves on rough permeable slopes. J. Coast. Res.149–166. (1990)
Wurjanto, A.; Kobayashi, N.: Irregular wave reflection and runup on permeable slopes. J. Waterw. Port Coast. Ocean Eng. 119(5), 537–557 (1993)Article Google Scholar
van Gent, M.R.: Numerical modelling of wave interaction with dynamically stable structures. In Coastal Engineering 1996. pp. 1930–1943. (1997)
Liu, P.L.F.; Wen, J.: Nonlinear diffusive surface waves in porous media. J. Fluid Mech. 347, 119–139 (1997)Article MathSciNet MATH Google Scholar
Troch, P.; De Rouck, J.: Development of two-dimensional numerical wave flume for wave interaction with rubble mound breakwaters. In Coastal Engineering. pp. 1638–1649. (1999)
Liu, P.L.F.; Lin, P.; Chang, K.A.; Sakakiyama, T.: Numerical modeling of wave interaction with porous structures. J. Waterw. Port Coast. Ocean Eng. 125(6), 322–330 (1999)Article Google Scholar
Abdolmaleki, K.; Thiagarajan, K.P.; Morris-Thomas, M.T.: Simulation of the dam break problem and impact flows using a Navier-Stokes solver. Simulation 13, 17 (2004)Google Scholar
Higuera, P.; Lara, J.L.; Losada, I.J.: Realistic wave generation and active wave absorption for Navier-Stokes models: application to OpenFOAM®. Coast. Eng. 71, 102–118 (2013)Article Google Scholar
Higuera, P.; Lara, J.L.; Losada, I.J.: Three-dimensional interaction of waves and porous coastal structures using OpenFOAM®. Part II: application. Coast. Eng. 83, 259–270 (2014)Article Google Scholar
Gui, Q.; Dong, P.; Shao, S.; Chen, Y.: Incompressible SPH simulation of wave interaction with porous structure. Ocean Eng. 110, 126–139 (2015)Article Google Scholar
Dentale, F.; Donnarumma, G.; Carratelli, E.P.; Reale, F.: A numerical method to analyze the interaction between sea waves and rubble mound emerged breakwaters. WSEAS Trans. Fluid Mech 10, 106–116 (2015)Google Scholar
Dentale, F.; Reale, F.; Di Leo, A.; Carratelli, E.P.: A CFD approach to rubble mound breakwater design. Int. J. Naval Archit. Ocean Eng. 10(5), 644–650 (2018)Article Google Scholar
Koley, S.: Wave transmission through multilayered porous breakwater under regular and irregular incident waves. Eng. Anal. Bound. Elem. 108, 393–401 (2019)Article MathSciNet MATH Google Scholar
Koley, S.; Panduranga, K.; Almashan, N.; Neelamani, S.; Al-Ragum, A.: Numerical and experimental modeling of water wave interaction with rubble mound offshore porous breakwaters. Ocean Eng. 218, 108218 (2020)Article Google Scholar
Pourteimouri, P.; Hejazi, K.: Development of an integrated numerical model for simulating wave interaction with permeable submerged breakwaters using extended Navier-Stokes equations. J. Mar. Sci. Eng. 8(2), 87 (2020)Article Google Scholar
Cao, D.; Yuan, J.; Chen, H.: Towards modelling wave-induced forces on an armour layer unit of rubble mound coastal revetments. Ocean Eng. 239, 109811 (2021)Article Google Scholar
Díaz-Carrasco, P.; Eldrup, M.R.; Andersen, T.L.: Advance in wave reflection estimation for rubble mound breakwaters: the importance of the relative water depth. Coast. Eng. 168, 103921 (2021)Article Google Scholar
Vieira, F.; Taveira-Pinto, F.; Rosa-Santos, P.: Damage evolution in single-layer cube armoured breakwaters with a regular placement pattern. Coast. Eng. 169, 103943 (2021)Article Google Scholar
Booshi, S.; Ketabdari, M.J.: Modeling of solitary wave interaction with emerged porous breakwater using PLIC-VOF method. Ocean Eng. 241, 110041 (2021)Article Google Scholar
Aristodemo, F.; Filianoti, P.; Tripepi, G.; Gurnari, L.; Ghaderi, A.: On the energy transmission by a submerged barrier interacting with a solitary wave. Appl. Ocean Res. 122, 103123 (2022)Article Google Scholar
Teixeira, P.R.; Didier, E.: Numerical analysis of performance of an oscillating water column wave energy converter inserted into a composite breakwater with rubble mound foundation. Ocean Eng. 278, 114421 (2023)Article Google Scholar
Burgan, H.I.: Numerical modeling of structural irregularities on unsymmetrical buildings. Tehnički vjesnik 28(3), 856–861 (2021)Google Scholar
Jones, I.P.: CFDS-Flow3D user guide. (1994)
Al Shaikhli, H.I.; Khassaf, S.I.: Stepped mound breakwater simulation by using flow 3D. Eurasian J. Eng. Technol. 6, 60–68 (2022)Google Scholar
Hirt, C.W.; Nichols, B.D.: Volume of fluid (VOF) method for the dynamics of free boundaries. J. Comput. Phys. 39(1), 201–225 (1981)Article MATH Google Scholar
Ghaderi, A.; Dasineh, M.; Aristodemo, F.; Aricò, C.: Numerical simulations of the flow field of a submerged hydraulic jump over triangular macroroughnesses. Water 13(5), 674 (2021)Article Google Scholar
Yakhot, V.; Orszag, S.A.; Thangam, S.; Gatski, T.B.; Speziale, C.G.: Development of turbulence models for shear flows by a double expansion technique. Phys. Fluids A 4(7), 1510–1520 (1992)Article MathSciNet MATH Google Scholar
Van der Meer, J.W.; Stam, C.J.M.: Wave runup on smooth and rock slopes of coastal structures. J. Waterw. Port Coast. Ocean Eng. 118(5), 534–550 (1992)Article Google Scholar
Goda, Y.; Suzuki, Y. Estimation of incident and reflected waves in random wave experiments. In: ASCE, Proceedings of 15th International Conference on Coastal Engineering, (Honolulu, Hawaii). vol. 1, pp. 828–845. (1976)
Zanuttigh, B.; Van der Meer, J.W.: Wave reflection from coastal structures. In: AA.VV., Proceedings of the XXX International Conference on Coastal Engineering, World Scientific, (San Diego, CA, USA, September 2006). pp. 4337–4349. (2006)
Seelig W.N.; Ahrens J.P.: Estimation of wave reflection and energy dissipation coefficients for beaches, revetments, and breakwaters. CERC, Technical Paper, Fort Belvoir. vol. 81, p. 41 (1981)
Mase, H.: Random wave runup height on gentle slope. J. Waterw. Port Coast. Ocean Eng. 115(5), 649–661 (1989)Article Google Scholar

원본 바로가기

FLOW-3D 2024R1 의 새로운 기능

FLOW-3D 2024R1 의 새로운 기능

FLOW-3D 2024R1은 버블 및 상변화 모델의 수정을 통해 제품 및 공정 개발 소프트웨어를 계속 개선하고 있으며, 이를 통해 특히 열 전달 또는 액체-증기 상변화 옵션을 사용할 때 일반적인 설정 오류를 피하면서 더 쉽게 사용할 수 있습니다. 사용자 인터페이스를 재구성하여 액체-증기 상변화 옵션을 고체-액체 상변화 옵션으로 그룹화합니다. 단열 버블 및 열 버블 모델을 통합된 이상 기체 상태 방정식으로 대체하고, 유체 특성 입력을 통합했으며, 상태 방정식을 정의하는 데 사용되는 매개 변수를 제어하는 옵션을 추가했습니다. 이 개발은 엔지니어링 오류의 가능성을 줄이고, 입력을 단순화하며, 상전이 모델에 대한 보다 자연스러운 그룹화를 제공합니다. 두 번째 개발은 새로운 EXODUS II 기반 출력 파일에서 유체-구조 상호작용 및 열 응력 진화 모델을 지원하여 후처리 성능을 크게 향상시킵니다.

FLOW-3D 2023R2 의 새로운 기능

새로운 결과 파일 형식

FLOW-3D POST 2023R2 는 EXODUS II 형식을 기반으로 하는 완전히 새로운 결과 파일 형식을 도입하여 더 빠른 후처리를 가능하게 합니다. 이 새로운 파일 형식은 크고 복잡한 시뮬레이션의 후처리 작업에 소요되는 시간을 크게 줄이는 동시에(평균 최대 5배!) 다른 시각화 도구와의 연결성을 향상시킵니다.

FLOW-3D POST 2023R2 에서 사용자는 이제 selected data를 flsgrf , EXODUS II 둘중 하나 또는 flsgrf 와 EXODUS II 둘다 파일 형식으로 쓸 수 있습니다 . 새로운 EXODUS II 파일 형식은 각 객체에 대해 유한 요소 메쉬를 활용하므로 사용자는 다른 호환 가능한 포스트 프로세서 및 FEA 코드를 사용하여 FLOW-3D 결과를 열 수도 있습니다. 새로운 워크플로우를 통해 사용자는 크고 복잡한 사례를 신속하게 시각화하고 임의 위치에서의 슬라이싱, 볼륨 렌더링 및 통계를 사용하여 추가 정보를 추출할 수 있습니다.

레이 트레이싱을 이용한 화장품 크림 충전 — ***FLOW-3D*** **POST** 의 새로운 EXODUS II 파일 형식으로 채워진 화장품 크림 모델의 향상된 광선 추적 기능의 예

새로운 결과 파일 형식은 솔버 엔진의 성능을 저하시키지 않으면서 flsgrf 에 비해 시각화 작업 흐름에서 놀라운 속도 향상을 자랑합니다. 이 흥미로운 새로운 개발은 결과 분석의 속도와 유연성이 향상되어 원활한 시뮬레이션 경험을 제공합니다.

FLOW-3D POST 의 새로운 시각화 기능 에 대해 자세히 알아보세요 .

난류 모델 개선

FLOW-3D 2023R2는 two-equation(RANS) 난류 모델에 대한 dynamic mixing length 계산을 크게 개선했습니다. 거의 층류 흐름 체계와 같은 특정 제한 사례에서는 이전 버전의 코드 계산 한계가 때때로 과도하게 예측되어 사용자가 특정 mixing length를 수동으로 입력해야 할 수 있습니다.

새로운 dynamic mixing length 계산은 이러한 상황에서 난류 길이와 시간 척도를 더 잘 설명합니다. 이제 사용자는 고정된(물리 기반) mixing length를 설정하는 대신 더 넓은 범위의 흐름에 동적 모델을 적용할 수 있습니다.

접촉식 탱크 혼합 시뮬레이션 — 적절한 고정 mixing length와 비교하여 접촉 탱크의 혼합 시뮬레이션을 위한 기존 동적 mixing length 모델과 새로운 동적 mixing length 모델 간의 비교

정수압 초기화

사용자가 미리 정의된 유체 영역에서 정수압을 초기화해야 하는 경우가 많습니다. 이전에는 대규모의 복잡한 시뮬레이션에서 정수압 솔버의 수렴 속도가 느려지는 경우가 있었습니다. FLOW-3D 2023R2는 정수압 솔버의 성능을 크게 향상시켜 전처리 단계에서 최대 6배 빠르게 수렴할 수 있도록 해줍니다.

압축성 흐름 솔버 성능

FLOW-3D 2023R2는 최적화된 압력 솔버를 도입하여 압축성 흐름 문제에 대해 상당한 성능 향상을 제공합니다. 압축성 제트 흐름의 예에서 2023R2 솔버는 2023R1 버전보다 최대 4배 빠릅니다.

FLOW-3D 2023R2 의 새로운 기능

FLOW-3D 소프트웨어 제품군의 모든 제품은 2023R2에서 IT 관련 개선 사항을 받았습니다. FLOW-3D 2023R2은 이제 Windows 11 및 RHEL 8을 지원합니다. Linux 설치 프로그램은 누락된 종속성을 보고하도록 개선되었으며 더 이상 루트 수준 권한이 필요하지 않으므로 설치가 더 쉽고 안전해집니다. 그리고 워크플로우를 자동화한 분들을 위해 입력 파일 변환기에 명령줄 인터페이스를 추가하여 스크립트 환경에서도 워크플로우가 업데이트된 입력 파일로 작동하는지 확인할 수 있습니다.

확장된 PQ ² 분석

제조에 사용되는 유압 시스템은 PQ ² 곡선을 사용하여 모델링할 수 있습니다. 장치의 세부 사항을 건너뛰고 흐름에 미치는 영향을 포함하기 위해 질량 운동량 소스 또는 속도 경계 조건을 사용하여 유압 시스템을 근사화하는 것이 편리하도록 단순화하는 경우가 많습니다. 우리는 기존 PQ ² 분석 모델을 확장하여 이러한 유형의 기하학적 단순화를 허용하면서도 현실적인 결과를 제공했습니다. 이로써 시뮬레이션 시간을 줄이고 모델 복잡성의 감소시킬 수 있습니다.

FLOW-3D 2022R2 의 새로운 기능

FLOW-3D 2022R2 제품군 출시로 Flow Science는 FLOW-3D 의 워크스테이션과 HPC 버전을 통합하여 노드 병렬 고성능 컴퓨팅 실행할 수 있도록 단일 노드 CPU 구성에서 다중 노드에 이르기까지 모든 유형의 하드웨어 아키텍처를 활용할 수 있는 단일 솔버 엔진을 제공합니다. 추가 개발에는 점탄성 흐름을 위한 새로운 로그 형태 텐서 방법, 지속적인 솔버 속도 성능 개선, 고급 냉각 채널 및 팬텀 구성요소 제어, entrained air 기능이 개선되었습니다.

통합 솔버

FLOW-3D 제품을 단일 통합 솔버로 마이그레이션하여 로컬 워크스테이션이나 고성능 컴퓨팅 하드웨어 환경에서 원활하게 실행할 수 있습니다.

많은 사용자가 노트북이나 로컬 워크스테이션에서 모델을 실행하지만, 고성능 컴퓨팅 클러스터에서 더 큰 모델을 실행합니다. 2022R2 릴리스에서는 통합 솔버를 통해 사용자가 HPC 솔루션의 Open MP/MPI 하이브리드 병렬화와 동일한 이점을 활용하여 워크스테이션과 노트북에서 실행할 수 있습니다.

메쉬 분해의 예 — Open MP/MPI 하이브리드 병렬화를 위한 메시 분해의 예

솔버 성능 개선

멀티 소켓 워크스테이션

다중 소켓 워크스테이션은 이제 매우 일반적이며 대규모 시뮬레이션을 실행할 수 있습니다. 새로운 통합 솔버를 사용하면 이러한 유형의 하드웨어를 사용하는 사용자는 일반적으로 HPC 클러스터 구성에서만 사용할 수 있었던 OpenMP/MPI 하이브리드 병렬화를 활용하여 모델을 실행할 수 있어 성능이 향상되는 것을 확인할 수 있습니다.

낮은 수준의 루틴으로 향상된 벡터화 및 메모리 액세스

대부분의 테스트 사례에서 10~20% 정도의 성능 향상이 관찰되었으며 일부 사례에서는 20%를 초과하는 런타임 이점이 나타났습니다.

정제된 체적 대류 안정성 한계

Time step 안정성 한계는 모델 런타임의 주요 요인이며, 2022R2에서는 새로운 time step 안정성 한계인 3D 대류 안정성 한계를 Numerics 탭에서 사용할 수 있습니다. 실행 중이고 대류가 제한된(cx, cy 또는 cz 제한) 모델의 경우 새 옵션은 일반적인 속도 향상을 30% 정도 보여줍니다.

압력 솔버 프리컨디셔너

경우에 따라 까다로운 유동 해석의 경우 과도한 압력 솔버 반복으로 인해 실행 시간이 길어질 수 있습니다. 이러한 어려운 경우 2022R2에서는 모델이 너무 많이 반복되면 FLOW-3D가 자동으로 새로운 프리컨디셔너 기능을 활성화하여 압력 수렴을 돕습니다. 런타임이 1.9~335배 더 빨라졌습니다!

점탄성 유체에 대한 로그 형태 텐서 방법

점탄성 유체에 대한 새로운 솔버 옵션을 사용자가 사용할 수 있으며 특히 높은 Weissenberg 수에 효과적입니다.

점탄성 흐름을 위한 개선된 솔루션 — 로그 구조 텐서 솔루션을 사용하여 점탄성 흐름에 대한 높은 Weissenberg 수의 개선된 솔루션의 예입니다. 제공: MF Tome 외, J. Non-Newton. Fluid. Mech. 175-176 (2012) 44–54

활성 시뮬레이션 제어 확장

Active simulation 제어 기능이 확장되어 연속 주조 및 적층 제조 응용 분야에 일반적으로 사용되는 팬텀 개체는 물론 주조 및 기타 여러 열 관리 응용 분야에 사용되는 냉각 채널에도 사용됩니다.

팬텀 물체 속도 제어의 예 — 연속 주조 응용 분야에 대한 가상 물체 속도 제어의 예

향상된 공기 동반 기능

디퓨저 및 이와 유사한 산업용 기포 흐름 응용 분야의 경우 이제 질량 공급원을 사용하여 물기둥에 공기를 유입할 수 있습니다. 또한, 동반된 공기 및 용존 산소의 난류 확산에 대한 기본값이 업데이트되었으며 매우 낮은 공기 농도에 대한 모델 정확도가 향상되었습니다.

디퓨저 모델의 예: 이제 질량 소스를 사용하여 물기둥에 공기를 유입할 수 있습니다.

FLOW-3D 아카이브 의 새로운 기능

FLOW-3D 2022R1 의 새로운 기능

FLOW-3D v12.0 의 새로운 기능

Figure 11. Sketch of scour mechanism around USAF under random waves.

Scour Characteristics and Equilibrium Scour Depth Prediction around Umbrella Suction Anchor Foundation under Random Waves

by Ruigeng Hu¹,Hongjun Liu²,Hao Leng¹,Peng Yu³ andXiuhai Wang^1,2,*

¹College of Environmental Science and Engineering, Ocean University of China, Qingdao 266000, China

²Key Lab of Marine Environment and Ecology (Ocean University of China), Ministry of Education, Qingdao 266000, China

³Qingdao Geo-Engineering Survering Institute, Qingdao 266100, China

^*Author to whom correspondence should be addressed.

J. Mar. Sci. Eng. 2021, 9(8), 886; https://doi.org/10.3390/jmse9080886

Received: 6 July 2021 / Revised: 8 August 2021 / Accepted: 13 August 2021 / Published: 17 August 2021

(This article belongs to the Section Ocean Engineering)

Download

Abstract

A series of numerical simulation were conducted to study the local scour around umbrella suction anchor foundation (USAF) under random waves. In this study, the validation was carried out firstly to verify the accuracy of the present model. Furthermore, the scour evolution and scour mechanism were analyzed respectively. In addition, two revised models were proposed to predict the equilibrium scour depth S_eq around USAF. At last, a parametric study was carried out to study the effects of the Froude number F_r and Euler number E_u for the S_eq. The results indicate that the present numerical model is accurate and reasonable for depicting the scour morphology under random waves. The revised Raaijmakers’s model shows good agreement with the simulating results of the present study when KC_s,p < 8. The predicting results of the revised stochastic model are the most favorable for n = 10 when KC_rms,a < 4. The higher F_r and E_u both lead to the more intensive horseshoe vortex and larger S_eq.

Keywords:

scour; numerical investigation; random waves; equilibrium scour depth; KC number

1. Introduction

The rapid expansion of cities tends to cause social and economic problems, such as environmental pollution and traffic jam. As a kind of clean energy, offshore wind power has developed rapidly in recent years. The foundation of offshore wind turbine (OWT) supports the upper tower, and suffers the cyclic loading induced by waves, tides and winds, which exerts a vital influence on the OWT system. The types of OWT foundation include the fixed and floating foundation, and the fixed foundation was used usually for nearshore wind turbine. After the construction of fixed foundation, the hydrodynamic field changes in the vicinity of the foundation, leading to the horseshoe vortex formation and streamline compression at the upside and sides of foundation respectively [1,2,3,4]. As a result, the neighboring soil would be carried away by the shear stress induced by vortex, and the scour hole would emerge in the vicinity of foundation. The scour holes increase the cantilever length, and weaken the lateral bearing capacity of foundation [5,6,7,8,9]. Moreover, the natural frequency of OWT system increases with the increase of cantilever length, causing the resonance occurs when the system natural frequency equals the wave or wind frequency [10,11,12]. Given that, an innovative foundation called umbrella suction anchor foundation (USAF) has been designed for nearshore wind power. The previous studies indicated the USAF was characterized by the favorable lateral bearing capacity with the low cost [6,13,14]. The close-up of USAF is shown in Figure 1, and it includes six parts: 1-interal buckets, 2-external skirt, 3-anchor ring, 4-anchor branch, 5-supporting rod, 6-telescopic hook. The detailed description and application method of USAF can be found in reference [13].

Figure 1. The close-up of umbrella suction anchor foundation (USAF).

Numerical and experimental investigations of scour around OWT foundation under steady currents and waves have been extensively studied by many researchers [1,2,15,16,17,18,19,20,21,22,23,24]. The seabed scour can be classified as two types according to Shields parameter θ, i.e., clear bed scour (θ < θ_cr) or live bed scour (θ > θ_cr). Due to the set of foundation, the adverse hydraulic pressure gradient exists at upstream foundation edges, resulting in the streamline separation between boundary layer flow and seabed. The separating boundary layer ascended at upstream anchor edges and developed into the horseshoe vortex. Then, the horseshoe vortex moved downstream gradually along the periphery of the anchor, and the vortex shed off continually at the lee-side of the anchor, i.e., wake vortex. The core of wake vortex is a negative pressure center, liking a vacuum cleaner. Hence, the soil particles were swirled into the negative pressure core and carried away by wake vortexes. At the same time, the onset of scour at rear side occurred. Finally, the wake vortex became downflow when the turbulence energy could not support the survival of wake vortex. According to Tavouktsoglou et al. [25], the scale of pile wall boundary layer is proportional to 1/ln(R_d) (R_d is pile Reynolds), which means the turbulence intensity induced by the flow-structure interaction would decrease with R_d increases, but the effects of R_d can be neglected only if the flow around the foundation is fully turbulent [26]. According to previous studies [1,15,27,28,29,30,31,32], the scour development around pile foundation under waves was significantly influenced by Shields parameter θ and KC number simultaneously (calculated by Equation (1)). Sand ripples widely existed around pile under waves in the case of live bed scour, and the scour morphology is related with θ and KC. Compared with θ, KC has a greater influence on the scour morphology [21,27,28]. The influence mechanism of KC on the scour around the pile is reflected in two aspects: the horseshoe vortex at upstream and wake vortex shedding at downstream.

KC=UwmTD��=�wm��(1)

where, U_wm is the maximum velocity of the undisturbed wave-induced oscillatory flow at the sea bottom above the wave boundary layer, T is wave period, and D is pile diameter.

There are two prerequisites to satisfy the formation of horseshoe vortex at upstream pile edges: (1) the incoming flow boundary layer with sufficient thickness and (2) the magnitude of upstream adverse pressure gradient making the boundary layer separating [1,15,16,18,20]. The smaller KC results the lower adverse pressure gradient, and the boundary layer cannot separate, herein, there is almost no horseshoe vortex emerging at upside of pile. Sumer et al. [1,15] carried out several sets of wave flume experiments under regular and irregular waves respectively, and the experiment results show that there is no horseshoe vortex when KC is less than 6. While the scale and lifespan of horseshoe vortex increase evidently with the increase of KC when KC is larger than 6. Moreover, the wake vortex contributes to the scour at lee-side of pile. Similar with the case of horseshoe vortex, there is no wake vortex when KC is less than 6. The wake vortex is mainly responsible for scour around pile when KC is greater than 6 and less than O(100), while horseshoe vortex controls scour nearly when KC is greater than O(100).

Sumer et al. [1] found that the equilibrium scour depth was nil around pile when KC was less than 6 under regular waves for live bed scour, while the equilibrium scour depth increased with the increase of KC. Based on that, Sumer proposed an equilibrium scour depth predicting equation (Equation (2)). Carreiras et al. [33] revised Sumer’s equation with m = 0.06 for nonlinear waves. Different with the findings of Sumer et al. [1] and Carreiras et al. [33], Corvaro et al. [21] found the scour still occurred for KC ≈ 4, and proposed the revised equilibrium scour depth predicting equation (Equation (3)) for KC > 4.

Rudolph and Bos [2] conducted a series of wave flume experiments to investigate the scour depth around monopile under waves only, waves and currents combined respectively, indicting KC was one of key parameters in influencing equilibrium scour depth, and proposed the equilibrium scour depth predicting equation (Equation (4)) for low KC (1 < KC < 10). Through analyzing the extensive data from published literatures, Raaijmakers and Rudolph [34] developed the equilibrium scour depth predicting equation (Equation (5)) for low KC, which was suitable for waves only, waves and currents combined. Khalfin [35] carried out several sets of wave flume experiments to study scour development around monopile, and proposed the equilibrium scour depth predicting equation (Equation (6)) for low KC (0.1 < KC < 3.5). Different with above equations, the Khalfin’s equation considers the Shields parameter θ and KC number simultaneously in predicting equilibrium scour depth. The flow reversal occurred under through in one wave period, so sand particles would be carried away from lee-side of pile to upside, resulting in sand particles backfilled into the upstream scour hole [20,29]. Considering the backfilling effects, Zanke et al. [36] proposed the equilibrium scour depth predicting equation (Equation (7)) around pile by theoretical analysis, and the equation is suitable for the whole range of KC number under regular waves and currents combined.

S/D=1.3(1−exp([−m(KC−6)])�/�=1.3(1−exp(−�(��−6))(2)

where, m = 0.03 for linear waves.

S/D=1.3(1−exp([−0.02(KC−4)])�/�=1.3(1−exp(−0.02(��−4))(3)

S/D=1.3γKwaveKhw�/�=1.3��wave�ℎw(4)

where, γ is safety factor, depending on design process, typically γ = 1.5, K_wave is correction factor considering wave action, K_hw is correction factor considering water depth.

S/D=1.5[tanh(hwD)]KwaveKhw�/�=1.5tanh(ℎw�)�wave�ℎw(5)

where, h_w is water depth.

S/D=0.0753(θθcr−−−√−0.5)0.69KC0.68�/�=0.0753(��cr−0.5)0.69��0.68(6)

where, θ is shields parameter, θ_cr is critical shields parameter.

S/D=2.5(1−0.5u/uc)xrelxrel=xeff/(1+xeff)xeff=0.03(1−0.35ucr/u)(KC−6)⎫⎭⎬⎪⎪�/�=2.5(1−0.5�/��)��=��/(1+��)��=0.03(1−0.35�cr/�)(��−6)(7)

where, u is near-bed orbital velocity amplitude, u_c is critical velocity corresponding the onset of sediment motion.

S/D=1.3{1−exp[−0.03(KC2lnn+36)1/2−6]}�/�=1.31−exp−0.03(��2ln�+36)1/2−6(8)

where, n is the 1/n’th highest wave for random waves

For predicting equilibrium scour depth under irregular waves, i.e., random waves, Sumer and Fredsøe [16] found it’s suitable to take Equation (2) to predict equilibrium scour depth around pile under random waves with the root-mean-square (RMS) value of near-bed orbital velocity amplitude U_m and peak wave period T_P to calculate KC. Khalfin [35] recommended the RMS wave height H_rms and peak wave period T_P were used to calculate KC for Equation (6). References [37,38,39,40] developed a series of stochastic theoretical models to predict equilibrium scour depth around pile under random waves, nonlinear random waves plus currents respectively. The stochastic approach thought the 1/n’th highest wave were responsible for scour in vicinity of pile under random waves, and the KC was calculated in Equation (8) with U_m and mean zero-crossing wave period T_z. The results calculated by Equation (8) agree well with experimental values of Sumer and Fredsøe [16] if the 1/10′th highest wave was used. To author’s knowledge, the stochastic approach proposed by Myrhaug and Rue [37] is the only theoretical model to predict equilibrium scour depth around pile under random waves for the whole range of KC number in published documents. Other methods of predicting scour depth under random waves are mainly originated from the equation for regular waves-only, waves and currents combined, which are limited to the large KC number, such as KC > 6 for Equation (2) and KC > 4 for Equation (3) respectively. However, situations with relatively low KC number (KC < 4) often occur in reality, for example, monopile or suction anchor for OWT foundations in ocean environment. Moreover, local scour around OWT foundations under random waves has not yet been investigated fully. Therefore, further study are still needed in the aspect of scour around OWT foundations with low KC number under random waves. Given that, this study presents the scour sediment model around umbrella suction anchor foundation (USAF) under random waves. In this study, a comparison of equilibrium scour depth around USAF between this present numerical models and the previous theoretical models and experimental results was presented firstly. Then, this study gave a comprehensive analysis for the scour mechanisms around USAF. After that, two revised models were proposed according to the model of Raaijmakers and Rudolph [34] and the stochastic model developed by Myrhaug and Rue [37] respectively to predict the equilibrium scour depth. Finally, a parametric study was conducted to study the effects of the Froude number (F_r) and Euler number (E_u) to equilibrium scour depth respectively.

2. Numerical Method

2.1. Governing Equations of Flow

The following equations adopted in present model are already available in Flow 3D software. The authors used these theoretical equations to simulate scour in random waves without modification. The incompressible viscous fluid motion satisfies the Reynolds-averaged Navier-Stokes (RANS) equation, so the present numerical model solves RANS equations:

∂u∂t+1VF(uAx∂u∂x+vAy∂u∂y+wAz∂u∂z)=−1ρf∂p∂x+Gx+fx∂�∂�+1��(��∂�∂�+��∂�∂�+��∂�∂�)=−1�f∂�∂�+��+��(9)

∂v∂t+1VF(uAx∂v∂x+vAy∂v∂y+wAz∂v∂z)=−1ρf∂p∂y+Gy+fy∂�∂�+1��(��∂�∂�+��∂�∂�+��∂�∂�)=−1�f∂�∂�+��+��(10)

∂w∂t+1VF(uAx∂w∂x+vAy∂w∂y+wAz∂w∂z)=−1ρf∂p∂z+Gz+fz∂�∂�+1��(��∂�∂�+��∂�∂�+��∂�∂�)=−1�f∂�∂�+��+��(11)

where, V_F is the volume fraction; u, v, and w are the velocity components in x, y, z direction respectively with Cartesian coordinates; A_i is the area fraction; ρ_f is the fluid density, f_i is the viscous fluid acceleration, G_i is the fluid body acceleration (i = x, y, z).

2.2. Turbulent Model

The turbulence closure is available by the turbulent model, such as one-equation, the one-equation k-ε model, the standard k-ε model, RNG k-ε turbulent model and large eddy simulation (LES) model. The LES model requires very fine mesh grid, so the computational time is large, which hinders the LES model application in engineering. The RNG k-ε model can reduce computational time greatly with high accuracy in the near-wall region. Furthermore, the RNG k-ε model computes the maximum turbulent mixing length dynamically in simulating sediment scour model. Therefore, the RNG k-ε model was adopted to study the scour around anchor under random waves [41,42].

∂kT∂T+1VF(uAx∂kT∂x+vAy∂kT∂y+wAz∂kT∂z)=PT+GT+DiffkT−εkT∂��∂�+1��(��∂��∂�+��∂��∂�+��∂��∂�)=��+��+��−��(12)

∂εT∂T+1VF(uAx∂εT∂x+vAy∂εT∂y+wAz∂εT∂z)=CDIS1εTkT(PT+CDIS3GT)+Diffε−CDIS2ε2TkT∂��∂�+1��(��∂��∂�+��∂��∂�+��∂��∂�)=��1��(��+��3��)+��−��2��2��(13)

where, k_T is specific kinetic energy involved with turbulent velocity, G_T is the turbulent energy generated by buoyancy; ε_T is the turbulent energy dissipating rate, P_T is the turbulent energy, Diff_ε and Diff_kT are diffusion terms associated with V_F, A_i; CDIS1, CDIS2 and CDIS3 are dimensionless parameters, and CDIS1, CDIS3 have default values of 1.42, 0.2 respectively. CDIS2 can be obtained from P_T and k_T.

2.3. Sediment Scour Model

The sand particles may suffer four processes under waves, i.e., entrainment, bed load transport, suspended load transport, and deposition, so the sediment scour model should depict the above processes efficiently. In present numerical simulation, the sediment scour model includes the following aspects:

2.3.1. Entrainment and Deposition

The combination of entrainment and deposition determines the net scour rate of seabed in present sediment scour model. The entrainment lift velocity of sand particles was calculated as [43]:

ulift,i=αinsd0.3∗(θ−θcr)1.5∥g∥di(ρi−ρf)ρf−−−−−−−−−−−−√�lift,i=��*0.3(�−�cr)1.5��(��−�f)�f(14)

where, α_i is the entrainment parameter, n_s is the outward point perpendicular to the seabed, d_* is the dimensionless diameter of sand particles, which was calculated by Equation (15), θ_cr is the critical Shields parameter, g is the gravity acceleration, d_i is the diameter of sand particles, ρ_i is the density of seabed species.

d∗=di(∥g∥ρf(ρi−ρf)μ2f)1/3�*=��(��f(��−�f)�f2)1/3(15)

where μ_f is the fluid dynamic viscosity.

In Equation (14), the entrainment parameter α_i confirms the rate at which sediment erodes when the given shear stress is larger than the critical shear stress, and the recommended value 0.018 was adopted according to the experimental data of Mastbergen and Von den Berg [43]. n_s is the outward pointing normal to the seabed interface, and n_s = (0,0,1) according to the Cartesian coordinates used in present numerical model.

The shields parameter was obtained from the following equation:

θ=U2f,m(ρi/ρf−1)gd50�=�f,m2(��/�f−1)��50(16)

where, U_f,m is the maximum value of the near-bed friction velocity; d₅₀ is the median diameter of sand particles. The detailed calculation procedure of θ was available in Soulsby [44].

The critical shields parameter θ_cr was obtained from the Equation (17) [44]

θcr=0.31+1.2d∗+0.055[1−exp(−0.02d∗)]�cr=0.31+1.2�*+0.0551−exp(−0.02�*)(17)

The sand particles begin to deposit on seabed when the turbulence energy weaken and cann’t support the particles suspending. The setting velocity of the particles was calculated from the following equation [44]:

usettling,i=νfdi[(10.362+1.049d3∗)0.5−10.36]�settling,�=�f��(10.362+1.049�*3)0.5−10.36(18)

where ν_f is the fluid kinematic viscosity.

2.3.2. Bed Load Transport

This is called bed load transport when the sand particles roll or bounce over the seabed and always have contact with seabed. The bed load transport velocity was computed by [45]:

ubedload,i=qb,iδicb,ifb�bedload,�=�b,��b,��b(19)

where, q_b,i is the bed load transport rate, which was obtained from Equation (20), δ_i is the bed load thickness, which was calculated by Equation (21), c_b,i is the volume fraction of sand i in the multiple species, f_b is the critical packing fraction of the seabed.

qb,i=8[∥g∥(ρi−ρfρf)d3i]1/2�b,�=8�(��−�f�f)��31/2(20)

δi=0.3d0.7∗(θθcr−1)0.5di��=0.3�*0.7(��cr−1)0.5��(21)

2.3.3. Suspended Load Transport

Through the following transport equation, the suspended sediment concentration could be acquired.

∂Cs,i∂t+∇(us,iCs,i)=∇∇(DfCs,i)∂�s,�∂�+∇(�s,��s,�)=∇∇(�f�s,�)(22)

where, C_s,i is the suspended sand particles mass concentration of sand i in the multiple species, u_s,i is the sand particles velocity of sand i, D_f is the diffusivity.

The velocity of sand i in the multiple species could be obtained from the following equation:

us,i=u¯¯+usettling,ics,i�s,�=�¯+�settling,��s,�(23)

where, u¯�¯ is the velocity of mixed fluid-particles, which can be calculated by the RANS equation with turbulence model, c_s,i is the suspended sand particles volume concentration, which was computed from Equation (24).

cs,i=Cs,iρi�s,�=�s,��(24)

3. Model Setup

The seabed-USAF-wave three-dimensional scour numerical model was built using Flow-3D software. As shown in Figure 2, the model includes sandy seabed, USAF model, sea water, two baffles and porous media. The dimensions of USAF are shown in Table 1. The sandy bed (210 m in length, 30 m in width and 11 m in height) is made up of uniform fine sand with median diameter d₅₀ = 0.041 cm. The USAF model includes upper steel tube with the length of 20 m, which was installed in the middle of seabed. The location of USAF is positioned at 140 m from the upstream inflow boundary and 70 m from the downstream outflow boundary. Two baffles were installed at two ends of seabed. In order to eliminate the wave reflection basically, the porous media was set at the outflow side on the seabed.

Figure 2. (a) The sketch of seabed-USAF-wave three-dimensional model; (b) boundary condation:Wv-wave boundary, S-symmetric boundary, O-outflow boundary; (c) USAF model.

Table 1. Numerical simulating cases.

3.1. Mesh Geometric Dimensions

In the simulation of the scour under the random waves, the model includes the umbrella suction anchor foundation, seabed and fluid. As shown in Figure 3, the model mesh includes global mesh grid and nested mesh grid, and the total number of grids is 1,812,000. The basic procedure for building mesh grid consists of two steps. Step 1: Divide the global mesh using regular hexahedron with size of 0.6 × 0.6. The global mesh area is cubic box, embracing the seabed and whole fluid volume, and the dimensions are 210 m in length, 30 m in width and 32 m in height. The details of determining the grid size can see the following mesh sensitivity section. Step 2: Set nested fine mesh grid in vicinity of the USAF with size of 0.3 × 0.3 so as to shorten the computation cost and improve the calculation accuracy. The encryption range is −15 m to 15 m in x direction, −15 m to 15 m in y direction and 0 m to 32 m in z direction, respectively. In order to accurately capture the free-surface dynamics, such as the fluid-air interface, the volume of fluid (VOF) method was adopted for tracking the free water surface. One specific algorithm called FAVORTM (Fractional Area/Volume Obstacle Representation) was used to define the fractional face areas and fractional volumes of the cells which are open to fluid flow.

Figure 3. The sketch of mesh grid.

3.2. Boundary Conditions

As shown in Figure 2, the initial fluid length is 210 m as long as seabed. A wave boundary was specified at the upstream offshore end. The details of determining the random wave spectrum can see the following wave parameters section. The outflow boundary was set at the downstream onshore end. The symmetry boundary was used at the top and two sides of the model. The symmetric boundaries were the better strategy to improve the computation efficiency and save the calculation cost [46]. At the seabed bottom, the wall boundary was adopted, which means the u = v = w= 0. Besides, the upper steel tube of USAF was set as no-slip condition.

3.3. Wave Parameters

The random waves with JONSWAP wave spectrum were used for all simulations as realistic representation of offshore conditions. The unidirectional JONSWAP frequency spectrum was described as [47]:

S(ω)=αg2ω5exp[−54(ωpω)4]γexp[−(ω−ωp)22σ2ω2p]�(�)=��2�5exp−54(�p�)4�exp−(�−�p)22�2�p2(25)

where, α is wave energy scale parameter, which is calculated by Equation (26), ω is frequency, ω_p is wave spectrum peak frequency, which can be obtained from Equation (27). γ is wave spectrum peak enhancement factor, in this study γ = 3.3. σ is spectral width factor, σ equals 0.07 for ω ≤ ω_p and 0.09 for ω > ω_p respectively.

α=0.0076(gXU2)−0.22�=0.0076(��2)−0.22(26)

ωp=22(gU)(gXU2)−0.33�p=22(��)(��2)−0.33(27)

where, X is fetch length, U is average wind velocity at 10 m height from mean sea level.

In present numerical model, the input key parameters include X and U for wave boundary with JONSWAP wave spectrum. The objective wave height and period are available by different combinations of X and U. In this study, we designed 9 cases with different wave heights, periods and water depths for simulating scour around USAF under random waves (see Table 2). For random waves, the wave steepness ε and Ursell number U_r were acquired form Equations (28) and (29) respectively

ε=2πgHsT2a�=2��s�a2(28)

Ur=Hsk2h3w�r=�s�2ℎw3(29)

where, H_s is significant wave height, T_a is average wave period, k is wave number, h_w is water depth. The Shield parameter θ satisfies θ _> θ_cr for all simulations in current study, indicating the live bed scour prevails.

Table 2. Numerical simulating cases.

3.4. Mesh Sensitivity

In this section, a mesh sensitivity analysis was conducted to investigate the influence of mesh grid size to results and make sure the calculation is mesh size independent and converged. Three mesh grid size were chosen: Mesh 1—global mesh grid size of 0.75 × 0.75, nested fine mesh grid size of 0.4 × 0.4, and total number of grids 1,724,000, Mesh 2—global mesh grid size of 0.6 × 0.6, nested fine mesh grid size of 0.3 × 0.3, and total number of grids 1,812,000, Mesh 3—global mesh grid size of 0.4 × 0.4, nested fine mesh grid size of 0.2 × 0.2, and total number of grids 1,932,000. The near-bed shear velocity U_* is an important factor for influencing scour process [1,15], so U_* at the position of (4,0,11.12) was evaluated under three mesh sizes. As the Figure 4 shown, the maximum error of shear velocity ∆U_*1,2 is about 39.8% between the mesh 1 and mesh 2, and 4.8% between the mesh 2 and mesh 3. According to the mesh sensitivity criterion adopted by Pang et al. [48], it’s reasonable to think the results are mesh size independent and converged with mesh 2. Additionally, the present model was built according to prototype size, and the mesh size used in present model is larger than the mesh size adopted by Higueira et al. [49] and Corvaro et al. [50]. If we choose the smallest cell size, it will take too much time. For example, the simulation with Mesh3 required about 260 h by using a computer with Intel Xeon Scalable Gold 4214 CPU @24 Cores, 2.2 GHz and 64.00 GB RAM. Therefore, in this case, considering calculation accuracy and computation efficiency, the mesh 2 was chosen for all the simulation in this study.

Figure 4. Comparison of near-bed shear velocity U_* with different mesh grid size.

The nested mesh block was adopted for seabed in vicinity of the USAF, which was overlapped with the global mesh block. When two mesh blocks overlap each other, the governing equations are by default solved on the mesh block with smaller average cell size (i.e., higher grid resolution). It is should be noted that the Flow 3D software used the moving mesh captures the scour evolution and automatically adjusts the time step size to be as large as possible without exceeding any of the stability limits, affecting accuracy, or unduly increasing the effort required to enforce the continuity condition [51].

3.5. Model Validation

In order to verify the reliability of the present model, the results of present study were compared with the experimental data of Khosronejad et al. [52]. The experiment was conducted in an open channel with a slender vertical pile under unidirectional currents. The comparison of scour development between the present results and the experimental results is shown in Figure 5. The Figure 5 reveals that the present results agree well with the experimental data of Khosronejad et al. [52]. In the first stage, the scour depth increases rapidly. After that, the scour depth achieves a maximum value gradually. The equilibrium scour depth calculated by the present model is basically corresponding with the experimental results of Khosronejad et al. [52], although scour depth in the present model is slightly larger than the experimental results at initial stage.

Figure 5. Comparison of time evolution of scour between the present study and Khosronejad et al. [52], Petersen et al. [17].

Secondly, another comparison was further conducted between the results of present study and the experimental data of Petersen et al. [17]. The experiment was carried out in a flume with a circular vertical pile in combined waves and current. Figure 4 shows a comparison of time evolution of scour depth between the simulating and the experimental results. As Figure 5 indicates, the scour depth in this study has good overall agreement with the experimental results proposed in Petersen et al. [17]. The equilibrium scour depth calculated by the present model is 0.399 m, which equals to the experimental value basically. Overall, the above verifications prove the present model is accurate and capable in dealing with sediment scour under waves.

In addition, in order to calibrate and validate the present model for hydrodynamic parameters, the comparison of water surface elevation was carried out with laboratory experiments conducted by Stahlmann [53] for wave gauge No. 3. The Figure 6 depicts the surface wave profiles between experiments and numerical model results. The comparison indicates that there is a good agreement between the model results and experimental values, especially the locations of wave crest and trough. Comparison of the surface elevation instructs the present model has an acceptable relative error, and the model is a calibrated in terms of the hydrodynamic parameters.

Figure 6. Comparison of surface elevation between the present study and Stahlmann [53].

Finally, another comparison was conducted for equilibrium scour depth or maximum scour depth under random waves with the experimental data of Sumer and Fredsøe [16] and Schendel et al. [22]. The Figure 7 shows the comparison between the numerical results and experimental data of Run01, Run05, Run21 and Run22 in Sumer and Fredsøe [16] and test A05 and A09 in Schendel et al. [22]. As shown in Figure 7, the equilibrium scour depth or maximum scour depth distributed within the ±30 error lines basically, meaning the reliability and accuracy of present model for predicting equilibrium scour depth around foundation in random waves. However, compared with the experimental values, the present model overestimated the equilibrium scour depth generally. Given that, a calibration for scour depth was carried out by multiplying the mean reduced coefficient 0.85 in following section.

Figure 7. Comparison of equilibrium (or maximum) scour depth between the present study and Sumer and Fredsøe [16], Schendel et al. [22].

Through the various examination for hydrodynamic and morphology parameters, it can be concluded that the present model is a validated and calibrated model for scour under random waves. Thus, the present numerical model would be utilized for scour simulation around foundation under random waves.

4. Numerical Results and Discussions

4.1. Scour Evolution

Figure 8 displays the scour evolution for case 1–9. As shown in Figure 8a, the scour depth increased rapidly at the initial stage, and then slowed down at the transition stage, which attributes to the backfilling occurred in scour holes under live bed scour condition, resulting in the net scour decreasing. Finally, the scour reached the equilibrium state when the amount of sediment backfilling equaled to that of scouring in the scour holes, i.e., the net scour transport rate was nil. Sumer and Fredsøe [16] proposed the following formula for the scour development under waves

St=Seq(1−exp(−t/Tc))�t=�eq(1−exp(−�/�c))(30)

where T_c is time scale of scour process.

Figure 8. Time evolution of scour for case 1–9: (a) Case 1–5; (b) Case 6–9.

The computing time is 3600 s and the scour development curves in Figure 8 kept fluctuating, meaning it’s still not in equilibrium scour stage in these cases. According to Sumer and Fredsøe [16], the equilibrium scour depth can be acquired by fitting the data with Equation (30). From Figure 8, it can be seen that the scour evolution obtained from Equation (30) is consistent with the present study basically at initial stage, but the scour depth predicted by Equation (30) developed slightly faster than the simulating results and the Equation (30) overestimated the scour depth to some extent. Overall, the whole tendency of the results calculated by Equation (30) agrees well with the simulating results of the present study, which means the Equation (30) is applicable to depict the scour evolution around USAF under random waves.

4.2. Scour Mechanism under Random Waves

The scour morphology and scour evolution around USAF are similar under random waves in case 1~9. Taking case 7 as an example, the scour morphology is shown in Figure 9.

Figure 9. Scour morphology under different times for case 7.

From Figure 9, at the initial stage (t < 1200 s), the scour occurred at upstream foundation edges between neighboring anchor branches. The maximum scour depth appeared at the lee-side of the USAF. Correspondingly, the sediments deposited at the periphery of the USAF, and the location of the maximum accretion depth was positioned at an angle of about 45° symmetrically with respect to the wave propagating direction in the lee-side of the USAF. After that, when t > 2400 s, the location of the maximum scour depth shifted to the upside of the USAF at an angle of about 45° with respect to the wave propagating direction.

According to previous studies [1,15,16,19,30,31], the horseshoe vortex, streamline compression and wake vortex shedding were responsible for scour around foundation. The Figure 10 displays the distribution of flow velocity in vicinity of foundation, which reflects the evolving processes of horseshoe vertex.

Figure 10. Velocity profile around USAF: (a) Flow runup and down stream at upstream anchor edges; (b) Horseshoe vortex at upstream anchor edges; (c) Flow reversal during wave through stage at lee side.

As shown in Figure 10, the inflow tripped to the upstream edges of the USAF and it was blocked by the upper tube of USAF. Then, the downflow formed the horizontal axis clockwise vortex and rolled on the seabed bypassing the tube, that is, the horseshoe vortex (Figure 11). The Figure 12 displays the turbulence intensity around the tube on the seabed. From Figure 12, it can be seen that the turbulence intensity was high-intensity with respect to the region of horseshoe vortex. This phenomenon occurred because of drastic water flow momentum exchanging in the horseshoe vortex. As a result, it created the prominent shear stress on the seabed, causing the local scour at the upstream edges of USAF. Besides, the horseshoe vortex moved downstream gradually along the periphery of the tube and the wake vortex shed off continually at the lee-side of the USAF, i.e., wake vortex.

Figure 11. Sketch of scour mechanism around USAF under random waves.

Figure 12. Turbulence intensity: (a) Turbulence intensity of horseshoe vortex; (b) Turbulence intensity of wake vortex; (c) Turbulence intensity of accretion area.

The core of wake vortex is a negative pressure center, liking a vacuum cleaner [11,42]. Hence, the soil particles were swirled into the negative pressure core and carried away by wake vortex. At the same time, the onset of scour at rear side occurred. Finally, the wake vortex became downflow at the downside of USAF. As is shown in Figure 12, the turbulence intensity was low where the downflow occurred at lee-side, which means the turbulence energy may not be able to support the survival of wake vortex, leading to accretion happening. As mentioned in previous section, the formation of horseshoe vortex was dependent with adverse pressure gradient at upside of foundation. As shown in Figure 13, the evaluated range of pressure distribution is −15 m to 15 m in x direction. The t = 450 s and t = 1800 s indicate that the wave crest and trough arrived at the upside and lee-side of the foundation respectively, and the t = 350 s was neither the wave crest nor trough. The adverse gradient pressure reached the maximum value at t = 450 s corresponding to the wave crest phase. In this case, it’s helpful for the wave boundary separating fully from seabed, which leads to the formation of horseshoe vortex with high turbulence intensity. Therefore, the horseshoe vortex is responsible for the local scour between neighboring anchor branches at upside of USAF. What’s more, due to the combination of the horseshoe vortex and streamline compression, the maximum scour depth occurred at the upside of the USAF with an angle of about 45° corresponding to the wave propagating direction. This is consistent with the findings of Pang et al. [48] and Sumer et al. [1,15] in case of regular waves. At the wave trough phase (t = 1800 s), the pressure gradient became positive at upstream USAF edges, which hindered the separating of wave boundary from seabed. In the meantime, the flow reversal occurred (Figure 10) and the adverse gradient pressure appeared at downstream USAF edges, but the magnitude of adverse gradient pressure at lee-side was lower than the upstream gradient pressure under wave crest. In this way, the intensity of horseshoe vortex behind the USAF under wave trough was low, which explains the difference of scour depth at upstream and downstream, i.e., the scour asymmetry. In other words, the scour asymmetry at upside and downside of USAF was attributed to wave asymmetry for random waves, and the phenomenon became more evident for nonlinear waves [21]. Briefly speaking, the vortex system at wave crest phase was mainly related to the scour process around USAF under random waves.

Figure 13. Pressure distribution around USAF.

4.3. Equilibrium Scour Depth

The KC number is a key parameter for horseshoe vortex emerging and evolving under waves. According to Equation (1), when pile diameter D is fixed, the KC depends on the maximum near-bed velocity U_wm and wave period T. For random waves, the U_wm can be denoted by the root-mean-square (RMS) value of near-bed velocity amplitude U_wm,rms or the significant value of near-bed velocity amplitude U_wm,s. The U_wm,rms and U_wm,s for all simulating cases of the present study are listed in Table 3 and Table 4. The T can be denoted by the mean up zero-crossing wave period T_a, peak wave period T_p, significant wave period T_s, the maximum wave period T_m, 1/10′th highest wave period T_{n = 1/10} and 1/5′th highest wave period T_{n = 1/5} for random waves, so the different combinations of U_wm and T will acquire different KC. The Table 3 and Table 4 list 12 types of KC, for example, the KC_rms,s was calculated by U_wm,rms and T_s. Sumer and Fredsøe [16] conducted a series of wave flume experiments to investigate the scour depth around monopile under random waves, and found the equilibrium scour depth predicting equation (Equation (2)) for regular waves was applicable for random waves with KC_rms,p. It should be noted that the Equation (2) is only suitable for KC > 6 under regular waves or KC_rms,p > 6 under random waves.

Table 3. U_wm,rms and KC for case 1~9.

Table 4. U_wm,s and KC for case 1~9.

Raaijmakers and Rudolph [34] proposed the equilibrium scour depth predicting model (Equation (5)) around pile under waves, which is suitable for low KC. The format of Equation (5) is similar with the formula proposed by Breusers [54], which can predict the equilibrium scour depth around pile at different scour stages. In order to verify the applicability of Raaijmakers’s model for predicting the equilibrium scour depth around USAF under random waves, a validation of the equilibrium scour depth S_eq between the present study and Raaijmakers’s equation was conducted. The position where the scour depth S_eq was evaluated is the location of the maximum scour depth, and it was depicted in Figure 14. The Figure 15 displays the comparison of S_eq with different KC between the present study and Raaijmakers’s model.

Figure 14. Sketch of the position where the S_eq was evaluated.

Figure 15. Comparison of the equilibrium scour depth between the present model and the model of Raaijmakers and Rudolph [34]: (a) KC_rms,s, KC_rms,a; (b) KC_rms,p, KC_rms,m; (c) KC_{rms,n = 1/10}, KC_{rms,n = 1/5}; (d) KC_s,s, KC_s,a; (e) KC_s,p, KC_s,m; (f) KC_{s,n = 1/10}, KC_{s,n = 1/5}.

As shown in Figure 15, there is an error in predicting S_eq between the present study and Raaijmakers’s model, and Raaijmakers’s model underestimates the results generally. Although the error exists, the varying trend of S_eq with KC obtained from Raaijmakers’s model is consistent with the present study basically. What’s more, the error is minimum and the Raaijmakers’s model is of relatively high accuracy for predicting scour around USAF under random waves by using KC_s,p. Based on this, a further revision was made to eliminate the error as much as possible, i.e., add the deviation value ∆S/D in the Raaijmakers’s model. The revised equilibrium scour depth predicting equation based on Raaijmakers’s model can be written as

S′eq/D=1.95[tanh(hD)](1−exp(−0.012KCs,p))+ΔS/D�eq′/�=1.95tanh(ℎ�)(1−exp(−0.012��s,p))+∆�/�(31)

As the Figure 16 shown, through trial-calculation, when ∆S/D = 0.05, the results calculated by Equation (31) show good agreement with the simulating results of the present study. The maximum error is about 18.2% and the engineering requirements have been met basically. In order to further verify the accuracy of the revised model for large KC (KC_s,p > 4) under random waves, a validation between the revised model and the previous experimental results [21]. The experiment was conducted in a flume (50 m in length, 1.0 m in width and 1.3 m in height) with a slender vertical pile (D = 0.1 m) under random waves. The seabed is composed of 0.13 m deep layer of sand with d₅₀ = 0.6 mm and the water depth is 0.5 m for all tests. The significant wave height is 0.12~0.21 m and the KC_s,p is 5.52~11.38. The comparison between the predicting results by Equation (31) and the experimental results of Corvaro et al. [21] is shown in Figure 17. From Figure 17, the experimental data evenly distributes around the predicted results and the prediction accuracy is favorable when KC_s,p < 8. However, the gap between the predicting results and experimental data becomes large and the Equation (31) overestimates the equilibrium scour depth to some extent when KC_s,p > 8.

Figure 16. Comparison of S_eq between the simulating results and the predicting values by Equation (31).

Figure 17. Comparison of S_eq/D between the Experimental results of Corvaro et al. [21] and the predicting values by Equation (31).

In ocean environment, the waves are composed of a train of sinusoidal waves with different frequencies and amplitudes. The energy of constituent waves with very large and very small frequencies is relatively low, and the energy of waves is mainly concentrated in a certain range of moderate frequencies. Myrhaug and Rue [37] thought the 1/n’th highest wave was responsible for scour and proposed the stochastic model to predict the equilibrium scour depth around pile under random waves for full range of KC. Noteworthy is that the KC was denoted by KC_rms,a in the stochastic model. To verify the application of the stochastic model for predicting scour depth around USAF, a validation between the simulating results of present study and predicting results by the stochastic model with n = 2,3,5,10,20,500 was carried out respectively.

As shown in Figure 18, compared with the simulating results, the stochastic model underestimates the equilibrium scour depth around USAF generally. Although the error exists, the varying trend of S_eq with KC_rms,a obtained from the stochastic model is consistent with the present study basically. What’s more, the gap between the predicting values by stochastic model and the simulating results decreases with the increase of n, but for large n, for example n = 500, the varying trend diverges between the predicting values and simulating results, meaning it’s not feasible only by increasing n in stochastic model to predict the equilibrium scour depth around USAF.

Figure 18. Comparison of S_eq between the simulating results and the predicting values by Equation (8).

The Figure 19 lists the deviation value ∆S_eq/D′ between the predicting values and simulating results with different KC_rms,a and n. Then, fitted the relationship between the ∆S′and n under different KC_rms,a, and the fitting curve can be written by Equation (32). The revised stochastic model (Equation (33)) can be acquired by adding ∆S_eq/D′ to Equation (8).

ΔSeq/D=0.052*exp(−n/6.566)+0.068∆�eq/�=0.052*exp(−�/6.566)+0.068(32)

S′eq¯/D=S′eq/D+0.052*exp(−n/6.566)+0.068�eq′¯/�=�eq′/�+0.052*exp(−�/6.566)+0.068(33)

Figure 19. The fitting line between ∆S′and n.

The comparison between the predicting results by Equation (33) and the simulating results of present study is shown in Figure 20. According to the Figure 20, the varying trend of S_eq with KC_rms,a obtained from the stochastic model is consistent with the present study basically. Compared with predicting results by the stochastic model, the results calculated by Equation (33) is favorable. Moreover, comparison with simulating results indicates that the predicting results are the most favorable for n = 10, which is consistent with the findings of Myrhaug and Rue [37] for equilibrium scour depth predicting around slender pile in case of random waves.

Figure 20. Comparison of S_eq between the simulating results and the predicting values by Equation (33).

In order to further verify the accuracy of the Equation (33) for large KC (KC_rms,a > 4) under random waves, a validation was conducted between the Equation (33) and the previous experimental results of Sumer and Fredsøe [16] and Corvaro et al. [21]. The details of experiments conducted by Corvaro et al. [21] were described in above section. Sumer and Fredsøe [16] investigated the local scour around pile under random waves. The experiments were conducted in a wave basin with a slender vertical pile (D = 0.032, 0.055 m). The seabed is composed of 0.14 m deep layer of sand with d₅₀ = 0.2 mm and the water depth was maintained at 0.5 m. The JONSWAP wave spectrum was used and the KC_rms,a was 5.29~16.95. The comparison between the predicting results by Equation (33) and the experimental results of Sumer and Fredsøe [16] and Corvaro et al. [21] are shown in Figure 21. From Figure 21, contrary to the case of low KC_rms,a (KC_rms,a < 4), the error between the predicting values and experimental results increases with decreasing of n for KC_rms,a > 4. Therefore, the predicting results are the most favorable for n = 2 when KC_rms,a > 4.

Figure 21. Comparison of S_eq between the experimental results of Sumer and Fredsøe [16] and Corvaro et al. [21] and the predicting values by Equation (33).

Noteworthy is that the present model was built according to prototype size, so the errors between the numerical results and experimental data of References [16,21] may be attribute to the scale effects. In laboratory experiments on scouring process, it is typically impossible to ensure a rigorous similarity of all physical parameters between the model and prototype structure, leading to the scale effects in the laboratory experiments. To avoid a cohesive behaviour, the bed material was not scaled geometrically according to model scale. As a consequence, the relatively large-scaled sediments sizes may result in the overestimation of bed load transport and underestimation of suspended load transport compared with field conditions. What’s more, the disproportional scaled sediment presumably lead to the difference of bed roughness between the model and prototype, and thus large influences for wave boundary layer on the seabed and scour process. Besides, according to Corvaro et al. [21] and Schendel et al. [55], the pile Reynolds numbers and Froude numbers both affect the scour depth for the condition of non fully developed turbulent flow in laboratory experiments.

4.4. Parametric Study

4.4.1. Influence of Froude Number

As described above, the set of foundation leads to the adverse pressure gradient appearing at upstream, leading to the wave boundary layer separating from seabed, then horseshoe vortex formatting and the horseshoe vortex are mainly responsible for scour around foundation (see Figure 22). The Froude number F_r is the key parameter to influence the scale and intensity of horseshoe vortex. The F_r under waves can be calculated by the following formula [42]

Fr=UwgD−−−√�r=�w��(34)

where U_w is the mean water particle velocity during 1/4 cycle of wave oscillation, obtained from the following formula. Noteworthy is that the root-mean-square (RMS) value of near-bed velocity amplitude U_wm,rms is used for calculating U_wm.

Uw=1T/4∫0T/4Uwmsin(t/T)dt=2πUwm�w=1�/4∫0�/4�wmsin(�/�)��=2��wm(35)

Figure 22. Sketch of flow field at upstream USAF edges.

Tavouktsoglou et al. [25] proposed the following formula between F_r and the vertical location of the stagnation y

yh∝Fer�ℎ∝�r�(36)

where e is constant.

The Figure 23 displays the relationship between S_eq/D and F_r of the present study. In order to compare with the simulating results, the experimental data of Corvaro et al. [21] was also depicted in Figure 23. As shown in Figure 23, the equilibrium scour depth appears a logarithmic increase as F_r increases and approaches the mathematical asymptotic value, which is also consistent with the experimental results of Corvaro et al. [21]. According to Figure 24, the adverse pressure gradient pressure at upstream USAF edges increases with the increase of F_r, which is benefit for the wave boundary layer separating from seabed, resulting in the high-intensity horseshoe vortex, hence, causing intensive scour around USAF. Based on the previous study of Tavouktsoglou et al. [25] for scour around pile under currents, the high F_r leads to the stagnation point is closer to the mean sea level for shallow water, causing the stronger downflow kinetic energy. As mentioned in previous section, the energy of downflow at upstream makes up the energy of the subsequent horseshoe vortex, so the stronger downflow kinetic energy results in the more intensive horseshoe vortex. Therefore, the higher F_r leads to the more intensive horseshoe vortex by influencing the position of stagnation point y presumably. Qi and Gao [19] carried out a series of flume tests to investigate the scour around pile under regular waves, and proposed the fitting formula between S_eq/D and F_r as following

lg(Seq/D)=Aexp(B/Fr)+Clg(�eq/�)=�exp(�/�r)+�(37)

where A, B and C are constant.

Figure 23. The fitting curve between S_eq/D and F_r.

Figure 24. Sketch of adverse pressure gradient at upstream USAF edges.

Took the Equation (37) to fit the simulating results with A = −0.002, B = 0.686 and C = −0.808, and the results are shown in Figure 23. From Figure 23, the simulating results evenly distribute around the Equation (37) and the varying trend of S_eq/D and F_r in present study is consistent with Equation (37) basically, meaning the Equation (37) is applicable to express the relationship of S_eq/D with F_r around USAF under random waves.

4.4.2. Influence of Euler Number

The Euler number E_u is the influencing factor for the hydrodynamic field around foundation. The E_u under waves can be calculated by the following formula. The E_u can be represented by the Equation (38) for uniform cylinders [25]. The root-mean-square (RMS) value of near-bed velocity amplitude U_m,rms is used for calculating U_m.

Eu=U2mgD�u=�m2��(38)

where U_m is depth-averaged flow velocity.

The Figure 25 displays the relationship between S_eq/D and E_u of the present study. In order to compare with the simulating results, the experimental data of Sumer and Fredsøe [16] and Corvaro et al. [21] were also plotted in Figure 25. As shown in Figure 25, similar with the varying trend of S_eq/D and F_r, the equilibrium scour depth appears a logarithmic increase as E_u increases and approaches the mathematical asymptotic value, which is also consistent with the experimental results of Sumer and Fredsøe [16] and Corvaro et al. [21]. According to Figure 24, the adverse pressure gradient pressure at upstream USAF edges increases with the increasing of E_u, which is benefit for the wave boundary layer separating from seabed, inducing the high-intensity horseshoe vortex, hence, causing intensive scour around USAF.

Figure 25. The fitting curve between S_eq/D and E_u.

Therefore, the variation of F_r and E_u reflect the magnitude of adverse pressure gradient pressure at upstream. Given that, the Equation (37) also was used to fit the simulating results with A = 8.875, B = 0.078 and C = −9.601, and the results are shown in Figure 25. From Figure 25, the simulating results evenly distribute around the Equation (37) and the varying trend of S_eq/D and E_u in present study is consistent with Equation (37) basically, meaning the Equation (37) is also applicable to express the relationship of S_eq/D with E_u around USAF under random waves. Additionally, according to the above description of F_r, it can be inferred that the higher F_r and E_u both lead to the more intensive horseshoe vortex by influencing the position of stagnation point y presumably.

5. Conclusions

A series of numerical models were established to investigate the local scour around umbrella suction anchor foundation (USAF) under random waves. The numerical model was validated for hydrodynamic and morphology parameters by comparing with the experimental data of Khosronejad et al. [52], Petersen et al. [17], Sumer and Fredsøe [16] and Schendel et al. [22]. Based on the simulating results, the scour evolution and scour mechanisms around USAF under random waves were analyzed respectively. Two revised models were proposed according to the model of Raaijmakers and Rudolph [34] and the stochastic model developed by Myrhaug and Rue [37] to predict the equilibrium scour depth around USAF under random waves. Finally, a parametric study was carried out with the present model to study the effects of the Froude number F_r and Euler number E_u to the equilibrium scour depth around USAF under random waves. The main conclusions can be described as follows.(1)

The packed sediment scour model and the RNG k−ε turbulence model were used to simulate the sand particles transport processes and the flow field around UASF respectively. The scour evolution obtained by the present model agrees well with the experimental results of Khosronejad et al. [52], Petersen et al. [17], Sumer and Fredsøe [16] and Schendel et al. [22], which indicates that the present model is accurate and reasonable for depicting the scour morphology around UASF under random waves.(2)

The vortex system at wave crest phase is mainly related to the scour process around USAF under random waves. The maximum scour depth appeared at the lee-side of the USAF at the initial stage (t < 1200 s). Subsequently, when t > 2400 s, the location of the maximum scour depth shifted to the upside of the USAF at an angle of about 45° with respect to the wave propagating direction.(3)

The error is negligible and the Raaijmakers’s model is of relatively high accuracy for predicting scour around USAF under random waves when KC is calculated by KC_s,p. Given that, a further revision model (Equation (31)) was proposed according to Raaijmakers’s model to predict the equilibrium scour depth around USAF under random waves and it shows good agreement with the simulating results of the present study when KC_s,p < 8.(4)

Another further revision model (Equation (33)) was proposed according to the stochastic model established by Myrhaug and Rue [37] to predict the equilibrium scour depth around USAF under random waves, and the predicting results are the most favorable for n = 10 when KC_rms,a < 4. However, contrary to the case of low KC_rms,a, the predicting results are the most favorable for n = 2 when KC_rms,a > 4 by the comparison with experimental results of Sumer and Fredsøe [16] and Corvaro et al. [21].(5)

The same formula (Equation (37)) is applicable to express the relationship of S_eq/D with E_u or F_r, and it can be inferred that the higher F_r and E_u both lead to the more intensive horseshoe vortex and larger S_eq.

Author Contributions

Conceptualization, H.L. (Hongjun Liu); Data curation, R.H. and P.Y.; Formal analysis, X.W. and H.L. (Hao Leng); Funding acquisition, X.W.; Writing—original draft, R.H. and P.Y.; Writing—review & editing, X.W. and H.L. (Hao Leng); The final manuscript has been approved by all the authors. All authors have read and agreed to the published version of the manuscript.

Funding

This research was funded by the Fundamental Research Funds for the Central Universities (grant number 202061027) and the National Natural Science Foundation of China (grant number 41572247).

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

The data presented in this study are available on request from the corresponding author.

Conflicts of Interest

The authors declare no conflict of interest.

References

Sumer, B.M.; Fredsøe, J.; Christiansen, N. Scour Around Vertical Pile in Waves. J. Waterw. Port. Coast. Ocean Eng. 1992, 118, 15–31. [Google Scholar] [CrossRef]
Rudolph, D.; Bos, K. Scour around a monopile under combined wave-current conditions and low KC-numbers. In Proceedings of the 6th International Conference on Scour and Erosion, Amsterdam, The Netherlands, 1–3 November 2006; pp. 582–588. [Google Scholar]
Nielsen, A.W.; Liu, X.; Sumer, B.M.; Fredsøe, J. Flow and bed shear stresses in scour protections around a pile in a current. Coast. Eng. 2013, 72, 20–38. [Google Scholar] [CrossRef]
Ahmad, N.; Bihs, H.; Myrhaug, D.; Kamath, A.; Arntsen, Ø.A. Three-dimensional numerical modelling of wave-induced scour around piles in a side-by-side arrangement. Coast. Eng. 2018, 138, 132–151. [Google Scholar] [CrossRef]
Li, H.; Ong, M.C.; Leira, B.J.; Myrhaug, D. Effects of Soil Profile Variation and Scour on Structural Response of an Offshore Monopile Wind Turbine. J. Offshore Mech. Arct. Eng. 2018, 140, 042001. [Google Scholar] [CrossRef]
Li, H.; Liu, H.; Liu, S. Dynamic analysis of umbrella suction anchor foundation embedded in seabed for offshore wind turbines. Géoméch. Energy Environ. 2017, 10, 12–20. [Google Scholar] [CrossRef]
Fazeres-Ferradosa, T.; Rosa-Santos, P.; Taveira-Pinto, F.; Vanem, E.; Carvalho, H.; Correia, J.A.F.D.O. Editorial: Advanced research on offshore structures and foundation design: Part 1. Proc. Inst. Civ. Eng. Marit. Eng. 2019, 172, 118–123. [Google Scholar] [CrossRef]
Chavez, C.E.A.; Stratigaki, V.; Wu, M.; Troch, P.; Schendel, A.; Welzel, M.; Villanueva, R.; Schlurmann, T.; De Vos, L.; Kisacik, D.; et al. Large-Scale Experiments to Improve Monopile Scour Protection Design Adapted to Climate Change—The PROTEUS Project. Energies 2019, 12, 1709. [Google Scholar] [CrossRef][Green Version]
Wu, M.; De Vos, L.; Chavez, C.E.A.; Stratigaki, V.; Fazeres-Ferradosa, T.; Rosa-Santos, P.; Taveira-Pinto, F.; Troch, P. Large Scale Experimental Study of the Scour Protection Damage Around a Monopile Foundation Under Combined Wave and Current Conditions. J. Mar. Sci. Eng. 2020, 8, 417. [Google Scholar] [CrossRef]
Sørensen, S.P.H.; Ibsen, L.B. Assessment of foundation design for offshore monopiles unprotected against scour. Ocean Eng. 2013, 63, 17–25. [Google Scholar] [CrossRef]
Prendergast, L.; Gavin, K.; Doherty, P. An investigation into the effect of scour on the natural frequency of an offshore wind turbine. Ocean Eng. 2015, 101, 1–11. [Google Scholar] [CrossRef][Green Version]
Fazeres-Ferradosa, T.; Chambel, J.; Taveira-Pinto, F.; Rosa-Santos, P.; Taveira-Pinto, F.; Giannini, G.; Haerens, P. Scour Protections for Offshore Foundations of Marine Energy Harvesting Technologies: A Review. J. Mar. Sci. Eng. 2021, 9, 297. [Google Scholar] [CrossRef]
Yang, Q.; Yu, P.; Liu, Y.; Liu, H.; Zhang, P.; Wang, Q. Scour characteristics of an offshore umbrella suction anchor foundation under the combined actions of waves and currents. Ocean Eng. 2020, 202, 106701. [Google Scholar] [CrossRef]
Yu, P.; Hu, R.; Yang, J.; Liu, H. Numerical investigation of local scour around USAF with different hydraulic conditions under currents and waves. Ocean Eng. 2020, 213, 107696. [Google Scholar] [CrossRef]
Sumer, B.M.; Christiansen, N.; Fredsøe, J. The horseshoe vortex and vortex shedding around a vertical wall-mounted cylinder exposed to waves. J. Fluid Mech. 1997, 332, 41–70. [Google Scholar] [CrossRef]
Sumer, B.M.; Fredsøe, J. Scour around Pile in Combined Waves and Current. J. Hydraul. Eng. 2001, 127, 403–411. [Google Scholar] [CrossRef]
Petersen, T.U.; Sumer, B.M.; Fredsøe, J. Time scale of scour around a pile in combined waves and current. In Proceedings of the 6th International Conference on Scour and Erosion, Paris, France, 27–31 August 2012. [Google Scholar]
Petersen, T.U.; Sumer, B.M.; Fredsøe, J.; Raaijmakers, T.C.; Schouten, J.-J. Edge scour at scour protections around piles in the marine environment—Laboratory and field investigation. Coast. Eng. 2015, 106, 42–72. [Google Scholar] [CrossRef]
Qi, W.; Gao, F. Equilibrium scour depth at offshore monopile foundation in combined waves and current. Sci. China Ser. E Technol. Sci. 2014, 57, 1030–1039. [Google Scholar] [CrossRef][Green Version]
Larsen, B.E.; Fuhrman, D.R.; Baykal, C.; Sumer, B.M. Tsunami-induced scour around monopile foundations. Coast. Eng. 2017, 129, 36–49. [Google Scholar] [CrossRef][Green Version]
Corvaro, S.; Marini, F.; Mancinelli, A.; Lorenzoni, C.; Brocchini, M. Hydro- and Morpho-dynamics Induced by a Vertical Slender Pile under Regular and Random Waves. J. Waterw. Port. Coast. Ocean Eng. 2018, 144, 04018018. [Google Scholar] [CrossRef]
Schendel, A.; Welzel, M.; Schlurmann, T.; Hsu, T.-W. Scour around a monopile induced by directionally spread irregular waves in combination with oblique currents. Coast. Eng. 2020, 161, 103751. [Google Scholar] [CrossRef]
Fazeres-Ferradosa, T.; Taveira-Pinto, F.; Romão, X.; Reis, M.; das Neves, L. Reliability assessment of offshore dynamic scour protections using copulas. Wind. Eng. 2018, 43, 506–538. [Google Scholar] [CrossRef]
Fazeres-Ferradosa, T.; Welzel, M.; Schendel, A.; Baelus, L.; Santos, P.R.; Pinto, F.T. Extended characterization of damage in rubble mound scour protections. Coast. Eng. 2020, 158, 103671. [Google Scholar] [CrossRef]
Tavouktsoglou, N.S.; Harris, J.M.; Simons, R.R.; Whitehouse, R.J.S. Equilibrium Scour-Depth Prediction around Cylindrical Structures. J. Waterw. Port. Coast. Ocean Eng. 2017, 143, 04017017. [Google Scholar] [CrossRef][Green Version]
Ettema, R.; Melville, B.; Barkdoll, B. Scale Effect in Pier-Scour Experiments. J. Hydraul. Eng. 1998, 124, 639–642. [Google Scholar] [CrossRef]
Umeda, S. Scour Regime and Scour Depth around a Pile in Waves. J. Coast. Res. Spec. Issue 2011, 64, 845–849. [Google Scholar]
Umeda, S. Scour process around monopiles during various phases of sea storms. J. Coast. Res. 2013, 165, 1599–1604. [Google Scholar] [CrossRef]
Baykal, C.; Sumer, B.; Fuhrman, D.R.; Jacobsen, N.; Fredsøe, J. Numerical simulation of scour and backfilling processes around a circular pile in waves. Coast. Eng. 2017, 122, 87–107. [Google Scholar] [CrossRef][Green Version]
Miles, J.; Martin, T.; Goddard, L. Current and wave effects around windfarm monopile foundations. Coast. Eng. 2017, 121, 167–178. [Google Scholar] [CrossRef][Green Version]
Miozzi, M.; Corvaro, S.; Pereira, F.A.; Brocchini, M. Wave-induced morphodynamics and sediment transport around a slender vertical cylinder. Adv. Water Resour. 2019, 129, 263–280. [Google Scholar] [CrossRef]
Yu, T.; Zhang, Y.; Zhang, S.; Shi, Z.; Chen, X.; Xu, Y.; Tang, Y. Experimental study on scour around a composite bucket foundation due to waves and current. Ocean Eng. 2019, 189, 106302. [Google Scholar] [CrossRef]
Carreiras, J.; Larroudé, P.; Seabra-Santos, F.; Mory, M. Wave Scour Around Piles. In Proceedings of the Coastal Engineering 2000, American Society of Civil Engineers (ASCE), Sydney, Australia, 16–21 July 2000; pp. 1860–1870. [Google Scholar]
Raaijmakers, T.; Rudolph, D. Time-dependent scour development under combined current and waves conditions—Laboratory experiments with online monitoring technique. In Proceedings of the 4th International Conference on Scour and Erosion, Tokyo, Japan, 5–7 November 2008; pp. 152–161. [Google Scholar]
Khalfin, I.S. Modeling and calculation of bed score around large-diameter vertical cylinder under wave action. Water Resour. 2007, 34, 357. [Google Scholar] [CrossRef][Green Version]
Zanke, U.C.; Hsu, T.-W.; Roland, A.; Link, O.; Diab, R. Equilibrium scour depths around piles in noncohesive sediments under currents and waves. Coast. Eng. 2011, 58, 986–991. [Google Scholar] [CrossRef]
Myrhaug, D.; Rue, H. Scour below pipelines and around vertical piles in random waves. Coast. Eng. 2003, 48, 227–242. [Google Scholar] [CrossRef]
Myrhaug, D.; Ong, M.C.; Føien, H.; Gjengedal, C.; Leira, B.J. Scour below pipelines and around vertical piles due to second-order random waves plus a current. Ocean Eng. 2009, 36, 605–616. [Google Scholar] [CrossRef]
Myrhaug, D.; Ong, M.C. Random wave-induced onshore scour characteristics around submerged breakwaters using a stochastic method. Ocean Eng. 2010, 37, 1233–1238. [Google Scholar] [CrossRef]
Ong, M.C.; Myrhaug, D.; Hesten, P. Scour around vertical piles due to long-crested and short-crested nonlinear random waves plus a current. Coast. Eng. 2013, 73, 106–114. [Google Scholar] [CrossRef]
Yakhot, V.; Orszag, S.A. Renormalization group analysis of turbulence. I. Basic theory. J. Sci. Comput. 1986, 1, 3–51. [Google Scholar] [CrossRef]
Yakhot, V.; Smith, L.M. The renormalization group, the e-expansion and derivation of turbulence models. J. Sci. Comput. 1992, 7, 35–61. [Google Scholar] [CrossRef]
Mastbergen, D.R.; Berg, J.V.D. Breaching in fine sands and the generation of sustained turbidity currents in submarine canyons. Sedimentology 2003, 50, 625–637. [Google Scholar] [CrossRef]
Soulsby, R. Dynamics of Marine Sands; Thomas Telford Ltd.: London, UK, 1998. [Google Scholar] [CrossRef]
Van Rijn, L.C. Sediment Transport, Part I: Bed Load Transport. J. Hydraul. Eng. 1984, 110, 1431–1456. [Google Scholar] [CrossRef][Green Version]
Zhang, Q.; Zhou, X.-L.; Wang, J.-H. Numerical investigation of local scour around three adjacent piles with different arrangements under current. Ocean Eng. 2017, 142, 625–638. [Google Scholar] [CrossRef]
Yu, Y.X.; Liu, S.X. Random Wave and Its Applications to Engineering, 4th ed.; Dalian University of Technology Press: Dalian, China, 2011. [Google Scholar]
Pang, A.; Skote, M.; Lim, S.; Gullman-Strand, J.; Morgan, N. A numerical approach for determining equilibrium scour depth around a mono-pile due to steady currents. Appl. Ocean Res. 2016, 57, 114–124. [Google Scholar] [CrossRef]
Higuera, P.; Lara, J.L.; Losada, I.J. Three-dimensional interaction of waves and porous coastal structures using Open-FOAM®. Part I: Formulation and validation. Coast. Eng. 2014, 83, 243–258. [Google Scholar] [CrossRef]
Corvaro, S.; Crivellini, A.; Marini, F.; Cimarelli, A.; Capitanelli, L.; Mancinelli, A. Experimental and Numerical Analysis of the Hydrodynamics around a Vertical Cylinder in Waves. J. Mar. Sci. Eng. 2019, 7, 453. [Google Scholar] [CrossRef][Green Version]
Flow3D User Manual, version 11.0.3; Flow Science, Inc.: Santa Fe, NM, USA, 2013.
Khosronejad, A.; Kang, S.; Sotiropoulos, F. Experimental and computational investigation of local scour around bridge piers. Adv. Water Resour. 2012, 37, 73–85. [Google Scholar] [CrossRef]
Stahlmann, A. Experimental and Numerical Modeling of Scour at Foundation Structures for Offshore Wind Turbines. Ph.D. Thesis, Franzius-Institute for Hydraulic, Estuarine and Coastal Engineering, Leibniz Universität Hannover, Hannover, Germany, 2013. [Google Scholar]
Breusers, H.N.C.; Nicollet, G.; Shen, H. Local Scour Around Cylindrical Piers. J. Hydraul. Res. 1977, 15, 211–252. [Google Scholar] [CrossRef]
Schendel, A.; Hildebrandt, A.; Goseberg, N.; Schlurmann, T. Processes and evolution of scour around a monopile induced by tidal currents. Coast. Eng. 2018, 139, 65–84. [Google Scholar] [CrossRef]

Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

© 2021 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https://creativecommons.org/licenses/by/4.0/).

Share and Cite

MDPI and ACS Style

Hu, R.; Liu, H.; Leng, H.; Yu, P.; Wang, X. Scour Characteristics and Equilibrium Scour Depth Prediction around Umbrella Suction Anchor Foundation under Random Waves. J. Mar. Sci. Eng. 2021, 9, 886. https://doi.org/10.3390/jmse9080886

AMA Style

Hu R, Liu H, Leng H, Yu P, Wang X. Scour Characteristics and Equilibrium Scour Depth Prediction around Umbrella Suction Anchor Foundation under Random Waves. Journal of Marine Science and Engineering. 2021; 9(8):886. https://doi.org/10.3390/jmse9080886Chicago/Turabian Style

Hu, Ruigeng, Hongjun Liu, Hao Leng, Peng Yu, and Xiuhai Wang. 2021. “Scour Characteristics and Equilibrium Scour Depth Prediction around Umbrella Suction Anchor Foundation under Random Waves” Journal of Marine Science and Engineering 9, no. 8: 886. https://doi.org/10.3390/jmse9080886

Find Other Styles

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

For more information on the journal statistics, click here.

Multiple requests from the same IP address are counted as one view.

원본 다운로드

Strain rate magnitude at the free surface, illustrating Kelvin-Helmoltz (KH) shear instabilities.

On the reef scale hydrodynamics at Sodwana Bay, South Africa

Environmental Fluid Mechanics (2022)Cite this article

277 Accesses
Metricsdetails

Abstract

The hydrodynamics of coral reefs strongly influences their biological functioning, impacting processes such as nutrient availability and uptake, recruitment success and bleaching. For example, coral reefs located in oligotrophic regions depend on upwelling for nutrient supply. Coral reefs at Sodwana Bay, located on the east coast of South Africa, are an example of high latitude marginal reefs. These reefs are subjected to complex hydrodynamic forcings due to the interaction between the strong Agulhas current and the highly variable topography of the region. In this study, we explore the reef scale hydrodynamics resulting from the bathymetry for two steady current scenarios at Two-Mile Reef (TMR) using a combination of field data and numerical simulations. The influence of tides or waves was not considered for this study as well as reef-scale roughness. Tilt current meters with onboard temperature sensors were deployed at selected locations within TMR. We used field observations to identify the dominant flow conditions on the reef for numerical simulations that focused on the hydrodynamics driven by mean currents. During the field campaign, southerly currents were the predominant flow feature with occasional flow reversals to the north. Northerly currents were associated with greater variability towards the southern end of TMR. Numerical simulations showed that Jesser Point was central to the development of flow features for both the northerly and southerly current scenarios. High current variability in the south of TMR during reverse currents is related to the formation of Kelvin-Helmholtz type shear instabilities along the outer edge of an eddy formed north of Jesser Point. Furthermore, downward vertical velocities were computed along the offshore shelf at TMR during southerly currents. Current reversals caused a change in vertical velocities to an upward direction due to the orientation of the bathymetry relative to flow directions.

Highlights

A predominant southerly current was measured at Two-Mile Reef with occasional reversals towards the north.
Field observations indicated that northerly currents are spatially varied along Two-Mile Reef.
Simulation of reverse currents show the formation of a separated flow due to interaction with Jesser Point with Kelvin–Helmholtz type shear instabilities along the seaward edge.

지금까지 Sodwana Bay에서 자세한 암초 규모 유체 역학을 모델링하려는 시도는 없었습니다. 이러한 모델의 결과는 규모가 있는 산호초 사이의 흐름이 산호초 건강에 어떤 영향을 미치는지 탐색하는 데 사용할 수 있습니다. 이 연구에서는 Sodwana Bay의 유체역학을 탐색하는 데 사용할 수 있는 LES 모델을 개발하기 위한 단계별 접근 방식을 구현합니다. 여기서 우리는 이 초기 단계에서 파도와 조수의 영향을 배제하면서 Agulhas 해류의 유체역학에 초점을 맞춥니다. 이 접근법은 흐름의 첫 번째 LES를 제시하고 Sodwana Bay의 산호초에서 혼합함으로써 향후 연구의 기초를 제공합니다.

This is a preview of subscription content, access via your institution.

References

Anarde K, Myres H, Figlus J (2016) Tilt current meter field validation in the surf zone. In: AGU fall meeting abstracts, vol 2016, pp EP23A—-0950
Blocken B (2018) LES over RANS in building simulation for outdoor and indoor applications: A foregone conclusion? Build Simul 11(5):821–870. https://doi.org/10.1007/s12273-018-0459-3 Article Google Scholar
Booij N, Ris RC, Holthuijsen LH (1999) A third-generation wave model for coastal regions: 1. Model description and validation. J Geophys Res Ocean 104(C4):7649–7666. https://doi.org/10.1029/98JC02622 Article Google Scholar
Bouffanais R (2010) Advances and challenges of applied large-eddy simulation. Comput Fluids 39:735–738. https://doi.org/10.1016/j.compfluid.2009.12.003 Article Google Scholar
Celliers L, Schleyer MH (2002) Coral bleaching on high-latitude marginal reefs at Sodwana Bay, South Africa. Mar Pollut Bull 44:1380–1387Article Google Scholar
Celliers L, Schleyer MH (2008) Coral community structure and risk assessment of high-latitude reefs at Sodwana Bay, South Africa. Biodivers Conserv 17(13):3097–3117. https://doi.org/10.1007/s10531-007-9271-6 Article Google Scholar
Chen SC (2018) Performance assessment of FLOW-3D and XFlow in the numerical modelling of fish-bone type fishway hydraulics https://doi.org/10.15142/T3HH1J
Corbella S, Pringle J, Stretch DD (2015) Assimilation of ocean wave spectra and atmospheric circulation patterns to improve wave modelling. Coast Eng 100:1–10. https://doi.org/10.1016/j.coastaleng.2015.03.003 Article Google Scholar
Davis KA, Pawlak G, Monismith SG (2021) Turbulence and coral reefs. Ann Rev Mar Sci. https://doi.org/10.1146/annurev-marine-042120-071823 Article Google Scholar
Flow Science Inc (2018) FLOW-3D, Version 12.0 Users Manual. Santa Fe, NM, https://www.flow3d.com/
Flow Science Inc (2019) FLOW-3D, Version 12.0 [Computer Software]. Santa Fe, NM, https://www.flow3d.com/
Franco A, Moernaut J, Schneider-Muntau B, Strasser M, Gems B (2020) The 1958 Lituya Bay tsunami – pre-event bathymetry reconstruction and 3D numerical modelling utilising the computational fluid dynamics software Flow-3D. Nat Hazards Earth Syst Sci 20(8):2255–2279Article Google Scholar
Fringer OB, Gerritsen M, Street RL (2006) An unstructured-grid, finite-volume, nonhydrostatic, parallel coastal ocean simulator. Ocean Model 14(3):139–173Article Google Scholar
Fringer OB, Dawson CN, He R, Ralston DK, Zhang YJ (2019) The future of coastal and estuarine modeling: findings from a workshop. Ocean Model 143(September):101458. https://doi.org/10.1016/j.ocemod.2019.101458 Article Google Scholar
Glassom D, Celliers L, Schleyer MH (2006) Coral recruitment patterns at Sodwana Bay, South Africa. Coral Reefs 25(3):485–492. https://doi.org/10.1007/s00338-006-0117-6 Article Google Scholar
Gomes A, Pinho JLS, Valente T, do Carmo JS, Hegde VA (2020) Performance assessment of a semi-circular breakwater through CFD modelling. J Mar Sci Eng. https://doi.org/10.3390/jmse8030226 Article Google Scholar
Green RH, Lowe RJ, Buckley ML (2018) Hydrodynamics of a tidally forced coral reef atoll. J Geophys Res Oceans 123(10):7084–7101. https://doi.org/10.1029/2018JC013946 Article Google Scholar
Hansen AB, Carstensen S, Christensen DF, Aagaard T (2017) Performance of a tilt current meter in the surf zone. Coastal dynamics
Hench JL, Rosman JH (2013) Observations of spatial flow patterns at the coral colony scale on a shallow reef flat. J Geophys Res Ocean 118(3):1142–1156. https://doi.org/10.1002/jgrc.20105 Article Google Scholar
Hirt CW (1993) Volume-fraction techniques: powerful tools for wind engineering. J Wind Eng Ind Aerodyn 46–47:327–338. https://doi.org/10.1016/0167-6105(93)90298-3 Article Google Scholar
Hirt CW, Sicilian JM (1985) A porosity technique for the definition of obstacles in rectangular cell meshes. In: Proceedings of 4th International Conference on Ship Hydrodynamics https://ci.nii.ac.jp/naid/10009570543/en/
Hocker LO, Hruska MA (2004) Interleaving synchronous data and asynchronous data in a single data storage file
Hossain MM, Staples AE (2020) Effects of coral colony morphology on turbulent flow dynamics. PLoS ONE 15(10):e0225676. https://doi.org/10.1371/journal.pone.0225676 Article Google Scholar
Jacob B, Stanev EV (2021) Understanding the impact of bathymetric changes in the german bight on coastal hydrodynamics: one step toward realistic morphodynamic modeling. Front Mar Sci. https://doi.org/10.3389/fmars.2021.640214 Article Google Scholar
Koehl MAR, Hadfield MG (2010) Hydrodynamics of larval settlement from a larva’s point of view. Integr Comp Biol 50(4):539–551. https://doi.org/10.1093/icb/icq101 Article Google Scholar
Lim A, Wheeler AJ, Price DM, O’Reilly L, Harris K, Conti L (2020) Influence of benthic currents on cold-water coral habitats: a combined benthic monitoring and 3D photogrammetric investigation. Sci Rep 10(1):19433. https://doi.org/10.1038/s41598-020-76446-y Article Google Scholar
Limer BD, Bloomberg J, Holstein DM (2020) The influence of eddies on coral larval retention in the flower garden banks. Front Mar Sci 7:372. https://doi.org/10.3389/fmars.2020.00372 Article Google Scholar
Monismith SG (2007) Hydrodynamics of coral reefs. Annu Rev Fluid Mech 39(1):37–55. https://doi.org/10.1146/annurev.fluid.38.050304.092125 Article Google Scholar
Morris T (2009) Physical oceanography of Sodwana Bay and its effect on larval transport and coral bleaching. PhD thesis, Cape Peninsula University of Technology
Morris T, Lamont T, Roberts MJ (2013) Effects of deep-sea eddies on the northern KwaZulu-Natal shelf, South Africa. Afr J Mar Sci 35(3):343–350. https://doi.org/10.2989/1814232X.2013.827991 Article Google Scholar
Perry C, Larcombe P (2003) Marginal and non-reef-building coral environments. Coral Reefs 22:427–432. https://doi.org/10.1007/s00338-003-0330-5 Article Google Scholar
Pope SB (2001) Turbulent flows. Cambridge University Press, CambridgeGoogle Scholar
Porter SN (2009) Biogeography and potential factors regulating shallow subtidal reef communities in the Western Indian Ocean. PhD thesis, University of Cape Town
Porter SN, Schleyer MH (2017) Long-term dynamics of a high-latitude coral reef community at Sodwana Bay, South Africa. Coral Reefs 36(2):369–382. https://doi.org/10.1007/s00338-016-1531-z Article Google Scholar
Porter SN, Schleyer MH (2019) Environmental variation and how its spatial structure influences the cross-shelf distribution of high-latitude coral communities in South Africa. Diversity. https://doi.org/10.3390/d11040057 Article Google Scholar
Ramsay PJ (1994) Marine geology of the Sodwana Bay shelf, southeast Africa. Mar Geol 120(3–4):225–247. https://doi.org/10.1016/0025-3227(94)90060-4 Article Google Scholar
Ramsay PJ, Mason TR (1990) Development of a type zoning model for Zululand coral reefs, Sodwana Bay, South Africa. J Coastal Res 6(4):829–852Google Scholar
Reguero BG, Beck MW, Agostini VN, Kramer P, Hancock B (2018) Coral reefs for coastal protection: a new methodological approach and engineering case study in Grenada. J Environ Manag 210:146–161. https://doi.org/10.1016/j.jenvman.2018.01.024 Article Google Scholar
Reidenbach M, Stocking J, Szczyrba L, Wendelken C (2021) Hydrodynamic interactions with coral topography and its impact on larval settlement. Coral Reefs 40:1–15. https://doi.org/10.1007/s00338-021-02069-y Article Google Scholar
Reidenbach MA, Koseff JR, Koehl MAR (2009) Hydrodynamic forces on larvae affect their settlement on coral reefs in turbulent, wave-driven flow. Limnol Oceanogr 54(1):318–330. https://doi.org/10.4319/lo.2009.54.1.0318 Article Google Scholar
Roberts H, Richardson J, Lagumbay R, Meselhe E, Ma Y (2013) Hydrodynamic and sediment transport modeling using FLOW-3D for siting and optimization of the LCA medium diversion at white ditch hydrodynamic and sediment transport modeling using FLOW-3D for siting and optimization of the LCA medium diversion at white D (December)
Roberts MJ, Ribbink AJ, Morris T, Berg MAVD, Engelbrecht DC, Harding RT (2006) Oceanographic environment of the Sodwana Bay coelacanths (Latimeria chalumnae), South Africa: coelacanth research. South Afr J Sci 102(9):435–443Google Scholar
Rogers JS, Monismith SG, Feddersen F, Storlazzi CD (2013) Hydrodynamics of spur and groove formations on a coral reef. J Geophys Res Ocean 118(6):3059–3073. https://doi.org/10.1002/jgrc.20225 Article Google Scholar
Rogers JS, Monismith SG, Koweek DA, Torres WI, Dunbar RB (2016) Thermodynamics and hydrodynamics in an atoll reef system and their influence on coral cover. Limnol Oceanogr 61(6):2191–2206. https://doi.org/10.1002/lno.10365 Article Google Scholar
Schleyer MH, Celliers L (2003) Coral dominance at the reef-sediment interface in marginal coral communities at Sodwana Bay, South Africa. Mar Freshw Res 54(8):967–972. https://doi.org/10.1071/MF02049 Article Google Scholar
Schleyer MH, Porter SN (2018) Chapter One – drivers of soft and stony coral community distribution on the high-latitude coral reefs of South Africa. advances in marine biology, vol 80, Academic Press, pp 1–55, https://doi.org/10.1016/bs.amb.2018.09.001
Scott F, Antolinez JAA, McCall R, Storlazzi C, Reniers A, Pearson S (2020) Hydro-morphological characterization of coral reefs for wave runup prediction. Front Mar Sci 7:361. https://doi.org/10.3389/fmars.2020.00361 Article Google Scholar
Sebens KP, Grace SP, Helmuth B, Maney EJ Jr, Miles JS (1998) Water flow and prey capture by three scleractinian corals, Madracis mirabilis, Montastrea cavernosa and Porites porites, in a field enclosure. Mar Biol 131(2):347–360Article Google Scholar
Smagorinsky J (1963) General circulation experiments with the primitive equations. Mon Weather Rev 91(3):99–164Article Google Scholar
Stocking J, Laforsch C, Sigl R, Reidenbach M (2018) The role of turbulent hydrodynamics and surface morphology on heat and mass transfer in corals. J R Soc Interface 15:20180448. https://doi.org/10.1098/rsif.2018.0448 Article Google Scholar
Van Leer B (1977) Towards the ultimate conservative difference scheme III. Upstream-centered finite-difference schemes for ideal compressible flow. J Comput Phys 23(3):263–275. https://doi.org/10.1016/0021-9991(77)90094-8 Article Google Scholar
Wells C, Pringle J, Stretch D (2021) Cold water temperature anomalies on the Sodwana reefs and their driving mechanisms. South Afr J Sci. https://doi.org/10.17159/sajs.2021/9304 Article Google Scholar
Wyatt ASJ, Lowe RJ, Humphries S, Waite AM (2010) Particulate nutrient fluxes over a fringing coral reef: relevant scales of phytoplankton production and mechanisms of supply. Mar Ecol Prog Ser 405:113–130Article Google Scholar
Yao Y, He T, Deng Z, Chen L, Guo H (2019) Large eddy simulation modeling of tsunami-like solitary wave processes over fringing reefs. Nat Hazards Earth Syst Sci 19(6):1281–1295. https://doi.org/10.5194/nhess-19-1281-2019 Article Google Scholar
Zhao Q, Tanimoto K (1998) Numerical simulation of breaking waves by large eddy simulation and vof method. Coastal Engineering Proceedings 1(26), 10.9753/icce.v26.%p, https://journals.tdl.org/icce/index.php/icce/article/view/5656

Text and image taken from Deoraj, et al. (2022), On the reef scale hydrodynamics at Sodwana Bay, South Africa. Preprint courtesy the authors.

Figure 2: 3D (left) and 2D (right) views of wave elevation using case C

CFD 접근법을 사용하여 파도에서 하이드로포일의 SEAKEEPING 성능

SYAFIQ ZIKRYAND FITRIADHY*
Faculty of Ocean Engineering Technology and Informatics, Universiti Malaysia Terengganu, 21030 Kuala
Terengganu, Terengganu, Malaysia
*
Corresponding author: naoe.afit@gmail.com http://doi.org/10.46754/umtjur.2021.07.017

Abstract

수중익선은 일반적으로 열악한 환경 조건으로 인해 승객의 편안함에 영향을 미칠 수 있는 높은 저항과 과도한 수직 운동(히브 및 피치)을 경험합니다. 따라서 복잡한 유체역학적 현상이 존재하기 때문에 파랑에서 수중익선의 내항성능을 규명할 필요가 있다.

이를 위해 수중익선 운동에 대한 CFD(Computational Fluid Dynamic) 해석을 제안한다. Froude Number 및 포일 받음각과 같은 여러 매개변수가 고려되었습니다.

그 결과 Froude Number의 후속 증가는 히브 및 피치 운동에 반비례한다는 것이 밝혀졌습니다. 본질적으로 이것은 높은 응답 진폭 연산자(RAO)의 형태로 제공되는 수중익선 항해 성능의 업그레이드로 이어졌습니다.

또한 포일 선수의 증가하는 각도는 히브 운동에 비례하는 반면, 포일 선미는 7.5o에서 낮은 히브 운동을 보였고, 그 다음으로 5o, 10o 순으로 나타났다. 피치모션의 경우 포일 보우의 증가는 5o에서 더 낮았고, 그 다음이 10o, 7.5o 순이었다. 포일 선미의 증가는 수중익선에 의한 피치 모션 경험에 비례했습니다.

일반적으로 이 CFD 시뮬레이션은 앞서 언급한 설계 매개변수와 관련하여 공해 상태에서 수중익선 설계의 운영 효율성을 보장하는 데 매우 유용합니다.

Keywords

CFD, hydrofoil, foil angle of attack, heave, pitch.

Figure 1: Overall mesh block being used in simulation

References

Djavareshkian, M. H., & Esmaeili, A. (2014). Heuristic optimization of submerged hydrofoil
using ANFIS–PSO. Ocean Engineering, 92, 55-63.
Fitriadhy, A., & Adam, N. A. (2017). Heave and pitch motions performance of a monotricat ship in
head-seas. International Journal of Automotive and Mechanical Engineering, 14, 4243-4258.
Islam, M., Jahra, F., & Hiscock, S. (2016). Data analysis methodologies for hydrodynamic
experiments in waves. Journal of Naval Architecture and Marine Engineering, 13(1),
1-15.
Koutsourakis, N., Bartzis, J. G., & Markatos, N. C. (2012). Evaluation of Reynolds stress, k-ε and
RNG k-ε turbulence models in street canyon flows using various experimental datasets.
Environmental fluid mechanics, 1-25.
Manual, F. D. U. (2011). Flow3D User Manual, v9. 4.2, Flow Science. Inc., Santa Fe, NM. Matveev, K., & Duncan, R. (2005). Development
of the tool for predicting hydrofoil system performance and simulating motion of hydrofoil-assisted boats. Paper presented at the High Speed and High Performance Ship and Craft Symposium, Everett/WA: ASNE, USA.
Seif, M., Mehdigholi, H., & Najafi, A. (2014). Experimental and numerical modeling of the
high speed planing vessel motion. Journal of Marine Engineering & Technology, 13(2), 62-
72.
Sun, X., Yao, C., Xiong, Y., & Ye, Q. (2017). Numerical and experimental study on
seakeeping performance of a swath vehicle in head waves. Applied Ocean Research, 68, 262-
275.
Vakilabadi, K. A., Khedmati, M. R., & Seif, M.S. (2014). Experimental study on heave and
pitch motion characteristics of a wave-piercing trimaran. Transactions of FAMENA, 38(3), 13-
26.
Yakhot, A., Rakib, S., & Flannery, W. (1994). LowReynolds number approximation for turbulent
eddy viscosity. Journal of scientific computing, 9(3), 283-292.
Yakhot, V., & Orszag, S. A. (1986). Renormalization group analysis of turbulence. I. Basic theory.
Journal of scientific computing, 1(1), 3-51.

PDF 원문보기

Fig. 2. Design of the grate inlet types studied: (a) R1, (b) R2, (c) R3, (d) R4, (e) R5, (f) R6, (g) R7 (source: based on geometries of Chaparro Andrade and Abaunza Tabares, 2021)

Three-dimensional Numerical Evaluation of Hydraulic Efficiency and Discharge Coefficient in Grate Inlets

쇠창살 격자 유입구의 수리효율 및 배출계수에 대한 3차원 수치적 평가

Melquisedec Cortés Zambrano*, Helmer Edgardo Monroy González,
Wilson Enrique Amaya Tequia
Faculty of Civil Engineering, Santo Tomas Tunja University. Address Av. Universitaria No. 45-202.
Tunja – Boyacá – Colombia

Abstract

홍수는 지반이동 및 이동의 원인 중 하나이며, 급속한 도시화 및 도시화로 인해 이전보다 빈번하게 발생할 수 있다. 도시 배수 시스템의 특성은 집수 요소가 결정적인 역할을 하는 범람의 발생 및 범위를 정의할 수 있습니다. 이 문서는 7가지 유형의 화격자 유입구의 수력 유입 효율 및 배출 계수에 대한 수치 조사를 제시합니다. FLOW-3D® 시뮬레이터는 Q = 24, 34.1, 44, 100, 200 및 300 L/s의 유속에서 풀 스케일로 격자를 테스트하는 데 사용되며 종방향 기울기가 1.0인 실험 프로토타입의 구성을 유지합니다. %, 1.5% 및 2.0% 및 고정 횡단 경사, 총 126개 모델. 그 결과를 바탕으로 종류별 및 종단경사 조건에 따른 수력유입구 효율곡선과 토출계수를 구성하였다. 결과는 다른 조사에서 제안된 경험적 공식으로 조정되어 프로토타입의 물리적 테스트 결과를 검증하는 역할을 합니다.

Floods are one of the causes of ground movement and displacement, and due to rapid urbanization and urban growth may occur more frequently than before. The characteristics of an urban drainage system can define the occurrence and extent of flooding, where catchment elements have a determining role. This document presents the numerical investigation of the hydraulic inlet efficiency and the discharge coefficient of seven types of grate inlets. The FLOW-3D® simulator is used to test the gratings at a full scale, under flow rates of Q = 24, 34.1, 44, 100, 200 and 300 L/s, preserving the configuration of the experimental prototype with longitudinal slopes of 1.0%, 1.5% and 2.0% and a fixed cross slope, for a total of 126 models. Based on the results, hydraulic inlet efficiency curves and discharge coefficients are constructed for each type and a longitudinal slope condition. The results are adjusted with empirical formulations proposed in other investigations, serving to verify the results of physical testing of prototypes.

Keywords

grate inlet, inlet efficiency, discharge coefficient, computational fluid dynamic, 3D modelling.

Fig. 1. Physical model of the experimental campaign (source: Chaparro Andrade and Abaunza Tabares, 2021)

Fig. 4. Comparison between the results obtained during physical experimentation in prototype 7 and simulation results with FLOW-3D® (source:
made with FlowSight® and photographic record by Chaparro Andrade and Abaunza Tabares, 2021) — Fig. 4. Comparison between the results obtained during physical experimentation in prototype 7 and simulation results with FLOW-3D® (source: made with FlowSight® and photographic record by Chaparro Andrade and Abaunza Tabares, 2021)

Fig. 6. Example of the results of flow depth and velocity vectors in the xy plane, for a stable flow condition in a grate inlet type and free surface
configuration and flow regime, of some grating types (source: produced with FlowSight®) — Fig. 6. Example of the results of flow depth and velocity vectors in the xy plane, for a stable flow condition in a grate inlet type and free surface configuration and flow regime, of some grating types (source: produced with FlowSight®)

References

Alia Md., S., and Sabtu, N. (2020). Comparison of Different Methodologies for Determining the Efficiency of Gully Inlets. In F. M.
Nazri (Ed.), Proceedings of AICCE‘19: Transforming the Nation
for a Sustainable Tomorrow (Vol. 53, pp. 1275-1284). Springer
Nature Switzerland AG. https://doi.org/10.1007/978-3-030-
32816-0_99
Antunes do Carmo, J. S. (2020). Physical Modelling vs. Numerical Modelling: Complementarity and Learning. July. https://doi.
org/10.20944/preprints202007.0753.v1
Aragón-Hernández, J. L. (2013). Modelación numérica integrada de los procesos hidráulicos en el drenaje urbano [Universidad Politécnica de Cataluña]. In Doctoral Tesis. https://
upcommons.upc.edu/handle/2117/95059?locale-attribute=es
Argue, J. R., and Pezzaniti, D. (1996). How reliable are inlet
(hydraulic) models at representing stormwater flow? Science
of the Total Environment, 189-190, 355-359. https://doi.org/10.1016/0048-9697(96)05231-X
Banco Mundial, O. (2019). Agua: Panorama general. https://
www.bancomundial.org/es/topic/water/overview
Cárdenas-Quintero, M., Carvajal-Serna, L. F., and Marbello-Pérez, R. (2018). Evaluación numérica tridimensional de un
sumidero de reja de fondo (Three-Dimensional Numerical Assessment of Grate Inlet). SSRN Electronic Journal, November.
https://doi.org/10.2139/ssrn.3112980
Carvalho, R. F., Lopes, P., Leandro, J., and David, L. M. (2019).
Numerical Research of Flows into Gullies with Different Outlet Locations. Water, 11(2), 794. https://doi.org/10.3390/
w11040794
Chaparro Andrade, F. G., and Abaunza Tabares, K. V. (2021). Importancia de los sumideros, su funcionamiento y diseño en redes de alcantarillado caso de estudio sector nororiental Tunja.
Universidad Santo Tomás.
Cortés Zambrano, M., Amaya Tequia, W. E., and Gamba Fernández, D. S. (2020). Implementation of the hydraulic modelling of
urban drainage in the northeast sector, Tunja, Boyacá. Revista
Facultad de Ingeniería Universidad de Antioquia. https://doi.
org/10.17533/udea.redin.20200578
Cosco, C., Gómez, M., Russo, B., Tellez-Alvarez, J., Macchione, F., Costabile, P., and Costanzo, C. (2020). Discharge coefficients for specific grated inlets. Influence of the Froude
number. Urban Water Journal, 17(7), 656-668. https://doi.org/10.1080/1573062X.2020.1811881
Despotovic, J., Plavsic, J., Stefanovic, N., and Pavlovic, D. (2005).
Inefficiency of storm water inlets as a source of urban floods.
Water Science and Technology, 51(2), 139-145. https://doi.
org/10.2166/wst.2005.0041
Ellis, J. B., and Marsalek, J. (1996). Overview of urban drainage:
Environmental impacts and concerns, means of mitigation and
implementation policies. Journal of Hydraulic Research, 34(6),
723-732. https://doi.org/10.1080/00221689609498446
Fang, X., Jiang, S., and Alam, S. R. (2010). Numerical simulations of efficiency of curb-opening inlets. Journal of Hydraulic
Engineering, 136(1), 62-66. https://doi.org/10.1061/(ASCE)
HY.1943-7900.0000131
Faram, M. G., and Harwood, R. (2000). CFD for the Water Industry; The Role of CFD as a Tool for the Development of Wastewater Treatment Systems. Hydro International, 21-22.
Faram, M. G., and Harwood, R. (2002). Assessment of the
effectiveness of stormwater treatment chambers using
computational fluid dynamics. Global Solutions for Urban Drainage, 40644(September 2002), 1-14. https://doi.
org/10.1061/40644(2002)7
Flow Science, I. (2018). FLOW-3D® Version 12.0 Users Manual.
In FLOW-3D [Computer software]. https://www.flow3d.com
Flow Science, I. (2019). FLOW-3D® Version 12.0 [Computer software] (No. 12). https://www.flow3d.com
Ghanbari, R., and Heidarnejad, M. (2020). Experimental and numerical analysis of flow hydraulics in triangular and rectangular
piano key weirs. Water Science, 00(00), 1-7. https://doi.org/10.
1080/11104929.2020.1724649

Gómez, M., and Russo, B. (2005a). Comparative study of methodologies to determine inlet efficiency from test data. HEC-12
methodology vs UPC method. Water Resources Management,
Algarve, Portugal., 80(October 2014), 623-632. https://doi.
org/10.2495/WRM050621
Gómez, M., and Russo, B. (2005b). Comparative study among
different methodologies to determine storm sewer inlet efficiency from test data. 10th International Conference on Urban
Drainage, August, 21-26. https://www.researchgate.net/publication/255602448_Comparative_study_among_different_methodologies_to_determine_storm_sewer_inlet_efficiency_
from_test_data
Gómez, M., Recasens, J., Russo, B., and Martínez-Gomariz, E.
(2016). Assessment of inlet efficiency through a 3D simulation: Numerical and experimental comparison. Water Science
and Technology, 74(8), 1926-1935. https://doi.org/10.2166/
wst.2016.326
Gómez, M., and Russo, B. (2011). Methodology to estimate hydraulic efficiency of drain inlets. Proceedings of the Institution of
Civil Engineers: Water Management, 164(2), 81-90. https://doi.
org/10.1680/wama.900070
Gómez Valentin, M. (2007). Hidrología urbana. In Hidrología Urbana (pp. 135-147). Instituto Flumen.
Jakeman, A. J., Letcher, R. A., and Norton, J. P. (2006). Ten iterative steps in development and evaluation of environmental
models. Environmental Modelling and Software, 21, 602-614.
https://doi.org/10.1016/j.envsoft.2006.01.004
Jang, J. H., Hsieh, C. T., and Chang, T. H. (2019). The importance of gully flow modelling to urban flood simulation. Urban Water Journal, 16(5), 377-388. https://doi.org/10.1080/1573062X.2019.1669198
Kaushal, D. R., Thinglas, T., Tomita, Y., Kuchii, S., and Tsukamoto, H. (2012). Experimental investigation on optimization of
invert trap configuration for sewer solid management. Powder Technology, 215-216, 1-14. https://doi.org/10.1016/j.powtec.2011.08.029
Khazaee, I., and Mohammadiun, M. (2010). Effects of flow field
on open channel flow properties using numerical investigation
and experimental comparison. International Journal of Energy
and Environment, 1(6), 1083-1096. https://doi.org/10.1016/
S0031-9384(10)00122-8
Kleidorfer, M., Tscheikner-Gratl, F., Vonach, T., and Rauch, W.
(2018). What can we learn from a 500-year event? Experiences
from urban drainage in Austria. Water Science and Technology,
77(8), 2146-2154. https://doi.org/10.2166/wst.2018.138
Leitão, J. P., Simões, N. E., Pina, R. D., Ochoa-Rodriguez, S.,
Onof, C., and Sá Marques, A. (2017). Stochastic evaluation of
the impact of sewer inlets‘ hydraulic capacity on urban pluvial
flooding. Stochastic Environmental Research and Risk Assessment, 31(8), 1907-1922. https://doi.org/10.1007/s00477-016-
1283-x
Lopes, P., Leandro, J., Carvalho, R. F., Russo, B., and Gómez, M.
(2016). Assessment of the ability of a volume of fluid model to
reproduce the efficiency of a continuous transverse gully with
grate. Journal of Irrigation and Drainage Engineering, 142(10),
1-9. https://doi.org/10.1061/(ASCE)IR.1943-4774.0001058
Mohsin, M., and Kaushal, D. R. (2016). 3D CFD validation of invert trap efficiency for sewer solid management using VOF model. Water Science and Engineering, 9(2), 106-114. https://doi.
org/10.1016/j.wse.2016.06.006
Palla, A., Colli, M., Candela, A., Aronica, G. T., and Lanza, L.
G. (2018). Pluvial flooding in urban areas: the role of surface
drainage efficiency. Journal of Flood Risk Management, 11,
S663-S676. https://doi.org/10.1111/jfr3.12246
Russo, B. (2010). Design of surface drainage systems according
to hazard criteria related to flooding of urban areas [Universitat
Politècnica de Catalunya]. https://dialnet.unirioja.es/servlet/
tesis?codigo=258828
Sedano-Cruz, K., Carvajal-Escoar, Y., and Ávila Díaz, A. J. (2013).
ANÁLISIS DE ASPECTOS QUE INCREMENTAN EL RIESGO
DE INUNDACIONES EN COLOMBIA. Luna Azul, 37, 219-218.
https://www.redalyc.org/articulo.oa?id=321729206014
Spaliviero, F., May, R. W. P., Escarameia, M. (2000). Spacing of road gullies. Hydraulic performance of BS EN 124 gully gratings. HR Walingford, 44(0). https://doi.org/10.13140/
RG.2.1.1344.0889
Téllez-Álvarez, J., Gómez, M., and Russo, B. (2020). Quantification of energy loss in two grated inlets under pressure. Water
(Switzerland), 12(6). https://doi.org/10.3390/w12061601
Téllez Álvarez, J., Gómez, V., Russo, B., and Redondo, J. M.
(2003). Performance assessment of numerical modelling
for hydraulic efficiency of a grated inlet. 1, 6-8. https://doi.org/10.16309/j.cnki.issn.1007-1776.2003.03.004
Téllez Álvarez, J., Gómez Valentin, M., Paindelli, A., and Russo,
B. (2017). ACTIVIDAD EXPERIMENTAL DE I+D+i EN INGENIERÍA
HIDRÁULICA EN ESPAÑA. In L. J. Balairón Pérez and D. López
Gómez (Eds.), Seminario 2017, Comunicaciones de las líneas prioritarias (pp. 41-43). Universitat Politècnica de València.
https://doi.org/10.1017/CBO9781107415324.004
Téllez Álvarez, J., Gómez Valentin, M., and Russo, B. (2019).
Modelling of Surcharge Flow Through Grated Inlet. In P. Gourbesville and G. Caignaert (Eds.), Advances in Hydroinformati-

cs. Springer, Singapore. https://doi.org/10.1007/978-981-
4451-42-0
UNDRR, I., and CRED, I. (2018). Pérdidas económicas, pobreza y
Desastres 1998 – 2017 (Vol. 6, Issue 1). https://doi.org/10.12962/
j23373520.v6i1.22451
Vyzikas, T., and Greaves, D. (2018). Numerial Modelling.
In D. Greaves and G. Iglesias (Eds.), Wave and Tidal Energy (pp. 289-363). John Wiley and Sons Ltd. https://doi.
org/10.1002/9781119014492
Yakhot, V., and Orszag, S. A. (1986). Renormalization Group Analysis of Turbulence. I . Basic Theory. Journal of Scientific Computing, 1(1), 3-51. https://doi.org/10.1007/BF01061452
Yakhot, V., and Smith, L. M. (1992). The renormalization group,
the ɛ-expansion and derivation of turbulence models. Journal
of Scientific Computing, 7(l), 35-61. https://doi.org/10.1007/
BF01060210
Yazdanfar, Z., and Sharma, A. (2015). Urban drainage system
planning and design – Challenges with climate change and urbanization: A review. Water Science and Technology, 72(2), 165-https://doi.org/10.2166/wst.2015.207

원문 PDF 보기

저수지 측수로형 여수로 불완전월류 정밀안전진단 수리 해석 ( 3차원 전산 수치해석 )

불완전 월류 조건의 저수지 측수로형 여수로에 대한 3차원 전산 추치해석

한국농어촌공사는 수리시설안전진단사업을 통하여 노후 및 기능 저하된 농업생산기반시설물에 대하여 정밀안전진단을 실시하여 사전에 재해, 재난을 대비하고 있습니다.

측수로형 여수로는 계획 홍수량이하의 홍수량이 유입시엔 안정적으로 방류가 일어나나 계획 홍수량 이상의 홍수량이 유입되면 물넘이에서 불완전 월류가 발생하며 방류량이 충분하지 않게 됩니다. 측수로형 여수로는 설계당시의 홍수량에 비해 늘어난 현재에 맞게 변경된 홍수량이 유입할 경우 물넘이에서 불완전월류가 발생하는지를 확인하게 됩니다.

현재 농어촌공사와 농어촌연구원, 수자원공사, 학계 등에서는 전 세계에서 오랜 기간 학계의 연구활동을 통한 수많은 논문 검증과 현장 사용을 통해 검증된 FLOW-3D 수치해석 프로그램을 이용하고 있습니다.

특히, 농어촌공사의 정밀안전진단을 실시할 때 설계홍수량의 저수지 유입 시 물넘이에서 불완전월류가 발생하는지를 확인하고, 불완전월류 발생 시 수위 상승 영향을 분석해 안전성 검토 후 문제가 발견되면 보수·보강 방안을 제시할 수 있는 대표적인 3차원 수치해석 프로그램인 FLOW-3D 를 사용하고 있습니다.

한국농어촌공사 재난안전진단본부 FLOW-3D 수치해석 교육 장면 — 2024년 한국농어촌공사 안전진단본부 여수로 불완전월류 정밀안전진단 FLOW-3D 수치 해석교육 장면

농어촌공사 정밀안전진단 업무 수행시 수치해석이 필요하십니까? 수치해석에 대해 궁금하신 사항이나 용역 의뢰가 필요하시면 언제든지 아래 연락처로 연락 주시기 바랍니다.

<담당자 연락처>

전화 : 02-2026-0455
Email : flow3d@stikorea.co.kr

당사에는 20년 이상 수치해석 수처리 분야의 수치해석 연구에 전념하고 있는 전문 연구인력과 다양한 기술적 경험과 전문 수치해석 용역 서비스를 제공하는 숙련된 기술팀이 준비되어 있습니다.

저수지 정밀안전진단 수치해석 과업 예시

과업의 범위

3차원 수치해석을 통한 OO저수지의 측수로부 수면 검토
측수로 불완전 월류 발생 여부 및 제방 여유고 검토

수치해석 과업 세부내용

가능최대홍수량과 200년, 100년 빈도의 홍수량에 대해 각각의 측수로부 3차원 수치해석

경계조건

가. 수위

만수위
홍수위
– 100년 빈도
– 200년 빈도
– 가능최대홍수량(PMF)

나. 홍수량

100년 빈도의 홍수량
200년 빈도의 홍수량
가능최대홍수량(PMF)

저수지 수위별 방류량 검토 및 제방 여유고 검토

경계조건에 대해 측수로부 물넘이 수면 형상 검토
수위별 방류량을 제공된 수리계산값과 수치해석 결과값을 비교하여 방류 능력 검토
수위에 따른 물넘이 수위를 검토하여 제방 여유고 검토

※ 수위별 수리계산값은 발주처에서 제공

성과물

100년빈도, 200년빈도 및 가능최대홍수량(PMF) 유입에 따른 측수로부 불완전 월류 여부로 인한 제방 여유고 안정성 검토
가능최대홍수량(PMF)을 고려 할 경우 검증된 3차원 수치해석 모델 Data 구축
과업보고서, 보고서 원본 파일 및 PDF 파일, 수치해석 원본 입력 파일 및 결과 파일
기타
※ 모든 성과물은 CD 및 이동저장장치에 별도 저장하여 납품

<수치해석 용역 문의 담당자 연락처>

전화 : 02-2026-0455
Email : flow3d@stikorea.co.kr

Figure 3: Wave pattern at sea surface at 20 knots (10.29 m/s) for mesh 1

Flow-3D에서 CFD 시뮬레이션을 사용한 선박 저항 분석

Ship resistance analysis using CFD simulations in Flow-3D

Author

Deshpande, Sujay; Sundsbø, Per-Arne; Das, Subhashis

Abstract

선박의 동력 요구 사항을 설계할 때 고려해야 할 가장 중요한 요소는 선박 저항 또는 선박에 작용하는 항력입니다. 항력을 극복하는 데 필요한 동력이 추진 시스템의 ‘손실’에 기여하기 때문에 추진 시스템을 설계하는 동안 선박 저항을 추정하는 것이 중요합니다. 선박 저항을 계산하는 세 가지 주요 방법이 있습니다:

Holtrop-Mennen(HM) 방법과 같은 통계적 방법, 수치 분석 또는 CFD(전산 유체 역학) 시뮬레이션 및 모델 테스트, 즉 예인 탱크에서 축소된 모델 테스트. 설계 단계 초기에는 기본 선박 매개변수만 사용할 수 있을 때 HM 방법과 같은 통계 모델만 사용할 수 있습니다.

수치 해석/CFD 시뮬레이션 및 모델 테스트는 선박의 완전한 3D 설계가 완료된 경우에만 수행할 수 있습니다. 본 논문은 Flow-3D 소프트웨어 패키지를 사용하여 CFD 시뮬레이션을 사용하여 잔잔한 수상 선박 저항을 예측하는 것을 목표로 합니다.

롤온/롤오프 승객(RoPax) 페리에 대한 사례 연구를 조사했습니다. 선박 저항은 다양한 선박 속도에서 계산되었습니다. 메쉬는 모든 CFD 시뮬레이션의 결과에 영향을 미치기 때문에 메쉬 민감도를 확인하기 위해 여러 개의 메쉬가 사용되었습니다. 시뮬레이션의 결과를 HM 방법의 추정치와 비교했습니다.

시뮬레이션 결과는 낮은 선박 속도에 대한 HM 방법과 잘 일치했습니다. 더 높은 선속을 위한 HM 방법에 비해 결과의 차이가 상당히 컸다. 선박 저항 분석을 수행하는 Flow-3D의 기능이 시연되었습니다.

While designing the power requirements of a ship, the most important factor to be considered is the ship resistance, or the sea drag forces acting on the ship. It is important to have an estimate of the ship resistance while designing the propulsion system since the power required to overcome the sea drag forces contribute to ‘losses’ in the propulsion system. There are three main methods to calculate ship resistance: Statistical methods like the Holtrop-Mennen (HM) method, numerical analysis or CFD (Computational Fluid Dynamics) simulations, and model testing, i.e. scaled model tests in towing tanks. At the start of the design stage, when only basic ship parameters are available, only statistical models like the HM method can be used. Numerical analysis/ CFD simulations and model tests can be performed only when the complete 3D design of the ship is completed. The present paper aims at predicting the calm water ship resistance using CFD simulations, using the Flow-3D software package. A case study of a roll-on/roll-off passenger (RoPax) ferry was investigated. Ship resistance was calculated at various ship speeds. Since the mesh affects the results in any CFD simulation, multiple meshes were used to check the mesh sensitivity. The results from the simulations were compared with the estimate from the HM method. The results from simulations agreed well with the HM method for low ship speeds. The difference in the results was considerably high compared to the HM method for higher ship speeds. The capability of Flow-3D to perform ship resistance analysis was demonstrated.

Figure 4: Ship Resistance (kN) vs Ship Speed (knots)

Publisher

International Society of Multiphysics

Citation

Deshpande SR, Sundsbø P, Das S. Ship resistance analysis using CFD simulations in Flow-3D. The International Journal of Multiphysics. 2020;14(3):227-236

REFERENCES

[1] K. Min and S. Kang, “Study on the form factor and full-scale ship resistance prediction
method,” Journal of Marine Science and Technology, vol. 15, pp. 108-118, June 2010.
[2] A. Molland, S. Turnock and D. Hudson, “Ship Resistance and Propulsion” Second
Edition. In Ship Resistance and Propulsion: Practical Estimation of Ship Propulsive
Power (pp. 12-69), August 2017, Cambridge University Press.
[3] K. Niklas and H. Pruszko, “Full-scale CFD simulations for the determination of ship
resistance as a rational, alternative method to towing tank experiments,” Ocean
Engineering, vol. 190, October 2019.
[4] A. Elkafas, M. Elgohary and A. Zeid, “Numerical study on the hydrodynamic drag force
of a container ship model,” Alexandria Engineering Journal, vol. 58, no. 3, pp. 849-859,
September 2019.
[5] J. Holtrop and G. Mennen, “An approximate power prediction method,” International
Shipbuilding Progress, vol. 29, no. 335, pp. 166-170, July 1982.
[6] E. Bøckmann and S. Steen, “Model test and simulation of a ship with wavefoils,” Applied
Ocean research, vol. 57, pp. 8-18, April 2016.
[7] K. Atreyapurapu, B. Tallapragada and K. Voonna, “Simulation of a Free Surface Flow
over a Container Vessel Using CFD,” International Journal of Engineering Trends and
Technology (IJETT), vol. 18, no. 7, pp. 334-339, December 2014.
[8] J. Petersen, D. Jacobsen and O. Winther, “Statistical modelling for ship propulsion
efficiency,” Journal of Marine Science and Technology, vol. 17, pp. 30-39, December
2011.
[9] H. Versteeg and W. Malalasekera, An introduction to computational fluid dynamics: the
finite volume method (second edition), Harlow, England: Pearson Education Ltd, 2007.
[10]C. Hirth and B. Nichols, “Volume of fluid (VOF) method for the dynamics of free
boundaries,” Journal of Computational Physics, vol. 39, no. 1, pp. 201-225, January 1981.
[11] A. Nordli and H. Khawaja, “Comparison of Explicit Method of Solution for CFD Euler
Problems using MATLAB® and FORTRAN 77,” International Journal of Multiphysics,
vol. 13, no. 2, 2019.
[12] FLOW-3D® Version 12.0 User’s Manual (2018). FLOW-3D [Computer software]. Santa
Fe, NM: Flow Science, Inc. https://www.flow3d.com.
[13] D. McCluskey and A. Holdø, “Optimizing the hydrocyclone for ballast water treatment
using computational fluid dynamics,” International Journal of Multiphysics, vol. 3, no. 3,
2009.
[14]M. Breuer, D. Lakehal and W. Rodi, “Flow around a Surface Mounted Cubical Obstacle:
Comparison of Les and Rans-Results,” Computation of Three-Dimensional Complex
Flows. Notes on Numerical Fluid Mechanics, vol. 49, p. 1996.
[15] G. Wei, “A Fixed-Mesh Method for General Moving Objects in Fluid Flow”, Modern
Physics Letters B, vol. 19, no. 28, pp. 1719-1722, 2005.
[16]J. Michell, “The wave-resistance of a ship,” The London, Edinburgh, and Dublin
Philosophical Magazine and Journal of Science, Vols. 45, 1898, no. 272, pp. 106-123,
May 2009.

원문보기

Figure 5. Schematic view of flap and support structure [32]

Design Optimization of Ocean Renewable Energy Converter Using a Combined Bi-level Metaheuristic Approach

결합된 Bi-level 메타휴리스틱 접근법을 사용한 해양 재생 에너지 변환기의 설계 최적화

Erfan Amini ^a1, Mahdieh Nasiri ^b1, Navid Salami Pargoo ^a, Zahra Mozhgani ^c, Danial Golbaz ^d, Mehrdad Baniesmaeil ^e, Meysam Majidi Nezhad ^f, Mehdi Neshat ^gj, Davide Astiaso Garcia ^h, Georgios Sylaios ⁱ

Abstract

In recent years, there has been an increasing interest in renewable energies in view of the fact that fossil fuels are the leading cause of catastrophic environmental consequences. Ocean wave energy is a renewable energy source that is particularly prevalent in coastal areas. Since many countries have tremendous potential to extract this type of energy, a number of researchers have sought to determine certain effective factors on wave converters’ performance, with a primary emphasis on ambient factors. In this study, we used metaheuristic optimization methods to investigate the effects of geometric factors on the performance of an Oscillating Surge Wave Energy Converter (OSWEC), in addition to the effects of hydrodynamic parameters. To do so, we used CATIA software to model different geometries which were then inserted into a numerical model developed in Flow3D software. A Ribed-surface design of the converter’s flap is also introduced in this study to maximize wave-converter interaction. Besides, a Bi-level Hill Climbing Multi-Verse Optimization (HCMVO) method was also developed for this application. The results showed that the converter performs better with greater wave heights, flap freeboard heights, and shorter wave periods. Additionally, the added ribs led to more wave-converter interaction and better performance, while the distance between the flap and flume bed negatively impacted the performance. Finally, tracking the changes in the five-dimensional objective function revealed the optimum value for each parameter in all scenarios. This is achieved by the newly developed optimization algorithm, which is much faster than other existing cutting-edge metaheuristic approaches.

Keywords

Wave Energy Converter

OSWEC

Hydrodynamic Effects

Geometric Design

Metaheuristic Optimization

Multi-Verse Optimizer

¹. ^Introduction

The increase in energy demand, the limitations of fossil fuels, as well as environmental crises, such as air pollution and global warming, are the leading causes of calling more attention to harvesting renewable energy recently [1], [2], [3]. While still in its infancy, ocean wave energy has neither reached commercial maturity nor technological convergence. In recent decades, remarkable progress has been made in the marine energy domain, which is still in the early stage of development, to improve the technology performance level (TPL) [4], [5]and technology readiness level (TRL) of wave energy converters (WECs). This has been achieved using novel modeling techniques [6], [7], [8], [9], [10], [11], [12], [13], [14] to gain the following advantages [15]: (i) As a source of sustainable energy, it contributes to the mix of energy resources that leads to greater diversity and attractiveness for coastal cities and suppliers. [16] (ii) Since wave energy can be exploited offshore and does not require any land, in-land site selection would be less expensive and undesirable visual effects would be reduced. [17] (iii) When the best layout and location of offshore site are taken into account, permanent generation of energy will be feasible (as opposed to using solar energy, for example, which is time-dependent) [18].

In general, the energy conversion process can be divided into three stages in a WEC device, including primary, secondary, and tertiary stages [19], [20]. In the first stage of energy conversion, which is the subject of this study, the wave power is converted to mechanical power by wave-structure interaction (WSI) between ocean waves and structures. Moreover, the mechanical power is transferred into electricity in the second stage, in which mechanical structures are coupled with power take-off systems (PTO). At this stage, optimal control strategies are useful to tune the system dynamics to maximize power output [10], [13], [12]. Furthermore, the tertiary energy conversion stage revolves around transferring the non-standard AC power into direct current (DC) power for energy storage or standard AC power for grid integration [21], [22]. We discuss only the first stage regardless of the secondary and tertiary stages. While Page 1 of 16 WECs include several categories and technologies such as terminators, point absorbers, and attenuators [15], [23], we focus on oscillating surge wave energy converters (OSWECs) in this paper due to its high capacity for industrialization [24].

Over the past two decades, a number of studies have been conducted to understand how OSWECs’ structures and interactions between ocean waves and flaps affect converters performance. Henry et al.’s experiment on oscillating surge wave energy converters is considered as one of the most influential pieces of research [25], which demonstrated how the performance of oscillating surge wave energy converters (OSWECs) is affected by seven different factors, including wave period, wave power, flap’s relative density, water depth, free-board of the flap, the gap between the tubes, gap underneath the flap, and flap width. These parameters were assessed in their two models in order to estimate the absorbed energy from incoming waves [26], [27]. In addition, Folly et al. investigated the impact of water depth on the OSWECs performance analytically, numerically, and experimentally. According to this and further similar studies, the average annual incident wave power is significantly reduced by water depth. Based on the experimental results, both the surge wave force and the power capture of OSWECs increase in shallow water [28], [29]. Following this, Sarkar et al. found that under such circumstances, the device that is located near the coast performs much better than those in the open ocean [30]. On the other hand, other studies are showing that the size of the converter, including height and width, is relatively independent of the location (within similar depth) [31]. Subsequently, Schmitt et al. studied OSWECs numerically and experimentally. In fact, for the simulation of OSWEC, OpenFOAM was used to test the applicability of Reynolds-averaged Navier-Stokes (RANS) solvers. Then, the experimental model reproduced the numerical results with satisfying accuracy [32]. In another influential study, Wang et al. numerically assessed the effect of OSWEC’s width on their performance. According to their findings, as converter width increases, its efficiency decreases in short wave periods while increases in long wave periods [33]. One of the main challenges in the analysis of the OSWEC is the coupled effect of hydrodynamic and geometric variables. As a result, numerous cutting-edge geometry studies have been performed in recent years in order to find the optimal structure that maximizes power output and minimizes costs. Garcia et al. reviewed hull geometry optimization studies in the literature in [19]. In addition, Guo and Ringwood surveyed geometric optimization methods to improve the hydrodynamic performance of OSWECs at the primary stage [14]. Besides, they classified the hull geometry of OSWECs based on Figure 1. Subsequently, Whittaker et al. proposed a different design of OSWEC called Oyster2. There have been three examples of different geometries of oysters with different water depths. Based on its water depth, they determined the width and height of the converter. They also found that in the constant wave period the less the converter’s width, the less power captures the converter has [34]. Afterward, O’Boyle et al. investigated a type of OSWEC called Oyster 800. They compared the experimental and numerical models with the prototype model. In order to precisely reproduce the shape, mass distribution, and buoyancy properties of the prototype, a 40th-scale experimental model has been designed. Overall, all the models were fairly accurate according to the results [35].

Inclusive analysis of recent research avenues in the area of flap geometry has revealed that the interaction-based designs of such converters are emerging as a novel approach. An initiative workflow is designed in the current study to maximizing the wave energy extrication by such systems. To begin with, a sensitivity analysis plays its role of determining the best hydrodynamic values for installing the converter’s flap. Then, all flap dimensions and characteristics come into play to finalize the primary model. Following, interactive designs is proposed to increase the influence of incident waves on the body by adding ribs on both sides of the flap as a novel design. Finally, a new bi-level metaheuristic method is proposed to consider the effects of simultaneous changes in ribs properties and other design parameters. We hope this novel approach will be utilized to make big-scale projects less costly and justifiable. The efficiency of the method is also compared with four well known metaheuristic algorithms and out weight them for this application.

This paper is organized as follows. First, the research methodology is introduced by providing details about the numerical model implementation. To that end, we first introduced the primary model’s geometry and software details. That primary model is later verified with a benchmark study with regard to the flap angle of rotation and water surface elevation. Then, governing equations and performance criteria are presented. In the third part of the paper, we discuss the model’s sensitivity to lower and upper parts width (we proposed a two cross-sectional design for the flap), bottom elevation, and freeboard. Finally, the novel optimization approach is introduced in the final part and compared with four recent metaheuristic algorithms.

2. Numerical Methods

In this section, after a brief introduction of the numerical software, Flow3D, boundary conditions are defined. Afterwards, the numerical model implementation, along with primary model properties are described. Finally, governing equations, as part of numerical process, are discussed.

^2.1. ^{Model Setup}

FLOW-3D is a powerful and comprehensive CFD simulation platform for studying fluid dynamics. This software has several modules to solve many complex engineering problems. In addition, modeling complex flows is simple and effective using FLOW-3D’s robust meshing capabilities [36]. Interaction between fluid and moving objects might alter the computational range. Dynamic meshes are used in our modeling to take these changes into account. At each time step, the computational node positions change in order to adapt the meshing area to the moving object. In addition, to choose mesh dimensions, some factors are taken into account such as computational accuracy, computational time, and stability. The final grid size is selected based on the detailed procedure provided in [37]. To that end, we performed grid-independence testing on a CFD model using three different mesh grid sizes of 0.01, 0.015, and 0.02 meters. The problem geometry and boundary conditions were defined the same, and simulations were run on all three grids under the same conditions. The predicted values of the relevant variable, such as velocity, was compared between the grids. The convergence behavior of the numerical solution was analyzed by calculating the relative L2 norm error between two consecutive grids. Based on the results obtained, it was found that the grid size of 0.02 meters showed the least error, indicating that it provided the most accurate and reliable solution among the three grids. Therefore, the grid size of 0.02 meters was selected as the optimal spatial resolution for the mesh grid.

In this work, the flume dimensions are 10 meters long, 0.1 meters wide, and 2.2 meters high, which are shown in figure2. In addition, input waves with linear characteristics have a height of 0.1 meters and a period of 1.4 seconds. Among the linear wave methods included in this software, RNGk-ε and k- ε are appropriate for turbulence model. The research of Lopez et al. shows that RNGk- ε provides the most accurate simulation of turbulence in OSWECs [21]. We use CATIA software to create the flap primary model and other innovative designs for this project. The flap measures 0.1 m x 0.65 m x 0.360 m in x, y and z directions, respectively. In Figure 3, the primary model of flap and its dimensions are shown. In this simulation, five boundaries have been defined, including 1. Inlet, 2. Outlet, 3. Converter flap, 4. Bed flume, and 5. Water surface, which are shown in figure 2. Besides, to avoid wave reflection in inlet and outlet zones, Flow3D is capable of defining some areas as damping zones, the length of which has to be one to one and a half times the wavelength. Therefore, in the model, this length is considered equal to 2 meters. Furthermore, there is no slip in all the boundaries. In other words, at every single time step, the fluid velocity is zero on the bed flume, while it is equal to the flap velocity on the converter flap. According to the wave theory defined in the software, at the inlet boundary, the water velocity is called from the wave speed to be fed into the model.

^2.2. ^Verification

In the current study, we utilize the Schmitt experimental model as a benchmark for verification, which was developed at the Queen’s University of Belfast. The experiments were conducted on the flap of the converter, its rotation, and its interaction with the water surface. Thus, the details of the experiments are presented below based up on the experimental setup’s description [38]. In the experiment, the laboratory flume has a length of 20m and a width of 4.58m. Besides, in order to avoid incident wave reflection, a wave absorption source is devised at the end of the left flume. The flume bed, also, includes two parts with different slops. The flap position and dimensions of the flume can be seen in Figure4. In addition, a wave-maker with 6 paddles is installed at one end. At the opposite end, there is a beach with wire meshes. Additionally, there are 6 indicators to extract the water level elevation. In the flap model, there are three components: the fixed support structure, the hinge, and the flap. The flap measures 0.1m x 0.65m x 0.341m in x, y and z directions, respectively. In Figure5, the details are given [32]. The support structure consists of a 15 mm thick stainless steel base plate measuring 1m by 1.4m, which is screwed onto the bottom of the tank. The hinge is supported by three bearing blocks. There is a foam centerpiece on the front and back of the flap which is sandwiched between two PVC plates. Enabling changes of the flap, three metal fittings link the flap to the hinge. Moreover, in this experiment, the selected wave is generated based on sea wave data at scale 1:40. The wave height and the wave period are equal to 0.038 (m) and 2.0625 (s), respectively, which are tantamount to a wave with a period of 13 (s) and a height of 1.5 (m).

Two distinct graphs illustrate the numerical and experi-mental study results. Figure6 and Figure7 are denoting the angle of rotation of flap and surface elevation in computational and experimental models, respectively. The two figures roughly represent that the numerical and experimental models are a good match. However, for the purpose of verifying the match, we calculated the correlation coefficient (C) and root mean square error (RMSE). According to Figure6, correlation coefficient and RMSE are 0.998 and 0.003, respectively, and in Figure7 correlation coefficient and RMSE are respectively 0.999 and 0.001. Accordingly, there is a good match between the numerical and empirical models. It is worth mentioning that the small differences between the numerical and experimental outputs may be due to the error of the measuring devices and the calibration of the data collection devices.

Including continuity equation and momentum conserva- tion for incompressible fluid are given as [32], [39]:(1)

where P represents the pressure, g denotes gravitational acceleration, u represents fluid velocity, and Di is damping coefficient. Likewise, the model uses the same equation. to calculate the fluid velocity in other directions as well. Considering the turbulence, we use the two-equation model of RNGK- ε. These equations are:

(3)��t(��)+��(��)=��[�eff��]+��-��and(4)��(��)+��(��)=��[�eff��]+�1�∗��-��2��2�Where �2� and �1� are constants. In addition, �� and �� represent the turbulent Prandtl number of � and k, respectively.

�� also denote the production of turbulent kinetic energy of k under the effect of velocity gradient, which is calculated as follows:(5)��=�eff[��+��]��(6)�eff=�+��(7)�eff=�+��where � is molecular viscosity,�� represents turbulence viscosity, k denotes kinetic energy, and ∊∊ is energy dissipation rate. The values of constant coefficients in the two-equation RNGK ∊-∊ model is as shown in the Table 1 [40].Table 2.

Table 1. Constant coefficients in RNGK-∊ model

Factors	�	�0	�1	�2	��	��	��
Quantity	0.012	4.38	1.42	1.68	1.39	1.39	0.084

Table 2. Flap properties

Joint height (m)	0.476
Height of the center of mass (m)	0.53
Weight (Kg)	10.77

It is worth mentioning that the volume of fluid method is used to separate water and air phases in this software [41]. Below is the equation of this method [40].(8)��+��(��)=0where α and 1 − _α are portion of water phase and air phase, respectively. As a weighting factor, each fluid phase portion is used to determine the mixture properties. Finally, using the following equations, we calculate the efficiency of converters [42], [34], [43]:(9)�=14|�|2�+�2+(�+�a)2(�n2-�2)2where �� represents natural frequency, I denotes the inertia of OSWEC, Ia is the added inertia, F is the complex wave force, and B denotes the hydrodynamic damping coefficient. Afterward, the capture factor of the converter is calculated by [44]:(10)��=�1/2��2��gw where �� represents the capture factor, which is the total efficiency of device per unit length of the wave crest at each time step [15], �� represent the dimensional amplitude of the incident wave, w is the flap’s width, and Cg is the group velocity of the incident wave, as below:(11)��=��0·121+2�0ℎsinh2�0ℎwhere �0 denotes the wave number, h is water depth, and H is the height of incident waves.

According to previous sections ∊,��-∊ modeling is used for all models simulated in this section. For this purpose, the empty boundary condition is used for flume walls. In order to preventing wave reflection at the inlet and outlet of the flume, the length of wave absorption is set to be at least one incident wavelength. In addition, the structured mesh is chosen, and the mesh dimensions are selected in two distinct directions. In each model, all grids have a length of 2 (cm) and a height of 1 (cm). Afterwards, as an input of the software for all of the models, we define the time step as 0.001 (s). Moreover, the run time of every simulation is 30 (s). As mentioned before, our primary model is Schmitt model, and the flap properties is given in table2. For all simulations, the flume measures 15 meters in length and 0.65 meters in width, and water depth is equal to 0.335 (m). The flap is also located 7 meters from the flume’s inlet.

Finally, in order to compare the results, the capture factor is calculated for each simulation and compared to the primary model. It is worth mentioning that capture factor refers to the ratio of absorbed wave energy to the input wave energy.

According to primary model simulation and due to the decreasing horizontal velocity with depth, the wave crest has the highest velocity. Considering the fact that the wave’s orbital velocity causes the flap to move, the contact between the upper edge of the flap and the incident wave can enhance its performance. Additionally, the numerical model shows that the dynamic pressure decreases as depth increases, and the hydrostatic pressure increases as depth increases.

To determine the OSWEC design, it is imperative to understand the correlation between the capture factor, wave period, and wave height. Therefore, as it is shown in Figure8, we plot the change in capture factor over the variations in wave period and wave height in 3D and 2D. In this diagram, the first axis features changes in wave period, the second axis displays changes in wave height, and the third axis depicts changes in capture factor. According to our wave properties in the numerical model, the wave period and wave height range from 2 to 14 seconds and 2 to 8 meters, respectively. This is due to the fact that the flap does not oscillate if the wave height is less than 2 (m), and it does not reverse if the wave height is more than 8 (m). In addition, with wave periods more than 14 (s), the wavelength would be so long that it would violate the deep-water conditions, and with wave periods less than 2 (s), the flap would not oscillate properly due to the shortness of wavelength. The results of simulation are shown in Figure 8. As it can be perceived from Figure 8, in a constant wave period, the capture factor is in direct proportion to the wave height. It is because of the fact that waves with more height have more energy to rotate the flap. Besides, in a constant wave height, the capture factor increases when the wave period increases, until a given wave period value. However, the capture factor falls after this point. These results are expected since the flap’s angular displacement is not high in lower wave periods, while the oscillating motion of that is not fast enough to activate the power take-off system in very high wave periods.

As is shown in Figure 9, we plot the change in capture factor over the variations in wave period (s) and water depth (m) in 3D. As it can be seen in this diagram, the first axis features changes in water depth (m), the second axis depicts the wave period (s), and the third axis displays OSWEC’s capture factor. The wave period ranges from 0 to 10 seconds based on our wave properties, which have been adopted from Schmitt’s model, while water depth ranges from 0 to 0.5 meters according to the flume and flap dimensions and laboratory limitations. According to Figure9, for any specific water depth, the capture factor increases in a varying rate when the wave period increases, until a given wave period value. However, the capture factor falls steadily after this point. In fact, the maximum capture factor occurs when the wave period is around 6 seconds. This trend is expected since, in a specific water depth, the flap cannot oscillate properly when the wavelength is too short. As the wave period increases, the flap can oscillate more easily, and consequently its capture factor increases. However, the capture factor drops in higher wave periods because the wavelength is too large to move the flap. Furthermore, in a constant wave period, by changing the water depth, the capture factor does not alter. In other words, the capture factor does not depend on the water depth when it is around its maximum value.

³. ^{Sensitivity Analysis}

Based on previous studies, in addition to the flap design, the location of the flap relative to the water surface (freeboard) and its elevation relative to the flume bed (flap bottom elevation) play a significant role in extracting energy from the wave energy converter. This study measures the sensitivity of the model to various parameters related to the flap design including upper part width of the flap, lower part width of the flap, the freeboard, and the flap bottom elevation. Moreover, as a novel idea, we propose that the flap widths differ in the lower and upper parts. In Figure10, as an example, a flap with an upper thickness of 100 (mm) and a lower thickness of 50 (mm) and a flap with an upper thickness of 50 (mm) and a lower thickness of 100 (mm) are shown. The influence of such discrepancy between the widths of the upper and lower parts on the interaction between the wave and the flap, or in other words on the capture factor, is evaluated. To do so, other parameters are remained constant, such as the freeboard, the distance between the flap and the flume bed, and the wave properties.

In Figure11, models are simulated with distinct upper and lower widths. As it is clear in this figure, the first axis depicts the lower part width of the flap, the second axis indicates the upper part width of the flap, and the colors represent the capture factor values. Additionally, in order to consider a sufficient range of change, the flap thickness varies from half to double the value of the primary model for each part.

According to this study, the greater the discrepancy in these two parts, the lower the capture factor. It is on account of the fact that when the lower part of the flap is thicker than the upper part, and this thickness difference in these two parts is extremely conspicuous, the inertia against the motion is significant at zero degrees of rotation. Consequently, it is difficult to move the flap, which results in a low capture factor. Similarly, when the upper part of the flap is thicker than the lower part, and this thickness difference in these two parts is exceedingly noticeable, the inertia is so great that the flap can not reverse at the maximum degree of rotation. As the results indicate, the discrepancy can enhance the performance of the converter if the difference between these two parts is around 20%. As it is depicted in the Figure11, the capture factor reaches its own maximum amount, when the lower part thickness is from 5 to 6 (cm), and the upper part thickness is between 6 and 7 (cm). Consequently, as a result of this discrepancy, less material will be used, and therefore there will be less cost.

As illustrated in Figure12, this study examines the effects of freeboard (level difference between the flap top and water surface) and the flap bottom elevation (the distance between the flume bed and flap bottom) on the converter performance. In this diagram, the first axis demonstrates the freeboard and the second axis on the left side displays the flap bottom elevation, while the colors indicate the capture factor. In addition, the feasible range of freeboard is between -15 to 15 (cm) due to the limitation of the numerical model, so that we can take the wave slamming and the overtopping into consideration. Additionally, based on the Schmitt model and its scaled model of 1:40 of the base height, the flap bottom should be at least 9 (cm) high. Since the effect of surface waves is distributed over the depth of the flume, it is imperative to maintain a reasonable flap height exposed to incoming waves. Thus, the maximum flap bottom elevation is limited to 19 (cm). As the Figure12 pictures, at constant negative values of the freeboard, the capture factor is in inverse proportion with the flap bottom elevation, although slightly.

Furthermore, at constant positive values of the freeboard, the capture factor fluctuates as the flap bottom elevation decreases while it maintains an overall increasing trend. This is on account of the fact that increasing the flap bottom elevation creates turbulence flow behind the flap, which encumbers its rotation, as well as the fact that the flap surface has less interaction with the incoming waves. Furthermore, while keeping the flap bottom elevation constant, the capture factor increases by raising the freeboard. This is due to the fact that there is overtopping with adverse impacts on the converter performance when the freeboard is negative and the flap is under the water surface. Besides, increasing the freeboard makes the wave slam more vigorously, which improves the converter performance.

Adding ribs to the flap surface, as shown in Figure13, is a novel idea that is investigated in the next section. To achieve an optimized design for the proposed geometry of the flap, we determine the optimal number and dimensions of ribs based on the flap properties as our decision variables in the optimization process. As an example, Figure13 illustrates a flap with 3 ribs on each side with specific dimensions.

Figure14 shows the flow velocity field around the flap jointed to the flume bed. During the oscillation of the flap, the pressure on the upper and lower surfaces of the flap changes dynamically due to the changing angle of attack and the resulting change in the direction of fluid flow. As the flap moves upwards, the pressure on the upper surface decreases, and the pressure on the lower surface increases. Conversely, as the flap moves downwards, the pressure on the upper surface increases, and the pressure on the lower surface decreases. This results in a cyclic pressure variation around the flap. Under certain conditions, the pressure field around the flap can exhibit significant variations in magnitude and direction, forming vortices and other flow structures. These flow structures can affect the performance of the OSWEC by altering the lift and drag forces acting on the flap.

⁴. ^{Design Optimization}

We consider optimizing the design parameters of the flap of converter using a nature-based swarm optimization method, that fall in the category of metaheuristic algorithms [45]. Accordingly, we choose four state-of-the-art algorithms to perform an optimization study. Then, based on their performances to achieve the highest capture factor, one of them will be chosen to be combined with the Hill Climb algorithm to carry out a local search. Therefore, in the remainder of this section, we discuss the search process of each algorithm and visualize their performance and convergence curve as they try to find the best values for decision variables.

4.1. Metaheuristic Approaches

As the first considered algorithm, the Gray Wolf Optimizer (GWO) algorithm simulates the natural leadership and hunting performance of gray wolves which tend to live in colonies. Hunters must obey the alpha wolf, the leader, who is responsible for hunting. Then, the beta wolf is at the second level of the gray wolf hierarchy. A subordinate of alpha wolf, beta stands under the command of the alpha. At the next level in this hierarchy, there are the delta wolves. They are subordinate to the alpha and beta wolves. This category of wolves includes scouts, sentinels, elders, hunters, and caretakers. In this ranking, omega wolves are at the bottom, having the lowest level and obeying all other wolves. They are also allowed to eat the prey just after others have eaten. Despite the fact that they seem less important than others, they are really central to the pack survival. Since, it has been shown that without omega wolves, the entire pack would experience some problems like fighting, violence, and frustration. In this simulation, there are three primary steps of hunting including searching, surrounding, and finally attacking the prey. Mathematically model of gray wolves’ hunting technique and their social hierarchy are applied in determined by optimization. this study. As mentioned before, gray wolves can locate their prey and surround them. The alpha wolf also leads the hunt. Assuming that the alpha, beta, and delta have more knowledge about prey locations, we can mathematically simulate gray wolf hunting behavior. Hence, in addition to saving the top three best solutions obtained so far, we compel the rest of the search agents (also the omegas) to adjust their positions based on the best search agent. Encircling behavior can be mathematically modeled by the following equations: [46].(12)�→=|�→·��→(�)-�→(�)|(13)�→(�+1)=��→(�)-�→·�→(14)�→=2.�2→(15)�→=2�→·�1→-�→Where �→indicates the position vector of gray wolf, ��→ defines the vector of prey, t indicates the current iteration, and �→and �→are coefficient vectors. To force the search agent to diverge from the prey, we use �→ with random values greater than 1 or less than -1. In addition, C→ contains random values in the range [0,2], and �→ ₁ and �2→ are random vectors in [0,1]. The second considered technique is the Moth Flame Optimizer (MFO) algorithm. This method revolves around the moths’ navigation mechanism, which is realized by positioning themselves and maintaining a fixed angle relative to the moon while flying. This effective mechanism helps moths to fly in a straight path. However, when the source of light is artificial, maintaining an angle with the light leads to a spiral flying path towards the source that causes the moth’s death [47]. In MFO algorithm, moths and flames are both solutions. The moths are actual search agents that fly in hyper-dimensional space by changing their position vectors, and the flames are considered pins that moths drop when searching the search space [48]. The problem’s variables are the position of moths in the space. Each moth searches around a flame and updates it in case of finding a better solution. The fitness value is the return value of each moth’s fitness (objective) function. The position vector of each moth is passed to the fitness function, and the output of the fitness function is assigned to the corresponding moth. With this mechanism, a moth never loses its best solution [49]. Some attributes of this algorithm are as follows:

•It takes different values to converge moth in any point around the flame.
•Distance to the flame is lowered to be eventually minimized.
•When the position gets closer to the flame, the updated positions around the flame become more frequent.

As another method, the Multi-Verse Optimizer is based on a multiverse theory which proposes there are other universes besides the one in which we all live. According to this theory, there are more than one big bang in the universe, and each big bang leads to the birth of a new universe [50]. Multi-Verse Optimizer (MVO) is mainly inspired by three phenomena in cosmology: white holes, black holes, and wormholes. A white hole has never been observed in our universe, but physicists believe the big bang could be considered a white hole [51]. Black holes, which behave completely in contrast to white holes, attract everything including light beams with their extremely high gravitational force [52]. In the multiverse theory, wormholes are time and space tunnels that allow objects to move instantly between any two corners of a universe (or even simultaneously from one universe to another) [53]. Based on these three concepts, mathematical models are designed to perform exploration, exploitation, and local search, respectively. The concept of white and black holes is implied as an exploration phase, while the concept of wormholes is considered as an exploitation phase by MVO. Additionally, each solution is analogous to a universe, and each variable in the solution represents an object in that universe. Furthermore, each solution is assigned an inflation rate, and the time is used instead of iterations. Following are the universe rules in MVO:

•The possibility of having white hole increases with the inflation rate.
•The possibility of having black hole decreases with the inflation rate.
•Objects tend to pass through black holes more frequently in universes with lower inflation rates.
•Regardless of inflation rate, wormholes may cause objects in universes to move randomly towards the best universe. [54]

Modeling the white/black hole tunnels and exchanging objects of universes mathematically was accomplished by using the roulette wheel mechanism. With every iteration, the universes are sorted according to their inflation rates, then, based on the roulette wheel, the one with the white hole is selected as the local extremum solution. This is accomplished through the following steps:

Assume that

(16)��=��1<��(��)��1≥��(��)

Where �� represents the jth parameter of the ith universe, Ui indicates the ith universe, NI(Ui) is normalized inflation rate of the ith universe, r1 is a random number in [0,1], and j x_k shows the jth parameter of the kth universe selected by a roulette wheel selection mechanism [54]. It is assumed that wormhole tunnels always exist between a universe and the best universe formed so far. This mechanism is as follows:(17)��=if�2<��:��+��×((��-��)×�4+��)�3<0.5��-��×((��-��)×�4+��)�3≥0.5��:��where Xj indicates the jth parameter of the best universe formed so far, TDR and WEP are coefficients, where Xj indicates the jth parameter of the best universelbjshows the lower bound of the jth variable, ubj is the upper bound of the jth variable, and r2, r3, and r4 are random numbers in [1], [54].

Finally, one of the newest optimization algorithms is WOA. The WOA algorithm simulates the movement of prey and the whale’s discipline when looking for their prey. Among several species, Humpback whales have a specific method of hunting [55]. Humpback whales can recognize the location of prey and encircle it before hunting. The optimal design position in the search space is not known a priori, and the WOA algorithm assumes that the best candidate solution is either the target prey or close to the optimum. This foraging behavior is called the bubble-net feeding method. Two maneuvers are associated with bubbles: upward spirals and double loops. A unique behavior exhibited only by humpback whales is bubble-net feeding. In fact, The WOA algorithm starts with a set of random solutions. At each iteration, search agents update their positions for either a randomly chosen search agent or the best solution obtained so far [56], [55]. When the best search agent is determined, the other search agents will attempt to update their positions toward that agent. It is important to note that humpback whales swim around their prey simultaneously in a circular, shrinking circle and along a spiral-shaped path. By using a mathematical model, the spiral bubble-net feeding maneuver is optimized. The following equation represents this behavior:(18)�→(�+1)=�′→·�bl·cos(2��)+�∗→(�)

Where:(19)�′→=|�∗→(�)-�→(�)|

X→_(t_{+ 1)} indicates the distance of the it h whale to the prey (best solution obtained so far),� is a constant for defining the shape of the logarithmic spiral, l is a random number in [−1_, _1], and dot (.) is an element-by-element multiplication [55].

Comparing the four above-mentioned methods, simulations are run with 10 search agents for 400 iterations. In Figure 15, there are 20 plots the optimal values of different parameters in optimization algorithms. The five parameters of this study are freeboard, bottom elevations, number of ribs on the converter, rib thickness, and rib Height. The optimal value for each was found by optimization algorithms, naming WOA, MVO, MFO, and GWO. By looking through the first row, the freeboard parameter converges to its maximum possible value in the optimization process of GWO after 300 iterations. Similarly, MFO finds the same result as GWO. In contrast, the freeboard converges to its minimum possible value in MVO optimizing process, which indicates positioning the converter under the water. Furthermore, WOA found the optimal value of freeboard as around 0.02 after almost 200 iterations. In the second row, the bottom elevation is found at almost 0.11 (m) in all algorithms; however, the curves follow different trends in each algorithm. The third row shows the number of ribs, where results immediately reveal that it should be over 4. All algorithms coincide at 5 ribs as the optimal number in this process. The fourth row displays the trends of algorithms to find optimal rib thickness. MFO finds the optimal value early and sets it to around 0.022, while others find the same value in higher iterations. Finally, regarding the rib height, MVO, MFO, and GWO state that the optimal value is 0.06 meters, but WOA did not find a higher value than 0.039.

4.2. HCMVO Bi-level Approach

Despite several strong search characteristics of MVO and its high performance in various optimization problems, it suffers from a few deficiencies in local and global search mechanisms. For instance, it is trapped in the local optimum when wormholes stochastically generate many solutions near the best universe achieved throughout iterations, especially in solving complex multimodal problems with high dimensions [57]. Furthermore, MVO needs to be modified by an escaping strategy from the local optima to enhance the global search abilities. To address these shortages, we propose a fast and effective meta-algorithm (HCMVO) to combine MVO with a Random-restart hill-climbing local search. This meta-algorithm uses MVO on the upper level to develop global tracking and provide a range of feasible and proper solutions. The hill-climbing algorithm is designed to develop a comprehensive neighborhood search around the best-found solution proposed by the upper-level (MVO) when MVO is faced with a stagnation issue or falling into a local optimum. The performance threshold is formulated as follows.(20)Δ��THD=∑�=1��TH��-��TH��-1�where BestTHDis the best-found solution per generation, andM is related to the domain of iterations to compute the average performance of MVO. If the proposed best solution by the local search is better than the initial one, the global best of MVO will be updated. HCMVO iteratively runs hill climbing when the performance of MVO goes down, each time with an initial condition to prepare for escaping such undesirable situations. In order to get a better balance between exploration and exploitation, the search step size linearly decreases as follows:(21)��=��-��Ma�iter��+1where iter and Maxiter are the current iteration and maximum number of evaluation, respectively. �� stands for the step size of the neighborhood search. Meanwhile, this strategy can improve the convergence rate of MVO compared with other algorithms.

Algorithm 1 shows the technical details of the proposed optimization method (HCMVO). The initial solution includes freeboard (�), bottom elevation (�), number of ribs (Nr), rib thickness (�), and rib height(�).

5. Conclusion

The high trend of diminishing worldwide energy resources has entailed a great crisis upon vulnerable societies. To withstand this effect, developing renewable energy technologies can open doors to a more reliable means, among which the wave energy converters will help the coastal residents and infrastructure. This paper set out to determine the optimized design for such devices that leads to the highest possible power output. The main goal of this research was to demonstrate the best design for an oscillating surge wave energy converter using a novel metaheuristic optimization algorithm. In this regard, the methodology was devised such that it argued the effects of influential parameters, including wave characteristics, WEC design, and interaction criteria.

To begin with, a numerical model was developed in Flow 3D software to simulate the response of the flap of a wave energy converter to incoming waves, followed by a validation study based upon a well-reputed experimental study to verify the accuracy of the model. Secondly, the hydrodynamics of the flap was investigated by incorporating the turbulence. The effect of depth, wave height, and wave period are also investigated in this part. The influence of two novel ideas on increasing the wave-converter interaction was then assessed: i) designing a flap with different widths in the upper and lower part, and ii) adding ribs on the surface of the flap. Finally, four trending single-objective metaheuristic optimization methods

Empty Cell	Algorithm 1: Hill Climb Multiverse Optimization
01:	procedure HCMVO
02:	�=30,�=5▹��
03:	�=〈F1,B1,N,R,H1〉,…〈FN,B2,N,R,HN〉⇒lb1N⩽�⩽ubN
04:	Initialize parameters�ER,�DR,�EP,Best�,��ite��▹Wormhole existence probability (WEP)
05:	��=��(��)
06:	��=Normalize the inflation rate��
07:	for iter in[1,⋯,��iter]do
08:	for�in[1,⋯,�]do
09:	Update�EP,�DR,Black��Index=�
10:	for��[1,⋯,�]��
11:	�1=��()
12:	if�1≤��(��)then
13:	White HoleIndex=Roulette�heelSelection(-��)
14:	�(Black HoleIndex,�)=��(White HoleIndex,�)
15:	*end if*
16:	�2=��([0,�])
17:	if�2≤�EPthen
18:	�3=��(),�4=��()
19:	if�3<0.5then
20:	�1=((��(�)-��(�))×�4+��(�))
21:	�(�,�)=Best�(�)+�DR×�
22:	*else*
23:	�(�,�)=Best�(�)-�DR×�
24:	*end if*
25:	*end if*
26:	*end for*
27:	*end for*
28:	�HD=��([�1,�2,⋯,�Np])
29:	Bes�TH�itr=��HD
30:	ΔBestTHD=∑�=1�BestTII��-BestTII��-1�
31:	ifΔBestTHD<��then▹Perform hill climbing local search
32:	BestTHD=��-�lim��THD
33:	*end if*
34:	*end for*
35:	return�,BestTHD▹Final configuration
36:	end procedure

The implementation details of the hill-climbing algorithm applied in HCMPA can be seen in Algorithm 2. One of the critical parameters isg, which denotes the resolution of the neighborhood search around the proposed global best by MVO. If we set a small step size for hill-climbing, the convergence speed will be decreased. On the other hand, a large step size reinforces the exploration ability. Still, it may reduce the exploitation ability and in return increase the act of jumping from a global optimum or surfaces with high-potential solutions. Per each decision variable, the neighborhood search evaluates two different direct searches, incremental or decremental. After assessing the generated solutions, the best candidate will be selected to iterate the search algorithm. It is noted that the hill-climbing algorithm should not be applied in the initial iteration of the optimization process due to the immense tendency for converging to local optima. Meanwhile, for optimizing largescale problems, hill-climbing is not an appropriate selection. In order to improve understanding of the proposed hybrid optimization algorithm’s steps, the flowchart of HCMVO is designed and can be seen in Figure 16.

Figure 17 shows the observed capture factor (which is the absorbed energy with respect to the available energy) by each optimization algorithm from iterations 1 to 400. The algorithms use ten search agents in their modified codes to find the optimal solutions. While GWO and MFO remain roughly constant after iterations 54 and 40, the other three algorithms keep improving the capture factor. In this case, HCMVO and MVO worked very well in the optimizing process with a capture factor obtained by the former as 0.594 and by the latter as 0.593. MFO almost found its highest value before the iteration 50, which means the exploration part of the algorithm works out well. Similarly, HCMVO does the same. However, it keeps finding the better solution during the optimization process until the last iteration, indicating the strong exploitation part of the algorithm. GWO reveals a weakness in exploration and exploitation because not only does it evoke the least capture factor value, but also the curve remains almost unchanged throughout 350 iterations.

Figure 18 illustrates complex interactions between the five optimization parameters and the capture factor for HCMVO (a), MPA (b), and MFO (c) algorithms. The first interesting observation is that there is a high level of nonlinear relationships among the setting parameters that can make a multi-modal search space. The dark blue lines represent the best-found configuration throughout the optimisation process. Based on both HCMVO (a) and MVO (b), we can infer that the dark blue lines concentrate in a specific range, showing the high convergence ability of both HCMVO and MVO. However, MFO (c) could not find the exact optimal range of the decision variables, and the best-found solutions per generation distribute mostly all around the search space.

Empty Cell	Algorithm 1: Hill Climb Multiverse Optimization
01:	procedure HCMVO
02:	Initialization
03:	Initialize the constraints��1�,��1�
04:	�1�=Mi�1�+��1�/�▹Compute the step size,�is search resolution
05:	So�1=〈�,�,�,�,�〉▹��
06:	��1=��So�1▹��ℎ��
07:	Main loop
08:	for iter≤��ita=do
09:	��=��±��
10:	while�≤��(Sol1)do
11:	��=��+�,▹��ℎ��ℎ��ℎ
12:	fitness��iter=��
13:	t = t+1
14:	end while
15:	〈��,��max〉=��
16:	��itev=��Inde�max▹��ℎ��ℎ��
17:	��=��-��Max��+1▹��
18:	end for
19:	return��iter,��
20:	end procedure

were utilized to illuminate the optimum values of the design parameters, and the best method was chosen to develop a new algorithm that performs both local and global search methods.

The correlation between hydrodynamic parameters and the capture factor of the converter was supported by the results. For any given water depth, the capture factor increases as the wave period increases, until a certain wave period value (6 seconds) is reached, after which the capture factor gradually decreases. It is expected since the flap cannot oscillate effectively when the wavelength is too short for a certain water depth. Conversely, when the wavelength is too long, the capture factor decreases. Furthermore, under a constant wave period, increasing the water depth does not affect the capture factor. Regarding the sensitivity analysis, the study found that increasing the flap bottom elevation causes turbulence flow behind the flap and limitation of rotation, which leads to less interaction with the incoming waves. Furthermore, while keeping the flap bottom elevation constant, increasing the freeboard improves the capture factor. Overtopping happens when the freeboard is negative and the flap is below the water surface, which has a detrimental influence on converter performance. Furthermore, raising the freeboard causes the wave impact to become more violent, which increases converter performance.

In the last part, we discussed the search process of each algorithm and visualized their performance and convergence curves as they try to find the best values for decision variables. Among the four selected metaheuristic algorithms, the Multi-verse Optimizer proved to be the most effective in achieving the best answer in terms of the WEC capture factor. However, the MVO needed modifications regarding its escape approach from the local optima in order to improve its global search capabilities. To overcome these constraints, we presented a fast and efficient meta-algorithm (HCMVO) that combines MVO with a Random-restart hill-climbing local search. On a higher level, this meta-algorithm employed MVO to generate global tracking and present a range of possible and appropriate solutions. Taken together, the results demonstrated that there is a significant degree of nonlinearity among the setup parameters that might result in a multimodal search space. Since MVO was faced with a stagnation issue or fell into a local optimum, we constructed a complete neighborhood search around the best-found solution offered by the upper level. In sum, the newly-developed algorithm proved to be highly effective for the problem compared to other similar optimization methods. The strength of the current findings may encourage future investigation on design optimization of wave energy converters using developed geometry as well as the novel approach.

CRediT authorship contribution statement

Erfan Amini: Conceptualization, Methodology, Validation, Data curation, Writing – original draft, Writing – review & editing, Visualization. Mahdieh Nasiri: Conceptualization, Methodology, Validation, Data curation, Writing – original draft, Writing – review & editing, Visualization. Navid Salami Pargoo: Writing – original draft, Writing – review & editing. Zahra Mozhgani: Conceptualization, Methodology. Danial Golbaz: Writing – original draft. Mehrdad Baniesmaeil: Writing – original draft. Meysam Majidi Nezhad: . Mehdi Neshat: Supervision, Conceptualization, Writing – original draft, Writing – review & editing, Visualization. Davide Astiaso Garcia: Supervision. Georgios Sylaios: Supervision.

Declaration of Competing Interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Acknowledgement

This research has been carried out within ILIAD (Inte-grated Digital Framework for Comprehensive Maritime Data and Information Services) project that received funding from the European Union’s H2020 programme.

Data availability

Data will be made available on request.

References

[1]F. Ahmadpoor, P. SharmaA perspective on the statistical mechanics of 2d materialsExtreme Mechanics Letters, 14 (2017), pp. 38-43View PDF View article View in Scopus Google Scholar
[2]A. Morteza, M. SadipourReza Saadati Fard, Saman Taheri, and Amirhossein AhmadiA dagging-based deep learning framework for transmission line flexibility assessment, IET Renewable Power Generation (2022)Google Scholar
[3]A. Morteza, M. Ilbeigi, J. SchwedA blockchain information management framework for construction safety. Comput-ingCivil Engineering (2021, 2022.)Google Scholar
[4]Jochem Weber, Ronan Costello, and John Ringwood. Wec technology performance levels (tpls)-metric for successful development of economic wec technology. Proceedings EWTEC 2013, 2013.Google Scholar
[5]K. Rahgooy, A. Bahmanpour, M. Derakhshandi, A.a. Bagherzadeh-KhalkhaliDistribution of elastoplastic modulus of subgrade reaction for analysis of raft foundationsGeomechanics and Engineering, 28 (1) (2022), pp. 89-105View in Scopus Google Scholar
[6]EugenRusuandFlorinOnea. Estimationofthewaveenergyconversion efficiencyintheatlanticoceanclosetotheeuropeanislands. Renewable Energy, 85:687–703, 2016.Google Scholar
[7]M. Penalba, G. Giorgi, J.V. RingwoodMathematical modelling of wave energy converters: A review of nonlinear approachesRenewable and Sustainable Energy Reviews, 78 (2017), pp. 1188-1207View PDF View article View in Scopus Google Scholar
[8]C. Windt, J. Davidson, J.V. RingwoodHigh-fidelity numerical modelling of ocean wave energy systems: A review of computational fluid dynamics-based numerical wave tanksRenewable and Sustainable Energy Reviews, 93 (2018), pp. 610-630View PDF View article View in Scopus Google Scholar
[9]W. ShengWave energy conversion and hydrodynamics modelling technologies: A reviewRenewable and Sustainable Energy Reviews, 109 (2019), pp. 482-498View PDF View article View in Scopus Google Scholar
[10]J.V. Ringwood, G. Bacelli, F. FuscoEnergymaximizing control of wave-energy converters: The development of control system technology to optimize their operationIEEE control systems magazine, 34 (5) (2014), pp. 30-55View article CrossRef View in Scopus Google Scholar
[11]N. Faedo, S. Olaya, J.V. RingwoodOptimal control, mpc and mpc-like algorithms for wave energy systems: An overviewIFAC Journal of Systems and Control, 1 (2017), pp. 37-56View PDF View article View in Scopus Google Scholar
[12]L. Wang, J. Isberg, E. TedeschiReview of control strategies for wave energy conversion systems and their validation: the wave-to-wire approachRenewable and Sustainable Energy Reviews, 81 (2018), pp. 366-379View PDF View article View in Scopus Google Scholar
[13]J.V. RingwoodWave energy control: status and perspectives 2020IFAC-PapersOnLine, 53 (2) (2020), pp. 12271-12282View PDF View article View in Scopus Google Scholar
[14]B. Guo, J.V. RingwoodGeometric optimisation of wave energy conversion devices: A surveyApplied Energy, 297 (2021), Article 117100View PDF View article View in Scopus Google Scholar
[15]E. Amini, D. Golbaz, R. Asadi, M. Nasiri, O. Ceylan, M.M. Nezhad, et al.A comparative study of metaheuristic algorithms for wave energy converter power take-off optimisation: A case study for eastern australiaJournal of Marine Science and Engineering, 9(5):490 (2021)Google Scholar
[16]Arthur Pecher and Jens Peter KofoedHandbook of ocean wave energySpringer Nature (2017)Google Scholar
[17]G. Chang, C.A. Jones, J.D. Roberts, V.S. NearyA comprehensive evaluation of factors affecting the levelized cost of wave energy conversion projectsRenewable Energy, 127 (2018), pp. 344-354View PDF View article View in Scopus Google Scholar
[18]E. Amini, H. Mehdipour, E. Faraggiana, D. Golbaz, S. Mozaffari, G. Bracco, et al.Optimization of hydraulic power take-off system settings for point absorber wave energy converterRenewable Energy, 194 (2022), pp. 938-954View PDF View article View in Scopus Google Scholar
[19]A. Garcia-Teruel, D.I.M. ForehandA review of geometry optimisa-tion of wave energy convertersRenewable and Sustainable Energy Reviews, 139 (2021), Article 110593View PDF View article View in Scopus Google Scholar
[20]M.M. Nezhad, A. Heydari, M. Neshat, F. Keynia, G. Piras, D.A. GarciaA mediterranean sea offshore wind classification using merra-2 and machine learning modelsRenewable Energy, 190 (2022), pp. 156-166Google Scholar
[21]I. López, J. Andreu, S. Ceballos, I.M.D. Alegría, I. KortabarriaReview of wave energy technologies and the necessary power-equipmentRenewable and sustainable energy reviews, 27 (2013), pp. 413-434View PDF View article View in Scopus Google Scholar
[22]R. Ekström, B. Ekergård, M. LeijonElectrical damping of linear generators for wave energy converters—a reviewRenewable and Sustainable Energy Reviews, 42 (2015), pp. 116-128View PDF View article Google Scholar
[23]Danial Golbaz, Rojin Asadi, Erfan Amini, Hossein Mehdipour, Mahdieh Nasiri, Meysam Majidi Nezhad, Seyed Taghi Omid Naeeni, and Mehdi Neshat. Ocean wave energy converters optimization: A comprehensive review on research directions. arXiv preprint arXiv:2105.07180, 2021.Google Scholar
[24]Michael Choiniere, Jacob Davis, Nhu Nguyen, Nathan Tom, Matthew Fowler, and Krish Thiagarajan Sharman. Hydrodynamics and load shedding behavior of a variable geometry oscillating surge wave energy converter (oswec). Available at SSRN 3900951, 2022.Google Scholar
[25]Alan Henry, Olivier Kimmoun, Jonathan Nicholson, Guillaume Dupont, Yanji Wei, andFrederic Dias. A two dimensional experimental investigation of slamming of an oscillating wave surge converter. In The Twenty-fourth International Ocean and Polar Engineering Conference. OnePetro, 2014.Google Scholar
[26]S. Doyle, G.A. AggidisDevelopment of multioscillating water columns as wave energy convertersRenewable and Sustainable Energy Reviews, 107 (2019), pp. 75-86View PDF View article View in Scopus Google Scholar
[27]S. Doyle, G.A. AggidisExperimental investigation and performance comparison of a 1 single owc, array and m-owcRenewable Energy, 168 (2021), pp. 365-374View PDF View article View in Scopus Google Scholar
[28]Matthew Folley, TJT Whittaker, and Alan Henry. The effect of water depth on the performance of a small surging wave energy converter. Ocean Engineering, 34(8-9):1265–1274, 2007.Google Scholar
[29]Q. Li, J. Mi, X. Li, S. Chen, B. Jiang, L. ZuoA self-floating oscillating surge wave energy converterEnergy, 230 (2021), Article 120668View PDF View article View in Scopus Google Scholar
[30]D. Sarkar, E. Renzi, F. DiasEffect of a straight coast on the hydrodynamics and performance of the oscillating wave surge converterOcean Engineering, 105 (2015), pp. 25-32View PDF View article View in Scopus Google Scholar
[31]Adrian de Andres, Jéromine Maillet, Jørgen Hals Todalshaug, Patrik Möller, and Henry Jeffrey. On the optimum sizing of a real wec from a techno-economic perspective. In International Conference on Offshore Mechanics and Arctic Engineering, volume 49972, page V006T09A013. American Society of Mechanical Engineers, 2016.Google Scholar
[32]P. Schmitt, B. ElsaesserOn the use of openfoam to model oscillating wave surge convertersOcean Engineering, 108 (2015), pp. 98-104View PDF View article View in Scopus Google Scholar
[33]D.-J. Wang, S.-Q. Qiu, J.-W. YeWidth effects on hydrodynamics of pendulum wave energy converterApplied Mathematics and Mechanics, 35 (9) (2014), pp. 1167-1176View article CrossRef View in Scopus Google Scholar
[34]T. Whittaker, M. FolleyNearshore oscillating wave surge converters and the development of oysterPhilosophical Transactions Sciences of the Royal Society A: Mathematical, Physical and Engineering, 370 (1959) (2012), pp. 345-364View article CrossRef View in Scopus Google Scholar
[35]Louise O’Boyle, Kenneth Doherty, Jos van’t Hoff, and Jessica Skelton. The value of full scale prototype data-testing oyster 800 at emec, orkney. In Proceedings of the 11th European wave and tidal energy conference (_{EWTEC), Nantes, France}, pages 6–11, 2015.Google Scholar
[36]Why flow-3d? https://www.flow3d.com/products/flow-3d_/ why-flow-3d_/, 2022.Google Scholar
[37]Ishmail B Celik, Urmila Ghia, Patrick J Roache, and Christopher J Freitas. Procedure for estimation and reporting of uncertainty due to discretization in cfd applications. Journal of fluids EngineeringTransactions of the ASME, 130(7), 2008.Google Scholar
[38]Pal Schmitt, K Doherty, Darragh Clabby, and T Whittaker. The opportunities and limitations of using cfd in the development of wave energyconverters. Marine&OffshoreRenewableEnergy, pages 89–97, 2012.Google Scholar
[39]M. Choiniere, J. Davis, N.u. Nguyen, N. Tom, M. Fowler, K. ThiagarajanHydrodynamics and load shedding behavior of a variable-geometry oscillating surge wave energy converter (oswec)Renewable Energy (2022)Google Scholar
[40]Yong Li and Mian Lin. Regular and irregular wave impacts on floating body. Ocean Engineering, 42:93–101, 2012. Pal Manuel Schmitt. Investigation of the near flow field of bottom hinged flap type wave energy converters. PhD thesis, Queen’s University Belfast, 2014.Google Scholar
[41]Alan John Henry. The hydrodynamics of small seabed mounted bottom hinged wave energy conerverters in shallow water. PhD thesis, Queen’s University Belfast, 2009.Google Scholar
[42]N. Ghorbani, A. Korzeniowski, et al.Adaptive risk hedging for call options under cox-ingersoll-ross interest ratesJournal of Mathematical Finance, 10 (04) (2020), p. 697 View PDF CrossRef View in Scopus Google Scholar
[43]E. Renzi, F. DiasRelations for a periodic array of flap-type wave energy convertersApplied Ocean Research, 39 (2013), pp. 31-39View PDF View article View in Scopus Google Scholar
[44]M. Abdel-Basset, L. Abdel-Fatah, A.K. SangaiahChapter 10metaheuristic algorithms: a comprehensive reviewcomputational intelligence for multimedia big data on the cloud with engineering applications (2018)Google Scholar
[45]S. Mirjalili, S.M. Mirjalili, A. LewisGrey wolf optimizerAdvances in engineering software, 69 (2014), pp. 46-61View PDF View article Google Scholar
[46]S. MirjaliliMoth-flame optimization algorithm: A novel natureinspired heuristic paradigmKnowledge-based systems, 89 (2015), pp. 228-249View PDF View article View in Scopus Google Scholar
[47]Mohammad Shehab, Laith Abualigah, Husam Al Hamad, Hamzeh Alabool, Mohammad Alshinwan, and Ahmad M Khasawneh. Moth– flame optimization algorithm: variants and applications. Neural Computing and Applications, 32(14):9859–9884, 2020.Google Scholar
[48]Betül Sultan Yıldız and Ali Rıza YıldızMoth-flame optimization algorithm to determine optimal machining parameters in manufacturing processesMaterials Testing, 59 (5) (2017), pp. 425-429Google Scholar
[49]M Tegmark. Barrow, jd davies, pc harper, cl, jr eds. Science and Ultimate Reality Cambridge University Press Cambridge, 2004.Google Scholar
[50]DouglasMEardley. Deathofwhiteholesintheearlyuniverse. Physical Review Letters, 33(7):442, 1974.Google Scholar
[51]P.CW. DaviesThermodynamics of black holesReportsonProgress in Physics, 41 (8) (1978), p. 1313View article CrossRef View in Scopus Google Scholar
[52]M.S. Morris, K.S. ThorneWormholes in spacetime and their use for interstellar travel: A tool for teaching general relativityAmerican Journal of Physics, 56 (5) (1988), pp. 395-412View article CrossRef View in Scopus Google Scholar
[53]S. Mirjalili, S.M. Mirjalili, A. HatamlouMulti-verse optimizer: a nature-inspired algorithm for global optimizationNeural Computing and Applications, 27 (2) (2016), pp. 495-513View article CrossRef View in Scopus Google Scholar
[54]S. Mirjalili, A. LewisThe whale optimizationalgorithmAdvances in engineering software, 95 (2016), pp. 51-67View PDF View article View in Scopus Google Scholar
[55]Farhad Soleimanian Gharehchopogh and Hojjat GholizadehA comprehensive survey: Whale optimization algorithm and its applicationsSwarm and Evolutionary Computation, 48 (2019), pp. 1-24Google Scholar
[56]L. AbualigahMulti-verse optimizer algorithm: a comprehensive survey of its results, variants, and applicationsNeural Computing and Applications, 32 (16) (2020), pp. 12381-12401View article CrossRef View in Scopus Google Scholar

원문보기

PDF View

Figure 15. Velocity distribution of impinging jet on a wall under different Reynolds numbers.

Hydraulic Characteristics of Continuous Submerged Jet Impinging on a Wall by Using Numerical Simulation and PIV Experiment

by Hongbo Mi^1,2, Chuan Wang^1,3, Xuanwen Jia^3,*, Bo Hu², Hongliang Wang⁴, Hui Wang³ and Yong Zhu⁵

¹College of Mechatronics Engineering, Hainan Vocational University of Science and Technology, Haikou 571126, China

²Department of Energy and Power Engineering, Tsinghua University, Beijing 100084, China

³College of Hydraulic Science and Engineering, Yangzhou University, Yangzhou 225009, China

⁴School of Aerospace and Mechanical Engineering/Flight College, Changzhou Institute of Technology, Changzhou 213032, China

⁵National Research Center of Pumps, Jiangsu University, Zhenjiang 212013, China

^*Author to whom correspondence should be addressed.Sustainability2023, 15(6), 5159; https://doi.org/10.3390/su15065159

Received: 30 January 2023 / Revised: 4 March 2023 / Accepted: 10 March 2023 / Published: 14 March 2023(This article belongs to the Special Issue Advanced Technologies of Renewable Energy and Water Management for Sustainable Environment

Abstract

Due to their high efficiency, low heat loss and associated sustainability advantages, impinging jets have been used extensively in marine engineering, geotechnical engineering and other engineering practices. In this paper, the flow structure and impact characteristics of impinging jets with different Reynolds numbers and impact distances are systematically studied by Flow-3D based on PIV experiments. In the study, the relevant state parameters of the jets are dimensionlessly treated, obtaining not only the linear relationship between the length of the potential nucleation zone and the impinging distance, but also the linear relationship between the axial velocity and the axial distance in the impinging zone. In addition, after the jet impinges on the flat plate, the vortex action range caused by the wall-attached flow of the jet gradually decreases inward with the increase of the impinging distance. By examining the effect of Reynolds number Re on the hydraulic characteristics of the submerged impact jet, it can be found that the structure of the continuous submerged impact jet is relatively independent of the Reynolds number. At the same time, the final simulation results demonstrate the applicability of the linear relationship between the length of the potential core region and the impact distance. This study provides methodological guidance and theoretical support for relevant engineering practice and subsequent research on impinging jets, which has strong theoretical and practical significance.

Keywords:

PIV; Flow-3D; impinging jet; hydraulic characteristics; impinging distance

원문 보기

Figure 1. Geometric model.

Figure 2. Model grid schematic.

Figure 3. (a) Schematic diagram of the experimental setup; (b) PIV images of vertical impinging jets with velocity fields.

Figure 4. (a) Velocity distribution verification at the outlet of the jet pipe; (b) Distribution of flow angle in the mid-axis of the jet [39].

Figure 5. Along-range distribution of the dimensionless axial velocity of the jet at different impact distances.Figure 6 shows the variation of H

Figure 6. Relationship between the distribution of potential core region and the impact height H/D.

Figure 7. The relationship between the potential core length

Figure 8. Along-range distribution of the flow angle φ of the jet at different impact distances.

Figure 9. Velocity distribution along the axis of the jet at different impinging regions.

Figure 10. The absolute value distribution of slope under different impact distances.

Figure 11. Velocity distribution of impinging jet on wall under different impinging distances.

Figure 12. Along-range distribution of the dimensionless axial velocity of the jet at different Reynolds numbers.

Figure 13. Along-range distribution of the flow angle φ of the jet at different Reynolds numbers.

Figure 14. Velocity distribution along the jet axis at different Reynolds numbers.

Figure 15. Velocity distribution of impinging jet on a wall under different Reynolds numbers.

References

Zhang, J.; Li, Y.; Zhang, Y.; Yang, F.; Liang, C.; Tan, S. Using a high-pressure water jet-assisted tunnel boring machine to break rock. Adv. Mech. Eng. 2020, 12, 1687814020962290. [Google Scholar] [CrossRef]
Shi, X.; Zhang, G.; Xu, G.; Ma, Y.; Wu, X. Inactivating Microorganism on Medical Instrument Using Plasma Jet. High Volt. Eng. 2009, 35, 632–635. [Google Scholar]
Gao, Y.; Han, P.; Wang, F.; Cao, J.; Zhang, S. Study on the Characteristics of Water Jet Breaking Coal Rock in a Drilling Hole. Sustainability 2022, 14, 8258. [Google Scholar] [CrossRef]
Xu, W.; Wang, C.; Zhang, L.; Ge, J.; Zhang, D.; Gao, Z. Numerical study of continuous jet impinging on a rotating wall based on Wray—Agarwal turbulence model. J. Braz. Soc. Mech. Sci. Eng. 2022, 44, 433. [Google Scholar] [CrossRef]
Hu, B.; Wang, C.; Wang, H.; Yu, Q.; Liu, J.; Zhu, Y.; Ge, J.; Chen, X.; Yang, Y. Numerical Simulation Study of the Horizontal Submerged Jet Based on the Wray—Agarwal Turbulence Model. J. Mar. Sci. Eng. 2022, 10, 1217. [Google Scholar] [CrossRef]
Dahiya, A.K.; Bhuyan, B.K.; Kumar, S. Perspective study of abrasive water jet machining of composites—A review. J. Mech. Sci. Technol. 2022, 36, 213–224. [Google Scholar] [CrossRef]
Abushanab, W.S.; Moustafa, E.B.; Harish, M.; Shanmugan, S.; Elsheikh, A.H. Experimental investigation on surface characteristics of Ti6Al4V alloy during abrasive water jet machining process. Alex. Eng. J. 2022, 61, 7529–7539. [Google Scholar] [CrossRef]
Hu, B.; Wang, H.; Liu, J.; Zhu, Y.; Wang, C.; Ge, J.; Zhang, Y. A numerical study of a submerged water jet impinging on a stationary wall. J. Mar. Sci. Eng. 2022, 10, 228. [Google Scholar] [CrossRef]
Peng, J.; Shen, H.; Xie, W.; Zhai, S.; Xi, G. Influence of flow fluctuation characteristics on flow and heat transfer in different regions. J. Drain. Irrig. Mach. Eng. 2022, 40, 826–833. [Google Scholar]
Zhai, S.; Xie, F.; Yin, G.; Xi, G. Effect of gap ratio on vortex-induced vibration characteristics of different blunt bodies near-wall. J. Drain. Irrig. Mach. Eng. 2021, 39, 1132–1138. [Google Scholar]
Lin, W.; Zhou, Y.; Wang, L.; Tao, L. PIV experiment and numerical simulation of trailing vortex structure of improved INTER-MIG impeller. J. Drain. Irrig. Mach. Eng. 2021, 39, 158–164. [Google Scholar]
Han, B.; Yao, Z.; Tang, R.; Xu, H. On the supersonic impinging jet by laser Doppler velocimetry. Exp. Meas. Fluid Mech. 2002, 16, 99–103. [Google Scholar]
Darisse, A.; Lemay, J.; Benaissa, A. LDV measurements of well converged third order moments in the far field of a free turbulent round jet. Exp. Therm. Fluid Sci. 2013, 44, 825–833. [Google Scholar] [CrossRef]
Kumar, S.; Kumar, A. Effect of initial conditions on mean flow characteristics of a three dimensional turbulent wall jet. Proc. Inst. Mech. Eng. Part C J. Mech. Eng. Sci. 2021, 235, 6177–6190. [Google Scholar] [CrossRef]
Tao, D.; Zhang, R.; Ying, C. Development and application of the pollutant diffusion testing apparatus based on the image analysis. J. Saf. Environ. 2016, 16, 247–251. [Google Scholar]
Seo, H.; Kim, K.C. Experimental study on flow and turbulence characteristics of bubbly jet with low void fraction. Int. J. Multiph. Flow 2021, 142, 103738. [Google Scholar] [CrossRef]
Wen, Q.; Sha, J.; Liu, Y. TR-PIV measurement of the turbulent submerged jet and POB analysis of the dynamic structure. J. Exp. Fluid Mech. 2014, 4, 16–24. [Google Scholar]
Yang, Y.; Zhou, L.; Shi, W.; He, Z.; Han, Y.; Xiao, Y. Interstage difference of pressure pulsation in a three-stage electrical submersible pump. J. Petrol. Sci. Eng. 2021, 196, 107653. [Google Scholar] [CrossRef]
Tang, S.; Zhu, Y.; Yuan, S. An improved convolutional neural network with an adaptable learning rate towards multi-signal fault diagnosis of hydraulic piston pump. Adv. Eng. Inform. 2021, 50, 101406. [Google Scholar] [CrossRef]
Han, Y.; Song, X.; Li, K.; Yan, X. Hybrid modeling for submergence depth of the pumping well using stochastic configuration networks with random sampling. J. Petrol. Sci. Eng. 2022, 208, 109423. [Google Scholar] [CrossRef]
Tang, S.; Zhu, Y.; Yuan, S. A novel adaptive convolutional neural network for fault diagnosis of hydraulic piston pump with acoustic images. Adv. Eng. Inform. 2022, 52, 101554. [Google Scholar] [CrossRef]
Long, J.; Song, X.; Shi, J.; Chen, J. Optimization and CFD Analysis on Nozzle Exit Position of Two-phase Ejector. J. Refrig. 2022, 43, 39–45. [Google Scholar]
Ni, Q.; Ruan, W. Optimization design of desilting jet pump parameters based on response surface model. J. Ship Mech. 2022, 26, 365–374. [Google Scholar]
Zhang, K.; Zhu, X.; Ren, X.; Qiu, Q.; Shen, S. Numerical investigation on the effect of nozzle position for design of high performance ejector. Appl. Therm. Eng. 2017, 126, 594–601. [Google Scholar] [CrossRef]
Fu, W.; Liu, Z.; Li, Y.; Wu, H.; Tang, Y. Numerical study for the influences of primary steam nozzle distance and mixing chamber throat diameter on steam ejector performance. Int. J. Therm. Sci. 2018, 132, 509–516. [Google Scholar] [CrossRef]
Lucas, C.; Rusche, H.; Schroeder, A.; Koehler, J. Numerical investigation of a two-phase CO₂ ejector. Int. J. Refrigeration 2014, 43, 154–166. [Google Scholar] [CrossRef]
Ma, X.; Zhu, T.; Fu, Y.; Yan, Y.; Chen, W. Numerical simulation of rock breaking by abrasive water jet. J. Coast. Res. 2019, 93, 274–283. [Google Scholar] [CrossRef]
He, L.; Liu, Y.; Shen, K.; Yang, X.; Ba, Q.; Xiong, W. Numerical research on the dynamic rock-breaking process of impact drilling with multi-nozzle water jets. J. Pet. Sci. Eng. 2021, 207, 109145. [Google Scholar] [CrossRef]
Yu, Z.; Wang, Z.; Lei, C.; Zhou, Y.; Qiu, X. Numerical Simulation on Internal Flow Field of a Self-excited Oscillation Pulsed Jet Nozzle with Back-flow. Mech. Sci. Technol. Aerosp. Eng. 2022, 41, 998–1002. [Google Scholar]
Huang, J.; Ni, F.; Gu, L. Numerical method of FLOW-3D for sediment erosion simulation. China Harb. Eng. 2019, 39, 6–11. [Google Scholar]
Al Shaikhli, H.I.; Khassaf, S.I. Using of flow 3d as CFD materials approach in waves generation. Mater. Today Proc. 2022, 49, 2907–2911. [Google Scholar] [CrossRef]
Kosaj, R.; Alboresha, R.S.; Sulaiman, S.O. Comparison Between Numerical Flow3d Software and Laboratory Data, For Sediment Incipient Motion. IOP Conf. Ser. Earth Environ. Sci. 2022, 961, 012031. [Google Scholar] [CrossRef]
Du, C.; Liu, X.; Zhang, J.; Wang, B.; Chen, X.; Yu, X. Long-distance water hammer protection of pipeline after pump being first lowered and then rasied. J. Drain. Irrig. Mach. Eng. 2022, 40, 1248–1253, 1267. [Google Scholar]
Gao, F.; Li, X.; Gao, Q. Experiment and numerical simulation on hydraulic characteristics of novel trapezoidal measuring weir. J. Drain. Irrig. Mach. Eng. 2022, 40, 1104–1111. [Google Scholar]
Tu, A.; Nie, X.; Li, Y.; Li, H. Experimental and simulation study on water infiltration characteristics of layered red soil. J. Drain. Irrig. Mach. Eng. 2021, 39, 1243–1249. [Google Scholar]
Chen, J.; Zeng, B.; Liu, L.; Tao, K.; Zhao, H.; Zhang, C.; Zhang, J.; Li, D. Investigating the anchorage performance of full-grouted anchor bolts with a modified numerical simulation method. Eng. Fail. Anal. 2022, 141, 106640. [Google Scholar] [CrossRef]
Hu, B.; Yao, Y.; Wang, M.; Wang, C.; Liu, Y. Flow and Performance of the Disk Cavity of a Marine Gas Turbine at Varying Nozzle Pressure and Low Rotation Speeds: A Numerical Investigation. Machines 2023, 11, 68. [Google Scholar] [CrossRef]
Yao, J.; Wang, X.; Zhang, S.; Xu, S.; Jin, B.; Ding, S. Orthogonal test of important parameters affecting hydraulic performance of negative pressure feedback jet sprinkler. J. Drain. Irrig. Mach. Eng. 2021, 39, 966–972. [Google Scholar]
Wang, C.; Wang, X.; Shi, W.; Lu, W.; Tan, S.K.; Zhou, L. Experimental investigation on impingement of a submerged circular water jet at varying impinging angles and Reynolds numbers. Exp. Therm. Fluid Sci. 2017, 89, 189–198. [Google Scholar] [CrossRef]
Speziale, C.G.; Thangam, S. Analysis of an RNG based turbulence model for separated flows. Int. J. Eng. Sci. 1992, 30, 1379–1388. [Google Scholar] [CrossRef]
El Hassan, M.; Assoum, H.H.; Sobolik, V.; Vétel, J.; Abed-Meraim, K.; Garon, A.; Sakout, A. Experimental investigation of the wall shear stress and the vortex dynamics in a circular impinging jet. Exp. Fluids 2012, 52, 1475–1489. [Google Scholar] [CrossRef]
Fairweather, M.; Hargrave, G. Experimental investigation of an axisymmetric, impinging turbulent jet. 1. Velocity field. Exp. Fluids 2002, 33, 464–471. [Google Scholar] [CrossRef]
Ashforth-Frost, S.; Jambunathan, K. Effect of nozzle geometry and semi-confinement on the potential core of a turbulent axisymmetric free jet. Int. Commun. Heat Mass Transf. 1996, 23, 155–162. [Google Scholar] [CrossRef]
Chen, M.; Huang, H.; Wang, D.; Lv, S.; Chen, Y. PIV tests for flow characteristics of impinging jet in a semi-closed circular pipe. J. Vib. Shock 2021, 40, 90–97, 113. [Google Scholar]

Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

© 2023 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https://creativecommons.org/licenses/by/4.0/).

Share and Cite

MDPI and ACS Style

Mi, H.; Wang, C.; Jia, X.; Hu, B.; Wang, H.; Wang, H.; Zhu, Y. Hydraulic Characteristics of Continuous Submerged Jet Impinging on a Wall by Using Numerical Simulation and PIV Experiment. Sustainability 2023, 15, 5159. https://doi.org/10.3390/su15065159

AMA Style

Mi H, Wang C, Jia X, Hu B, Wang H, Wang H, Zhu Y. Hydraulic Characteristics of Continuous Submerged Jet Impinging on a Wall by Using Numerical Simulation and PIV Experiment. Sustainability. 2023; 15(6):5159. https://doi.org/10.3390/su15065159Chicago/Turabian Style

Mi, Hongbo, Chuan Wang, Xuanwen Jia, Bo Hu, Hongliang Wang, Hui Wang, and Yong Zhu. 2023. “Hydraulic Characteristics of Continuous Submerged Jet Impinging on a Wall by Using Numerical Simulation and PIV Experiment” Sustainability 15, no. 6: 5159. https://doi.org/10.3390/su15065159

Figure 5 A schematic of the water model of reactor URO 200.

Physical and Numerical Modeling of the Impeller Construction Impact on the Aluminum Degassing Process

알루미늄 탈기 공정에 미치는 임펠러 구성의 물리적 및 수치적 모델링

Kamil Kuglin,¹ Michał Szucki,² Jacek Pieprzyca,³ Simon Genthe,² Tomasz Merder,³ and Dorota Kalisz^1,^*

Mikael Ersson, Academic Editor

Author information Article notes Copyright and License information Disclaimer

Associated Data

Data Availability Statement

Abstract

This paper presents the results of tests on the suitability of designed heads (impellers) for aluminum refining. The research was carried out on a physical model of the URO-200, followed by numerical simulations in the FLOW 3D program. Four design variants of impellers were used in the study. The degree of dispersion of the gas phase in the model liquid was used as a criterion for evaluating the performance of each solution using different process parameters, i.e., gas flow rate and impeller speed. Afterward, numerical simulations in Flow 3D software were conducted for the best solution. These simulations confirmed the results obtained with the water model and verified them.

Keywords: aluminum, impeller construction, degassing process, numerical modeling, physical modeling

1. Introduction

Constantly increasing requirements concerning metallurgical purity in terms of hydrogen content and nonmetallic inclusions make casting manufacturers use effective refining techniques. The answer to this demand is the implementation of the aluminum refining technique making use of a rotor with an original design guaranteeing efficient refining [1,2,3,4]. The main task of the impeller (rotor) is to reduce the contamination of liquid metal (primary and recycled aluminum) with hydrogen and nonmetallic inclusions. An inert gas, mainly argon or a mixture of gases, is introduced through the rotor into the liquid metal to bring both hydrogen and nonmetallic inclusions to the metal surface through the flotation process. Appropriately and uniformly distributed gas bubbles in the liquid metal guarantee achieving the assumed level of contaminant removal economically. A very important factor in deciding about the obtained degassing effect is the optimal rotor design [5,6,7,8]. Thanks to the appropriate geometry of the rotor, gas bubbles introduced into the liquid metal are split into smaller ones, and the spinning movement of the rotor distributes them throughout the volume of the liquid metal bath. In this solution impurities in the liquid metal are removed both in the volume and from the upper surface of the metal. With a well-designed impeller, the costs of refining aluminum and its alloys can be lowered thanks to the reduced inert gas and energy consumption (optimal selection of rotor rotational speed). Shorter processing time and a high degree of dehydrogenation decrease the formation of dross on the metal surface (waste). A bigger produced dross leads to bigger process losses. Consequently, this means that the choice of rotor geometry has an indirect impact on the degree to which the generated waste is reduced [9,10].

Another equally important factor is the selection of process parameters such as gas flow rate and rotor speed [11,12]. A well-designed gas injection system for liquid metal meets two key requirements; it causes rapid mixing of the liquid metal to maintain a uniform temperature throughout the volume and during the entire process, to produce a chemically homogeneous metal composition. This solution ensures effective degassing of the metal bath. Therefore, the shape of the rotor, the arrangement of the nozzles, and their number are significant design parameters that guarantee the optimum course of the refining process. It is equally important to complete the mixing of the metal bath in a relatively short time, as this considerably shortens the refining process and, consequently, reduces the process costs. Another important criterion conditioning the implementation of the developed rotor is the generation of fine diffused gas bubbles which are distributed throughout the metal volume, and whose residence time will be sufficient for the bubbles to collide and adsorb the contaminants. The process of bubble formation by the spinning rotors differs from that in the nozzles or porous molders. In the case of a spinning rotor, the shear force generated by the rotor motion splits the bubbles into smaller ones. Here, the rotational speed, mixing force, surface tension, and fluid density have a key effect on the bubble size. The velocity of the bubbles, which depends mainly on their size and shape, determines their residence time in the reactor and is, therefore, very important for the refining process, especially since gas bubbles in liquid aluminum may remain steady only below a certain size [13,14,15].

The impeller designs presented in the article were developed to improve the efficiency of the process and reduce its costs. The impellers used so far have a complicated structure and are very pricey. The success of the conducted research will allow small companies to become independent of external supplies through the possibility of making simple and effective impellers on their own. The developed structures were tested on the water model. The results of this study can be considered as pilot.

2. Materials and Methods

Rotors were realized with the SolidWorks computer design technique and a 3D printer. The developed designs were tested on a water model. Afterward, the solution with the most advantageous refining parameters was selected and subjected to calculations with the Flow3D package. As a result, an impeller was designed for aluminum refining. Its principal lies in an even distribution of gas bubbles in the entire volume of liquid metal, with the largest possible participation of the bubble surface, without disturbing the metal surface. This procedure guarantees the removal of gaseous, as well as metallic and nonmetallic, impurities.

2.1. Rotor Designs

The developed impeller constructions, shown in Figure 1, Figure 2, Figure 3 and Figure 4, were printed on a 3D printer using the PLA (polylactide) material. The impeller design models differ in their shape and the number of holes through which the inert gas flows. Figure 1, Figure 2 and Figure 3 show the same impeller model but with a different number of gas outlets. The arrangement of four, eight, and 12 outlet holes was adopted in the developed design. A triangle-shaped structure equipped with three gas outlet holes is presented in Figure 4.

An external file that holds a picture, illustration, etc.
Object name is materials-15-05273-g001.jpg

Figure 1

A 3D model—impeller with four holes—variant B4.

An external file that holds a picture, illustration, etc.
Object name is materials-15-05273-g002.jpg

Figure 2

A 3D model—impeller with eight holes—variant B8.

An external file that holds a picture, illustration, etc.
Object name is materials-15-05273-g003.jpg

Figure 3

A 3D model—impeller with twelve holes—variant B12.

An external file that holds a picture, illustration, etc.
Object name is materials-15-05273-g004.jpg

Figure 4

A 3D model—‘red triangle’ impeller with three holes—variant RT3.

2.2. Physical Models

Investigations were carried out on a water model of the URO 200 reactor of the barbotage refining process (see Figure 5).

An external file that holds a picture, illustration, etc.
Object name is materials-15-05273-g005.jpg

Figure 5

A schematic of the water model of reactor URO 200.

The URO 200 reactor can be classified as a cyclic reactor. The main element of the device is a rotor, which ends the impeller. The whole system is attached to a shaft via which the refining gas is supplied. Then, the shaft with the rotor is immersed in the liquid metal in the melting pot or the furnace chamber. In URO 200 reactors, the refining process lasts 600 s (10 min), the gas flow rate that can be obtained ranges from 5 to 20 dm³·min⁻¹, and the speed at which the rotor can move is 0 to 400 rpm. The permissible quantity of liquid metal for barbotage refining is 300 kg or 700 kg [8,16,17]. The URO 200 has several design solutions which improve operation and can be adapted to the existing equipment in the foundry. These solutions include the following [8,16]:

URO-200XR—used for small crucible furnaces, the capacity of which does not exceed 250 kg, with no control system and no control of the refining process.
URO-200SA—used to service several crucible furnaces of capacity from 250 kg to 700 kg, fully automated and equipped with a mechanical rotor lift.
URO-200KA—used for refining processes in crucible furnaces and allows refining in a ladle. The process is fully automated, with a hydraulic rotor lift.
URO-200KX—a combination of the XR and KA models, designed for the ladle refining process. Additionally, refining in heated crucibles is possible. The unit is equipped with a manual hydraulic rotor lift.
URO-200PA—designed to cooperate with induction or crucible furnaces or intermediate chambers, the capacity of which does not exceed one ton. This unit is an integral part of the furnace. The rotor lift is equipped with a screw drive.

Studies making use of a physical model can be associated with the observation of the flow and circulation of gas bubbles. They require meeting several criteria regarding the similarity of the process and the object characteristics. The similarity conditions mainly include geometric, mechanical, chemical, thermal, and kinetic parameters. During simulation of aluminum refining with inert gas, it is necessary to maintain the geometric similarity between the model and the real object, as well as the similarity related to the flow of liquid metal and gas (hydrodynamic similarity). These quantities are characterized by the Reynolds, Weber, and Froude numbers. The Froude number is the most important parameter characterizing the process, its magnitude is the same for the physical model and the real object. Water was used as the medium in the physical modeling. The factors influencing the choice of water are its availability, relatively low cost, and kinematic viscosity at room temperature, which is very close to that of liquid aluminum.

The physical model studies focused on the flow of inert gas in the form of gas bubbles with varying degrees of dispersion, particularly with respect to some flow patterns such as flow in columns and geysers, as well as disturbance of the metal surface. The most important refining parameters are gas flow rate and rotor speed. The barbotage refining studies for the developed impeller (variants B4, B8, B12, and RT3) designs were conducted for the following process parameters:

Rotor speed: 200, 300, 400, and 500 rpm,
Ideal gas flow: 10, 20, and 30 dm³·min⁻¹,
Temperature: 293 K (20 °C).

These studies were aimed at determining the most favorable variants of impellers, which were then verified using the numerical modeling methods in the Flow-3D program.

2.3. Numerical Simulations with Flow-3D Program

Testing different rotor impellers using a physical model allows for observing the phenomena taking place while refining. This is a very important step when testing new design solutions without using expensive industrial trials. Another solution is modeling by means of commercial simulation programs such as ANSYS Fluent or Flow-3D [18,19]. Unlike studies on a physical model, in a computer program, the parameters of the refining process and the object itself, including the impeller design, can be easily modified. The simulations were performed with the Flow-3D program version 12.03.02. A three-dimensional system with the same dimensions as in the physical modeling was used in the calculations. The isothermal flow of liquid–gas bubbles was analyzed. As in the physical model, three speeds were adopted in the numerical tests: 200, 300, and 500 rpm. During the initial phase of the simulations, the velocity field around the rotor generated an appropriate direction of motion for the newly produced bubbles. When the required speed was reached, the generation of randomly distributed bubbles around the rotor was started at a rate of 2000 per second. Table 1 lists the most important simulation parameters.

Table 1

Values of parameters used in the calculations.

Parameter	Value	Unit
Maximum number of gas particles	1,000,000	–
Rate of particle generation	2000	1·s⁻¹
Specific gas constant	287.058	J·kg⁻¹·K⁻¹
Atmospheric pressure	1.013 × 10⁵	Pa
Water density	1000	kg·m⁻³
Water viscosity	0.001	kg·m⁻¹·s⁻¹
Boundary condition on the walls	No-slip	–
Size of computational cell	0.0034	m

In the case of the CFD analysis, the numerical solutions require great care when generating the computational mesh. Therefore, computational mesh tests were performed prior to the CFD calculations. The effect of mesh density was evaluated by taking into account the velocity of water in the tested object on the measurement line A (height of 0.065 m from the bottom) in a characteristic cross-section passing through the object axis (see Figure 6). The mesh contained 3,207,600, 6,311,981, 7,889,512, 11,569,230, and 14,115,049 cells.

An external file that holds a picture, illustration, etc.
Object name is materials-15-05273-g006.jpg

Figure 6

The velocity of the water depending on the size of the computational grid.

The quality of the generated computational meshes was checked using the criterion skewness angle Q_EAS [18]. This criterion is described by the following relationship:

QEAS=max{βmax−βeq180−βeq,βeq−βminβeq},

(1)

where β_max, β_min are the maximal and minimal angles (in degrees) between the edges of the cell, and β_eq is the angle corresponding to an ideal cell, which for cubic cells is 90°.

Normalized in the interval [0;1], the value of Q_EAS should not exceed 0.75, which identifies the permissible skewness angle of the generated mesh. For the computed meshes, this value was equal to 0.55–0.65.

Moreover, when generating the computational grids in the studied facility, they were compacted in the areas of the highest gradients of the calculated values, where higher turbulence is to be expected (near the impeller). The obtained results of water velocity in the studied object at constant gas flow rate are shown in Figure 6.

The analysis of the obtained water velocity distributions (see Figure 6) along the line inside the object revealed that, with the density of the grid of nodal points, the velocity changed and its changes for the test cases of 7,889,512, 11,569,230, and 14,115,049 were insignificant. Therefore, it was assumed that a grid containing not less than 7,900,000 (7,889,512) cells would not affect the result of CFD calculations.

A single-block mesh of regular cells with a size of 0.0034 m was used in the numerical calculations. The total number of cells was approximately 7,900,000 (7,889,512). This grid resolution (see Figure 7) allowed the geometry of the system to be properly represented, maintaining acceptable computation time (about 3 days on a workstation with 2× CPU and 12 computing cores).

An external file that holds a picture, illustration, etc.
Object name is materials-15-05273-g007.jpg

Figure 7

Structured equidistant mesh used in numerical calculations: (a) mesh with smoothed, surface cells (the so-called FAVOR method) used in Flow-3D; (b) visualization of the applied mesh resolution.

The calculations were conducted with an explicit scheme. The timestep was selected by the program automatically and controlled by stability and convergence. From the moment of the initial velocity field generation (start of particle generation), it was 0.0001 s.

When modeling the degassing process, three fluids are present in the system: water, gas supplied through the rotor head (impeller), and the surrounding air. Modeling such a multiphase flow is a numerically very complex issue. The necessity to overcome the liquid backpressure by the gas flowing out from the impeller leads to the formation of numerical instabilities in the volume of fluid (VOF)-based approach used by Flow-3D software. Therefore, a mixed description of the analyzed flow was used here. In this case, water was treated as a continuous medium, while, in the case of gas bubbles, the discrete phase model (DPM) model was applied. The way in which the air surrounding the system was taken into account is later described in detail.

The following additional assumptions were made in the modeling:

—The liquid phase was considered as an incompressible Newtonian fluid.
—The effect of chemical reactions during the refining process was neglected.
—The composition of each phase (gas and liquid) was considered homogeneous; therefore, the viscosity and surface tension were set as constants.
—Only full turbulence existed in the liquid, and the effect of molecular viscosity was neglected.
—The gas bubbles were shaped as perfect spheres.
—The mutual interaction between gas bubbles (particles) was neglected.

2.3.1. Modeling of Liquid Flow

The motion of the real fluid (continuous medium) is described by the Navier–Stokes Equation [20].

dudt=−1ρ∇p+ν∇2u+13ν∇(∇⋅ u)+F,

(2)

where du/dt is the time derivative, u is the velocity vector, t is the time, and F is the term accounting for external forces including gravity (unit components denoted by X, Y, Z).

In the simulations, the fluid flow was assumed to be incompressible, in which case the following equation is applicable:

∂u∂t+(u⋅∇)u=−1ρ∇p+ν∇2u+F.

(3)

Due to the large range of liquid velocities during flows, the turbulence formation process was included in the modeling. For this purpose, the k–ε model turbulence kinetic energy k and turbulence dissipation ε were the target parameters, as expressed by the following equations [21]:

∂(ρk)∂t+∂(ρkvi)∂xi=∂∂xj[(μ+μtσk)⋅∂k∂xi]+Gk+Gb−ρε−Ym+Sk,

(4)

∂(ρε)∂t+∂(ρεui)∂xi=∂∂xj[(μ+μtσε)⋅∂k∂xi]+C1εεk(Gk+G3εGb)+C2ερε2k+Sε,

(5)

where ρ is the gas density, σ_κ and σ_ε are the Prandtl turbulence numbers, k and ε are constants of 1.0 and 1.3, and G_k and G_b are the kinetic energy of turbulence generated by the average velocity and buoyancy, respectively.

As mentioned earlier, there are two gas phases in the considered problem. In addition to the gas bubbles, which are treated here as particles, there is also air, which surrounds the system. The boundary of phase separation is in this case the free surface of the water. The shape of the free surface can change as a result of the forming velocity field in the liquid. Therefore, it is necessary to use an appropriate approach to free surface tracking. The most commonly used concept in liquid–gas flow modeling is the volume of fluid (VOF) method [22,23], and Flow-3D uses a modified version of this method called TrueVOF. It introduces the concept of the volume fraction of the liquid phase f_l. This parameter can be used for classifying the cells of a discrete grid into areas filled with liquid phase (f_l = 1), gaseous phase, or empty cells (f_l = 0) and those through which the phase separation boundary (f_l ∈ (0, 1)) passes (free surface). To determine the local variations of the liquid phase fraction, it is necessary to solve the following continuity equation:

dfldt=0.

(6)

Then, the fluid parameters in the region of coexistence of the two phases (the so-called interface) depend on the volume fraction of each phase.

ρ=flρl+(1−fl)ρg,

(7)

ν=flνl+(1−fl)νg,

(8)

where indices l and g refer to the liquid and gaseous phases, respectively.

The parameter of fluid velocity in cells containing both phases is also determined in the same way.

u=flul+(1−fl)ug.

(9)

Since the processes taking place in the surrounding air can be omitted, to speed up the calculations, a single-phase, free-surface model was used. This means that no calculations were performed in the gas cells (they were treated as empty cells). The liquid could fill them freely, and the air surrounding the system was considered by the atmospheric pressure exerted on the free surface. This approach is often used in modeling foundry and metallurgical processes [24].

2.3.2. Modeling of Gas Bubble Flow

As stated, a particle model was used to model bubble flow. Spherical particles (gas bubbles) of a given size were randomly generated in the area marked with green in Figure 7b. In the simulations, the gas bubbles were assumed to have diameters of 0.016 and 0.02 m corresponding to the gas flow rates of 10 and 30 dm³·min⁻¹, respectively.

Experimental studies have shown that, as a result of turbulent fluid motion, some of the bubbles may burst, leading to the formation of smaller bubbles, although merging of bubbles into larger groupings may also occur. Therefore, to be able to observe the behavior of bubbles of different sizes (diameter), the calculations generated two additional particle types with diameters twice smaller and twice larger, respectively. The proportion of each species in the system was set to 33.33% (Table 2).

Table 2

Data assumed for calculations.

No	Rotor Speed (Rotational Speed) rpm	Bubbles Diameter m	Corresponding Gas Flow Rate dm³·min⁻¹	No	Rotor Speed (Rotational Speed) rpm	Bubbles Diameter m	Corresponding Gas Flow Rate dm³·min⁻¹
A	200	0.016	10	D	200	0.02	30
0.008	0.01
0.032	0.04
B	300	0.016	10	E	300	0.02	30
0.008	0.01
0.032	0.04
C	500	0.016	10	F	500	0.02	30
0.008	0.01
0.032	0.04

The velocity of the particle results from the generated velocity field (calculated from Equation (3) in the liquid ul around it and its velocity resulting from the buoyancy force ub. The effect of particle radius r on the terminal velocity associated with buoyancy force can be determined according to Stokes’ law.

ub=29 (ρg−ρl)μlgr2,

(10)

where g is the acceleration (9.81).

The DPM model was used for modeling the two-phase (water–air) flow. In this model, the fluid (water) is treated as a continuous phase and described by the Navier–Stokes equation, while gas bubbles are particles flowing in the model fluid (discrete phase). The trajectories of each bubble in the DPM system are calculated at each timestep taking into account the mass forces acting on it. Table 3 characterizes the DPM model used in our own research [18].

Table 3

Characteristic of the DPM model.

Method	Equations
Euler–Lagrange	Balance equation: dugdt=FD(u−ug)+g(ϱg−ϱ)ϱg+F. F_D (u − u_p) denotes the drag forces per mass unit of a bubble, and the expression for the drag coefficient F_D is of the form FD=18μCDReϱ⋅gd2g24. The relative Reynolds number has the form Re≡ρdg\|ug−u\|μ. On the other hand, the force resulting from the additional acceleration of the model fluid has the form F=12dρdtρg(u−ug), where u_g is the gas bubble velocity, u is the liquid velocity, d_g is the bubble diameter, and C_D is the drag coefficient.

Method

Equations

Euler–Lagrange

Balance equation:
dugdt=FD(u−ug)+g(ϱg−ϱ)ϱg+F.
F_D (u − u_p) denotes the drag forces per mass unit of a bubble, and the expression for the drag coefficient F_D is of the form
FD=18μCDReϱ⋅gd2g24.
The relative Reynolds number has the form
Re≡ρdg|ug−u|μ.
On the other hand, the force resulting from the additional acceleration of the model fluid has the form
F=12dρdtρg(u−ug),
where u_g is the gas bubble velocity, u is the liquid velocity, d_g is the bubble diameter, and C_D is the drag coefficient.

3. Results and Discussion

3.1. Calculations of Power and Mixing Time by the Flowing Gas Bubbles

One of the most important parameters of refining with a rotor is the mixing power induced by the spinning rotor and the outflowing gas bubbles (via impeller). The mixing power of liquid metal in a ladle of height (h) by gas injection can be determined from the following relation [15]:

pgVm=ρ⋅g⋅uB,

(11)

where p_g is the mixing power, V_m is the volume of liquid metal in the reactor, ρ is the density of liquid aluminum, and u_B is the average speed of bubbles, given below.

uB=n⋅R⋅TAc⋅Pm⋅t,

(12)

where n is the number of gas moles, R is the gas constant (8.314), A_c is the cross-sectional area of the reactor vessel, T is the temperature of liquid aluminum in the reactor, and P_m is the pressure at the middle tank level. The pressure at the middle level of the tank is calculated by a function of the mean logarithmic difference.

Pm=(Pa+ρ⋅g⋅h)−Paln(Pa+ρ⋅g⋅h)Pa,

(13)

where P_a is the atmospheric pressure, and h is the the height of metal in the reactor.

Themelis and Goyal [25] developed a model for calculating mixing power delivered by gas injection.

pg=2Q⋅R⋅T⋅ln(1+m⋅ρ⋅g⋅hP),

(14)

where Q is the gas flow, and m is the mass of liquid metal.

Zhang [26] proposed a model taking into account the temperature difference between gas and alloy (metal).

pg=QRTgVm[ln(1+ρ⋅g⋅hPa)+(1−TTg)],

(15)

where T_g is the gas temperature at the entry point.

Data for calculating the mixing power resulting from inert gas injection into liquid aluminum are given below in Table 4. The design parameters were adopted for the model, the parameters of which are shown in Figure 5.

Table 4

Data for calculating mixing power introduced by an inert gas.

Parameter	Value	Unit
Height of metal column	0.7	m
Density of aluminum	2375	kg·m⁻³
Process duration	20	s
Gas temperature at the injection site	940	K
Cross-sectional area of ladle	0.448	m²
Mass of liquid aluminum	546.25	kg
Volume of ladle	0.23	M³
Temperature of liquid aluminum	941.15	K

Table 5 presents the results of mixing power calculations according to the models of Themelis and Goyal and of Zhang for inert gas flows of 10, 20, and 30 dm³·min⁻¹. The obtained calculation results significantly differed from each other. The difference was an order of magnitude, which indicates that the model is highly inaccurate without considering the temperature of the injected gas. Moreover, the calculations apply to the case when the mixing was performed only by the flowing gas bubbles, without using a rotor, which is a great simplification of the phenomenon.

Table 5

Mixing power calculated from mathematical models.

Mathematical Model	Mixing Power (W·t⁻¹) for a Given Inert Gas Flow (dm³·min⁻¹)
	10	20	30
Themelis and Goyal	11.49	23.33	35.03
Zhang	0.82	1.66	2.49

The mixing time is defined as the time required to achieve 95% complete mixing of liquid metal in the ladle [27,28,29,30]. Table 6 groups together equations for the mixing time according to the models.

Table 6

Models for calculating mixing time.

Authors	Model	Remarks
Szekely [31]	τ=800ε−0.4	ε—W·t⁻¹
Chiti and Paglianti [27]	τ=CVQl	V—volume of reactor, m³ Ql—flow intensity, m³·s⁻¹
Iguchi and Nakamura [32]	τ=1200⋅Q−0.4D1.97h−1.0υ0.47	υ—kinematic viscosity, m²·s⁻¹ D—diameter of ladle, m h—height of metal column, m Q—liquid flow intensity, m³·s⁻¹

Figure 8 and Figure 9 show the mixing time as a function of gas flow rate for various heights of the liquid column in the ladle and mixing power values.

An external file that holds a picture, illustration, etc.
Object name is materials-15-05273-g008.jpg

Figure 8

Mixing time as a function of gas flow rate for various heights of the metal column (Iguchi and Nakamura model).

An external file that holds a picture, illustration, etc.
Object name is materials-15-05273-g009.jpg

Figure 9

Mixing time as a function of mixing power (Szekly model).

3.2. Determining the Bubble Size

The mechanisms controlling bubble size and mass transfer in an alloy undergoing refining are complex. Strong mixing conditions in the reactor promote impurity mass transfer. In the case of a spinning rotor, the shear force generated by the rotor motion separates the bubbles into smaller bubbles. Rotational speed, mixing force, surface tension, and liquid density have a strong influence on the bubble size. To characterize the kinetic state of the refining process, parameters k and A were introduced. Parameters k, A, and u_B can be calculated using the below equations [33].

k=2D⋅uBdB⋅π−−−−−−√,

(16)

A=6Q⋅hdB⋅uB,

(17)

uB=1.02g⋅dB,−−−−−√

(18)

where D is the diffusion coefficient, and d_B is the bubble diameter.

After substituting appropriate values, we get

dB=3.03×104(πD)−2/5g−1/5h4/5Q0.344N−1.48.

(19)

According to the last equation, the size of the gas bubble decreases with the increasing rotational speed (see Figure 10).

An external file that holds a picture, illustration, etc.
Object name is materials-15-05273-g010.jpg

Figure 10

Effect of rotational speed on the bubble diameter.

In a flow of given turbulence intensity, the diameter of the bubble does not exceed the maximum size d_max, which is inversely proportional to the rate of kinetic energy dissipation in a viscous flow ε. The size of the gas bubble diameter as a function of the mixing energy, also considering the Weber number and the mixing energy in the negative power, can be determined from the following equations [31,34]:

—Sevik and Park:

dBmax=We0.6kr⋅(σ⋅103ρ⋅10−3)0.6⋅(10⋅ε)−0.4⋅10−2.

(20)

—Evans:

dBmax=⎡⎣Wekr⋅σ⋅1032⋅(ρ⋅10−3)13⎤⎦35 ⋅(10⋅ε)−25⋅10−2.

(21)

The results of calculating the maximum diameter of the bubble d_Bmax determined from Equation (21) are given in Table 7.

Table 7

The results of calculating the maximum diameter of the bubble using Equation (21).

Model	Mixing Energy ĺ (m²·s⁻³)	Weber Number (We_kr)
0.59	1.0	1.2
Zhang and Taniguchi d_max	0.1	0.0167	0.0230	0.026
0.5	0.0088	0.0121	0.013
1.0	0.0067	0.0091	0.010
1.5	0.0057	0.0078	0.009
Sevik and Park d_B_max	0.1	0.265	0.36	0.41
0.5	0.139	0.19	0.21
1.0	0.106	0.14	0.16
1.5	0.090	0.12	0.14
Evans d_B_max	0.1	0.247	0.340	0.38
0.5	0.130	0.178	0.20
1.0	0.098	0.135	0.15
1.5	0.084	0.115	0.13

3.3. Physical Modeling

The first stage of experiments (using the URO-200 water model) included conducting experiments with impellers equipped with four, eight, and 12 gas outlets (variants B4, B8, B12). The tests were carried out for different process parameters. Selected results for these experiments are presented in Figure 11, Figure 12, Figure 13 and Figure 14.

An external file that holds a picture, illustration, etc.
Object name is materials-15-05273-g011.jpg

Figure 11

Impeller variant B4—gas bubbles dispersion registered for a gas flow rate of 10 dm³·min⁻¹ and rotor speed of (a) 200, (b) 300, (c) 400, and (d) 500 rpm.

An external file that holds a picture, illustration, etc.
Object name is materials-15-05273-g012.jpg

Figure 12

Impeller variant B8—gas bubbles dispersion registered for a gas flow rate of 10 dm³·min⁻¹ and rotor speed of (a) 200, (b) 300, (c) 400, and (d) 500 rpm.

An external file that holds a picture, illustration, etc.
Object name is materials-15-05273-g013.jpg

Figure 13

Gas bubble dispersion registered for different processing parameters (impeller variant B12).

An external file that holds a picture, illustration, etc.
Object name is materials-15-05273-g014.jpg

Figure 14

Gas bubble dispersion registered for different processing parameters (impeller variant RT3).

The analysis of the refining variants presented in Figure 11, Figure 12, Figure 13 and Figure 14 reveals that the proposed impellers design model is not useful for the aluminum refining process. The number of gas outlet orifices, rotational speed, and flow did not affect the refining efficiency. In all the variants shown in the figures, very poor dispersion of gas bubbles was observed in the object. The gas bubble flow had a columnar character, and so-called dead zones, i.e., areas where no inert gas bubbles are present, were visible in the analyzed object. Such dead zones were located in the bottom and side zones of the ladle, while the flow of bubbles occurred near the turning rotor. Another negative phenomenon observed was a significant agitation of the water surface due to excessive (rotational) rotor speed and gas flow (see Figure 13, cases 20; 400, 30; 300, 30; 400, and 30; 500).

Research results for a ‘red triangle’ impeller equipped with three gas supply orifices (variant RT3) are presented in Figure 14.

In this impeller design, a uniform degree of bubble dispersion in the entire volume of the modeling fluid was achieved for most cases presented (see Figure 14). In all tested variants, single bubbles were observed in the area of the water surface in the vessel. For variants 20; 200, 30; 200, and 20; 300 shown in Figure 14, the bubble dispersion results were the worst as the so-called dead zones were identified in the area near the bottom and sidewalls of the vessel, which disqualifies these work parameters for further applications. Interestingly, areas where swirls and gas bubble chains formed were identified only for the inert gas flows of 20 and 30 dm³·min⁻¹ and 200 rpm in the analyzed model. This means that the presented model had the best performance in terms of dispersion of gas bubbles in the model liquid. Its design with sharp edges also differed from previously analyzed models, which is beneficial for gas bubble dispersion, but may interfere with its suitability in industrial conditions due to possible premature wear.

3.4. Qualitative Comparison of Research Results (CFD and Physical Model)

The analysis (physical modeling) revealed that the best mixing efficiency results were obtained with the RT3 impeller variant. Therefore, numerical calculations were carried out for the impeller model with three outlet orifices (variant RT3). The CFD results are presented in Figure 15 and Figure 16.

An external file that holds a picture, illustration, etc.
Object name is materials-15-05273-g015.jpg

Figure 15

Simulation results of the impeller RT3, for given flows and rotational speeds after a time of 1 s: simulation variants (a) A, (b) B, (c) C, (d) D, (e) E, and (f) F.

An external file that holds a picture, illustration, etc.
Object name is materials-15-05273-g016.jpg

Figure 16

Simulation results of the impeller RT3, for given flows and rotational speeds after a time of 5.4 s.: simulation variants (a) A, (b) B, (c) C, (d) D, (e) E, and (f) F.

CFD results are presented for all analyzed variants (impeller RT3) at two selected calculation timesteps of 1 and 5.40 s. They show the velocity field of the medium (water) and the dispersion of gas bubbles.

Figure 15 shows the initial refining phase after 1 s of the process. In this case, the gas bubble formation and flow were observed in an area close to contact with the rotor. Figure 16 shows the phase when the dispersion and flow of gas bubbles were advanced in the reactor area of the URO-200 model.

The quantitative evaluation of the obtained results of physical and numerical model tests was based on the comparison of the degree of gas dispersion in the model liquid. The degree of gas bubble dispersion in the volume of the model liquid and the areas of strong turbulent zones formation were evaluated during the analysis of the results of visualization and numerical simulations. These two effects sufficiently characterize the required course of the process from the physical point of view. The known scheme of the below description was adopted as a basic criterion for the evaluation of the degree of dispersion of gas bubbles in the model liquid.

Minimal dispersion—single bubbles ascending in the region of their formation along the ladle axis; lack of mixing in the whole bath volume.
Accurate dispersion—single and well-mixed bubbles ascending toward the bath mirror in the region of the ladle axis; no dispersion near the walls and in the lower part of the ladle.
Uniform dispersion—most desirable; very good mixing of fine bubbles with model liquid.
Excessive dispersion—bubbles join together to form chains; large turbulence zones; uneven flow of gas.

The numerical simulation results give a good agreement with the experiments performed with the physical model. For all studied variants (used process parameters), the single bubbles were observed in the area of water surface in the vessel. For variants presented in Figure 13 (200 rpm, gas flow 20 and dm³·min⁻¹) and relevant examples in numerical simulation Figure 16, the worst bubble dispersion results were obtained because the dead zones were identified in the area near the bottom and sidewalls of the vessel, which disqualifies these work parameters for further use. The areas where swirls and gas bubble chains formed were identified only for the inert gas flows of 20 and 30 dm³·min⁻¹ and 200 rpm in the analyzed model (physical model). This means that the presented impeller model had the best performance in terms of dispersion of gas bubbles in the model liquid. The worst bubble dispersion results were obtained because the dead zones were identified in the area near the bottom and side walls of the vessel, which disqualifies these work parameters for further use.

Figure 17 presents exemplary results of model tests (CFD and physical model) with marked gas bubble dispersion zones. All variants of tests were analogously compared, and this comparison allowed validating the numerical model.

An external file that holds a picture, illustration, etc.
Object name is materials-15-05273-g017.jpg

Figure 17

Compilations of model research results (CFD and physical): A—single gas bubbles formed on the surface of the modeling liquid, B—excessive formation of gas chains and swirls, C—uniform distribution of gas bubbles in the entire volume of the tank, and D—dead zones without gas bubbles, no dispersion. (a) Variant B; (b) variant F.

It should be mentioned here that, in numerical simulations, it is necessary to make certain assumptions and simplifications. The calculations assumed three particle size classes (Table 2), which represent the different gas bubbles that form due to different gas flow rates. The maximum number of particles/bubbles (Table 1) generated was assumed in advance and related to the computational capabilities of the computer. Too many particles can also make it difficult to visualize and analyze the results. The size of the particles, of course, affects their behavior during simulation, while, in the figures provided in the article, the bubbles are represented by spheres (visualization of the results) of the same size. Please note that, due to the adopted Lagrangian–Eulerian approach, the simulation did not take into account phenomena such as bubble collapse or fusion. However, the obtained results allow a comprehensive analysis of the behavior of gas bubbles in the system under consideration.

The comparative analysis of the visualization (quantitative) results obtained with the water model and CFD simulations (see Figure 17) generated a sufficient agreement from the point of view of the trends. A precise quantitative evaluation is difficult to perform because of the lack of a refraction compensating system in the water model. Furthermore, in numerical simulations, it is not possible to determine the geometry of the forming gas bubbles and their interaction with each other as opposed to the visualization in the water model. The use of both research methods is complementary. Thus, a direct comparison of images obtained by the two methods requires appropriate interpretation. However, such an assessment gives the possibility to qualitatively determine the types of the present gas bubble dispersion, thus ultimately validating the CFD results with the water model.

A summary of the visualization results for impellers RT3, i.e., analysis of the occurring gas bubble dispersion types, is presented in Table 8.

Table 8

Summary of visualization results (impeller RT3)—different types of gas bubble dispersion.

No Exp.	A	B	C	D	E	F
Gas flow rate, dm³·min⁻¹	10	30
Impeller speed, rpm	200	300	500	200	300	500
Type of dispersion	Accurate	Uniform	Uniform/excessive	Minimal	Excessive	Excessive

Tests carried out for impeller RT3 confirmed the high efficiency of gas bubble distribution in the volume of the tested object at a low inert gas flow rate of 10 dm³·min⁻¹. The most optimal variant was variant B (300 rpm, 10 dm³·min⁻¹). However, the other variants A and C (gas flow rate 10 dm³·min⁻¹) seemed to be favorable for this type of impeller and are recommended for further testing. The above process parameters will be analyzed in detail in a quantitative analysis to be performed on the basis of the obtained efficiency curves of the degassing process (oxygen removal). This analysis will give an unambiguous answer as to which process parameters are the most optimal for this type of impeller; the results are planned for publication in the next article.

It should also be noted here that the high agreement between the results of numerical calculations and physical modelling prompts a conclusion that the proposed approach to the simulation of a degassing process which consists of a single-phase flow model with a free surface and a particle flow model is appropriate. The simulation results enable us to understand how the velocity field in the fluid is formed and to analyze the distribution of gas bubbles in the system. The simulations in Flow-3D software can, therefore, be useful for both the design of the impeller geometry and the selection of process parameters.

4. Conclusions

The results of experiments carried out on the physical model of the device for the simulation of barbotage refining of aluminum revealed that the worst results in terms of distribution and dispersion of gas bubbles in the studied object were obtained for the black impellers variants B4, B8, and B12 (multi-orifice impellers—four, eight, and 12 outlet holes, respectively).

In this case, the control of flow, speed, and number of gas exit orifices did not improve the process efficiency, and the developed design did not meet the criteria for industrial tests. In the case of the ‘red triangle’ impeller (variant RT3), uniform gas bubble dispersion was achieved throughout the volume of the modeling fluid for most of the tested variants. The worst bubble dispersion results due to the occurrence of the so-called dead zones in the area near the bottom and sidewalls of the vessel were obtained for the flow variants of 20 dm³·min⁻¹ and 200 rpm and 30 dm³·min⁻¹ and 200 rpm. For the analyzed model, areas where swirls and gas bubble chains were formed were found only for the inert gas flow of 20 and 30 dm³·min⁻¹ and 200 rpm. The model impeller (variant RT3) had the best performance compared to the previously presented impellers in terms of dispersion of gas bubbles in the model liquid. Moreover, its design differed from previously presented models because of its sharp edges. This can be advantageous for gas bubble dispersion, but may negatively affect its suitability in industrial conditions due to premature wearing.

The CFD simulation results confirmed the results obtained from the experiments performed on the physical model. The numerical simulation of the operation of the ‘red triangle’ impeller model (using Flow-3D software) gave good agreement with the experiments performed on the physical model. This means that the presented model impeller, as compared to other (analyzed) designs, had the best performance in terms of gas bubble dispersion in the model liquid.

In further work, the developed numerical model is planned to be used for CFD simulations of the gas bubble distribution process taking into account physicochemical parameters of liquid aluminum based on industrial tests. Consequently, the obtained results may be implemented in production practice.

Funding Statement

This paper was created with the financial support grants from the AGH-UST, Faculty of Foundry Engineering, Poland (16.16.170.654 and 11/990/BK_22/0083) for the Faculty of Materials Engineering, Silesian University of Technology, Poland.

Author Contributions

Conceptualization, K.K. and D.K.; methodology, J.P. and T.M.; validation, M.S. and S.G.; formal analysis, D.K. and T.M.; investigation, J.P., K.K. and S.G.; resources, M.S., J.P. and K.K.; writing—original draft preparation, D.K. and T.M.; writing—review and editing, D.K. and T.M.; visualization, J.P., K.K. and S.G.; supervision, D.K.; funding acquisition, D.K. and T.M. All authors have read and agreed to the published version of the manuscript.

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

Data are contained within the article.

Conflicts of Interest

The authors declare no conflict of interest.

Footnotes

Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.