A Mapping Approach for Efﬁcient CFD Simulation of Low-Speed Large-Bore Marine Engine with Pre-Chamber and Dual-Fuel Operation

: AbstractA natural-gas-diesel dual-fuel marine engine with a pre-chamber is a promising solution for ocean transportation to meet the International Maritime Organization (IMO) emission regulations. This engine system employs a pre-chamber with direct injection of diesel to ignite premixed natural gas due to its higher ignition energy, which can enable lower lean limit and higher thermal efﬁciency. The dual-fuel pre-chamber marine engine presents complex multi-regime combustion characteristics in- and outside the pre-chamber, thus posing challenges in its numerical simulation in a cost-effective manner. Therefore, this paper presents a three-dimensional modeling study for the multi-regime combustion in a large-bore two-stroke marine dual-fuel engine, proposing a novel mapping approach, which couples the well-stirred reactor (WSR) model with the G-equation model to achieve high computational accuracy and efﬁciency simultaneously. In-depth analysis is performed using representative exothermic reaction (RXR) analysis and premixed turbulent combustion fundamentals to better understand the combustion process and to provide guidance in the selection of mapping timing. The results show that the use of mapping to switch from the WSR to the G-equation model can effectively reduce the runtime signiﬁcantly by 71.5%, meanwhile maintaining similar accuracies in predictions of in-cylinder pressure traces, HRR and NOx emissions, compared to using WSR all along. Additionally, the choice of mapping timing based on several parameters is preliminarily discussed.


Introduction
In recent decades, the emission regulations of the International Maritime Organization (IMO) are becoming drastically stricter [1,2], posing challenges for marine engines to work in a more environmentally friendly manner and more efficiently. For example, the third-stage standard (IMO Tier III) has been implemented since 2016 of which a limit on oxynitride (NO x ) emission in the emission control areas (ECAs) corresponds to an 80% reduction in IMO Tier I standard [3], and the global fuel sulfur limit of 0.5% also came into force in 2020. Meanwhile, in order to contribute to the climate change mitigation, IMO announced its initial strategy on greenhouse gas (GHG) reduction targeting at >50% reduction in GHG emissions in 2050 compared to 2008 [4]. Therefore, the maritime sector is in urgent need of clean fuels with a low carbon footprint alternative to the heavy oil fuel (HOF) currently used in marine engines. Liquefied natural gas (LNG) has been a promising solution, since it provides great potentials for the reduction in NO x emissions and has near-zero sulfur content [5][6][7]. In addition, the H/C mass ratio of NG is much higher compared to conventional fossil fuels, leading to an approximately 20-25% reduction in CO 2 emissions during combustion compared to HOF [8].

Engine Specifications and Operating Conditions
The test engine is a two-stroke pre-chamber marine engine manufactured by China Shipbuilding Power Engineering Institute, which is designed for flexible operation on both HOF single-fuel and diesel-LNG dual-fuel engines. The gas admission valves (GAVs) are located at the two sides of cylinder and supply the gaseous fuel into the cylinder while the piston moving upward when operating in the dual-fuel mode. Pilot diesel fuel is injected into the pre-chamber through a three-hole, 0.42 mm-orifice injector under 500 bar injection pressure. The detailed configuration of the test engine is shown in Figure 1a, and the close-up of the pre-chamber is shown in Figure 1b. The key parameters of the test engine are listed in Table 1, and the baseline operating conditions are listed in Table 2. In this work, the top dead center (TDC) was defined at the 0 • CA, and the engine simulation runs from the exhaust valve opening (EVO, seen in Table 2) and for a 360-• CA-period complete cycle, including the scavenging process to provide accurate in-cylinder flow and thermodynamic conditions. Data on the test engine were used to evaluate and validate numerical models.

