Numerical Investigation into the Effect of Damage Openings on Ship Hydrodynamics by the Overset Mesh Technique

Damage stability is difficult to assess due to the complex hydrodynamic phenomena regarding interactions between fluid and structures. Therefore, a detailed analysis of the flooding progression and motion responses is important for improving ship safety. In this paper, numerical simulations are performed on the damaged DTMB 5415 ship at zero speed. All calculation are carried out using CD Adapco Star CCM + software, investigating the effect of damage openings on ship hydrodynamics, including the side damage and the bottom damage. The computational domain is modelled by the overset mesh and solved using the unsteady Reynold-average Navier-Stokes (URANS) solver. An implicit solver is used to find the field of all hydrodynamics unknown quantities, in conjunction with an iterative solver to solve each time step. The Volume of Fluid (VOF) method is applied to visualize the flooding process and capture the complex hydrodynamics behaviors. The simulation results indicated that two damage locations produce the characteristic flooding processes, and the motion responses corresponding to the hydrodynamic behaviors are different. Through comparative analysis, due to the difference between the horizontal impact on the longitudinal bulkhead and the vertical impact on the bottom plate, the bottom damage scenario always has a larger heel angle than the side damage scenario in the same period. However, the pitch motions are basically consistent. Generally, the visualization of the flooding process is efficient to explain the causes of the motion responses. Also, when the damage occurs, regardless of the bottom damage or the side damage, the excessive heel angle due to asymmetric flooding is often a threat to ship survivability with respect to the pitch angle.


Introduction
Nowadays, ship safety is of high priority to the maritime industry. However, despite many efforts being to improve ship design in recent years, damage accidents continue to occur due to collision, grounding, or the unpredictable sea environment (wind, current and waves). The loss of hull integrity leading to damage flooding can be a severe risk to ship stability [1], even making the damaged ship sink or capsize. For a damaged ship, different damage scenarios correspond to the special opening locations. The resulting flooding processes and ship motion responses are also characteristic. In the flooding process, the ship motions impact the water flooding and sloshing in the flooded compartments. Simultaneously, the liquid loads acting on the compartments also influence the ship motions [2]. Therefore, the accurate prediction of the hydrodynamic behavior coupled with the ship motions is crucial to assess the remaining survivability of the damaged ship. Also, the complex hydrodynamic behaviors caused by the coupled motion has attracted significant attention at several recent International Towing Tank Conferences [3,4].
In order to enhance the understanding of the flooding process and motion responses of the damaged ship, a series of model experiments were performed while various numerical simulation methods were developed and implemented. An experimental campaign was carried out on a passenger ferry hull to underline the effects of the damage opening on the ship roll response. The damaged ship was placed in still water and beam regular waves at zero speed [5]. Lim et al. [6] used a course-keeping model ship to measure the advance speed and motion response of the damaged ship in head and following seas. Through the free-running tests, the motion characteristics under safe return to port (SRTP) regulations were identified. Siddiqui et al. [7] performed a detailed series of experiments in a wave flume on a thin walled prismatic hull form. The obtained results demonstrate the occurrence of sloshing and piston mode resonance in the tests and their influence on the hydrodynamics load of a damaged ship. The effect of air compressibility in the airtight compartment on local floodwater behavior was also investigated. Rodrigues, et al. [8] presented an experimental procedure to measure progressive flooding on a small-scale damaged floating body. Their focus was on the estimation of the discharge coefficient at different opening geometries in still water. Additionally, a numerical simulation based on the Reynolds averaged Navier-Stokes (RANS) solver was applied to validate its capability to reasonably reproduce the physical experiments. Korkut et al. [9] carried out six degrees of freedom motion responses tests in regular waves for intact and damaged conditions. The effect of the damage opening and waves with different wave heights and wave frequencies on the motion responses of the damaged model was explored. The obtained experimental results indicated that the damage opening has an adverse influence depending on the directionality of the waves and the applied wave frequency.
Generally, the applied experimental models above are created on the simplified assumption that the damaged compartment is empty, not considering the effect of permeability on the flooding and motion responses. In the real compartment layout, obstacles in the compartment and inner subdivision may affect the flooding path and quantity. Therefore, accounting for more realistic modeling of the damage flooding, Acanfora et al. [10] carried out an experimental investigation on the dynamic response of a damaged ship with a realistic arrangement of the flooded compartment. The results presented the effects of obstacles in the engine room compartment, such as decks and engine, on the roll responses. Similarly, Domesh et al. [11] used a damaged segmented ship model to study the effect of permeability and the internal arrangement of the damaged compartment on the pitch and heave responses. When the effect of the internal structures on the flooding process and motion responses is taken into account, the bearing capacity of different components will determine whether the secondary water ingress will occur. When the flooding water pressure exceeds the bearing limitation of components, the components will leak or collapse. Based on this aspect, Risto et al. [12] conducted unique full-scale tests to determine the leakage and collapse characteristics of various typical non-watertight structures, when subjected to the flooding water pressure. The obtained results can provide guideline values to determine when the structure may collapse with the accumulation of the flooding water. These well-designed model tests can accurately assess the damaged stability with complicated physical phenomena, establishing a database for the motion responses of the damaged ship.
Although the experimental tests can well investigate damaged ships, their flexibility and economic efficiency are very limited. In most cases, full-scale experiments are impossible, and model experiments are associated with problems of the scale effect. In the last decade, owing to the development of high-performance computers, there has been an increasing interest in the application of computational fluid dynamics (CFD) to investigate the multi-phenomena hydrodynamic problem of the damaged ship. A Navier-Stokes (NS) solver with a free surface capturing technique, i.e., the volume of fluid method, was developed to numerically simulate water flooding into a damaged vessel. The proposed method can be used to predict the dynamic behavior of the flooding water and its impact forces on the flooded compartment [13]. Sadat-Hosseini et al. [14] performed unsteady Reynold averaged Navier-Stokes (URANS) simulations for zero-speed damaged passenger ships in calm water and waves. The flooding procedure and roll decay in calm water were studied, and the motions in regular beam waves for various wavelength were analyzed. Even though the simulation demands a larger computation costs, the predicted results coincide better with the experiment results than those reported for potential flow solver. Santos et al. [15] described a mathematical model in the time domain of the motions and flooding of ships in a seaway. Different factors affecting the survivability of the damaged ships were assessed. Ming et al. [16] applied the weakly compressible smoothed particle hydrodynamics (SPH) method to explore the influences of transversal waves on the dynamic flooding process of a damaged compartment. The simulation results indicate that when the waves slam against the damaged ship, the relative position between the damage opening and the free surface will be changed. Further, different wave directions will result in different flooding processes. This method has the advantage of dealing with large deformation problems of free surface flow and fluid-structure interaction. In addition, Manderbacka et al. [17] and Acanfora et al. [18] presented a non-linear time domain simulation method for damaged ships. The flooding water motion is based on the lumped mass method with a moving free surface, and the ship's transient response to an abrupt flooding is simulated. It has been proven that these numerical methods have been an alternative approach to study the damage flooding.
In this paper, the URANS method combining the overset mesh technique in conjunction with the 6 Degree-of-Freedom (DOF) solver is applied to investigate the effect of side damage and bottom damage on the flooding process and motion responses. The paper is organized as follows. The overset mesh methodology is depicted in Section 2, including the definition of the overset mesh and the interpolation options. Section 3 introduces the utilized 5415 scale model. Section 4 underlies the whole simulation process and the critical settings, including the creation of the simulation domain, the boundary conditions, the choice of mesh types and the relevant solver settings. The captured flooding process and measured motion responses are presented and discussed in Section 5. Finally, conclusions and future research directions are taken in Section 6.

