Large-eddy simulation of split injection strategies in RCCI conditions

In this study, we investigate the effect of different split injection strategies on ignition delay time (IDT) and heat release rate (HRR) characteristics in Reactivity Controlled Compression Ignition conditions via large-eddy simulation and finite-rate chemistry. A diesel surrogate (n-dodecane) is injected into a domain with premixed methane and oxidiser in two separate injection pulses. Three different split injection strategies are investigated by fixing the amount of total fuel mass: varying the first injection timing, varying the second injection timing, and changing the fuel mass ratio between the two injections at a fixed injection timing. A compression heating mass source term approach is utilised to take compression heating into account. The main findings of the study are as follows: (1) In general, the IDT shifts towards the top-dead centre when the first injection is advanced or the second injection is retarded. The size and spatial pattern of the ignition kernels are shown to depend on the dwell time between the injections. (2) A precisely timed first injection offered the best control over ignition and HRR characteristics. However, advancing the first injection may lead to over-dilution downstream, preventing volumetric ignition and reducing the peak HRR value. (3) Approximately 21% decrease in the maximum HRR value, as well as a factor of 2.8 increase in combustion duration could be achieved by advancing the first injection timing. (4) As indicated by frozen-flow chemistry analysis, in the investigated configurations, the reactivity stratification is controlled by mixture stratification rather than temperature. The findings indicate that the first injection controls the downstream reactivity stratification, offering ignition and HRR control.


Introduction
Spray assisted ignition has been investigated extensively in the past in relation to lowtemperature combustion (LTC) technologies [1][2][3]. LTC technologies aim at emission reduction through lower burning temperatures and leaner mixture conditions to mitigate NOx and soot emissions, while targeting performance comparable to Conventional Diesel Combustion (CDC) [3,4]. The most significant LTC strategies are Homogeneous Charge Compression Ignition (HCCI), Premixed Charge Compression Ignition (PCCI) and Reactivity Controlled Compression Ignition (RCCI) [5]. In the present paper, fundamental numerical studies are carried out to gain insights to robust ignition control in RCCI combustion [5]. RCCI combustion aims at mitigating emissions while providing high efficiency through LTC for compression ignition (CI) engines. The working principle of RCCI engines is to create an in-cylinder mixture and reactivity stratification via single or multiple direct injection(s) of a high-reactivity fuel (HRF) (e.g. diesel) into a homogeneous mixture of a low-reactivity fuel (LRF) (e.g. methane) and oxidiser. The difference in the auto-ignition characteristics of the HRF and LRF along with different operational parameters provide control over the ignition timing, as well as the combustion phasing and the duration. The scope of the present work is to provide fundamental insights on the mechanisms that control ignition timing and combustion duration in a dual-fuel multiple injection system through large-eddy simulations (LES).
According to the literature, the use of multiple injections and injection timing control is reported to be the most effective ways of creating the desired reactivity stratification in RCCI [8]. Numerical investigations of RCCI engine with natural gas and diesel by Nieman et al. [10] highlighted the importance of SOI timing, where complete combustion with high unburned hydrocarbon (UHC) and NOx emissions was reported for very early SOI. These observations were supported by Li et al. [14], where advancing SOI reduced UHC and soot, although it had an adverse effect on NOx and CO emissions. According to Paykani et al. [16], retarded ignition timing by adjusting SOI is preferred for decreasing the residence time of the mixture at high ambient temperatures to avoid NOx formation, while advanced SOI can avoid fuel-rich regions and reduce soot emissions. In an experimental study using iso-octane and n-heptane, Kokjohn et al. [12] reported the highest peak of HRR for either early or late SOI's and the minimum peak at mid-range injections. It was explained that for early SOI's, the mixing time is long and the charge is well-mixed at lean mixture conditions, while for late SOI's, the richer diesel charge has similar ignition delays with short mixing time, both leading to volumetric ignition. For mid-range SOI's, however, ignition initiates in the downstream of the diesel jet in local auto-ignition pockets and propagates to the rest of the charge at a slower rate than in other SOI conditions.
The common observation in the literature is that adopting an advanced SOI strategy creates a leaner mixture and lower burning temperature at the time of ignition, hence, reducing the high NOx and soot emissions associated with rich mixtures and high burning temperatures. However, advanced SOI may also reduce the reactivity gradient through long mixing times and lead to volumetric ignition with rapid energy release, which is a common problem in HCCI combustion that prevents its operation at higher loads [12]. Our previous single injection study has provided evidence that although ignition timing control and leaner mixture conditions can be achieved by advancing SOI, the initial ignition kernels of very early injections may not result in strong ignition kernels, indicating slow combustion initiation. [21]. We have also shown consistent trends between our numerical results and experimental data in terms of ignition characteristics. In addition, the local auto-ignition regions are usually found at the tip of the injected spray, which could leave a large amount of the injected HRF unburned between the local auto-ignition kernel and the injector.
To mitigate the above-mentioned issues, multiple injection strategies have been proposed to enhance the ignition/combustion control [3,5]. In split injection strategies, an early injection (SOI 1 ) is utilised to create a lean charge, favouring low soot and NOx emissions. A second injection (SOI 2 ) closer to the top-dead centre (TDC) is used to create a higher reactivity region and a subsequent reactivity gradient, enabling the control of HRR and ignition timing. In addition, for shorter dwell times, the jet-jet interaction between SOI 1 and SOI 2 can be used to further control the ignition timing and overall mixture. Hadadpour et al. [22] investigated various mass split and dwell time configurations to study jet-jet interactions and reported that with longer dwell time, the entrainment time is increased and leaner mixtures can be obtained. They also reported that SOI 2 can be optimised to either push the lean near-nozzle residue gas away from the nozzle to reduce the UHC, or enrich this region and initiate a local auto-ignition zone.
Two main parameters used to control split injection strategies are the dwell time (time between SOI 1 and SOI 2 ) and the mass split ratio (ratio of injected mass in each injection). In a recent numerical study, Li et al. [20] applied various split injection strategies with different dwell time and mass split ratios to RCCI engines. They reported longer ignition times with a split injection strategy compared to single injection, regardless of the dwell time and mass split ratio. In addition, it was noted that extending the dwell time can gradually change the combustion mode from reactivity-controlled to reactivity and mixing controlled, offering potential control over the combustion duration. Increasing the mass fraction of SOI 1 was also reported to be beneficial for reducing soot and NOx emissions [16].
Despite the proven benefits of multiple injections, the effects of different injection strategies on ignition and combustion characteristics are still not completely predictable, which is mainly due to the high dimensionality of the test space and the limited number of relevant numerical studies. Based on our literature survey, a fine-grained understanding of the effects of changing the dwell time and mass split ratio systematically on the ignition mechanism in multi-injection dual-fuel (DF) sprays is still lacking. In particular, high-fidelity numerical simulations along with detailed chemistry can provide a better understanding on the interplay between mixing and reactivity in DF spray-assisted ignition with different dwell time and mass split ratios. Following the numerical framework set in our earlier single injection RCCI study [21], here, we investigate various split injection strategies for diesel surrogate sprays and the effects of these strategies on ignition and combustion characteristics using LES. While keeping the total injected mass constant, different mass split ratios and dwell time between SOI 1 and SOI 2 are simulated. The main objectives of the study are to: (1) Carry out three-dimensional (3D) LES to analyse effects of different split injection strategies on ignition delay time (IDT) and HRR characteristics under RCCI conditions. (2) Investigate the effect of injection strategy on mixture/reactivity distribution and formation of ignition kernels.
(3) Analyse HRR characteristics based on the injection strategy in both physical and mixture fraction spaces. (4) Perform a frozen-flow ignition index analysis to characterise the relative importance of two potential reactivity stratification sources (thermal and mixture), and draw conclusions from reactivity stratification to overall HRR characteristics.