Engine Specifications and Operating Conditions
The test engine is a two-stroke pre-chamber marine engine manufactured by China Shipbuilding Power Engineering Institute, which is designed for flexible operation on both HOF single-fuel and diesel-LNG dual-fuel engines. The gas admission valves (GAVs) are located at the two sides of cylinder and supply the gaseous fuel into the cylinder while the piston moving upward when operating in the dual-fuel mode. Pilot diesel fuel is injected into the pre-chamber through a three-hole, 0.42 mm-orifice injector under 500 bar injection pressure. The detailed configuration of the test engine is shown in Figure 1a, and the close-up of the pre-chamber is shown in Figure 1b. The key parameters of the test engine are listed in Table 1, and the baseline operating conditions are listed in Table 2. In this work, the top dead center (TDC) was defined at the 0 °CA, and the engine simulation runs from the exhaust valve opening (EVO, seen in Table 2) and for a 360-°CA-period complete cycle, including the scavenging process to provide accurate in-cylinder flow and thermodynamic conditions. Data on the test engine were used to evaluate and validate numerical models.

Mapping Method and Numerical Models
In this work, a three-dimensional (3D) computational fluid dynamic (CFD) model of the test engine is set up on CONVERGE V2.3 platform. Some of the sub-model setup can be found in Table 3. For each simulation case, a first simulation run is performed using the WSR model with a reduced mechanism to simulate the diesel spray auto-ignition and part of jet flame development, and it runs from EVO until the map timing (T m ). Then, a second simulation run is started at T m with the initial fields of temperature, pressure, species and flow mapped from the first run and employs the G-equation model for flame propagation, as illustrated in Figure 2. The choice of T m is usually around the early stage of jet flame emergence into the main chamber when the combustion regime transitions from the turbulent flame jet to flame propagation, as will be discussed in detail in Sections 3.3 and 4.2. For comparison, simulations using WSR alone and G-equation alone are also performed.

Mapping Method and Numerical Models
In this work, a three-dimensional (3D) computational fluid dynamic (CFD) model of the test engine is set up on CONVERGE V2.3 platform. Some of the sub-model setup can be found in Table 3. For each simulation case, a first simulation run is performed using the WSR model with a reduced mechanism to simulate the diesel spray auto-ignition and part of jet flame development, and it runs from EVO until the map timing (Tm). Then, a second simulation run is started at Tm with the initial fields of temperature, pressure, species and flow mapped from the first run and employs the G-equation model for flame propagation, as illustrated in Figure 2. The choice of Tm is usually around the early stage of jet flame emergence into the main chamber when the combustion regime transitions from the turbulent flame jet to flame propagation, as will be discussed in detail in Sections 3.3 and 4.2. For comparison, simulations using WSR alone and G-equation alone are also performed. The physical and chemical properties of diesel are approximated as those of n-dodecane and n-heptane, respectively, and the natural gas is approximated as methane. For WSR model, a reduced chemical kinetic mechanism of methane/n-heptane with 79 species and 262 reactions [24] was employed to simulate the natural-gas/diesel dual-fuel combustion. This chemical mechanism was validated for both methane and n-heptane combustion experiments and showed quite good agreement with the experiments in terms of laminar flame speed and ignition delay. In addition, the multi-zone model [25] is employed with WSR in order to accelerate the solution of detailed chemical kinetics, which groups the computational cells with similar composition and thermodynamic states in the chemistry calculation. The physical and chemical properties of diesel are approximated as those of ndodecane and n-heptane, respectively, and the natural gas is approximated as methane. For WSR model, a reduced chemical kinetic mechanism of methane/n-heptane with 79 species and 262 reactions [24] was employed to simulate the natural-gas/diesel dual-fuel combustion. This chemical mechanism was validated for both methane and n-heptane combustion experiments and showed quite good agreement with the experiments in terms of laminar flame speed and ignition delay. In addition, the multi-zone model [25] is employed with WSR in order to accelerate the solution of detailed chemical kinetics, which groups the computational cells with similar composition and thermodynamic states in the chemistry calculation. Table 3. CFD sub-model setup.