Definition of the Overset Mesh
With the development of computational fluid dynamics, time-domain simulation approaches based on the finite volume method (FOM) have been constantly evolving. However, for unstructured grids, the mesh-partitioning stage can be challenging due to the memory limitations of the massively parallel architectures. In this case, in order to run a simulation with increasingly fine grids and increasingly complex physics modelling, the high-performance computing represents a crucial capability to solve this problem [19]. Moreover, among the different motion-mesh techniques, the overset mesh (also known as Chimera or overlapping grids) has been considered as an efficient way to accurately describe the rotation motion of the damaged ship. Taking the simulation case as an example, as illustrated in Figure 1, when the overset mesh is activated, two individual regions are created: a background region and an overset region surrounding the damaged ship. The motion specification can be assigned to the overset region, rather than the background region. Because when the overset mesh is not applied, the motion specification can only be assigned to the background region. The background region will rotate relative to the stationary hull under the influence of the flooding water, monitoring the roll and pitch motion of the damaged hull. This will make the free surface out of the original mesh refinement block, resulting in poor simulation accuracy. However, the overset mesh can avoid this problem well. When the motion specification is assigned to the overset region, the background region is always stationary. The limited overset region will follow the movement of the damaged ship.
Even if the damaged ship has a large roll or pitch motion, good overset mesh quality will ensure the accuracy of the simulation results. In addition, two regions are meshed separately, and the overset interface is created between them. For implicitly coupling the background region and the overset region, the interface is set to the overset mesh boundary condition. In this case, as the damaged ship moves within the background region, the overset region will correspondingly change. Simultaneously, the data information between the regions is exchanged through the overlapping cells. As illustrated in Figure 2, once the overset mesh is performed, the hole-cutting process in STAR-CCM + automatically couples the overset region with the background region through the overset interface. Four types of cells from the hole-cutting process are created. Active cells (cyan and yellow): Discretizing governing equations are solved here. Passive cells (dark blue). Donor cells (green): These provide interpolation information to the mesh acceptor cells. Acceptor cells (red): The boundary cells that receive information from the donor cells [20]. meshed separately, and the overset interface is created between them. For implicitly coupling the background region and the overset region, the interface is set to the overset mesh boundary condition. In this case, as the damaged ship moves within the background region, the overset region will correspondingly change. Simultaneously, the data information between the regions is exchanged through the overlapping cells. As illustrated in Figure 2, once the overset mesh is performed, the hole-cutting process in STAR-CCM + automatically couples the overset region with the background region through the overset interface. Four types of cells from the hole-cutting process are created. Active cells (cyan and yellow): Discretizing governing equations are solved here. Passive cells (dark blue). Donor cells (green): These provide interpolation information to the mesh acceptor cells. Acceptor cells (red): The boundary cells that receive information from the donor cells [20].

Interpolation Option
For the overset mesh, the interpolation function determines the data transfer relationship between the acceptor cells and the donor cells, ensuring implicit coupling of the background region and the overset region. In this case, a solution is computed on all grids simultaneously, leading to improved robustness and convergence. Of the specified interpolation options, the alternative interpolation options include distance-weighted interpolation, linear interpolation, and least-squares interpolation [21]. For the distance-weight interpolation, the interpolation factors are inversely proportional to the distance from the acceptor to the donor cell center, resulting in the closest cell giving the largest contribution. If the simulation involves moving mesh without great motion, the linear interpolation will be a better choice as it can ensure that interpolation elements do not overlap. This choice is more accurate but also more expensive due to the computational effort required. Finally, the least-squares interpolation transfers data from source to target meshes by data mappers. This method is suitable where there is a large variation of the moving grid with respect to the background mesh, as indicated in CD Adapco User's Guide [21] and De Luca et al. [22].
Based on the features of the simulation cases, linear interpolation is adopted to obtain an accurate solution.

Model Description
Numerical simulations have been performed using the well-known benchmark US Navy Destroyer Hull DTMB 5415 with a corresponding scale ratio of 1:25. Figures 3 and 4  meshed separately, and the overset interface is created between them. For implicitly coupling the background region and the overset region, the interface is set to the overset mesh boundary condition. In this case, as the damaged ship moves within the background region, the overset region will correspondingly change. Simultaneously, the data information between the regions is exchanged through the overlapping cells. As illustrated in Figure 2, once the overset mesh is performed, the hole-cutting process in STAR-CCM + automatically couples the overset region with the background region through the overset interface. Four types of cells from the hole-cutting process are created.

Interpolation Option
For the overset mesh, the interpolation function determines the data transfer relationship between the acceptor cells and the donor cells, ensuring implicit coupling of the background region and the overset region. In this case, a solution is computed on all grids simultaneously, leading to improved robustness and convergence. Of the specified interpolation options, the alternative interpolation options include distance-weighted interpolation, linear interpolation, and least-squares interpolation [21]. For the distance-weight interpolation, the interpolation factors are inversely proportional to the distance from the acceptor to the donor cell center, resulting in the closest cell giving the largest contribution. If the simulation involves moving mesh without great motion, the linear interpolation will be a better choice as it can ensure that interpolation elements do not overlap. This choice is more accurate but also more expensive due to the computational effort required. Finally, the least-squares interpolation transfers data from source to target meshes by data mappers. This method is suitable where there is a large variation of the moving grid with respect to the background mesh, as indicated in CD Adapco User's Guide [21] and De Luca et al. [22].
Based on the features of the simulation cases, linear interpolation is adopted to obtain an accurate solution.

Model Description
Numerical simulations have been performed using the well-known benchmark US Navy Destroyer Hull DTMB 5415 with a corresponding scale ratio of 1:25. Figures 3 and 4

