A comparison of numerical methods for damage index based residual ultimate limit state assessment of grounded ship hulls

Considerable efforts have been devoted to developing rapid methodologies for predicting the residual strength of ship hull girders for a given damage scenario (e.g., R–D diagram). This task is usually challenged by the difficulty of a proper definition of damaged scenarios, which is a function of the damage location and extent. The concept of damage index (DI) was proposed to resolve this issue, and its application has been demonstrated previously through the incorporation of the Intelligent Supersize Finite Element Method (ISFEM) and Modified Paik–Mansour (P–M) approach. Alternatively, a Smith-type progressive collapse can be adopted. In fact, this may enable a more comprehensive application of the DI concept, as the Smith-type approach is codified in the Common Structural Rule (CSR) for assessing the longitudinal strength of the hull girder. In light of this, the damage index based assessment tools (i.e., R–D diagram) incorporated with a Smith-type progressive collapse method is presented in this paper. This also facilitates a comparison with the R–D diagram previously developed by ISFEM and Modified P–M approach. Discussion is given regarding the discrepancy between different methods, and recommendations for future research are outlined.

One of the key steps in the required action is to examine the remaining capacity (= residual ultimate longitudinal strength) of damaged ship hulls [40][41][42][43] so as to comply with the ultimate limit state (ULS) criterion [44][45][46][47][48][49][50][51]. Since a direct strength calculation requires dedicated expertise, marine structural engineers seek to devise various closedform formulas and/or design diagrams for efficient first-cut estimation. However, this is often challenged by the appropriate definition of a damaged scenario. Ship grounding, which is the focus of this work, is a complex accidental event that could damage the ship hull at random locations with different shapes. * Corresponding authors.
An innovative damage index (DI) concept was proposed by Paik et al. [52] to address these issues concerning the condition assessment for damaged structures, and they validated its applicability for grounding damaged oil tankers. In this damage index technique, a reliable number of damage scenarios are sampled from the probabilistic distributions of relevant random variables. The sampled variables, including damaged location and extent, define the damage index. The corresponding residual ultimate strength is then computed by analytical or numerical methods for all selected scenarios to develop the Residual strength versus Damage index (R-D) diagram. Based on this, the condition assessment of grounding damaged oil tankers was performed, and practical R-D diagrams were proposed, including its user guide [52]. The same concept was applied to container ships [53] and bulk carriers [54]. The secondary effects, i.e., the ageing deterioration [55] and the low-temperature condition [56], were investigated by continued studies. The use of the R-D diagram was also extended to the vessels damaged by collision accidents [57].
In the previous works abovementioned, the residual ultimate strength was calculated by Modified Paik-Mansour Method [52] and by Intelligent Supersize Finite Element Method (ISFEM) [53][54][55][56][57]. A comparison of the R-D diagrams developed by these two methods, i.e., modified P-M method and ISFEM, was also conducted by Kim et al. [58]. An alternative approach is the Smith-type progressive collapse method [59], which is a prevailing approach codified in the Harmonised Common Structural Rules (CSR-H), issued by IACS for bulk carriers and oil tankers [60], and container ships [61]. Many kinds of research are available in the literature on the development and application of the Smith method for a range of ultimate strength analyses of hull girders [62][63][64][65][66][67][68][69][70][71][72][73]. However, the integration between the Smith method and damage index concept-based residual strength assessment for grounded ships has not yet been completed. A more comprehensive application of the DI-based concept to deal with the residual ultimate strength assessment of ship hull girders can be realised by incorporating with the Smith-type progressive collapse method. In the meantime, whist all these approach (ISFEM, modified P-M and Smith method) are wellestablished for predicting the ship hull girder strength, quantifying their computational discrepancy is important. However, a thorough comparison between the ISFEM, modified P-M approach and Smith method seems to be lacking in the literature.
In the light of this, this paper contributes the R-D diagram developed from a Smith-type progressive collapse method for four classes of double-hull oil tankers, i.e. VLCC, Suezmax, Aframax, and Panamax. A comparison of the developed R-D diagrams based on the ISFEM, modified P-M and Smith method is presented. The remaining part of the present study is structured as follows. Section 2 reviews the fundamentals of three analysis methods, i.e., ISFEM, Modified P-M method and Smith method. Section 3 summarises the key steps in the damage index concept-based condition assessment procedure. In Section 4, the development of the R-D diagram based on the Smith method is presented, including the comparison results with other methods, i.e., ISFEM and Modified P-M method. Section 5 elaborates the insights developed from the comparison between the three methods. The conclusions drawn from the present study are documented in Section 6.