Methodology
In this section, information on the utilised numerical framework is presented. The Open-FOAM [23] framework is used and all in-house implementations are built on top of this framework. The proposed case setup is similar to the setup utilised in our previous studies [21,[24][25][26].

Governing equations and discretisation
The gas-phase fluid flow is solved using the compressible Navier-Stokes equations. The mass conservation, momentum, species transport, and enthalpy equations are solved together. The system of equations is closed by the ideal gas law and thermal equation of state. The LES formulation of these equations with Favre filtering is given as: whereρ, u i ,p, Y k ,h s , τ ij are filtered density, velocity, pressure, mass fraction of k th species, sensible enthalpy and viscous stress sensor, respectively. Overbar (¯) denotes the unweighted ensemble average and tilde (˜) denotes the density-weighted ensemble average. The variables c p and λ in the energy equation (Equation (4)) represents the heat capacity and thermal conductivity of the mixture. In reacting cases, the production rate of each specie is denoted byω k and heat release rate (HRR) is calculated asω h = k h o f ,kω k , where h o f ,k is the enthalpy of formation. A unity Lewis number is assumed with the diffusion coefficient D = λ/(ρc p ).
A 2 nd order time and space discretisation is utilised to discretise the governing equations using the finite volume method. The pressure-velocity coupling is achieved by using the PISO algorithm. The convection terms of the governing equations are discretised using the nonlinear and dissipative Gamma scheme [27], similar to our previous work [21,[24][25][26]28,29]. An implicit LES (ILES) approach is adopted with no explicit turbulence model, In ILES, the underlying assumption is that the dissipation introduced by the Gamma scheme is in the same order as the subgrid-scale dissipation. The Gamma model parameters were selected based on previous tests performed in Ref. [28]. Similar usage of ILES on free shear flows [30], supersonic jets [31], and spray injection problems [25,29,32] can be found in the literature. Further theoretical work comparing the implicit and explicit subgrid-scale models can be found in [33].

Spray and chemistry models
Lagrangian Particle Tracking (LPT) is used to inject n-dodecane liquid spray to a cylindrical tank. Frössling correlation [34] and Ranz and Marshall correlation [35,36] are used to model the evaporation between the two phases. Unlike our previous studies [24][25][26], an explicit droplet break-up model (e.g. KHRT [37,38]) is not utilised due to the dependency of the model constants on the ambient pressure and temperature, which are changing throughout the simulation. Alternatively, a constant droplet diameter (0.5 μ m) approach is utilised with the assumption that the secondary break-up process has already occurred. This modelling decision is based on the previous findings in the literature [39,40] and further investigated in our previous study in detail [21]. However, the salient details are reported here. To validate the applicability of this approach, four different experimental conditions from the ECN spray dataset [41] with a range of different pressures (50 and 150 Mpa), ambient densities (7.6, 15.2 and 22.8 kg/m 3 ) and ambient temperatures (900 and 1100 K) were investigated in our previous study. For all the investigated conditions we obtained good agreement with both liquid and vapour penetration experimental profiles. While our validation efforts are limited to four existing operating conditions in the ECN database, these conditions cover a wide range of injection pressure, ambient density and ambient temperature values with close relevance to the present work. It is also important the note that while validation efforts are carried out for the liquid injection process, the role of the Lagrangian particles is only to introduce the fuel into the numerical framework. The point of focus of our study is the gas phase chemical kinetics and the combustion process.
The chemical reactivity is introduced to the governing equations through filtered species reaction termsω k and HRRω h in Equations (3) and (4). An operator splitting method is used to split the chemical reaction rate calculation from gas-phase flow calculation, due to the large timescale difference between the two timescales. A skeletal chemical mechanism developed by Yao et al. [42] (54 species and 269 reactions) were used with a finite-rate chemistry model. This mechanism has been deemed suitable for n-dodecanemethane blends and utilised in our previous studies [21,[24][25][26], where it predicted the auto-ignition times of n-dodecane/air and methane/air mixtures successfully. Moreover, based on the experimental comparison analysis performed by Desantes et al. [43], Yao mechanism includes the essential low-temperature chemistry reaction pathways, although shows slight over-prediction on low-temperature heat release rate characteristics compared to other mechanisms they investigated.
The local thermochemical state of the fuel and the oxidiser are represented using a normalised mixture fraction (Z) based on conservation of the elemental mass fractions using Bilger's definition [44], where the value 0 corresponds to pure oxidiser (homogeneous methane/air mixture), and the value 1 corresponds to pure fuel. Turbulence-chemistry interaction (TCI) is utilised using a first-order closure approach with no subgrid-scale model for the chemical source terms in Equations (3) and (4). This decision is made due to the focus of this study being targeted to the spray auto-ignition, which is assumed to be less sensitive to micro-mixing issues [45], compared to e.g. quasi-steady flame lift-off length estimation in the Spray A context [46]. In addition, our comparison of the current nomodel approach with a Flamelet Generated Manifold (FGM) method using a presumed PDF approach indicated an ignition delay time (IDT) difference around 1% [25]. Moreover, in [25], the low-temperature combustion features were noted to be very similar for both ILES no-model and FGM methods, both agreeing with the experimental CH2O PLIF data.
The chemical Jacobian used by the finite-rate chemistry ODE solver is obtained by linking the open-source pyJac library [47] to our in-house OpenFOAM code. A linearly implicit extrapolation method (Seulex) [48] is used the solve the chemical ODE in each grid point. The computational load imbalance arising from the chemistry CPU load in different processors is solved by implementing a dynamic load balancing model using MPI routines. This model is called DLBFoam and the implementation details can be found in Refs. [49] and [50].

Compression heating model
In RCCI applications, the injection of the high-reactivity fuel (HRF) generally occurs well before TDC [21]. Therefore, it is important to take the increase in temperature and pressure during compression into account to accurately predict the chemical oxidation the HRF undergoes prior to ignition. While a moving mesh approach with topological changes is the most common way to include compression effects, it can also be modelled by a source term added to the governing flow equations. Following the mass source term approach introduced by Bhagatwala et al. [51,52] and further utilised by Luong et al. [53], a compression heating model is implemented to account for the increase in temperature and pressure throughout the simulation. Using the kinematic piston motion equations, the motored pressure inside the cylinder around TDC can be represented by: where P 0,m is the pressure at TDC, t c is the time it takes for one crank rotation, t 0 is the time at TDC and g and n are the engine geometry model constants.
Using Equation (5), a mass source termṁ is derived to be added to the right-hand side of all governing equations as:ṁ Figure 1 shows the validation of the source term compared to the motored pressure curve used in this study, obtained from Cantera using kinematic piston equations. It is also important to point out that the source term added to the governing equations only introduces the motored pressure for prescribed conditions, without the pressure rise originating from the combustion of the injected fuel. The pressure rise due to combustion is allowed to evolve on its own, following the approach adopted by Luong et al. [53]. This approach is valid since the ECN Spray A, the baseline of our numerical configuration, reports an overall increase in the constant volume chamber around 0.2 bar, corresponding roughly to a marginal ( < 1%) pressure rise [28]. Therefore, we chose the approach adopted by Ref. [53] and let the local pressure rise originating from combustion to evolve on its own.

Simulation setup
This study adopts a split HRF injection strategy to control the ignition timing, HRR distribution, and burning rate of the injected spray, while keeping the injected mass constant (∼ 1.1 mg). Different injection strategies include i) fixing the first injection timing (SOI 1 ) and varying the second injection timing (SOI 2 ), ii) fixing the SOI 2 and varying the SOI 1 , and iii) fixing both SOI 1 , SOI 2 and varying the injected mass ratio were investigated. A schematic representing the adopted strategy is presented in Figure 2. The LRF (methane) is assumed to be port injected earlier in the compression cycle and homogeneously mixed with the ambient oxidiser. HRF is direct-injected in two split injections. Using this approach, 9 different cases are investigated, the details of which can be found in Table 1.
Further details on the simulation setup, including engine compression and injection parameters as well as ambient conditions are provided in Table 2. The engine parameters and ambient conditions are identical to our previous single injection study [21]. A constant grid size of 62.5 μ m is used around the spray region and it is represented in Figure 3.