Interpolation Option
For the overset mesh, the interpolation function determines the data transfer relationship between the acceptor cells and the donor cells, ensuring implicit coupling of the background region and the overset region. In this case, a solution is computed on all grids simultaneously, leading to improved robustness and convergence. Of the specified interpolation options, the alternative interpolation options include distance-weighted interpolation, linear interpolation, and least-squares interpolation [21]. For the distance-weight interpolation, the interpolation factors are inversely proportional to the distance from the acceptor to the donor cell center, resulting in the closest cell giving the largest contribution. If the simulation involves moving mesh without great motion, the linear interpolation will be a better choice as it can ensure that interpolation elements do not overlap. This choice is more accurate but also more expensive due to the computational effort required. Finally, the least-squares interpolation transfers data from source to target meshes by data mappers. This method is suitable where there is a large variation of the moving grid with respect to the background mesh, as indicated in CD Adapco User's Guide [21] and De Luca et al. [22].
Based on the features of the simulation cases, linear interpolation is adopted to obtain an accurate solution.

Model Description
Numerical simulations have been performed using the well-known benchmark US Navy Destroyer Hull DTMB 5415 with a corresponding scale ratio of 1:25. Figures 3 and 4 respectively show the side view and body lines of the ship model, and Table 1 presents the principal dimensions of the ship model [23]. The created damaged compartment is located near the bow. The damaged compartment is assumed to be empty, not considering the influence of permeability and internal arrangements on the flooding water motion and damaged stability. However, in the real damage scenario, compartments could be full of equipment and obstacles that modify the flow of the water ingress. The simulation results will visualize the flooding process and measure the motion responses of the damaged ship in still water. The influence of the forward speed and external wave conditions on the coupled motion of the damaged ship and flooding water is not taken into account. the side view and body lines of the ship model, and Table 1 presents the principal dimensions of the ship model [23]. The created damaged compartment is located near the bow. The damaged compartment is assumed to be empty, not considering the influence of permeability and internal arrangements on the flooding water motion and damaged stability. However, in the real damage scenario, compartments could be full of equipment and obstacles that modify the flow of the water ingress. The simulation results will visualize the flooding process and measure the motion responses of the damaged ship in still water. The influence of the forward speed and external wave conditions on the coupled motion of the damaged ship and flooding water is not taken into account.   As shown in Figure 5, two comparative damage scenarios are modeled separately. It can be found that the scenarios are specific to the location of the damage opening. Side damage and bottom damage will cause different types of flooding, while the corresponding motion responses are different. In the simulation process, in order to eliminate the composite influence arising from the air compressibility, an appropriate ventilation hole is constructed on the upper deck. In [24][25][26][27], it has the side view and body lines of the ship model, and Table 1 presents the principal dimensions of the ship model [23]. The created damaged compartment is located near the bow. The damaged compartment is assumed to be empty, not considering the influence of permeability and internal arrangements on the flooding water motion and damaged stability. However, in the real damage scenario, compartments could be full of equipment and obstacles that modify the flow of the water ingress. The simulation results will visualize the flooding process and measure the motion responses of the damaged ship in still water. The influence of the forward speed and external wave conditions on the coupled motion of the damaged ship and flooding water is not taken into account.   As shown in Figure 5, two comparative damage scenarios are modeled separately. It can be found that the scenarios are specific to the location of the damage opening. Side damage and bottom damage will cause different types of flooding, while the corresponding motion responses are different. In the simulation process, in order to eliminate the composite influence arising from the air compressibility, an appropriate ventilation hole is constructed on the upper deck. In [24][25][26][27], it has  As shown in Figure 5, two comparative damage scenarios are modeled separately. It can be found that the scenarios are specific to the location of the damage opening. Side damage and bottom damage will cause different types of flooding, while the corresponding motion responses are different. In the simulation process, in order to eliminate the composite influence arising from the air compressibility, an appropriate ventilation hole is constructed on the upper deck. In [24][25][26][27], it has been elaborated that the air compression in the flooded compartments will delay the flooding process and affect the dynamic behaviors of the damaged ship. According to the much-simplified assumption in MSC.362 (92) [28], if the total ventilation hole sectional area is 10% or more of the damage opening, the air compression may be neglected and the flooded compartment can be considered to be fully ventilated. Therefore, in order to ensure the adequate ventilation condition, a relatively large ventilation hole is set in the damage scenarios. The detailed dimensions of the damage opening and ventilation hole are shown in Table 2. It is worth noting that, except for the location of the damage opening, the other properties of the hull in the two scenarios are completely consistent, including the size and location of the ventilation hole, the size of the damage opening, the weight, the center of the gravity, and the inertia moments. In this case, it is meaningful to investigate the effect of the damage location on the flooding process and motion responses of the damaged ship. For the characteristics of the hull, a preliminary analysis was conducted by the means of a fine CAD (Computer aided design) model using CATIA CAD software. The density of the low carbon steel is assigned to each part in the CAD environment. Subsequently, the total weight of the hull, the center of mass, and moments of inertia can be accurately calculated, as presented in Table 2. been elaborated that the air compression in the flooded compartments will delay the flooding process and affect the dynamic behaviors of the damaged ship. According to the much-simplified assumption in MSC.362 (92) [28], if the total ventilation hole sectional area is 10% or more of the damage opening, the air compression may be neglected and the flooded compartment can be considered to be fully ventilated. Therefore, in order to ensure the adequate ventilation condition, a relatively large ventilation hole is set in the damage scenarios. The detailed dimensions of the damage opening and ventilation hole are shown in Table 2. It is worth noting that, except for the location of the damage opening, the other properties of the hull in the two scenarios are completely consistent, including the size and location of the ventilation hole, the size of the damage opening, the weight, the center of the gravity, and the inertia moments. In this case, it is meaningful to investigate the effect of the damage location on the flooding process and motion responses of the damaged ship. For the characteristics of the hull, a preliminary analysis was conducted by the means of a fine CAD (Computer aided design) model using CATIA CAD software. The density of the low carbon steel is assigned to each part in the CAD environment. Subsequently, the total weight of the hull, the center of mass, and moments of inertia can be accurately calculated, as presented in Table 2.  After determining the characteristics of the damaged hull and ensuring that the damage location is the single variation, the actual draft of the damaged hull at the corresponding weight needs to be calculated. In the simulation setup, the height of the free surface needs to be consistent with the actual draft height calculated above. Specifying an accurate draft value will ensure that the damaged hull does not instantaneously heave due to the difference between weight and displacement at the moment of release, which also reproduces the actual physical process to some extent. Therefore, before carrying out the damage simulations, it is necessary to specify the weight of the damaged hull for the intact hull, calculating the actual draft in the case of the 665.64 kg displacement. As shown in Figure 6, the initial draft specified by the simulation is 0.3 m. Under the action of gravity and buoyancy, the final actual draft is stable at 0.269737 m. However, in the process of calculating the draft, only the vertical Z-axis motion is released, and other degrees of freedom are restrained. Such  After determining the characteristics of the damaged hull and ensuring that the damage location is the single variation, the actual draft of the damaged hull at the corresponding weight needs to be calculated. In the simulation setup, the height of the free surface needs to be consistent with the actual draft height calculated above. Specifying an accurate draft value will ensure that the damaged hull does not instantaneously heave due to the difference between weight and displacement at the moment of release, which also reproduces the actual physical process to some extent. Therefore, before carrying out the damage simulations, it is necessary to specify the weight of the damaged hull for the intact hull, calculating the actual draft in the case of the 665.64 kg displacement. As shown in Figure 6, the initial draft specified by the simulation is 0.3 m. Under the action of gravity and buoyancy, the final actual draft is stable at 0.269737 m. However, in the process of calculating the draft, only the vertical Z-axis motion is released, and other degrees of freedom are restrained. Such simplification will result in the transient pitch motion of the damaged hull due to the difference in bow and stern weight distribution. The specific simulation results will be analyzed in Section 5. simplification will result in the transient pitch motion of the damaged hull due to the difference in bow and stern weight distribution. The specific simulation results will be analyzed in Section 5.