Description Model
Turbulence RNG k-ε [26] Spray break-up KH-RT [27] Ignition/combustion Well-stirred reactor (WSR) model [24,28] G-equation with Gulder [29] correlation NO x emission Extended Zeldovic's [30] In the second run, the G-equation model is used to describe flame propagation, and the WSR model is only applied in the burnt zone for emission prediction. In the G-equation model, the flame front can be tracked by solving the transport equation of distance function G, as follows [31]: ∂ρG ∂t where S t is the turbulent flamespeed, ρ u is the unburned density, and k is the turbulent kinetic energy. The turbulent diffusivity terms are given by: where S c is the turbulent Schmidt number, and c s is a user-supplied constant. The turbulent flame speed S t is calculated in the RANS framework as follows [32]: where S l is the laminar flame speed, and the second term on the right hand side is the turbulence correction including three modeling constants a 4 , b 1 , b 3 and the non-dimensional Damkohler number (Da). Correspondingly, the turbulence correction term is positively correlated with the turbulence intensity (u ) and Da. The expression of S l is given by Gulder [29] as: where T u_ref is the reference unburned temperature, T u is the unburned temperature, P ref is the reference pressure, and γ _dil is the mass fraction of dilution species. γ and β are correction exponents for temperature and pressure, respectively. The typical values of the constants mentioned are summarized in Table 4 for simulating the CH 4 flame. In addition, for jet flame simulation in this work, the turbulence correction modeling constants b 1 and b 3 are adjusted in calibration with the experiments.

Meshing Details
For simulation of a two-stroke low-speed marine engine, the computational mesh needs to be carefully setup due to its huge size compared to an automobile engine. CON-VERGE provides a modified cut-cell Cartesian grid generation method to automatically generate the computational grid at runtime and allows for regional grid embedding and adaptive mesh refinement (AMR) according to the local condition [32]. In this work, in addition to the base grid, fixed embedding of level 2 is applied in the pre-chamber and AMR of level 2 with respect to temperature, and velocity is applied in the main chamber. A series of cases with a different choice of base grid size of 20 mm and 10 mm are set up with both two models as shown in Table 5 to investigate the effects of meshing on simulation results using different models.

Tm and Turbulent Flame Speed Model Constant
As mentioned in Section 3.1, the T m is a key parameter that determines when Gequation starts to work in the computational model. In order to shed insight into the choice of T m , three cases with different T m are set up in this section and compared to the baseline case using WSR alone (case 1) and the case using G-equation alone (case 3), as seen in Table 6. 351. 8 40 As introduced in Equation (3), the turbulent flame speed prediction is dependent on the turbulence correction term, which is affected by model constants a 4 , b 1 and b 3 . These constants are often used as tuning knob to accommodate the variation in turbulence-chemistry interaction and to optimize prediction accuracy. Hence, typical values as recommended by Peters [31], as well as 2 groups of constants using higher values for b 1 and b 3 , are adopted in this study, as detailed in Table 6. The use of large values for b 1 and b 3 is necessary to match the experiments, which is also observed in other studies [33]. Note that the mesh setup of case 1 is adopted for all the cases in Tables 6 and 7. Table 7. Modeling constants study cases.

ID Combustion Model
Another parameter, C, is defined to describe the combustion heat release progress of pilot fuel as follows: where in q is the integrated heat release, and Q is the calorific value of pilot fuel. Considering the heat release by both two fuels cannot be decoupled, C is inevitable larger than 1 when pilot fuel is completed burned.

Flame Structure Analysis
Da is a key parameter for elucidation of the flame structure and is defined as the ratio of mixing time scale (τ l ) and chemical time scale (τ F ): where s l is the laminar flame speed and can be estimated using the Gulder's model [29]. l F is the flame front thickness, and l t is the integral length scale. u is the turbulent fluctuating velocity from the RNG k-ε model [26]. When Da is less than 1, the combustion is mainly controlled by chemical kinetics, while the combustion process is characterized by fast chemistry when Da is larger than 1 [34]. In this work, flame front is defined as the 1400 K temperature isosurface [35,36], as shown by the red isosurface in Figure 3. In addition, the flame tip is defined as the cross of flame front and the axial centerline of channel orifice. Along the orifice direction, simulation data of Da, u , turbulent flame speed (St) in each cell and the corresponding distance from orifice exit are extracted for several time instants to illuminate the flame structure variation. Refer to Xu et al. [37] for the Da calculation method. When the G-equation model is applied, the turbulent flame speed in each cell can be output directly, while when the WSR model is applied, the turbulent flame speed is derived as the movement rate of 1400 K isosurface through post-processing. where in q is the integrated heat release, and Q is the calorific value of pilot fuel. Considering the heat release by both two fuels cannot be decoupled, C is inevitable larger than 1 when pilot fuel is completed burned.