Background
This section provides a brief summary of the three prevailing numerical methods to compute the ultimate strength of ship hull girders, namely the Intelligent Supersize Finite Element Method, Modified Paik-Mansour method and Smith-type progressive collapse method. Their fundamentals, common features and distinct differences are highlighted.

Intelligent Supersize Finite Element Method (ISFEM)
The Intelligent Supersize Finite Element Method (ISFEM) is developed to resolve the high computational requirement of the conventional finite element method (FEM) [2]. In contrast to the conventional FEM, the ISFEM is said to be intelligent because the highly nonlinear behaviour of the large elements is educated or formulated in advance, generating a high level of intelligence in terms of judging the failure status and modes of such a large element. The ISFEM overcomes the high computational requirement in conventional FEM but still retains the versatility to analyse the nonlinear inelastic large deflection behaviours of different steel and aluminium structures, such as ship hulls, bridges, cranes, and others [74,75]. Pioneering development of this approach, which was named the Idealised Structural Unit Method (ISUM), was reported by Ueda and Rashed [76,77]. The ISFEM is implemented into ALPS/HULL computer simulation code for the progressive collapse analysis of hull girders [78].
Two common elements of ISFEM are illustrated in Fig. 1, namely the beam-column unit and rectangular plate unit. The former is usually adopted to simulate the nonlinear behaviour of stiffener without attached plating, and the latter is generally employed for predicting the nonlinear response of locally attached plating. The beam-column unit is formulated with three degrees of freedom at each node: The rectangular plate unit is formulated with six degrees of freedom at each node as given in Box I. The general relationship between nodal displacement and nodal force is given by where [ ] is the tangent stiffness matrix that depends on the nodal displacement. A complete formulation of the stiffness matrix is documented in [2]. As the applied load increases, the ISFEM element will undergo buckling and/or yielding. Thus, the element tangent stiffness is updated considering these failure modes. In principle, the procedure is similar to the conventional finite element method. The strain increment will first be calculated for a given increment of displacement vector via the strain-displacement relationship. Subsequently, the membrane stress increment will be evaluated using a pre-defined and efficient stress-strain relationship that encapsulates both geometric and material nonlinearities. This is also a distinct feature of ISFEM as compared with conventional FEM. Finally, an updated Lagrangian approach is applied to formulate the tangent stiffness matrix. Similar to the conventional FEM, the element stiffness matrix [ ] of ISFEM will be transformed from local coordinate to global coordinate via the transformation matrix [ ]: The stiffness matrix of all elements in the global coordinate is assembled to formulate the overall stiffness matrix and the stiffness equation for the target structures. The structural response can be estimated by solving the stiffness equation of the entire structure with prescribed boundary conditions and loading.