Simulation Domain and Physical Models
The applied overset mesh technique requires two different regions, including the background region and the overset region, as shown in Figure 7. The overset region is obtained by Boolean operation (subtraction) between the cylinder block and the damaged hull. The overset region will rotate and translate with the movements of the damaged ship. The background region is stationary, only providing the external flow field information. The background region and the overset region are implicit through the interface, while the connectivity between them takes place through the interpolation scheme specified for the interface. The interface mentioned here is the surface of the cylinder block. According to the Mancini et al. [29], Handschel et al. [30], and the available ITTC recommended procedure and guidelines [31], the dimensions of the background region and the overset region are summarized in Table 3. The background region is usually designed in compliance with the "Practical Guidelines for Ship CFD Application" [31]. However, there is no clear specification for the dimension of the overset region, as indicated by Tezdogan et al. [32]. It is worth mentioning that in the process of generating grids, the two regions use their individual mesh continuum to generate grids respectively, and the parameter properties of the two mesh continuums are independent of each other. However, in order to ensure the consistency of the external physical field, one identical physical continuum is applied to define the physical models in two regions.

Simulation Domain and Physical Models
The applied overset mesh technique requires two different regions, including the background region and the overset region, as shown in Figure 7. The overset region is obtained by Boolean operation (subtraction) between the cylinder block and the damaged hull. The overset region will rotate and translate with the movements of the damaged ship. The background region is stationary, only providing the external flow field information. The background region and the overset region are implicit through the interface, while the connectivity between them takes place through the interpolation scheme specified for the interface. The interface mentioned here is the surface of the cylinder block. According to the Mancini et al. [29], Handschel et al. [30], and the available ITTC recommended procedure and guidelines [31], the dimensions of the background region and the overset region are summarized in Table 3. The background region is usually designed in compliance with the "Practical Guidelines for Ship CFD Application" [31]. However, there is no clear specification for the dimension of the overset region, as indicated by Tezdogan et al. [32]. It is worth mentioning that in the process of generating grids, the two regions use their individual mesh continuum to generate grids respectively, and the parameter properties of the two mesh continuums are independent of each other. However, in order to ensure the consistency of the external physical field, one identical physical continuum is applied to define the physical models in two regions. with the "Practical Guidelines for Ship CFD Application" [31]. However, there is no clear specification for the dimension of the overset region, as indicated by Tezdogan et al. [32]. It is worth mentioning that in the process of generating grids, the two regions use their individual mesh continuum to generate grids respectively, and the parameter properties of the two mesh continuums are independent of each other. However, in order to ensure the consistency of the external physical field, one identical physical continuum is applied to define the physical models in two regions.  In the simulations, the behaviors of two fluids (liquid and air) are modelled in the same physical continuum by the volume of fluid (VOF) approach. Due to the presence of two fluids, the Euler multiphase flow model is activated, and the gravity model is used to consider the gravitational effects of two fluids. The liquid phase is modeled with constant density water. In order to reproduce the real physical process, although the air compression is not considered, the air phase is still characterized by an ideal gas model. The necessary user defined field function (UDFF) model is needed to distribute the water and air [33]. Finally, a realizable k-ε two-Layer turbulence model is applied to solve the Reynolds stress problem, which can provide a good compromise between robustness, computational cost, and accuracy [34].

Boundary Conditions and Solver Settings
According to Zhang et al. [33] and Begovic et al. [35], the boundaries of the simulation domain are represented in Figure 8. The chosen boundary conditions and optimal solver settings are presented in Table 4. For a clear description of the internal arrangement of the simulation domain, the boundaries on both sides are not shown. It can be found that by creating isospheric surfaces, the entire simulation domain is divided into two parts by the free surface, including the air part above and the water part below. In order to obtain the sharp interfaces between the air and water, the second-order convection term is recommended. In this case, the high-resolution interface capturing (HRIC) scheme is designed to mimic the convective transport of immiscible fluid components, forming a scheme that is suited for tracking sharp interfaces.  In the simulations, the behaviors of two fluids (liquid and air) are modelled in the same physical continuum by the volume of fluid (VOF) approach. Due to the presence of two fluids, the Euler multiphase flow model is activated, and the gravity model is used to consider the gravitational effects of two fluids. The liquid phase is modeled with constant density water. In order to reproduce the real physical process, although the air compression is not considered, the air phase is still characterized by an ideal gas model. The necessary user defined field function (UDFF) model is needed to distribute the water and air [33]. Finally, a realizable k-ε two-Layer turbulence model is applied to solve the Reynolds stress problem, which can provide a good compromise between robustness, computational cost, and accuracy [34].

Boundary Conditions and Solver Settings
According to Zhang et al. [33] and Begovic et al. [35], the boundaries of the simulation domain are represented in Figure 8. The chosen boundary conditions and optimal solver settings are presented in Table 4. For a clear description of the internal arrangement of the simulation domain, the boundaries on both sides are not shown. It can be found that by creating isospheric surfaces, the entire simulation domain is divided into two parts by the free surface, including the air part above and the water part below. In order to obtain the sharp interfaces between the air and water, the second-order convection term is recommended. In this case, the high-resolution interface capturing (HRIC) scheme is designed to mimic the convective transport of immiscible fluid components, forming a scheme that is suited for tracking sharp interfaces.    The mass conservation equation and the momentum conservation equations including the turbulence model were calculated in the incompressible based unsteady state. For the coupling of velocity and pressure, a semi-implicit method for pressure linked equations (SIMPLE) method was used [36]. In order to increase the convergence performance of the linear algebraic equation, AMG (Algebraic Multi-Grid) method [37] was used and, using the Gauss-Seidel method, the simultaneous linear equation was solved.
The URANS (unsteady Reynolds-averaged Navier-Stokes) equations have been applied to control the update at each physical time for the calculation. In order to converge the solution for that given instant of time, each physical time is set to involve some number of inner iterations. Considering the compromise between the computational accuracy and the computational cost, the simulation program in this paper used a constant time step of 0.002 s, while the maximum number of the inner iteration steps is 10. The numbers determined are also consistent with the related recommendations of practical guidelines for ship CFD applications [31]. In addition, the second-order temporal discretization is activated to perform the transient calculations, which uses the current time level and the solutions from the previous two-time levels. Therefore, when the solver performs first-step calculation utilizing the second-order temporal discretization, the first-order temporal discretization is temporarily activated to provide the solutions of the first two inner iteration steps.