Flame Structure Analysis
Da is a key parameter for elucidation of the flame structure and is defined as the ratio of mixing time scale ( ) and chemical time scale ( ): where is the laminar flame speed and can be estimated using the Gulder's model [29]. is the flame front thickness, and is the integral length scale. is the turbulent fluctuating velocity from the RNG k-ε model [26]. When Da is less than 1, the combustion is mainly controlled by chemical kinetics, while the combustion process is characterized by fast chemistry when Da is larger than 1 [34]. In this work, flame front is defined as the 1400 K temperature isosurface [35,36], as shown by the red isosurface in Figure 3. In addition, the flame tip is defined as the cross of flame front and the axial centerline of channel orifice. Along the orifice direction, simulation data of Da, , turbulent flame speed (St) in each cell and the corresponding distance from orifice exit are extracted for several time instants to illuminate the flame structure variation. Refer to Xu et al. [37] for the Da calculation method. When the G-equation model is applied, the turbulent flame speed in each cell can be output directly, while when the WSR model is applied, the turbulent flame speed is derived as the movement rate of 1400 K isosurface through post-processing. In addition, the RXR analysis is applied here to further investigate the flame structure at the mapping timing to illustrate the chemical process. The RXR in each cell was identified as the reaction that makes the highest contribution to the local heat release using a post-processing code developed by Liu et al. [38]. In addition, the RXR analysis is applied here to further investigate the flame structure at the mapping timing to illustrate the chemical process. The RXR in each cell was identified as the reaction that makes the highest contribution to the local heat release using a postprocessing code developed by Liu et al. [38].

Comparative Analysis of WSR, G-Equation and the Mapping Approach
The model accuracy and grid sensitivity are evaluated by comparing the simulation results to the experimental data. Figure 4 compares the experimental and predicted incylinder pressure and HRR of case setups in Table 5, employing two combustion models on meshing setup with minimum grid sizes 5 and 2.5 mm, respectively. It is obvious that the WSR model is more sensitive to the grid size, wherein simulation results using a finer computational grid match the experimental data better compared to the coarser one. In contrast, negligible differences between two mesh setups are noticed when using the G-equation model. Furthermore, even under a 2.5 mm-scale computational grid, the predicted results using the WSR model match the experiment only for the early combustion stage (before 353 • CA) but show a much slower heat release rate from 353 to 365 • CA, which leads to the delayed HRR at around 375 • CA. Employing a much finer grid could be beneficial to the WSR simulation results; however, it is prohibitively expensive in computational cost, since the marine engine as investigated in this study is significantly larger than a typical automobile engine, and the speed is also much lower.