Modified Paik-Mansour method
The Modified Paik and Mansour method (= modified P-M method) is usually categorised as the presumed stress distribution-based method. Pioneering development was presented by Caldwell [79], which is also the first systematic approach to predict the ultimate longitudinal strength of ships. In this method, a stress distribution pattern at the ultimate limit state of the ship hull should be pre-defined. An example of typical bending stress distribution across the cross-section of a ship's hull at the ultimate limit state under a hogging bending moment is shown in Fig. 2.
The presumed stress distribution is usually divided into several regions: compressive buckling collapse region, tensile yielding region and linear elastic region. The heights of different regions are estimated by imposing zero axial force condition, i.e., The vertical coordinate of the neutral axis where the resultant stress is zero can be estimated by The ultimate bending moment is calculated by taking the first moment of the resultant stresses about the neutral axis as given by the following: In the original Paik-Mansour method [81], the only unknown is the height of the compressive buckling collapse region. Thus, Eq. (6) is sufficient to determine the unknown. The modified Paik-Mansour method was adapted from the original P-M approach. It improves the original approach by considering the expansion of the yielding zone  There are two unknowns in this modified method, i.e. the height of the compressive buckling region and the height of the tensile yielding region. Eq. (6) is insufficient for two unknowns, and an iterative process is required. In principle, the cross-section will be sub-divided into structural components, i.e., plate-stiffener combination (PSC) or platestiffener separation (PSS). The heights of the compressive and tensile failure regions are searched iteratively until the minimal difference between the compression and tensile axial forces is found. Once the height of the compressive and tensile failure regions are determined, the computation of the ultimate bending strength is completed by Eq. (8). Details may be referred to [82].

Smith-type progressive collapse method
The Smith-type progressive collapse method is a generalisation of the elementary beam theory to consider the influence of local compressive failure due to buckling. With the pioneering paper published by Smith [59], the Smith method was formalised and validated with smallscale testing by Dow et al. [83]. In the early days, the Smith method was developed in a descriptive format. Smith and Dow subsequently introduced a stiffness equation to deal with the bi-axial bending problem, which becomes the governing equation of this class of method as  presented in Eq. (9).
where and are the horizontal and vertical bending moment, respectively. and are the horizontal and vertical curvature, respectively. and are the vertical and horizontal bending stiffness, respectively. and are the interactive bending stiffness. Regarding the ultimate strength of damaged ship hulls, this bi-axial bending equation should be adopted to allow for the rotation of the neutral axis. Fujikubo et al. [84] derived the solution to this bi-axial bending for a damaged ship hull. With a vertical curvature increment ( ), the increment of vertical bending moment ( ) is given as In an unconstrained vertical bending case, horizontal bending might be induced due to the asymmetric damage in the cross-section.
To evaluate the four bending stiffness terms (i.e., , , and ), a sub-division of the ship hull girder is needed, as illustrated in Fig. 4. Plate-stiffener combination (PSC) or plate-stiffener separation (PSS) techniques can be utilised for the cross-section sub-division. In this paper, the PSC technique is applied. A load-shortening curve (LSC) characterising the average stress versus average strain response under in-plane loading will be assigned to each structural component. For all structural components, their tensile load-shortening behaviour is defined as elastic-perfectly plastic. The inelastic buckling behaviour will be considered for most structural components in terms of compressive behaviour, except the so-called hard corner elements. It is assumed that the hard corner element will not experience any buckling, and therefore its compressive behaviour follows the elastic-perfectly plastic. The pre-defined load-shortening curves evaluate the tangent stiffness of structural component at a given applied strain With the tangent stiffness of all components, the instantaneous neutral axis of the cross-section can be estimated as follows: where and are the vertical and horizontal coordinates of the centroid of the neutral axis, respectively. , , and are the elemental tangent stiffness, cross-sectional area, vertical coordinate and horizontal coordinate, respectively. With the elemental tangent stiffness and instantaneous neural axis of hull girder cross-section, the bending stiffness can be computed as As given by Eq. (10), the bending moment increment can be solved. To continue the iteration, the resultant strain of each structural component due to the bending moment increment, which drives the update of element tangent stiffness and instantaneous neutral axis, is calculated as follows: As can be seen in Eq. (10), the interactive stiffness terms are included in calculating the vertical bending moment increment for a given curvature increment. Additionally, a horizontal curvature is induced as a result of the applied vertical curvature (Eq. (11)), which will also be taken into account in the applied strain increment on each element (Eq. (18)). Meanwhile, caused by the horizontal bending (in an indirect way), the position of the horizontal neutral axis will be updated at each incremental step in a similar way as the vertical neutral axis. These effectively account for the effects of the neutral axial translation in both directions and its rotation around the principal vertical and horizontal axes. Note that an explicit calculation of the rotation angle is not required, as the bending curvature is still applied about the upright axes (i.e., principal axes) rather than the rotated one. However, it may be calculated as follows: As shown in Fig. 5, the LSC adopted in the present Smith method is evaluated by the adaptable algorithm proposed and systematically validated in [67], in which the LSC is described by three distinct responses: linear elastic response, arc-shape ultimate response and the asymptotic post-collapse response. This method was developed by idealising the observations in parametric nonlinear finite element analyses and offers an efficient and robust approach to be incorporated in the Smith-type method. The formulation of this method is given by Eqs. (20)- (22) in which the ultimate compressive strength of stiffened panels ( ) is estimated by the empirical formula introduced in [85].
Arc-shape ultimate response Asymptotic post-collapse response where is the normalised initial stiffness which can be taken as unity if the initial imperfection is not excessively large; is the normalised instantaneous stiffness which incrementally reduces from the initial stiffness to null at the ultimate limit state. Meanwhile, the following expressions define the arc radius of the nonlinear response close to the ultimate limit state (i.e., ), the linear strain limit which separates the linear elastic response and arc-shape ultimate response (i.e., ∕ ), the ultimate strain at collapse state (i.e., ∕ ) and the parameter related to the asymptotic convergence of the post-collapse strength (i.e., ): where is the plate slenderness ratio of stiffened panel element, and is the column slenderness ratio of stiffened panel element. A Smith method example is given in Fig. 6, which shows the difference between the analysis results with and without considering neutral axis rotation for the Suezmax tanker with damage scenario No. 35 (explained in Section 3). The two bending moment/curvature relationships are nearly identical, and the deviation is also minimal in terms of the progressive translation of the vertical neutral axis. There is a sudden change in the rotational angle of the neutral axis when the cross-section reaches the ultimate limit state. This is likely due to the local failure of compressed structural elements. But soon after that, the rotational angle recovers to relatively small, which is probably caused by the gross yielding failure of the elements under tension. In most damaged scenarios analysed in this paper, the effects of the neutral axis rotation seem to be relatively small. This is because only grounding damage is considered, i.e. only bottom damage, in which case a highly un-symmetrical cross-section does not present. Thus, the effect of the neutral axis is generally minor. It is acknowledged that different techniques were also proposed in several studies for the progressive collapse analysis of hull girder under unsymmetrical bending by a Smith-type approach considering the neutral axis rotation [86][87][88].