Mesh Type and Mesh Size
For the mesh type, the trimmed hexahedral type is used to generate the mesh. In Begovic et al. [35], detailed sensitivity analysis of the mesh types has been performed on the roll damping assessments of the damaged ship with two hybrid meshes (polyhedral and trimmed) and two trimmed meshes. The simulation results indicate that the hybrid meshes are prohibitive due to the high time consumption and poor simulation accuracy, while the trimmed meshes are recommended. Based on this conclusion, the generated mesh in this paper is shown in Figure 5. It can be found that the entire domain is divided into three regions with different mesh densities, including the background region, the cylindrical overset region, and the overlapping region. In order to optimize the discretization of the overset region, the mesh density of the overset region is denser than the other two regions. To minimize the errors that occur when interpolating variables between two meshes, the same mesh density order of magnitude is used in the overlapping region and the background region. Therefore, it can be seen from Figure 9 that the mesh density of the overlapping region is denser than the background region and coarser than the overset region. However, fundamentally, the overlapping region is still part of the background area. It is just a refinement block extracted from the background region. In addition, the meshes around the free surface are also locally refined, which can prevent the floating-point exception due to the free surface fluctuation or breakage. Finally, in order to avoid large computational costs, the free surface of the overset region is finer than that of the background region. The mesh sizes in different parts are summarized in Table 5.
recommendations of practical guidelines for ship CFD applications [31]. In addition, the second-order temporal discretization is activated to perform the transient calculations, which uses the current time level and the solutions from the previous two-time levels. Therefore, when the solver performs firststep calculation utilizing the second-order temporal discretization, the first-order temporal discretization is temporarily activated to provide the solutions of the first two inner iteration steps.

Mesh Type and Mesh Size
For the mesh type, the trimmed hexahedral type is used to generate the mesh. In Begovic et al. [35], detailed sensitivity analysis of the mesh types has been performed on the roll damping assessments of the damaged ship with two hybrid meshes (polyhedral and trimmed) and two trimmed meshes. The simulation results indicate that the hybrid meshes are prohibitive due to the high time consumption and poor simulation accuracy, while the trimmed meshes are recommended. Based on this conclusion, the generated mesh in this paper is shown in Figure 5. It can be found that the entire domain is divided into three regions with different mesh densities, including the background region, the cylindrical overset region, and the overlapping region. In order to optimize the discretization of the overset region, the mesh density of the overset region is denser than the other two regions. To minimize the errors that occur when interpolating variables between two meshes, Figure 9. Visualization of the hexahedral mesh. In addition to optimizing the mesh sizes of the different regions described above, the mesh quality of the damaged ship also determines the calculation accuracy. As shown in Table 5, the wrapper model, the remesh model and the trim model are activated in the mesh continuum. The appropriate wrapper and remesh sizes will restore the original geometry of the damaged ship. It can also be found that the wrapper model is not activated in the mesh continuum of the background region. This is because the damaged ship is located in the overset region, and only the wrapper characteristics in the mesh continuum of the overset region can be defined to the damaged ship. The setting of a mesh size requires a comprehensive combination of mesh quality and the mesh number. Too fine meshes will cause a longer computation time, while too coarse meshes will result in poor computation convergence. Taking the mesh size of the overlapping region as an example, when the trim size is set to 0.100 m, the mesh number will be much larger than when the mesh size is set to 0.110 m. Consequently, the time consumption is very high. So, in order to achieve the balance between the simulation accuracy and the time consumption, the final trim size is to set 0.110 m. Similarly, the mesh sizes of other parts are determined after repeated attempts. Especially for the damage opening, the ventilation hole, and the bulbous bow, since their outer contours display a large curvature, the wrapper and remesh sizes need to be modified repeatedly to ensure the real opening shape. For the flooded compartment, the trim size is locally refined. Because only when the mesh is fine enough, the complex hydrodynamics phenomena in the flooding process can be accurately captured. The generated surface of the volume mesh is shown in Figure 10. It can be seen from the figure that the original surface of the damaged ship is restored well, and the necessary block is correspondingly refined.  In addition to optimizing the mesh sizes of the different regions described above, the mesh quality of the damaged ship also determines the calculation accuracy. As shown in Table 5, the wrapper model, the remesh model and the trim model are activated in the mesh continuum. The appropriate wrapper and remesh sizes will restore the original geometry of the damaged ship. It can also be found that the wrapper model is not activated in the mesh continuum of the background region. This is because the damaged ship is located in the overset region, and only the wrapper characteristics in the mesh continuum of the overset region can be defined to the damaged ship. The setting of a mesh size requires a comprehensive combination of mesh quality and the mesh number. Too fine meshes will cause a longer computation time, while too coarse meshes will result in poor computation convergence. Taking the mesh size of the overlapping region as an example, when the trim size is set to 0.100 m, the mesh number will be much larger than when the mesh size is set to 0.110 m. Consequently, the time consumption is very high. So, in order to achieve the balance between the simulation accuracy and the time consumption, the final trim size is to set 0.110 m. Similarly, the mesh sizes of other parts are determined after repeated attempts. Especially for the damage opening, the ventilation hole, and the bulbous bow, since their outer contours display a large curvature, the wrapper and remesh sizes need to be modified repeatedly to ensure the real opening shape. For the flooded compartment, the trim size is locally refined. Because only when the mesh is fine enough, the complex hydrodynamics phenomena in the flooding process can be accurately

Near-Wall Treatment
The wall function approach is used for the near-wall treatment, in particular, the All wall Y + model. This approach is formulated to assure reasonable answers for meshes of intermediate resolution and is considered as the best compromise between description of the boundary layer with acceptable quality and the time required for the calculation [29]. The wall y + is a non-dimensional distance similar to the local Reynolds number, often used in CFD to describe how coarse or fine a mesh is for a particular flow [38]. As indicated in the User's Guide [21], values of y+ ≈ 30 are most desirable for wall functions, whereas values of y+ ≈ 1 are most desirable for near-wall modeling. The values of wall y + on the hull surface is shown in Figure 11. It can be found that the y + values on the hull are very close to 1. For this reason, the realizable k-ε two layers turbulence model is applied. This turbulence model represents an improved treatment of the near-wall region for turbulent flows at low Reynolds numbers. This model is characterized for the layer next to the wall, where the turbulent dissipation rate and the turbulent viscosity are specified as functions of wall distance. More details about this model are available in Rodi [39]. captured. The generated surface of the volume mesh is shown in Figure 10. It can be seen from the figure that the original surface of the damaged ship is restored well, and the necessary block is correspondingly refined.