Comparative Analysis of WSR, G-Equation and the Mapping Approach
The model accuracy and grid sensitivity are evaluated by comparing the simulation results to the experimental data. Figure 4 compares the experimental and predicted in-cylinder pressure and HRR of case setups in Table 5, employing two combustion models on meshing setup with minimum grid sizes 5 and 2.5 mm, respectively. It is obvious that the WSR model is more sensitive to the grid size, wherein simulation results using a finer computational grid match the experimental data better compared to the coarser one. In contrast, negligible differences between two mesh setups are noticed when using the G-equation model. Furthermore, even under a 2.5 mm-scale computational grid, the predicted results using the WSR model match the experiment only for the early combustion stage (before 353 °CA) but show a much slower heat release rate from 353 to 365 °CA, which leads to the delayed HRR at around 375 °CA. Employing a much finer grid could be beneficial to the WSR simulation results; however, it is prohibitively expensive in computational cost, since the marine engine as investigated in this study is significantly larger than a typical automobile engine, and the speed is also much lower.  Table 5.
On the other hand, cases using G-equation cost less run time (see Table 5), which is more favorable for large-bore engine R&D activity. In addition, the results of cases employing G-equation are closer to experimental data during the entire combustion process. However, the G-equation model also has its drawbacks for such a multi-stage combustion process. Figure 5b shows that the flame front in the pre-chamber develops spherical flame propagation using the G-equation model in case 3, and it is quite different compared to the results of the WSR model, as shown in Figure 5a. The diesel spray combustion process in premixed methane/air mixture is expected to occur in the pre-chamber, which is similar to the results in Figure 5a according to optical engine studies [38,39]. The flame development in case 3 is significantly slower than that in case 1, whereas the jet flame in case 1 has already entered the main chamber at 351.7 °CA but the flame in case 3 has not reached the exit yet at 352 °CA. Therefore, it is not appropriate to apply the G-equation model alone for the multi-regime combustion in the dual-fuel pre-chamber engine, as it is designed for turbulent premixed combustion rather than liquid fuel spray combustion. In addition, a slow combustion rate inside the pre-chamber leads to underpredicted NOx emission with the G-equation model (seen in Section 4.2 in details), since most NOx is generated inside the pre-chamber according to the previous work [40].  Table 5.
On the other hand, cases using G-equation cost less run time (see Table 5), which is more favorable for large-bore engine R&D activity. In addition, the results of cases employing G-equation are closer to experimental data during the entire combustion process. However, the G-equation model also has its drawbacks for such a multi-stage combustion process. Figure 5b shows that the flame front in the pre-chamber develops spherical flame propagation using the G-equation model in case 3, and it is quite different compared to the results of the WSR model, as shown in Figure 5a. The diesel spray combustion process in premixed methane/air mixture is expected to occur in the pre-chamber, which is similar to the results in Figure 5a according to optical engine studies [38,39]. The flame development in case 3 is significantly slower than that in case 1, whereas the jet flame in case 1 has already entered the main chamber at 351.7 • CA but the flame in case 3 has not reached the exit yet at 352 • CA. Therefore, it is not appropriate to apply the G-equation model alone for the multi-regime combustion in the dual-fuel pre-chamber engine, as it is designed for turbulent premixed combustion rather than liquid fuel spray combustion. In addition, a slow combustion rate inside the pre-chamber leads to underpredicted NOx emission with the G-equation model (seen in Section 4.2 in details), since most NOx is generated inside the pre-chamber according to the previous work [40].
Therefore, neither WSR nor G-equation used solely can accurately predict the complex combustion process in the engine using premixed natural gas fuel ignited by diesel spray inside the pre-chamber. More specifically, the WSR model is more appropriate to predict diesel spray combustion, while G-equation can simulate jet flame and the subsequent flame propagation process in the main combustion chamber. Hence, a mapping method combining both models for different combustion stages in such a dual fuel pre-chamber engine is set up. Figure 6 shows the in-cylinder pressure and HRR of experimental and simulation results with the single G-equation model, WSR model, and mapping method connecting both models. The mapping timing applied here is 353.8 • CA, and the choice of mapping timing will be discussed in the next section. In order to optimize simulation results, two groups of different model constants named G 1 and G 2 , respectively, are used for the G-  Table 7. From Figure 6, we can see that using the combined WSR and G-equation model through the mapping method (case 6) shows better prediction accuracy compared to the other two. Specifically, although it shows little difference in heat release rate in early combustion stage, the process inside the pre-chamber is better described by the mapping approach using the WSR model for this period compared to the G-equation model, as discussed above and illustrated in Figure 5. For the following combustion stage (from 353.8 to 360 • CA), the mapping approach shows faster heat release than employing WSR alone, because the mapping approach employs G-equation for this period to achieve less grid-dependency. Figure 7 summarizes the temporal S t at the flame tip (defined as shown in Figure 3) after the flame develops into the main chamber of the cases with different models. The starting points of lines correspond to the timing of flame emergence into the main chamber, and the earlier flame emergence using the mapping approach is consistent with the previous discussion. It is obvious that the flame development process inside the main chamber includes three stages according to the flame speed profile, which was also observed in direct numerical simulation results of pre-chamber jet flame [21]. This three-stage flame development exists for all three model setups as seen in Figure 7; although the specific timings and durations are different from the results of Figure 7, we can see that employing the G-equation model after map timing (353.8 • CA) can slow down the reduction rate of flame speed and achieve a better match against the experiment for the heat release prediction. Therefore, neither WSR nor G-equation used solely can accurately predict the complex combustion process in the engine using premixed natural gas fuel ignited by diesel spray inside the pre-chamber. More specifically, the WSR model is more appropriate to predict diesel spray combustion, while G-equation can simulate jet flame and the subsequent flame propagation process in the main combustion chamber. Hence, a mapping method combining both models for different combustion stages in such a dual fuel pre-chamber engine is set up. Figure 6 shows the in-cylinder pressure and HRR of experimental and simulation results with the single G-equation model, WSR model, and mapping method connecting both models. The mapping timing applied here is 353.8 °CA, and the choice of mapping timing will be discussed in the next section. In order to optimize simulation results, two groups of different model constants named G1 and G2, respectively, are used for the G-equation model, as listed in Table 7. From Figure 6, we can see that using the combined WSR and G-equation model through the mapping method (case 6) shows better prediction accuracy compared to the other two. Specifically, although it shows little difference in heat release rate in early combustion stage, the process inside the pre-chamber is better described by the mapping approach using the WSR model for this period compared to the G-equation model, as discussed above and illustrated in Figure 5. For the following combustion stage (from 353.8 to 360 °CA), the mapping approach shows faster heat release than employing WSR alone, because the mapping approach employs G-equation for this period to achieve less grid-dependency. Figure 7 summarizes the temporal St at the flame tip (defined as shown in Figure 3) after the flame develops into the main chamber of the cases with different models. The starting points of lines correspond to the timing of flame emergence into the main chamber, and the earlier flame emergence us-