Results and discussion 3.1. Overview of different injection strategies
The main focus of this section is to identify and analyse the different ignition modes attained by various split injection strategies. In particular, ignition modes may span between volumetric and stratified modes, corresponding to a narrow distribution and high peak HRR (typical in HCCI) and well-distributed HRR with lower peaks (typical in RCCI), respectively. Commonly, RCCI applications aim at a wider HRR distribution with a lower peak value without compromising the thermal efficiency. We remind that RCCI was originally introduced to mitigate the issues related to operational range and high peak pressure/HRR in HCCI applications [5]. According to Karimkashi et al. [54] the reactivity stratification level affects the mode of combustion which will be further discussed in Section 3.3.
In this section, we investigate the three different split injection strategies described in Figure 2, and examine their effects on HRR characteristics. These three injection strategies, introduced in Table 1, are respectively sweeps of SOI 1 timing (Case 1), SOI 2 timing (Case 2), and mass split ratio (Case 3). The main purpose behind exploring these three split injection strategies is to find the most efficient strategy for controlling the HRR distribution under the present simulation conditions. At the end of this section, we select the strategy that provides the best control over HRR for further analysis in the remainder of the study. Figure 4 presents a qualitative comparison of the 3 different split injection strategies (Cases 1-3), each sweeping over 3 different values of the selected parameters is presented in Table 1, at their corresponding IDT + 0.2 ms time instances. The presented visualisations for the 9 LES cases are based on the superimposing of hydrogen peroxide H2O2 and temperature to indicate the locations of low-temperature and high-temperature combustion (LTC, HTC), respectively.
For Case 1 (SOI 1 = 21.6, 28.8 and 36 CAD while SOI 2 is fixed at 14.4 CAD), provided in Figure 4(a), we observe different ignition regimes based on the timing of SOI 1 . For SOI 1 = 21.6 CAD, the short dwell time between the two injections increases the local fuel reactivity, creating a merged volumetric ignition kernel contributed by both injections. However, with advanced SOI 1 (28.8 CAD), the jet-jet interaction is restricted due to the increased dwell time between the two injections and the longer mixing time of the first injection. In particular, the earlier SOI 1 extends the spray volume along the spray axis and reduces the overall reactivity of the first injection through creating a leaner mixture downstream, resulting in spatially separated ignition kernels originating from different injections. Finally, with the earliest SOI 1 (36.6 CAD), the reactivity difference between the two injections is much more visible. It can be seen that the ignition kernel associated with the first injection is not yet present due to the very high mixing time and subsequent over-leaning and reduced reactivity.
For Case 2 (SOI 2 = 21.6, 14.4 and 7.2 CAD while SOI 1 is fixed at 28.8 CAD), Figure 4(b) shows that with retarding SOI 2 (from 21.6 to 14.4 CAD), the increased dwell time creates explicit ignition zones separated in space. With further retarded SOI 2 (from 14.4 to 7.2 CAD), the ignition timing is controlled by the first injection, while the ignition due to the second injection is not present due to the short mixing time.
For Case 3 (SOI 1 and SOI 2 fixed at 36 and 14 CAD, with varying mass ratios), Figure 4(c) shows that for the 60/40% mass split, there is a temporally stratified ignition structure with the second injection ignition first, followed by the first injection. Increasing the mass fraction of the first injection increases the associated local reactivity downstream, while the reactivity closer to the nozzle is reduced due to a shorter second injection. This behaviour creates multiple ignition kernels originating at various locations due to temporally less stratified reactivity (mass split ratio 80/20%).
The qualitative observations presented above are supported by HRR profiles provided in Figure 5 (left column) along with the quantification provided in Table 3. To distinguish between the ignition modes, IDT, maximum HRR, and peak HRR duration are evaluated and compared. Here, IDT is defined as the time instance when 2% of the cumulative heat release (CHR) is reached. To quantify the peak HRR duration in CAD, we define the peak HRR window as the time duration between 80% of the maximum HRR before and after maximum HRR. IDT (CAD), maximum HRR (kJ/s), and peak HRR window (CAD) are compared for each case in Table 3. In addition, the CHR in each case is presented in Figure 5 (right column). Since the total injected diesel amount is constant for all investigated cases, CHR value provides an insight on the amount of LRF (methane) entrainment within the HRF charge due to mixing. We note that four different characteristics time instances associated with the HRR are also marked on HRR profiles to be used in our later analysis. Also, in all HRR profiles, the initial smaller peaks correspond to the first stage ignition. In the following, we discuss the observations in Figure 5 along with the quantified values in Table 3.
In Figure 5(a) for Case 1, HRR is investigated when SOI 2 and mass split ratio are fixed and SOI 1 is varied. According to Table 3, with advancing SOI 1 , IDT is advanced and the peak HRR value reduces ≈ 21% and the peak HRR window increases by a factor of 2.8; i.e. ignition approaches the stratified mode. This is due to the increased dwell time between the first and second injection, which reduces jet-jet interactions between the injected sprays and lowers the local fuel concentration, providing reduced overall reactivity. In contrast, retarding SOI 1 towards TDC results in a retarded IDT, higher HRR peak, and a narrower peak HRR window, i.e. approaching volumetric ignition, due to the short mixing time and increased jet-jet interactions. Moreover, the CHR distribution shows that the early SOI 1 (36 CAD) case starts with a lower CHR value, but in the end, it releases ≈ 25% more total heat than SOI 1 = 21.6 CAD case. For SOI 1 = 36 CAD, since volumetric ignition is avoided, a smaller CHR initial value is attained. However, due to the long mixing time of early SOI 1 , more LRF is entrained, which increases the CHR despite the same amount of injected HRF.  Table 1. The markers indicate peak HRR at first stage ignition, 10% of CHR, peak HRR and 80% of CHR for each case. In particular, the effect of higher dwell time on lowering peak HRR can be observed in Case 1.
The effect of SOI 2 timing on HRR characteristics, with fixed SOI 1 and mass split ratio (Case 2), can be observed in Figure 5(b). In contrast to the SOI 1 sweep reported in Figure 5(a), it can be seen that although retarding SOI 2 advances the ignition timing, the peak HRR value is first decreased ≈ 20% and then increased ≈ 10% (c.f. Table 3 for quantification). Similarly, the peak HRR window prolongs first by a factor of 1.9 and then becomes shorter by a factor of 2.3. Therefore, fixing SOI 1 and sweeping SOI 2 does not lead to a systematic change from one ignition mode to another, like in Case 1. While for advanced SOI 2 the dwell time between the injection pulses is small and the local reactivity is higher, retarding the SOI 2 increases the dwell and reduces reactivity, resulting in retarded Table 3. Ignition delay time (IDT), peak HRR value, and peak HRR duration for SOI 1 (Case 1), SOI 2 (Case2) and mass split fraction (Case 3) sweeps. Peak HRR duration denotes the HRR window where HRR is higher than 80% of the peak HRR value.  Figure 5(c), we consider the HRR for Case 3 in which SOI 1 and SOI 2 are fixed and the mass split ratio between the injections is varied. A two-stage HRR structure can be observed for the 60-40% split between the first and second injection, indicating a temporally stratified ignition structure with the second injection igniting first, followed by the first injection, as seen in Figure 4(c). Increasing the mass fraction of the first injection increases the associated local reactivity downstream, while the reactivity closer to the nozzle is reduced due to a shorter second injection. This behaviour creates multiple ignition kernels originating at various locations due to temporally less stratified reactivity (c.f. Figure 4(c)), which leads to a higher peak HRR value ( ≈ 47%) and a shorter peak HRR window (by a factor of 1.9); i.e. approaching from stratified to volumetric ignition mode, which is unfavourable in the RCCI context. Additionally, while a decrease in CHR can be observed for the 60/40% case due to the shorter first injection and less LRF entrainment, the other two cases (70/30% and 80/20%) yield similar CHR values.