Overall framework and procedure
Paik et al. [52] proposed an innovative structural condition assessment technique by developing the residual strength versus damage index (R-D) diagram technique, and the general procedure is presented in Fig. 7. For the application to the grounding damage ship structures in this study, the definition of ship structure characteristics in principle refers to the modelling for strength analysis in accord with the chosen method, including model discretisation and local behaviour definition etc. The characterisation of the damage parameters aims to select the most critical influence parameters for the rational definition of a grounding scenario. Once this is completed, the selected parameters are sampled from their probability distributions with an efficient sampling scheme (e.g. Latin Hypercube Sampling, LHS). A damaged index is defined, representing the severity of a damage scenario. The residual strength versus damage index diagram can be established by computing the target vessel's residual ultimate strength in all sampled damaged scenarios.

Definition of grounding scenarios
As illustrated by Fig. 8, a ship grounding scenario may be featured by (1) location and extent of the seabed obstacle in the longitudinal direction; (2) location and extent of the seabed obstacle in the transverse direction; (3) location and extent of the seabed obstacle in the vertical direction; (4) shape of the seabed obstacle [52]. The statistics of ship grounding may refer to the studies reported in [90,91]. Concerning Fig. 9, the following four parameters are considered in the DI-based approach for grounded ship hull girders: • Transverse location of grounding damage affected zone ( 1 ) • Height of grounding damage affected zone ( 2 ) • Breadth of grounding damage affected zone ( 3 ) • Angle of the affected zone ( 4 ) Fig. 7. Flowchart of the damage index (DI) concept-based condition assessment of damaged structures [52,89].
In some scenarios, it is not possible to configure the height of the damage affected zone because its angle is too large. In these cases, the following relationship will be used instead: A grounding damage index (GDI) is calculated to indicate the damage severity for each grounding damage scenario. As given by Eq. (19), the GDI is a function of the ratio between the damaged area of the inner bottom ( ,damaged ) and the total area of the original intact inner bottom ( ,original ), the ratio between the damaged area of the outer bottom ( ,damaged ) and the total area of the original intact outer bottom ( ,original ), and a correction factor ( ) to account for the contribution difference between the inner and outer bottom to the strength of the hull girders.