Flame Structure Analysis and Choice of Tm
The mapping timing Tm is determined as 353.8 °CA for this studied condition mainly based on calibration rather than in a predictive manner, and it is necessary to further discuss the choice of Tm to improve the predictivity of this approach. In this section, the flame structure evolution based on parameters of , Da and the Borghi-Peters diagram is analyzed in detail.
Local temperature, and Da in each simulation cell along the orifice centerline (seen in Figure 3) using WSR alone or the G-equation model (G1) alone are provided in Figures 8-10, respectively. The temperature profiles show that flame development speeds predicted by the two models are remarkably different. Specifically, at 352 °CA, jet flame length is a bit longer using the WSR model than using G-equation; however, flame length shows the opposite result after 353 °CA, which indicates that the simulated flame develops faster with the WSR model at the early stage of the jet flame development process, while later on, the simulated flame develops more slowly using WSR model. Figure  9 shows local results at different times using the two models. It is obvious that the of the flame is profoundly larger at 352 °CA than later on. Results using two models show very similar trends, while only the flame lengths are different, as earlier mentioned. The

Flame Structure Analysis and Choice of Tm
The mapping timing T m is determined as 353.8 • CA for this studied condition mainly based on calibration rather than in a predictive manner, and it is necessary to further discuss the choice of T m to improve the predictivity of this approach. In this section, the flame structure evolution based on parameters of u , Da and the Borghi-Peters diagram is analyzed in detail.
Local temperature, u and Da in each simulation cell along the orifice centerline (seen in Figure 3) using WSR alone or the G-equation model (G 1 ) alone are provided in Figures 8-10, respectively. The temperature profiles show that flame development speeds predicted by the two models are remarkably different. Specifically, at 352 • CA, jet flame length is a bit longer using the WSR model than using G-equation; however, flame length shows the opposite result after 353 • CA, which indicates that the simulated flame develops faster with the WSR model at the early stage of the jet flame development process, while later on, the simulated flame develops more slowly using WSR model. Figure 9 shows local u results at different times using the two models. It is obvious that the u of the flame is profoundly larger at 352 • CA than later on. Results using two models show very similar trends, while only the flame lengths are different, as earlier mentioned. The change in u indicates that the flame structure changed significantly after 352 • CA, which is captured by the simulations with either the WSR model or G-equation. In Figure 10, Da results show that large Da appears at the location nearby the flame tip, especially after 354 • CA. Although the magnitudes are quite different using the two models, which is due to the differences in predicted flame speed and l t /l f , the overall trends of Da changing with time for each model are consistent. The scatter points in the Borghi-Peters diagram of two models at different times as shown in Figure 11 can further illustrate the flame structure change. The jet flame is initially close to the broken reaction zones and is always in thin reaction zones. Then, it moves towards the corrugated flamelet zones, which is more profound when using the G-equation model. Furthermore, after 353 • CA, the scatter points cross the line of Da = 1 towards larger values, meaning that fast chemistry dominates the flame, and the flame front can be tracked with a level-set approach, i.e., the G-equation model, as suggested by Peters [31]. flame structure change. The jet flame is initially close to the broken reaction zones and is always in thin reaction zones. Then, it moves towards the corrugated flamelet zones, which is more profound when using the G-equation model. Furthermore, after 353 °CA, the scatter points cross the line of Da = 1 towards larger values, meaning that fast chemistry dominates the flame, and the flame front can be tracked with a level-set approach, i.e., the G-equation model, as suggested by Peters [31].   flame structure change. The jet flame is initially close to the broken reaction zones and is always in thin reaction zones. Then, it moves towards the corrugated flamelet zones, which is more profound when using the G-equation model. Furthermore, after 353 °CA, the scatter points cross the line of Da = 1 towards larger values, meaning that fast chemistry dominates the flame, and the flame front can be tracked with a level-set approach, i.e., the G-equation model, as suggested by Peters [31].    In order to further understand the flame structure, the RXR distribution on the clip plane along the jet direction for case 3 at 353.8 °CA (Tm in case 5) is shown in Figure 12, and the corresponding formation and production rate of methyl are shown in Figure 13. As can be seen, similar RXR occurs inside the pre-chamber and jet flame core near the exit, which suggests the burning mixture in the pre-chamber is pushed into the main chamber along with the jet, and the reaction status is not changed much during this process. In addition, the reaction at the periphery of the jet is mainly participated by methyl, which is mainly from the premixed CH4 reaction according to the analysis in Figure 13. This justifies the use of G-equation to describe flame propagation in premixed CH4 in the mapping approach. In order to further understand the flame structure, the RXR distribution on the clip plane along the jet direction for case 3 at 353.8 • CA (T m in case 5) is shown in Figure 12, and the corresponding formation and production rate of methyl are shown in Figure 13. As can be seen, similar RXR occurs inside the pre-chamber and jet flame core near the exit, which suggests the burning mixture in the pre-chamber is pushed into the main chamber along with the jet, and the reaction status is not changed much during this process. In addition, the reaction at the periphery of the jet is mainly participated by methyl, which is mainly from the premixed CH 4 reaction according to the analysis in Figure 13. This justifies the use of G-equation to describe flame propagation in premixed CH 4 in the mapping approach.  Using a key parameter as the formally simple approach to choose the T m is ideal in engineering application. In this section, the mean turbulence intensity u , the overall NO x emission level and the combustion schedule parameter C are considered tentatively to determine T m . Figure 14 shows the simulated mean turbulence intensity u results, defined as statistics data of the region inside the flame surface in main chamber. u declines as the flame development proceeds. In addition, the results of two models show a considerable difference when u reduces below 4 m/s, and the corresponding crank angle is at 353.8 • CA (T m in case 5). Thus, choosing the crank angle when u = 4 m/s as T m may be reasonable. Table 8 shows the IMEP and NO x emission prediction results of cases 3, 5, 6 and 7. As can be seen, using the G-equation model alone (case 3) or when T m is early (cases 6 and 7), which means WSR does not run long enough for flame development prediction, leads to NO x emissions a bit lower than in experimental data, while T m = 353.8 • CA (case 5) shows the best simulated result. On the other hand, although the IMEP result of case 5 shows a higher discrepancy against experimental data, it is still considered within the acceptable range compared to the errors seen for NO x prediction. As mentioned, the NO x emission is mainly produced in the pre-chamber and is dominated by the first combustion stage duration. Hence, choosing T m based on the calibration of NO x emissions may also be a reasonable option. Besides, C at T m is also listed in Table 8, as discussed previously, and the flame structure is changed significantly at 354 • CA; hence, C = 13.3 (case 5, T m = 353.8 • CA) can also be considered as an indicator for combustion mode transition. In addition, when T m is too small (case 7), C is less than 1, meaning pilot fuel has not completely burned yet. When using the crank angle at C = 4.8 as T m (case 6), the result shows underpredicted NO x emissions compared to experimental data. Therefore, given the limited experimental data, the choice of T m = 353.8 • CA provides the most accurate and reasonable results, and three indicators of u , NO x emissions and C can be potentially used to determine the mapping timing T m for engineering application.  Using a key parameter as the formally simple approach to choose the Tm is ideal in engineering application. In this section, the mean turbulence intensity , the overall NOx emission level and the combustion schedule parameter C are considered tentatively to determine Tm. Figure 14 shows the simulated mean turbulence intensity results, defined as statistics data of the region inside the flame surface in main chamber.
declines as the flame development proceeds. In addition, the results of two models show a considerable difference when reduces below 4 m/s, and the corresponding crank angle is at 353.8 °CA (Tm in case 5). Thus, choosing the crank angle when = 4 m/s as Tm may be reasonable. Table 8 shows the IMEP and NOx emission prediction results of cases 3, 5, 6 and 7. As can be seen, using the G-equation model alone (case 3) or when Tm is early (cases 6 and 7), which means WSR does not run long enough for flame development prediction, leads to NOx emissions a bit lower than in experimental data, while Tm = 353.8 °CA (case 5) shows the best simulated result. On the other hand, although the IMEP result of case 5 shows a higher discrepancy against experimental data, it is still considered within the acceptable range compared to the errors seen for NOx prediction. As mentioned, the NOx emission is mainly produced in the pre-chamber and is dominated by the first combustion stage duration. Hence, choosing Tm based on the calibration of NOx emissions may also be a reasonable option. Besides, C at Tm is also listed in Table 8, as discussed previously, and the flame structure is changed significantly at 354 °CA; hence, C = 13.3 (case 5, Tm = 353.8 °CA) can also be considered as an indicator for combustion mode transition. In Figure 13. Formation rate of methyl for major reactions that account for more than 1% at 353.8 • CA in the clip plane, same as Figure 12. addition, when Tm is too small (case 7), C is less than 1, meaning pilot fuel has not completely burned yet. When using the crank angle at C = 4.8 as Tm (case 6), the result shows underpredicted NOx emissions compared to experimental data. Therefore, given the limited experimental data, the choice of Tm = 353.8 °CA provides the most accurate and reasonable results, and three indicators of , NOx emissions and C can be potentially used to determine the mapping timing Tm for engineering application.