IDT [CAD BTDC] Peak HRR [kJ/s] Peak HRR window [CAD]
The results presented for Cases 1-3 demonstrate that depending on the injection strategy, a good amount of control can be obtained over the shape of the HRR profile and ignition mode. In particular, adjusting the SOI 1 timing shows clear potential for IDT and HRR control in the RCCI context. According to the results in Figure 5(a) and Table 3, fixing SOI 2 and mass split ratio while sweeping SOI 1 provides a more systematic control over maximum HRR and HRR peak window to switch from volumetric to stratified ignition mode. More obvious differences in CHR profiles were also observed in Case 1 compared to the other two cases. Therefore, in the remainder of this study, we focus only on the effect of SOI 1 (Case 1) on the ignition mode.

HRR analysis for different SOI 1 timings
As justified in the previous section, we next focus solely on the connection between SOI 1 and the ignition mechanism. To further investigate HRR characteristics in Case 1 (sweep of SOI 1 with fixed SOI 2 and mass split ratios), the spatial distribution of HRR for different SOI 1 timings is provided in Figure 6 at different time instances (t 1 -t 4 ) marked in Figure 5(a). The marked points correspond to the peak of the first stage ignition HRR (t 1 ), 10% of total HRR (t 2 ), maximum HRR (t 3 ), and 80% of total HRR (t 4 ), respectively. First, it can be observed that for SOI 1 = 21.6 CAD (left), the overall high reactivity due to short dwell time results in a volumetric ignition initiating approximately 40-60 mm from the injector (t 1 -t 2 ). The ignition then propagates within the spray charge, consuming the HRF (t 3 ). Finally, the HRR is noted to primarily occur at lower mixture fraction values and at the boundary of the spray, indicating HRR originating from the ambient methane-oxidiser charge (t 4 ), which is common in all cases.
For SOI 1 = 28.8 CAD (middle), the initial HRR originates from the rich downstream region associated with the first injection (t 1 ). Subsequently, two independent ignition locations spatially confined from one another can be observed (t 2 -t 3 ). It can be seen that the ignition locations correspond to rich pockets within the spray charge, demonstrating the direct effect of local fuel concentration (and reactivity) on the ignition characteristics. Advancing SOI 1 not only changes the ignition characteristics from a single volumetric ignition pocket to two independent ignition events in space, it also changes the HRR profile associated with the case. Consistently, Figure 5(a) showed an HRR profile with two lower peaks for SOI 1 = 28.8 case, compared to the single peak in SOI 1 = 21.6 case.
Finally, for SOI 1 = 36.6 CAD (right), the initial HRR can be observed in the very small rich pockets downstream, associated with the first injection (t 1 ). The smaller rich regions downstream are due to the early SOI 1 timing and subsequent long mixing time of the first injection. However, the initial ignition pocket is formed upstream closer to the nozzle, originating from the second injection (t 2 ). The downstream region associated with the first injection ignites later in small pockets at multiple locations (stratified ignition), which indicates the reactivity stratification caused by mixing of the HRF and LRF (t 3 ). Such ignition characteristics create the double-peaked HRR curve observed in Figure 5(a) with a lower peak HRR value and much wider distribution, often desired in RCCI applications.
In relation to Figure 6, Figure 7 displays scatter plots and conditional means of the HRR distribution in the mixture fraction space for three SOI 1 values. The data points are coloured with the axial distance from the nozzle. For t 1 , all three SOI 1 timings show similar trends for the HRR conditional mean. It can be seen from the scatter data that the downstream HRR contribution for SOI 1 = 21.6 CAD is closer to the injector and distributed over a broader spatial range, indicating the local high reactivity region formed due to the short dwell time. At time t 2 , SOI 1 = 21.6 CAD and 28.8 CAD present a bi-modal HRR in the mixture fraction space with distinct contributions on the rich and lean mixture fractions. However, for SOI 1 = 36 CAD, the distribution of HRR is more biased towards the lean mixture fraction values and the colour range of the data points indicates a single ignition pocket. At t 3 , the effect of reactivity stratification can be observed for the SOI 1 = 36 CAD case (right). The HRR has two distinct peaks: the first peak is associated with the second injection and closer to the injector with a narrow distribution around Z st , and a second peak below Z st with a wider distribution in Z. The second peak below Z st indicates stratified reactivity, which can also be observed in Figure 6 as multiple ignition pockets downstream. Finally, for t 4 , the conditional means depict more similarity between the cases for the latest combustion phase. In particular, the majority of the injected HRF is already consumed at that point, and the main source of HRR is the spray charge boundary propagating through the ambient methane-oxidiser charge.
The observations on the spatial and mixture fraction space distribution of HRR suggest that the injection timing of the first injection has a significant impact on the HRR mode. In particular, the different HRR modes are correlated to different reactivity stratification in the system. It has been discussed in the literature [4,5] that mixture and/or thermal stratification are the main factors that determine the level of reactivity stratification. In the following section, we provide a frozen-flow ignition index analysis for the three cases discussed in this section (sweep of SOI 1 with fixed SOI 2 and mass split ratios) to shed more light on the effect of injection on the reactivity stratification, as well as its dependence on thermal and mixture stratification.