Near-Wall Treatment
The wall function approach is used for the near-wall treatment, in particular, the All wall Y + model. This approach is formulated to assure reasonable answers for meshes of intermediate resolution and is considered as the best compromise between description of the boundary layer with acceptable quality and the time required for the calculation [29]. The wall y + is a non-dimensional distance similar to the local Reynolds number, often used in CFD to describe how coarse or fine a mesh is for a particular flow [38]. As indicated in the User's Guide [21], values of y 30 are most desirable for wall functions, whereas values of y 1 are most desirable for near-wall modeling.
The values of wall y + on the hull surface is shown in Figure 11. It can be found that the y + values on the hull are very close to 1. For this reason, the realizable k-ε two layers turbulence model is applied. This turbulence model represents an improved treatment of the near-wall region for turbulent flows at low Reynolds numbers. This model is characterized for the layer next to the wall, where the turbulent dissipation rate and the turbulent viscosity are specified as functions of wall distance. More details about this model are available in Rodi [39].

Simulation Results and Discussion
Based on the appropriate solver settings and the optimized mesh generation, the simulation results respectively analyze the effect of damage locations on the flooding process and the motion responses. The complex hydrodynamic behaviors in the flooding process are visualized. The corresponding motion responses are compared and discussed, including the roll and pitch motion. Figure 12 presents the distribution of the flooding water at different time points. In the early stage of the damage flooding, due to large pressure difference between inside and outside of the damaged opening, the water ingress flows into the flooded compartment with a jet form in the opening section. The water ingress impacts the bottom plate and the longitudinal bulkhead, resulting in the splashing of the flooding water. Since the flooding water is violent, the compressed air in the damaged compartment cannot smoothly escape from the ventilation hole, and the complex hydrodynamic behavior such as a bubble is formed in the flooding water. This stage is complex but short, which is often defined as the transient flooding stage. As shown in Figure 12, the period from 0 s to 4.9 s can be roughly referred to as the transient flooding stage. With the effect of asymmetric water ingress, the damaged ship heels towards the starboard side, causing the pressure at the opening section to become larger. From the four graphs corresponding to 1, 2.9, 4, and 4.9 s, it can be found that the increased pressure at the opening section makes the slamming point (1, 2, 3, 4) on the

Simulation Results and Discussion
Based on the appropriate solver settings and the optimized mesh generation, the simulation results respectively analyze the effect of damage locations on the flooding process and the motion responses. The complex hydrodynamic behaviors in the flooding process are visualized. The corresponding motion responses are compared and discussed, including the roll and pitch motion. Figure 12 presents the distribution of the flooding water at different time points. In the early stage of the damage flooding, due to large pressure difference between inside and outside of the damaged opening, the water ingress flows into the flooded compartment with a jet form in the opening section. The water ingress impacts the bottom plate and the longitudinal bulkhead, resulting in the splashing of the flooding water. Since the flooding water is violent, the compressed air in the damaged compartment cannot smoothly escape from the ventilation hole, and the complex hydrodynamic behavior such as a bubble is formed in the flooding water. This stage is complex but short, which is often defined as the transient flooding stage. As shown in Figure 12, the period from 0 s to 4.9 s can be roughly referred to as the transient flooding stage. With the effect of asymmetric water ingress, the damaged ship heels towards the starboard side, causing the pressure at the opening section to become larger. From the four graphs corresponding to 1, 2.9, 4, and 4.9 s, it can be found that the increased pressure at the opening section makes the slamming point (1, 2, 3, 4) on the longitude move upward. Such a slamming effect will produce a restoring moment that makes the damaged ship heel towards the port side. When the flooding water develops to a certain extent, the pressure difference between the inside and outside of the damaged compartment will gradually smaller, the flooding will become slower. This stage is often referred to as the progressive flooding stage, as shown in Figure 12 for the period from 4.9 s to 20 s. In this stage, the flooding water continuously flooded the damaged compartment. The free surface exhibits a wave-propagating form, producing a reflective behavior when it touches the longitudinal bulkhead. Although the flooding process is almost completed in about 20 s, from the capture of the flooding process at 30 s, the free surface is still sloshing due to the roll motion of the damaged ship. Conversely, the sloshing of the free surface also affects the motion response of the damaged ship. Finally, if the damaged ship can keep afloat, not capsizing or sinking due to the added flooding water, the final equilibrium state will be characterized. This stage is often referred to as the steady stage. Such detailed descriptions of the hydrodynamic behaviors in the flooding process can be applied to explain the causes of the damaged ship's motion responses. The specific and comprehensive explanation is elaborated in the motion response analysis in Section 5.3.1.