Sampling of grounding scenarios
The probability distribution of four grounding damage parameters is shown in Figs. 10(a) to 10(d), respectively. The distributions of 1 , 2 Fig. 8. Schematics of ship grounding scenario. Fig. 9. Schematic illustration of sharp and blunt rocks [52]. and 3 are consistent with the requirement by International Maritime Organisation (IMO) [92]. The probability distribution of 4 was developed in [52]. As the sampling statistics covers all grounding damaged scenarios recorded by IMO, the grounding damaged scenarios would include the sharp rock-type scenario and the blunt shoal-type scenario with a large contact surface. However, it is assumed that the damaged structures will lose their effectiveness completely, regardless of the grounding damage failure mode, i.e., tearing or denting. Nevertheless, this aligns with most studies on the post-accident residual strength assessment of grounded ship hull girders.
To sample the grounding scenarios, a Latin-Hypercube Sampling (LHS) technique is utilised [93]. In this sampling method, each dimension is divided into n intervals with equal probability. The LHS sampling is carried out by first randomly selecting an interval that has not been selected in the previous sample. Once an interval is selected, a random value within the corresponding interval will be generated. The main benefit of using LHS is the efficiency to simulate the desired probability distribution with a limited number of samples. Fifty grounding damage scenarios are sampled in this work, which is summarised in Appendix. A similar random sampling technique is also applied to the stochastic imperfection [94], plate selection [95], ship berthing [96], corroded structures [97,98], current profile [99] and offshore riser [100].

Case study models
In this work, the GDI-based residual strength assessment is performed for four different classes of double hull oil tankers, including VLCC class, Suezmax class, Aframax class and Panamax class. Their midship cross-sections are shown in Fig. 11, in which L = ship length, B = ship breadth, D = ship depth, b = double side shell width and h = double bottom height.
Consistent with the previous studies, it is assumed that the grounding takes place at the midship section. Therefore, the R-D diagram   will be developed in association with the residual ultimate strength performance of the midship cross-section.

Validation of numerical methods
Paik et al. [52] performed a comparison study on the ultimate longitudinal strength behaviour of the VLCC class and Suezmax class oil tankers suffered by a minor (Scenario No: 24) and major (Scenario No: 35) grounding damages (Fig. 12). Various simulations were performed by the Modified P-M method, ALPS/HULL (= ISFEM), ANSYS (= NLFEM), and CSR method, which are presented in Figs. 13 to 16. In addition, the predictions by the present Smith method are also plotted. Note that the CSR method is a variant of the Smith-type progressive collapse method. The distinct difference between the present Smith method and CSR is the formulation for load-shortening curves of structural components.
The calculation based on the present Smith method would implicitly include an average-level initial geometric imperfection and residual stress. In the ANSYS simulation, a single-bay finite element model is adopted to model the ship hull cross section subjected to longitudinal bending. The bending load is applied through reference points coupled with the cross-section boundary. The rotation-controlled technique is adopted. An average-level initial geometric imperfection is considered, including local plate distortion, column-type distortion and stiffener sideway distortion. Conversely, the effect of residual stress is ignored. The finite element analysis is executed by the arc-length solver.
It can be seen that the calculations by ISFEM, modified P-M method, and the present Smith method are generally in reasonable agreement with the reference prediction by ANSYS nonlinear finite element analysis and CSR method. As suggested by the benchmark studies regarding the ultimate longitudinal strength of ship hull girder [102][103][104], the present computational discrepancy should be acceptable, and hence the rationality of the three numerical methods can be verified.