Frozen flow ignition index analysis for different SOI 1 timings
As mentioned earlier, reactivity stratification in the investigated multi-injection sprays can influence HRR distribution in time and space thus affecting the combustion duration. Reactivity stratification in multi-injection sprays may be related to temperature and mixture stratification generated by different dwell time and mass split ratios. In a recent study, Karimkashi et al. [54] showed that in DF mixtures, different mixture stratification levels can change the mode of combustion; i.e. deflagrative versus volumetric. In particular, it was shown that a broader distribution of reactivity in the mixture fraction and physical space is favourable, leading to longer combustion duration and flame initiation.
In this section, a simple approach is utilised to elaborate on the role of the mixture and temperature stratification on reactivity distribution and ignition mechanism with SOI 1 timing sweep (Case 1). For this purpose, we isolate reactivity from mixing by convection and diffusion. A frozen-flow ignition index analysis, similar to the analysis performed in our previous work [21,24] is utilised. Such a frozen-flow approach uses the LES data at a time instance close to the second-stage ignition and solves only for the chemistry in each computational cell. Each cell is assumed as an independent 0D homogeneous reactor, which enables investigating the chemistry and reacting characteristics of the flow without the influence of mixing. The analysis is expected to provide further knowledge on the role of temperature versus mixture stratification on reactivity distribution and ignition mechanism.
The analysis solves only the species and enthalpy equations at a specific time instance taken from the LES data. For each case, the time instance of IDT-0.2 ms is used as a frozen-flow instance containing the thermochemical properties f (Y , p, T) such as mass fraction of species, temperature and pressure of each computational cell. Then, individual 0D homogeneous reactors are simulated in each computational cell for 6 ms (42 CAD), using these f (Y , p, T) conditions. Considering the selected time instance from the LES data, it is expected that ignition starts after ≈ 0.2 ms in the frozen-flow simulations. In order to distinguish between non-ignitive and prone to ignition cells, here 6 ms is considered which is a very long ignition time compared to 0.2 ms until the spray IDT. For each cell, an IDT 0D value is calculated by finding the maximum gradient of temperature (dT/dt) max as a function of time. The calculated IDT 0D is interpreted as an indicator of the reactivity in each cell. It should be noted that the turbulent strain is neglected and the obtained IDT 0D results only give information on the reactivity characteristics of individual computational cells. Figure 8 shows the PDF distribution of IDT 0D for each SOI 1 timing. We note that here only Case 1 (SOI 1 timing sweep) is analysed. The analysis is carried out by selecting all spray cloud data points (Z ≥ 1e −4 ) from the 3D IDT 0D fields. It can be seen that for short dwell (SOI 1 = 21.6 CAD), PDF of IDT 0D has a higher peak at lower IDT 0D , indicating an overall high reactivity within the spray charge. With advanced SOI 1 , it can be observed that the PDF peak gets smaller around lower IDT 0D values with a more uniform distribution, indicating stratified reactivity. Therefore, with advanced SOI 1 , a favourably longer combustion duration is expected as previously discussed in Figure 5.
While the PDF analysis presented in Figure 8 gives information on the overall reactivity within the spray charge, the effect of temperature and mixture fraction on reactivity as well as the spatial distribution of reactivity should be investigated to get a better understanding of the effect of chemical stratification on ignition characteristics. The distribution of IDT 0D in both T-Z space and physical space is presented in Figure 9 for Case 1 (SOI 1 timing sweep). The presented IDT 0D data is coloured with its corresponding values up to 4 ms. Data points with IDT 0D > 4 ms are considered 'non-ignitive' and marked with the black colour. The distribution of IDT in T-Z space in Figure 9 (left) shows that the reactivity stratification is mainly due to the changes in mixture fraction Z, rather than temperature. This observation demonstrates that the main factor in creating reactivity stratification is through creating different blending ratios of LRF (methane) and HRF (n-dodecane), while thermal stratification plays a secondary role. Our findings are consisted with our previous studies on HRF-LRF mixing in Refs. [21,54], wherein mixture fraction stratification is identified as the dominant factor in determining the reactivity distribution in DF mixtures relevant to RCCI engines. According to Figure 9 (left), while for Z > Z st the reactivity is relatively high (lower IDT 0D values), around and below Z st the IDT 0D starts to increase with leaner conditions, indicating lower reactivity. It is also important to point out that the effect of T and Z on IDT 0D is independent of the injection strategy: the injection strategy determines how the range of different reactivities is distributed in space, which further determines the reactivity stratification required for RCCI. Figure 9. Scatter plot of IDT 0D in temperature-mixture fraction space (left) and spatial distribution of IDT 0D (right) at IDT-0.2 ms time instance for different SOI 1 conditions in Case 1. It can be observed that while the reactivity stratification is very similar in mixture fraction space for all the cases, the distribution of reactivity in physical space determines the ignition characteristics.
The distribution of IDT 0D in spatial coordinates is also provided in Figure 9 (right) for 3 different SOI 1 timings. For SOI 1 = 21.6 CAD, the spray is observed to be mostly rich posing high reactivity throughout region (I). The low-reactivity regions with higher IDT 0D values are clustered around the periphery of the spray instead of being distributed in a stratified manner. The uniform distribution of high reactivity conditions throughout the spray explains the volumetric ignition and high peak HRR observed in Figures 5(a) and 6. With advancing the SOI 1 to 28.8 CAD, a stratified lower reactivity region downstream (II) is formed in between the high reactivity regions confined within the rich mixture conditions (III). The high reactivity rich pockets denoted with (III) correspond to the major HRR locations demonstrated earlier in Figure 6 for SOI 1 = 28.8 CAD, which further proves that the stratified reactivity within the two high reactivity pockets prevents a fully volumetric ignition and provides lower peak HRR. Finally, for SOI 1 = 36.6 CAD, it can be observed that while there is a high reactivity region near the nozzle associated with the second injection, the downstream spray mixture is heavily stratified (IV) due to the very early first injection and the subsequent over-leaning as a result of the long mixing time. Such a stratified reactivity distribution downstream is consistent with the small HRR pockets originating downstream for SOI 1 = 36 CAD case, introduced earlier in Figure 6. The lower reactivity and stratified reactivity distribution due to early SOI 1 timing results in longer combustion duration and non-volumetric ignition within small pockets, reducing the peak HRR value.
As a closing note, while the short dwell time in SOI 1 = 21.6 CAD leads to volumetric ignition, advancing the first injection in SOI 1 = 36.6 CAD provides a well-distributed spatial reactivity which is a pre-condition for flame initiation in the ambient charge according to the discussions in Ref. [54]. We note that in the present study, our simulations are carried out slightly after the ignition point and flame initiation in the ambient charge will be considered in our future research.