The Analysis of the Flooding Process in the Bottom Damage Scenario
After comparison, there are both similarities and differences between the side flooding process illustrated in Figure 12 and the bottom flooding process illustrated in Figure 13. The overall similarity is that the bottom flooding process also experiences three flooding stages, including the transient stage, the progressing stage, and the steady stage. In the transient flooding stage (0-2.0 s), due to the large pressure difference between the inside and outside of the damaged opening, the violent seawater flooded the damaged compartment in a short time. Then, with the accumulation of the flooding water, the pressure difference gradually reduces while the flooding water occupied the damaged compartment at a slow rate (2.0-15.2 s). Finally, under the coupled influence of the tank sloshing, the damaged ship tends to be stable in the roll decay motion (15.

The Analysis of the Flooding Process in the Bottom Damage Scenario
After comparison, there are both similarities and differences between the side flooding process illustrated in Figure 12 and the bottom flooding process illustrated in Figure 13. The overall similarity is that the bottom flooding process also experiences three flooding stages, including the transient stage, the progressing stage, and the steady stage. In the transient flooding stage (0-2.0 s), due to the large pressure difference between the inside and outside of the damaged opening, the violent seawater flooded the damaged compartment in a short time. Then, with the accumulation of the flooding water, the pressure difference gradually reduces while the flooding water occupied the damaged compartment at a slow rate (2.0-15.2 s). Finally, under the coupled influence of the tank sloshing, the damaged ship tends to be stable in the roll decay motion (15.2 s-30.0 s). At 24.6 s and 30.0 s, the internal wave propagation caused by the sloshing of the free surface can be clearly seen in the damaged compartment. At the same time, the flooding characteristics in the bottom damage scenario are also evident. Corresponding to the normal direction of the bottom opening section, the flooding water is sprayed from the ship bottom in the form of a water column. When the flooding water slams the longitudinal bulkhead and reaches the highest point at 0.3 s, the flooding waterfalls under the influence of gravity. Complex dynamic behaviors can be observed when the falling flooding water touches the ship bottom from 0.5 to 1.2 s, including splashing and bubble phenomena. In addition, in contrast to the side damage scenario, the damage opening in the bottom damage scenario is located deeper below the waterline. In this case, the hydrostatic pressure at the bottom opening section is much greater than that at the side bottom opening section. This also explains why under the premise of the same damage opening size, the bottom flooding process is completed in about 15 s, however, the side damage flooding takes about 20 s. Although the extra 5 s is relatively short, according to the Froude law, the converted flooding time for the real ship longer is than on the scale model. All simulation cases in this paper were carried out on the scale model (1/25). Therefore, the extra 5 s is about 25 s when converted to the full-scale ship. Once the damage accident occurs, 25 s can provide more rescue options. Therefore, accurately predicting the flooding time in different damage scenarios is meaningful for the emergency crew to take appropriate rescue managements. Generally, this detailed visualization of the hydrodynamic behaviors is helpful to enhance the understanding of the entire flooding process among the crew and ship designers.

Description of the Roll Motion Response
Based on the hydrodynamic behavior in Sections 5.1. and 5.2., the resulting motion responses are elaborated in this section. As shown in Figure 14, although the damage locations of the two damage scenarios are different, the asymmetric flooding occurs in both damage scenarios. The asymmetric moment generated by the flooding water causes the damaged ship to heel only in the starboard. There is no periodic reciprocating roll motion between the portside and starboard. According to the right-hand rule, the value of the heel angle is negative. The heel angle can reach 15 degrees or more, indicating that the asymmetric flooding in a damaged ship is a dangerous situation. Therefore, efficient and feasible cross-flooding arrangements are needed to provide the necessary equalization across the ship in order to decrease the heel angle [40]. After a separate analysis of the roll motion curves in the two damage scenarios, it can be found that the damaged ship does not always heel towards the starboard, but gradually heels during left and right shaking. On the one hand, this is due to the inherent restoring moment of the hull itself. On the other hand, the flooding water causes a leftward impact on the longitudinal bulkhead, which also causes the damaged ship to have a tendency to heel towards the portside. As can be seen from Figure 12, the side flooding water

Description of the Roll Motion Response
Based on the hydrodynamic behavior in Section 5.1. and Section 5.2, the resulting motion responses are elaborated in this section. As shown in Figure 14, although the damage locations of the two damage scenarios are different, the asymmetric flooding occurs in both damage scenarios. The asymmetric moment generated by the flooding water causes the damaged ship to heel only in the starboard. There is no periodic reciprocating roll motion between the portside and starboard. According to the right-hand rule, the value of the heel angle is negative. The heel angle can reach 15 degrees or more, indicating that the asymmetric flooding in a damaged ship is a dangerous situation. Therefore, efficient and feasible cross-flooding arrangements are needed to provide the necessary equalization across the ship in order to decrease the heel angle [40]. After a separate analysis of the roll motion curves in the two damage scenarios, it can be found that the damaged ship does not always heel towards the starboard, but gradually heels during left and right shaking. On the one hand, this is due to the inherent restoring moment of the hull itself. On the other hand, the flooding water causes a leftward impact on the longitudinal bulkhead, which also causes the damaged ship to have a tendency to heel towards the portside. As can be seen from Figure 12, the side flooding water continuously impacts the longitudinal bulkhead from the starboard. However, the horizontal impact effect due to side flooding is small compared to the vertical effect of the flooding water accumulation on the starboard bottom plate. In this case, the damaged ship only has a slight tendency to heel towards the portside, and never has a positive heel angle. At the same time, since the normal direction of the bottom opening section points to the longitudinal bulkhead, the water column in Figure 13 inevitably have an impact effect on the longitudinal bulkhead, so that the same shaking phenomenon as the side damage scenario occurs in the roll motion curve of the bottom damage scenario. In addition, the roll motion of the damaged ship and the flooding water affect each other. The roll motion makes the flooding water slosh in the flooded compartment. Conversely, the water sloshing has an impact effect on the internal bulkhead, including the longitudinal bulkhead and hull plate. Such coupled motion can present a risk to the ship survivability, even making the damaged ship capsize due to the parametric roll motion. This coupled analysis also explains why the damage scenarios still have a roll motion, even though the flooding water is no longer increased. However, due to the dissipation of energy, the entire roll motion will gradually decay. has an impact effect on the internal bulkhead, including the longitudinal bulkhead and hull plate. Such coupled motion can present a risk to the ship survivability, even making the damaged ship capsize due to the parametric roll motion. This coupled analysis also explains why the damage scenarios still have a roll motion, even though the flooding water is no longer increased. However, due to the dissipation of energy, the entire roll motion will gradually decay. Finally, by comprehensively comparing and analyzing the roll motion of the two damage scenarios, it can be found that different hydrodynamic behaviors produce distinct motion responses. The roll motion curves of the two damaged scenarios follow the similar periodic variation rule. However, there are certain differences in the peak and trough values for the same period. In Figure  14, the peaks and troughs in the same periods are connected by the black line segments and the yellow line segments. The magnitude of the peak value represents the extent to which the damaged ship heels towards the portside (the intact side), and the magnitude of the trough value represents the extent to which the damaged ship heels towards the starboard (the damaged side). It can be seen that the bottom damage flooding in the same period produces a larger heel angle with respect to the side damage flooding, regardless of the peak value or trough value. The reason for this difference is close to the hydrodynamic behaviors of the specified damage scenario. For the bottom damage flooding, the vertical effect of the upward flooding on the bottom plate is much larger than the horizontal effect of the water column on the longitudinal bulkhead. In this case, in comparison with the side damage scenario, the damaged ship with the bottom opening has a greater inclination towards the starboard. This also corresponds to the fact that the trough values of the bottom damage scenario are below those of the side damage scenario. For the side damage flooding, the flooding water strikes the longitudinal bulkhead from the starboard to the portside. The horizontal impact drives the damaged ship heel towards the portside. This explains why the peak values of the side damage scenario in the same period are above the bottom damage scenario. In general, the visualization of the flooding process is very helpful and meaningful for analyzing the causes of motion responses.

Description of the Pitch Motion Response
For analyzing the influence of the damage location on the pitch motion, Figure 15 presents the Finally, by comprehensively comparing and analyzing the roll motion of the two damage scenarios, it can be found that different hydrodynamic behaviors produce distinct motion responses. The roll motion curves of the two damaged scenarios follow the similar periodic variation rule. However, there are certain differences in the peak and trough values for the same period. In Figure 14, the peaks and troughs in the same periods are connected by the black line segments and the yellow line segments. The magnitude of the peak value represents the extent to which the damaged ship heels towards the portside (the intact side), and the magnitude of the trough value represents the extent to which the damaged ship heels towards the starboard (the damaged side). It can be seen that the bottom damage flooding in the same period produces a larger heel angle with respect to the side damage flooding, regardless of the peak value or trough value. The reason for this difference is close to the hydrodynamic behaviors of the specified damage scenario. For the bottom damage flooding, the vertical effect of the upward flooding on the bottom plate is much larger than the horizontal effect of the water column on the longitudinal bulkhead. In this case, in comparison with the side damage scenario, the damaged ship with the bottom opening has a greater inclination towards the starboard. This also corresponds to the fact that the trough values of the bottom damage scenario are below those of the side damage scenario. For the side damage flooding, the flooding water strikes the longitudinal bulkhead from the starboard to the portside. The horizontal impact drives the damaged ship heel towards the portside. This explains why the peak values of the side damage scenario in the same period are above the bottom damage scenario. In general, the visualization of the flooding process is very helpful and meaningful for analyzing the causes of motion responses.

Description of the Pitch Motion Response
For analyzing the influence of the damage location on the pitch motion, Figure 15 presents the pitch motion curves of the two damage scenarios. From the overall analysis, the maximum values of the pitch angle for the two damage scenarios are only about 2 degrees in the flooding process. Conversely, the maximum value of the heel angle can reach about 15 degrees. This validates a basic conclusion that the damaged ships rarely lose stability due to the excessive pitch angle. The damaged ships often capsized due to the excessive heel angle caused by the additional flooding water. Therefore, this also provides an empirical reference for the ship designers and the emergency personnel onboard. When the damage occurs, especially for asymmetric flooding, in order to reduce the risk of capsizing or sinking, appropriate countermeasures should be taken to equalize the heel angle. Although the pitch angle caused by the flooding water does not pose a threat to the safety of the damaged ship, its relevant effect deserves enough attention. Especially for the military ships, even if the ship is in a damaged situation, the ship must guarantee the corresponding operational missions and has to able to recover functionality following an incident (recoverability). Hence, the resulting 2 degree pitch angle is likely to affect the accuracy of the weapon strikes. So, the way to eliminate the extra pitch angle is of research significance. By comparing the pitch motion curves of the two damage scenarios, it can be found that the two pitch motion curves maintain the same variation rule, and the differences of the peak and trough values are small. This shows that the effect of the damage location on the pitch motion is small, even can be negligible. For analyzing the tendency of the pitch motion curves, it can be divided into two parts. The first part is that the damaged ship will have a transient head-pitching process at the beginning of the calculation. Because in the simulation settings, the damaged ship is placed horizontally by default. When the simulation runs, the damaged ship firstly have the initial pitch motion due to the uneven distribution of the fore and bow weights. Therefore, it can be seen from the Figure 15 that the pitch motion curves of the two damage scenarios are basically consistent at the beginning of the damage flooding, which is caused by the ship's own weight distribution described above and has little to do with the hydrodynamic behavior. The second part is that, as the flooding continues to develop, the flooding water gradually flooded the damaged compartment. Due to the strong nonlinear phenomenon of the flooding water, the flooding water spread irregularly in the damaged compartment. This is why the predicting pitch motion curves fluctuate up and down. However, the amplitude of the peak and trough values are very small. In this case, the crew onboard will not have the sloshing feeling caused by the pitch motion. Therefore, from the perspective of ensuring the survivability of the hull and the safety of human life, the threat posed by the roll motion is the problem that should be solved first.

Conclusion and Future Researches
The paper demonstrated the feasibility of CFD simulations to investigate complex flooding phenomena. The developed numerical approach in this paper can well capture the complex hydrodynamic behavior in the flooding process, including splash, jet, water column, and bubble. The Although the pitch angle caused by the flooding water does not pose a threat to the safety of the damaged ship, its relevant effect deserves enough attention. Especially for the military ships, even if the ship is in a damaged situation, the ship must guarantee the corresponding operational missions and has to able to recover functionality following an incident (recoverability). Hence, the resulting 2 degree pitch angle is likely to affect the accuracy of the weapon strikes. So, the way to eliminate the extra pitch angle is of research significance. By comparing the pitch motion curves of the two damage scenarios, it can be found that the two pitch motion curves maintain the same variation rule, and the differences of the peak and trough values are small. This shows that the effect of the damage location on the pitch motion is small, even can be negligible. For analyzing the tendency of the pitch motion curves, it can be divided into two parts. The first part is that the damaged ship will have a transient head-pitching process at the beginning of the calculation. Because in the simulation settings, the damaged ship is placed horizontally by default. When the simulation runs, the damaged ship firstly have the initial pitch motion due to the uneven distribution of the fore and bow weights. Therefore, it can be seen from the Figure 15 that the pitch motion curves of the two damage scenarios are basically consistent at the beginning of the damage flooding, which is caused by the ship's own weight distribution described above and has little to do with the hydrodynamic behavior. The second part is that, as the flooding continues to develop, the flooding water gradually flooded the damaged compartment. Due to the strong nonlinear phenomenon of the flooding water, the flooding water spread irregularly in the damaged compartment. This is why the predicting pitch motion curves fluctuate up and down. However, the amplitude of the peak and trough values are very small. In this case, the crew onboard will not have the sloshing feeling caused by the pitch motion. Therefore, from the perspective of ensuring the survivability of the hull and the safety of human life, the threat posed by the roll motion is the problem that should be solved first.

Conclusion and Future Researches
The paper demonstrated the feasibility of CFD simulations to investigate complex flooding phenomena. The developed numerical approach in this paper can well capture the complex hydrodynamic behavior in the flooding process, including splash, jet, water column, and bubble. The URANS solver involving the overset mesh technique is applied to monitor the motion responses of the damaged ship with different damage locations. Through analysis, the visualization of the flooding process can be efficient to explain the cause of the resulting motion response. Comparing the bottom damage and the side damage, the upward flooding in the bottom damage scenario causes a larger vertical impact on the bottom plate while the flooding water from starboard to portside in the side damage scenario causes a larger horizontal impact on the longitudinal bulkhead. This detailed visualization description explains why the bottom damage flooding in the same period produces a larger heel angle with respect to the side damage flooding, regardless of the peak value or trough value. In addition, due to the coupled motion between the damaged ship and the flooding water, the tank sloshing makes the damaged ship still roll even if the flooding water is not increased. And, the wave propagation in the flooded is seen clearly. After comparing the roll motion and pitch motion, it can be found that the symmetric flooding is a dangerous situation, even making the damaged ship capsize. Though the damage locations are different, the asymmetric flooding produces an excessive heel angle. However, asymmetric flooding has little effect on the pitch angle. Based on these summaries, the final suggestion is that when the damage flooding takes place, especially asymmetric flooding, appropriate counter measures should be taken first to equalize the heel angel. In this case, the damaged ship can keep a good floating state, which is helpful to improve the survivability of the damaged ship and ensure the safety of human life. From the perspective of the simulation validation, only numerical simulation results are introduced, without validation with experimental results. Subsequently, the specified model tests will be carried out to prove the reliability of the applied numerical simulation approach.
In the future, a more realistic compartment arrangement will be created, considering the effect of the permeability on the flooding process and motion responses. Also, more sea conditions will be included, including wind, wave, and forward speed. This gradual improvement process will also present new challenges to the current simulation approach.