Correction factors
The correction factor ( ) in Eq. (19) aims to account for the contribution difference between the inner bottom and otter bottom to the ultimate longitudinal strength of ship hull girders. To calibrate this correction factor, a series of computations are completed on the four oil tanker models with systematically varying damages in the inner bottom or outer bottom only. Linear regression curves are developed for the variation in ultimate longitudinal strength with respect to the damages in the inner bottom or outer bottom. As illustrated by Fig. 17, the correction factor is defined as the ratio between the slopes of the fitted linear curves for the results corresponding to the inner bottom ( ) and outer bottom ( ), respectively (Eq. (29)). The variation in ultimate longitudinal strength of the four oil tankers with respect to the damages in the inner bottom or outer bottom are plotted in Fig. 18.
As presented in Fig. 17(a) to (d), different correction factors are developed for different numerical methods. The computational discrepancy between the three numerical methods will be discussed in the later section. Nevertheless, it should be noted that there will be some differences in the GDI values between the three numerical methods due to the deviation in correction factors, even though the same damage parameters are defined.

Residual strength versus grounding damage index (R-D) relationships
The residual ultimate strength of four oil tanker models is calculated for the generated 50 grounding damaged scenarios. The developed R-D diagrams (= empirical formula) are summarised in Figs. 19(a) to 19(d). Note that the R-D formula is derived by a quadratic curve-fitting on the obtained results. Generally, it is found in all four vessels that the R-D diagram based upon the present Smith method is the most optimistic with regards to the residual ultimate strength of ship hull girders in sagging. When it comes to the residual strength of hull girders in Fig. 11. Cross-sections of the four different double-hull oil tankers models [78,101].     hogging, the R-D diagrams are close for GDI < 0.6. A larger discrepancy is observed in more severe grounding scenarios.
The developed R-D diagram provides efficient first-cut estimation for the residual ultimate strength of double hull oil tankers suffered from grounding damage, in which the approximate location and damage extent are known. Additionally, the R-D diagram can also be used to determine the acceptance criteria in relation to the post-grounding residual strength performance of ship hull girders. The IMO [105] proposed that the residual strength of all vessels should not be less than 90% of its intact ultimate strength (Fig. 20). In accord with this requirement, the upper limit of GDI of each vessel can be computed, as summarised in Table 1. These data are useful for developing acceptance criteria to assess the residual ultimate strength of grounded double-hull oil tankers based on different numerical methods.
When it comes to the assessment of double-hull oil tankers with different sizes, the relationship between the correction factor and the principal dimensions of ships (e.g. ship length, breadth, depth etc.) could provide a measure for size scaling, so that the damage index GDI of the target vessel can be calculated using Eq. (18). The relationship between ship's breadth and the correction factor is shown in Fig. 21, the best-fit one among different parameters. Note that the data shown in Fig. 18 corresponds to the analysis by the Fig. 17. Schematic illustration of the calibration of the correction factor . Smith-type progressive collapse method. The scaling relationships for ISFEM and the modified P-M method are given in Appendix (Fig. A.1 and Fig. A.2). The calculated GDI will be used in combination with a suitable R-D diagram/formula depending on the type of the vessels and the chosen strength analysis method. If the target vessel is out of the scope of the four different classes analysed in this work, a generic R-D diagram/formula is provided in Fig. 22, which are developed based on the calculation results of all four oil tankers. The generalised R-D diagram may support structural designers/engineers in predicting the residual strength capacity of the hull girders damaged by grounding. Application examples of the damage index-based assessment of grounded ship hull girders are given in Appendix.