Conclusions
An accurate and efficient numerical method to simulate combustion process for largebore dual-fuel marine engine with pre-chamber was proposed in this study. Firstly, a grid-dependency study was performed using the WSR model and G-equation with the RANS turbulence model to manifest that the WSR model is more sensitive to grid resolution compared to the G-equation model. Then, detailed RXR analysis and premixed turbulent combustion regime analysis were applied to better understand the characteristics of the pre-chamber flame jet, which also justifies the use of the mapping approach that couples the WSR and G-equation models to enable fast and accurate simulation. Furthermore, the choice of mapping timing (T m ) was discussed, and three potential approaches to determine T m were proposed to improve the predictivity of this mapping approach for engineering application. Major conclusions are as follows: 1.
The grid sensitivity study showed that a minimum mesh resolution of 2.5 mm cannot satisfactorily predict the combustion HRR and in-cylinder pressure using the RANS method coupled with the WSR model. However, it is quite unacceptable to use a further finer grid in engineering application. The G-equation model is not as grid sensitive as the WSR model, thus enabling faster simulation.

2.
Employing G-equation alone for the combustion prediction can significantly reduce simulation time cost. However, the prediction of flame development and NO x emission formation in the pre-chamber is unsatisfying using G-equation alone.

3.
The flame structure and RXR analysis shows that flame structure changes significantly at around 354 • CA. Besides, the flame propagation mainly occurs in the premixed CH4. Therefore, choosing T m around this crank angle shows the best calibrated result compared to the test data.

4.
Three potential approaches to determine T m for engineering application are proposed, viz., the crank angle when u = 4 m/s, when NO x emission is closest to experimental data or when the combustion progress parameter C = 13.3.
This work represents a first step towards a predictive and cost-effective simulation framework for a dual-fuel marine engine, which is very challenging, characterized by the complex combustion process and large computational domain. Further development and verification for the modelling approach based on more engine test data and optical diagnostic results will be conducted in a future work.