Conclusions
In this paper, LES and finite-rate chemistry were used to investigate different split injection strategies and their effect on ignition delay time and heat release rate characteristics in conditions emulating RCCI engine charge preparation process. A compression heating model was employed to create dynamically changing ambient pressure and temperature based on a given motored pressure profile. It was shown that the ignition timing and heat release rate characteristics can be controlled by utilising different split injection strategies. In particular, adjusting the first injection timing was found to be an effective way to control the reactivity stratification and peak heat release rate. The main findings of the study are: (1) In the studied split injection cases, increasing the dwell time between the injections causes lower local reactivity and delayed IDT. (2) Ignition and HRR characteristics can be controlled by adjusting the timing of the injections. In particular, adjusting the SOI 1 timing was shown to provide good control over HRR characteristics. (3) With advanced SOI 1 , the mixing time of the first injection is increased, leading to over-dilution downstream. This delays the downstream reactivity compared to the subsequent second injection, preventing a volumetric ignition event; i.e. a more stratified ignition. In Case 1, when shifting SOI 1 from 21.6 to 36 CAD BTDC, approximately 21% decrease in the maximum HRR is reported along with a factor of 2.8 increase in combustion duration. (4) The frozen-flow ignition index analysis shows that the reactivity stratification is mostly controlled by mixture fraction rather than temperature stratification. For advanced SOI 1 , the downstream mixture was found to be low reactivity and heavily stratified, consistent with the findings from the HRR analysis.
Although the numerical framework utilised in this study takes the dynamic thermophysical conditions due to piston motion into account, the wall effects and the in-cylinder mixing of the two fuels also play an important role in mixture formation in RCCI conditions. In addition, the results presented here provide an ignition structure for the particular operating conditions, and the in-cylinder mixing and ignition characteristics are expected to be different for different operating conditions (e.g. compression ratio, injection timing, injected mass). While the current study focuses on the interactions between the injections in a split injection scenario and creating reactivity stratification within the mixture, in future studies the wall effects should also be taken into account to be able to make more practical conclusions on real RCCI engines.
for this study were provided by CSC -Finnish IT Center for Science. The first author has been financially supported by the Merenkulun Säätiö.

Disclosure statement
The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

Funding
The present study has been financially and academically supported by the Wärtsilä Finland Oy and the Academy of Finland (grant numbers 297248,318024,332784)