Discussions
As shown in Section 4, there is considerable uncertainty in the residual strength prediction between ISFEM, Modified P-M method and Smith-type progressive collapse method. The computational discrepancy is mainly driven by (1) differences in their fundamental theories; (2) sub-division technique; (3) hard corner element.
As reviewed in Section 2, there is a distinct difference between these methods regarding their fundamental concepts. The ISFEM is developed similarly to the conventional finite element method. The modified P-M method is a presumed stress distribution-based approach, while the Smith method is a generalisation of the elementary beam theory. In addition, the idealisation of the local behaviour of structural components also substantially differs. Although using a similar concept of pre-defining the load-shortening curves of structural elements, the ISFEM and the present Smith method employ different formulations to generate these average stress/average strain relationships for characterising the in-plane response of the local component. This difference may be highlighted in the residual strength calculation of ship hull girders with un-symmetrical damages, as the problem is complicated by the neutral axis rotation. As for the modified P-M method, only the ultimate limit state of the local structural component is utilised to evaluate the overall hull girder bending strength. In other words, the post-collapse behaviour of local structures is ignored. In the meantime, the effect of neutral axis rotation is not accounted for in the modified P-M approach.
Sub-division of the hull girder cross-section into structural elements is an essential step in all three numerical methods. Nevertheless, different sub-division techniques are employed in these methods. The ISFEM adopts a plate-stiffener separation (PSS) technique where the local plating and stiffener are considered different elements. The present Smith method utilises a plate-stiffener combination (PSC) approach where the stiffener and its attached plate is viewed as one stiffened panel element. The PSS technique is adopted for the sagging strength in the modified P-M method, whereas the PSC technique is adopted for the hogging strength. The difference in the subdivision technique has several implications regarding the calculation of residual ultimate strength of ship hull girders. Given the same grounding damage scenario, more damaged elements may be identified when using the PSC technique. This may occur when the longitudinal stiffeners in the inner bottom are damaged. In such a case, the whole stiffened panel, including stiffener and attached plate, needs to be removed according to the PSC technique.
However, only the stiffeners should be deleted when using the PSS technique. The conservatism of the PSC technique may partially explain why the R-D diagram for hogging developed by ISFEM is slightly optimistic than the other two numerical methods. Aside from this, the conservatism of the PSC technique regarding damaged element identification also affects the vertical position of the neutral axis of a damaged cross-section. Since more elements should be removed from the bottom, the vertical position of the neutral axis will further translate upward comparing with the case where PSS could be applied. This implies that the applied in-plane compression on the deck panel resulting from global bending will be smaller, possibly leading to a larger ultimate sagging strength of the hull girder. Hence, this may partially explain the optimistic sagging strength reduction predicted by the present Smith method.
Hard corner element, with elastic-perfectly plastic behaviour, is considered in most Smith-type progressive collapse methods. However, it is usually not included in the ISFEM and modified P-M approach. It is apparent that neglecting the buckling effect in those elements would increase the ultimate strength of the hull girder, which therefore adds a discrepancy to the Smith method as compared with the other methods. It is usually argued that the intersections between the deck and side shell, or bottom and longitudinal bulkhead, in which the hard corner elements appear, would gain a strengthening effect because of the relatively strong boundary constrain. Therefore, the plating in these locations will not suffer from a significant buckling failure. However, this is still an empirical statement, and no research has been conducted to verify this conclusion.
It is assumed in all three numerical methods that the affected zone by grounding damage will lose completely its structural effectiveness (= load-carrying capacity) regardless of the grounding damage modes, e.g., tearing and denting. This also aligns with most studies on the postaccident residual strength assessment of grounded ship hull girders. However, this may be somewhat conservative, especially in the case of minor dented structures. A complete loss of structural effectiveness can be unrealistic.
Additionally, the damage affected zone is defined only as those in direct contact with the seabed obstacles, whereas the adjacent structural elements are assumed to be intact. In fact, minor denting or residual stress may also be induced in those adjacent areas. There have been many advanced numerical simulations of ship hull structures that suffered from grounding or collision. Nevertheless, the assumption of removing damaged elements for the simplified residual strength calculation is still in place. Future research is needed to develop a   more reasonable mapping between the damage simulation and simplified residual strength calculation models. Concerning the DI-based concept, the definition of grounding damage index may be updated to accommodate the difference in damage modes.

Conclusions
A comparison of the numerical methods for the DI-based residual strength assessment of grounded ship hull girders is presented in this paper. The Intelligent Supersize Finite Element Method, Modified Paik-Mansour method, and Smith-type progressive collapse are compared regarding the developed residual strength versus damage index (R-D) diagram. Discussion upon the computational discrepancy between the three numerical methods is given. From this study, the following conclusions are drawn: • With damage location and extent embedded in the formulation, the GDI concept is a practical approach to evaluate the residual ultimate strength of ship hull girders suffered from accidental grounding. • Different numerical methods for the ultimate strength prediction of ship hull girders can be integrated with the GDI-based assessment. The present work supplements the R-D diagram based on the Smith method. Thus, the GDI-based assessment has now been     • There is an appreciable computational discrepancy between three numerical methods, likely caused by the differences in fundamental theories, sub-division techniques of structural components, and the effects of hard corner elements.   • Plate-stiffener combination sub-division may lead to a conservative prediction of the hogging strength reduction. Conversely, the sagging strength reduction may be overestimated as compared with a plate-stiffener separation technique.
Based upon the present study, the following future researches are recommended: (1) application of the other types of existing empirical formulas [106][107][108][109][110][111][112] in predicting the ultimate limit state (ULS) of the plate-stiffener combination (PSC) model; (2) developments of plate elements and stiffener elements for Smith-type progressive collapse method to allow for the use of plate-stiffener separation (PSS) subdivision technique; (3) investigation of hard corner element behaviour to better evaluate their buckling capacity; (4) development of R-D diagram for bi-axial bending; (5) development of time-dependent R-D diagram based on Smith method; (6) application of the R-D diagram for a probabilistic reliability assessment; (7) an improved mapping between the grounding damage simulation models and the simplified residual strength calculation models. Lastly, this study assumed that the vessel's mid-ship section would suffer from raking type grounding damage. It means that the grounding damage in the longitudinal direction has not been considered, and thus future research should be conducted to resolve this issue.

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.
To demonstrate the application of the damage index-based for grounded ship hull girders, two examples are provided in the following. The first example corresponds to damage scenario No. 24 and the second example corresponds to damage scenario No. 35, as given in  24). In this example, the VLCC-class double hull oil tanker is analysed, and its principal dimension is given as follows.  By using the R-D diagram developed by the Smith method shown in Fig. 17(a), the residual strength of the damaged VLCC in hogging and  .9632 > 0.9 (Satisfies the IMO requirement) From the above calculation, it is clear that the residual strength of the VLCC suffered from damage scenario No. 24 satisfies the IMO requirement [105]. ,damaged = 0.43 m 2 As given in the previous application example, the correction factors for hogging and sagging are: In accord with Eq. (19), the GDI in hogging and sagging can be calculated. By using the R-D diagram developed by the Smith method shown in Fig. 17(a), the residual strength of the damaged VLCC in hogging and sagging can be evaluated as follows: From the above calculation, it is found that the residual sagging strength of the VLCC suffered by damage scenario No. 35 satisfies the IMO requirement [105]. However, the residual hogging strength is not acceptable. Similarly, the R-D diagram developed by another method, i.e., IS-FEM or Modified P-M presented, can also be applied to assess the safety of the hull girder damaged by grounding.