CFD Simulation of Pre-chamber Spark-ignition Engines- A Perspective Review

1. Introduction

Given the current global energy demand and shortages, it is imperative to develop technological solutions that enhance fuel conversion efficiency in internal combustion engines. In this regard, lean burn combustion is a pivotal technique for boosting thermal efficiency and reducing pumping losses [1,2]. Leaner fuel-air mixtures have been proven in several studies [1,2] to reduce fuel consumption and NOx emissions and improve thermal efficiency. However, there are several issues that must be overcome when applying the lean burn strategy to various fueled gasoline engines. One issue is that three-way catalytic converters cannot be used in gasoline engines. Another challenge is that, in order to comply with current NOx emission regulations, λ levels need to be increased to the range of 1.8-2.0 [2], which exceeds the lean stability limit (λ~1.4). Additionally, leaner sparse air-fuel mixtures are known to pose serious challenges in terms of combustion stability degradation, increased cycle-by-cycle variation, reduction in thermal efficiency, and increased levels of unburned hydrocarbon (UHC) emissions.

Another prevalent method for enhancing the efficiency of SI engines is to sustain a stoichiometric air-to-fuel ratio (λ = 1) while attenuating the oxygen concentration in the mixture through exhaust gas recirculation (EGR) [3,4]. This approach offers the clear advantage of utilizing a three-way catalytic converter, which is particularly effective in mitigating NOx emissions [4]. Furthermore, an important effect of utilizing EGR as diluents is not only reducing knock tendency at low-speed high-load conditions but also improving fuel economy by eliminating fuel enrichment [5]. Nonetheless, it is also accompanied by challenges such as increased cycle-to-cycle variability, reduced flame propagation speed, and the occurrence of misfires, issues that are typically associated with lean combustion absent of EGR.*

Previous literature [6,7] has identified that the main problems in lean blending operations, where either air or EGR is used as diluents, include the need for increased energy to initiate combustion and the low flame propagation velocity. Therefore, the most significant obstacle to implementing lean combustion technology in practical engine designs is developing the ignition system’s capacity to reliably spark the fuel-lean blend. Lean-burn internal combustion engines targeting reduced nitrogen oxide (NOx) emissions necessitate substantial ignition power and widespread ignition sources to spark and incinerate the lean premixed primary charge effectively. This strategy aims to enhance efficiency and diminish residual hydrocarbon emissions. To circumvent this limitation, alternative ignition approaches such as laser-induced plasma ignition [8,9], diesel pilot injection [10], and pre-chamber ignition systems have been suggested and thoroughly examined by researchers. Despite extensive research efforts in this domain, only a select few of these methods have successfully transitioned into commercial applications. Laser-induced plasma ignition technology offers advantages such as precise control over ignition timing and potential for high energy ignition. However, drawbacks include cost implications, complexity in implementation, and sensitivity to environmental conditions. While extensively researched and considered promising, widespread commercial adoption of laser-induced ignition technology is currently limited due to these challenges. The technology still remains in research phase. The pilot injection of diesel also presents several drawbacks. Firstly, introducing an additional fuel source increases the engine’s complexity. Secondly, since ignition relies on the auto-ignition of diesel, it requires elevated in-cylinder temperatures and lacks precise control. Consequently, the advantages of pilot injection may be constrained under specific operating conditions.

All TJI(Turbulent Jet Ignition) systems are equipped with a spark plug, a compact pre-chamber(<3% of clearance volume), and one or more orifices through which the pre-chamber’s reacting contents are injected into the main chamber. The spark plug ignites the mixture in the pre-chamber, creating a pressure differential that drives the flame through the nozzle. This results in a jet of intermediate combustion products containing active radicals and high-temperature burned materials, which rapidly ignite the charge in the main chamber, initiating fast and turbulent combustion.

Recently, pre-chamber spark ignition (PCSI) technology has been actively implemented as an ignition method not only in highly-downsized TGI engines [11,12] and heavy duty gas engines [13,14] but also in hybrid vehicle engines [15]. Very recently, PCSI system have been also applied carbon-neutral fueled passenger car engine for automotive application.

Summarizing the research results of the aforementioned PCSI (Pre-Chamber Spark Ignition) system, it can be observed that by supplying the higher amount of energy in the main chamber at the start of combustion, stable combustion can be achieved under very lean conditions(λ>2.5) for gasoline fueled engine [16,17] and λ of up to 2.6 for heavy duty natural gas engine with near zero NOx emission [18]. Recently, a pre-chamber was applied to a gasoline-fueled PFI engine, achieving not only 52.5% indicated thermal efficiency but also stable ultra-lean burn combustion at a lambda of 2.4 [19]. These make PCSI a more practical and economical lean combustion technique compared to other lean combustion methods.

1.1. Pre-Chamber Ignition Concept

The first pre-chamber combustion engine is often attributed to the Ricardo Dolphin engine, developed in the early 20th century [20]. Sir Harry Ricardo, a pioneering British engineer, designed this engine to improve the combustion efficiency and reduce knock in gasoline engines.

Subsequently, pre-chamber combustion technology evolved into the jet igniter [21], characterized by a much smaller orifice connecting the main chamber and the pre-chamber combustion cavities. The smaller orifice size causes the burning mixture to travel quickly through the orifice, which extinguishes the flame and allows reacting active radical species to reignite at a certain distance away from the pre-chamber. The concept of jet ignition was introduced in the late 1950s by Nikolai Semenov, famous for developing the general theory of chemical chain reactions. In 1981, the LAG system was implemented into the powertrain of the Volga passenger vehicle, which was equipped with a cam-actuated injector that introduced a rich (λ=0.5) mixture into the pre-chamber, subsequently igniting an ultra-lean (λ=2) mixture in the cylinder [22,23].

The detailed development process and characteristics of the PCSI (Pre-Chamber Spark Ignition) system are thoroughly described in previously published review papers [22,23]; therefore, this paper will omit those details for brevity. The most successful representative example of a modern PCSI engine with a divided chamber (pre-chamber) is Honda’s Compound Vortex Controlled Combustion (CVCC) system, developed between 1968 and 1972 [24] and first introduced in the 1975 Honda Civic. This technology was innovative at the time as it met the suddenly stringent emission regulations of the California Air Resources Board without the need for a catalytic converter.

However, the CVCC (Compound Vortex Controlled Combustion) system was not continuously applied to Honda cars due to several reasons such as advancements in Emission control technology, cost and complexity and shift to direct Fuel Injection. As emission regulations became stricter, newer and more effective technologies were developed. Catalytic converters and advanced electronic fuel injection systems became standard, offering better performance, efficiency, and lower emissions compared to the CVCC system. Namely, Direct fuel injection systems offered better control over the combustion process, leading to improved fuel efficiency and lower emissions. These systems became more prevalent in the automotive industry, replacing older technologies like CVCC. Additionally, emission standards continued to evolve, requiring more sophisticated and integrated approaches to meet the new requirements. The CVCC system, while effective in its time, could not keep up with the increasingly stringent regulations without significant modifications.

Recently, pre-chamber combustion techniques have been gaining renewed attention as a key technology for environmentally friendly, next-generation engine-based vehicles. This can be attributed not only to advancements in machining and production technologies and air-fuel ratio control techniques but also to improvements in turbulence flow control both inside and outside the pre-chamber. These improvements have enhanced combustion stability under lean burn conditions and advanced knock control techniques. Furthermore, the development of combustion visualization techniques and chemiluminescence using RCEM (Rapid Compression Expansion Machine), along with advancements in CFD analysis techniques, have significantly contributed to a more detailed analysis and understanding of turbulent jet flame development and combustion characteristics.

This section aims to present a concise overview of the combustion strategies and principles employed in the PCSI (Pre-Chamber Spark Ignition) engine, which has recently attracted significant interest.

This ignition strategy can be implemented in two main ways: active [25,26](or scavenged) and passive (or unscavenged) systems [27,28,29] as shown in Fig.

The active system features a dedicated fuel injector or a miniature check valve within the pre-chamber that does not allow the backflow of pressurized gases and combustion products from the pre-chamber to the fuel line, ensuring the air-to-fuel mixture remains optimal for combustion, near stoichiometric conditions [31]. Hence, in active, fueled PCSI system, an auxiliary fueling event occurs in the pre-chamber via an injector or check valve, enabling the formation of a stoichiometric or fuel-rich mixture near the spark plug. In contrast, the passive system simplifies the design by excluding this secondary fuel injector, which lowers assembly and packaging costs, making it a viable option for passenger vehicles. Thus, fuel injection is performed solely in the main chamber, utilizing either port fuel injection (PFI) or direct injection (DI), ensuring an identical air-fuel ratio in both combustion chambers. The air-fuel mixture then flows into the pre-chamber through interconnection orifices. Consequently, the geometry of the pre-chamber and connecting pipelines plays a crucial role in the combustion process. The importance of this ignition system has been confirmed by numerous studies over the past decade. Several studies examined not only the fundamental aspects of jet ignition using both experimental and simulation methods [32] but also the effects of pre-chamber configuration on jet dynamics and combustion physics [33,34].

In recent, Benejes et al. [35,36] assessed the synergies between pre-chamber ignition and other strategies aimed at enhancing engine thermal efficiency, such as lean burn and CO₂-free fuel use such as hydrogen [37,38,40] and ammonia [36].

Furthermore, numerous studies have been also conducted to optimize the geometry of the pre-chamber to maximize the effects of turbulence-chemistry interaction by the shear stress of the flow through the holes. For this goal, numerous previous studies [39,40,41,42] investigated how nozzle diameter and pre-chamber volume influence combustion performance.

1.3. Review Objectives

Review papers on the performance enhancement and emission reduction effects of PCSI combustion engines, as well as on EGR (Exhaust Gas Recirculation) and lean burn limits, have been consistently published from 2010 to 2023 [2,23,43,44,45]. This indicates ongoing research and technological advancements in this field. Numerous studies have focused on improving the accuracy of CFD analysis techniques to numerically model the complex physical phenomena of PCSI engines more accurately.

Recently, pre-chamber combustion techniques have been gaining renewed attention as a key technology for environmentally friendly, next-generation engine-based vehicles [43,45]. However, PCSI (Pre-Chamber Spark Ignition) engines are strongly influenced by the turbulence-chemistry interactions between the turbulent flame jet generated in the pre-chamber and the turbulent characteristics such as swirl and tumble within the main chamber [2,45]. Therefore, the optimal design of PCSI engines requires a thorough analysis and prediction of complex physical phenomena, including multi-mode combustion, optimization of the nozzle diameter and number, orifice length, shape, and volume of the pre-chamber, the stretching and quenching of the turbulent jet torch, the rapid changes in turbulence length scales, and the thermal mixing between the turbulent jet and lean mixture in the main chamber [2,44]. For this purpose, the utilization of three-dimensional computational fluid dynamics (CFD) capable of simulating physico-chemical processes in PCSI system is essential. This study aims to review and analyze the technical achievements and limitations of CFD analysis techniques in PCSI engine research over the past 20 years, and to discuss the prospects for future advancements in CFD technology.

This paper reviews the evolution and application of CFD in PCSI engines, with a detailed look at the numerical modeling of the complex physico-chemical processes involved in the PCSI ignition system and the contemporary issues regarding the limitations of commercial CFD codes. Finally, the prospects for future advancements in CFD technology is discussed. We performed an in-depth exploration of relevant literature by utilizing widely used databases such as Web of Science Core Collection and Google Scholar.

2. Major Achievements of CFD Applications on PCSI Engine design

Since the early 2000s, research on PCSI using CFD simulation has primarily focused on optimizing the geometric configuration of the pre-chamber, including the nozzle diameter and number, orifice length, shape, and volume. Recently, these research themes have continued with the application of upgraded and more accurate numerical models. These research topics utilizing CFD are schematically illustrated in Figure 2. The research topics using CFD can be summarized into two categories.

First, some articles [46,47] presents a study of the influence of pre-chamber volume and nozzle diameter on the resultant ignition characteristics by using CFD simulation. The influence of the prechamber’s nozzle Orifice Diameter and thereby the total cross-sectional area of the nozzle orifices has been investigated several times. The influence of the orientation and number of nozzle orifices connecting the prechamber and the main combustion chamber, as well as the prechamber’s internal volume and shape is evaluated.

Silva et al. [48] constructed a RANS-based CFD model for a methane-fueled engine, integrating a well-stirred reactor combustion model with a methane oxidation mechanism. They evaluated the impact of geometric parameters of the passive pre-chamber, including throat diameter, nozzle length, and nozzle diameter, on combustion performance. The model was validated by comparing pressure traces in both chambers under motoring conditions at 1200 rpm. The results demonstrated that throat diameter significantly affects pressure build-up and residence time within the pre-chamber, whereas nozzle diameter influences both peak pressure and residence time.

Distaso et al. [49] performed 3D CFD simulations on an active pre-chamber ignition system in a lean-burn methane engine, analyzing six phases of scavenging and combustion: filling & scavenging, mixing, flame propagation, ejection, reburning, and expulsion & extraction. Thelen et al. [50,51,52] conducted extensive CFD modeling of the TJI process. Their 3D CFD simulations, incorporating detailed combustion chemistry, investigated the effects of various orifice diameters (1.0 mm, 1.5 mm, 2.0 mm, and 3.0 mm). The findings indicated that a 1.5 mm orifice diameter provides the quickest ignition and overall combustion, based on pressure data, while a 1.0 mm orifice diameter results in higher jet velocity but a longer burn duration compared to larger diameters.

Secondly, using 3D CFD simulation, research has been conducted on the combustion characteristics within the pre-chamber, which cannot be measured experimentally, and the effects of the turbulent jet generated from the nozzle on the combustion characteristics and NOx emissions in the main chamber. Through these studies, various parametric studies have reported on the optimal geometry of the pre-chamber that can extend the lean limit of the main chamber. Recently, research has been published analyzing the turbulence-flame interaction caused by the hot turbulent jet and the turbulence (swirl and tumble) within the main chamber, along with the predictive limitations of related combustion models used in commercial CFD software. Additionally, studies using optically accessible engines and engine-like geometries aim to capture the turbulence and mixture characteristics related to the hot turbulent jets from pre-chamber, while achieving engine-relevant thermodynamic conditions [52]. Investigations in a rapid compression machine with generic pre-chamber geometries have been conducted, providing insights into combustion and the effects of nozzle size. The advantage of optical data compared to pressure-only data from a metal engine is the additional validation opportunities for model development. Hence, a comparison of characteristic data such as jet exit timing, jet penetration velocity, and cyclic variations can be used to validate simulations and provide additional understanding of the phenomena observed in the engine under similar conditions. Currently, from a simulation perspective, this includes 0D models for turbulence generation/dissipation, heat transfer, and combustion; 3-D RANS calculations using level-set combustion models for industrial CFD and design optimization; LES calculations with detailed combustion models for further understanding of mixing and combustion phenomena; and 2-D/3-D DNS calculation tools that provide unprecedented high-fidelity data for fundamental phenomena studies and model development. However, a physically accurate analytical model for PCSI combustion has not yet been developed [43,45,53]. Consequently, there remain uncertainties in the understanding of combustion phenomena, as well as in the accuracy of turbulence and combustion models, especially concerning turbulence-flame interactions [45,53]. These aspects will be discussed in detail in the following sections. Recently, an experimental study has been conducted to optimize all previously researched pre-chamber shapes using Taguchi’s method, followed by redesigning to enhance performance [60].

3. CFD Modeling PCSI Engines

3.1. CFD Software

The development of user-friendly and feature-rich CFD software has greatly expanded its usage across different research fields. In order to assess the reasons for the increased usage of CFD in the PCSI engine field in recent years, this section provides an overview of the major achievements of CFD in the PCSI engine over the past about 30 years.

I reviewed various CFD software used for modeling the PCSI engines and the RCEM, as shown in Figure 3. Out of 54 studies [52,53,54,55,56,57,58,59,60,61,62,63,64,65,66,67,68,69,70,71,72,73,74,75,76,77,78,79,80,81,82,83,84,85,86,87,88,89,90,91,92,93,94,95,96,97,98,99,100,101,102,103,104,105,106,107,108,109,110,111], 1 did not specify the CFD software employed by the researchers. Among the remaining articles, the majority used the commercial software CONVERGE™ [54] for simulations. Other programs like STAR-CD™ [55], FIRE™ [56] and VECTICS [57] were used much less frequently. CFD codes that appeared in only one study were grouped into an ‘Others’ category, which includes two programs such as KIVA-V3 [58]. Additionally, some studies used open-source software and custom-built codes: three used OpenFOAM [59].

The selection of CFD software is primarily influenced by its capabilities, accessibility, and user preferences. Commercial programs, which are typically updated annually with enhanced features, are adept at handling a wide spectrum of simulation tasks, as illustrated in Figure 3. Their user-friendly interfaces further enhance their popularity by simplifying the simulation workflow. However, the substantial cost of licenses and limited customization options have driven some users towards open-source alternatives or the development of custom CFD codes. While the utilization of these software tools has enhanced the accessibility of CFD across a wider users, it is imperative to acknowledge that mere accessibility is insufficient in addressing the challenges outlined in PCSI engines as detailed earlier. The following sections discusses four key focal points investigated to overcome these challenges. These include the distinctive features and comparative analysis of LES and Reynolds-averaged Navier-Stokes (RANS) models, the establishment of appropriate inflow and outflow boundary conditions, multi-mode combustion models, wall heat transfer, turbulence-chemistry interaction, rigorous verification and validation processes, and the development of best practice guidelines in the field.

3.1. Turbulence Models

In computational fluid dynamics (CFD), turbulence modeling can generally be divided into direct numerical simulation (DNS), large eddy simulation (LES), and the Reynolds-averaged Navier-Stokes (RANS) model. Due to the challenges associated with applying DNS to real-world engineering problems, especially in the field of automotive engineering, LES and the RANS model are more frequently utilized.

Turbulence models, especially simpler ones like RANS models, often struggle to accurately predict not only flow separation and reattachment, which are critical in engine intake and exhaust processes but also turbulence length scale which is crucial for determining not only the rate of energy dissipation and mixing efficiency in engine flows but also combustion processes.

Moreover, the flow fields encountered in PCSI engines are highly complex, involving phenomena such as impinging flows, jet-like flow, flow separation, strong swirl, highly variation of turbulent length scale, and vortex shedding (refer to Section 2). Accurately and universally modeling the all these turbulent flow characteristics remains a formidable challenge.

As shown in Figure 4, in the analysis of 52 publications reviewed, , the majority, 52.8% studies, relied solely on the RNG k-ε turbulence model, while 11.91% studies opted for LES alone. Moreover, 7.14% studies integrated both LES and RANS models in their research. Common RANS turbulence models cited were standard k-ε, renormalization group (RNG) k-ε, shear stress transport (SST) k-ω. Sub-grid scale models for LES included Smagorinsky only.

3.1.1. RANS Turbulence Models

The RANS-type turbulence models are simple to use, computationally inexpensive, and economical [113,114].

Hence, given its simplicity and shorter computation time, the RANS-based model, which represents the average behavior of important parameters, is typically chosen over the LES (Large Eddy Simulation) model. However, the RANS-type k- ε turbulence models are based on Boussinesq’s isotropic eddy viscosity assumptions, and it is well known that have several problems with deteriorated prediction performance in cases of swirling flow, separation and reattached flow and flows with large rapid extra strains. Therefore, numerous variants of the k- ε model have been researched, reported, and are currently in use to overcome these challenges [115,116,117]. The existence of numerous k- ε model variants is due to the ε equation with 3-4 model coefficients, which are empirically derived or determined through ad-hoc methods based on various turbulent flow patterns [*]. Additionally, turbulent flows can vary significantly in their characteristics depending on the geometry, flow velocity, pressure gradients, and other factors. To accurately capture these diverse turbulent flow patterns, modifications to the standard k-epsilon model, particularly in the ε equation, are necessary [115,116,117]. Each variant attempts to address specific deficiencies of the standard model in certain flow regimes, such as free shear flows, boundary layers, or rotating flows. Consequently, a universally accurate RANS-type turbulence model has not yet been achieved. Therefore, engineers face the challenge of selecting a turbulence model that is appropriate for the specific turbulent flow characteristics of the geometry of interest.

The RNG (Re-Normalization Group) k-ε model [118], which is the most frequently adopted in previous studies for PCSI system for various engine types [115,116,117], is one of the variations of the standard k-ε model used for turbulence modeling in various kinds of CFD simulations. It incorporates additional terms that accounts for the interaction between turbulence and mean strain rate, which is not present in the standard k-ε model and modifications to the standard k-ε model to improve accuracy for certain flow conditions, particularly for swirling and highly strained flows. To achieve this goal, this model contains a strain-dependent correction term in the constant C_1ε of the production term in the RNG k-ε model’s dissipation rate(ε) equation [115,116,117,119]. The RNG k-ε model is particularly useful for simulations involving complex flow features such as swirling flows, recirculating flows, and flows with high strain rates, making it suitable for jet-like flows and complex industrial applications in which the velocity gradients are significant causing intense mixing and variations in velocity. [54,115,116,117]. In PCSI internal combustion engines, the scavenging processes in pre-chamber, along with turbulent flame jets and combustion, create regions of high strain rates. Particularly, turbulent jets issuing from nozzles of pre-chamber are typical examples of flows with high strain rates. Therefore, for these reasons, applying RNG k-ε model to PCSI engines extensively is considered practical for engine simulation, as it provides a good balance between accuracy and computational efficiency. However, noteworthy demerit of RNG k- ε model is the fact that near-wall treatment of this model can struggle with accurately predicting flows close to walls, particularly in cases involving adverse pressure gradients or separation [54,115,116,121,122,123]. Hence it should be careful to use RNG k- ε model if the turbulent jet interacts significantly with cylinder walls or piston head surface. Additionally, it is clear that RNG k- ε model is less accurate for detailed Structures: and does not capture the detailed eddy structures as well as LES(Large Eddy Simulation) [115,116,117,119,120]

The k- ζ -f turbulence model [121,122,123] is the second most commonly used turbulence model in PCSI engine CFD simulation after the RNG model. The k-ζ-f model is an extension of the eddy-viscosity concept and includes three transport equations for turbulence quantities, namely turbulent kinetic energy(k), Turbulence Frequency (ζ) and Dissipation Rate (f). This is selected as the turbulence model in previous studies due to its high accuracy and convergence stability. This model, optimized from Durbin’s near-wall turbulence closure model [124,125], is a variant of RSM(Reynolds Stress Model) turbulence models [124,125,128]. The k-ζ-f model is particularly useful in complex flow simulations, including those with significant near-wall effects and flow separation [117,121]. It enhances the standard k-ε turbulence model by introducing the wall-normal velocity fluctuation v² and its source term f, which incorporate near-wall turbulence anisotropy and non-local pressure-strain effects. The careful introduction of these relaxation terms eliminates the need for damping functions. Additionally, this model improves numerical stability over the original v²-f model by solving a transport equation for the velocity scale ratio ζ=v²/k instead of the velocity scale v². Moreover, This model demonstrates superior predictive accuracy in heat transfer, surface friction at low Reynolds numbers, adverse pressure gradients, and recirculation regions compared to traditional k-ε turbulence models [126,127]. While it shares characteristics with low-Reynolds-number models, it uniquely eliminates the need for wall functions by being applicable near the wall. Instead, it introduces a transport mechanism for turbulent energy from the wall, effectively representing near-wall viscous damping effects through an elliptic relaxation equation [121,122,123]. Hence, this model may offer superior performance in capturing the detailed flow features and interactions in not only the intake and exhaust processes but also scavenging process inside pre-chamber, crucial for predicting flame propagation, heat transfer, and emissions. Its improved near-wall treatment is advantageous for accurately predicting heat transfer inside orifices of pre-chamber and flow separation around valves and ports.

3.1.2. LES Turbulence Model

Due to the substantial turbulent energy and significant influence on momentum transfer and turbulent mixing carried by large eddies, LES methodologies offer superior accuracy compared to RANS turbulence models. [129]

LES captures flow structures from the domain scale down to the filter scale, necessitating significant resolution of high-frequency turbulent fluctuations. This requires the use of either high-order numerical methods or fine grid resolution when employing lower-order numerical techniques. Therefore, the implementation of Large Eddy Simulation (LES) methods in automotive and mechanical engineering necessitates finely resolved grids with grid points positioned in close proximity to the boundary layers. This results in significantly heightened computational costs compared to Reynolds-Averaged Navier-Stokes (RANS) methods [ over serveral applications [130]. Consequently, previous studies have reported only a limited number of simulations utilizing LES, with the majority being conducted on Rapid Compression Machines (RCEMs) [102,105] rather than on full-scale metal engines [67].

Recently, 3D CFD analysis using LES turbulence models was performed to investigate the engine characteristics of active and passive PCSI engines fueled by natural gas for large ships under a single operating condition, both at stoichiometric and lean conditions (λ=2). Additionally, the heat release rate curves obtained from the analysis were compared with experimental results [67].

However, to perform accurate CFD simulations of PCSI engines using the LES turbulence model, it is necessary to understand the limitations and characteristics of various sub-models.

The most important noteworthy point is that the accuracy of LES heavily relies on the Sub-grid Scale (SGS) models used to represent the unresolved scales. This is because the computational grid limits the size of eddies that can be physically represented. Despite their presence in the flow field, these eddies cannot be resolved because the CFD mesh lacks the resolution to accurately capture and depict them in CFD simulations. In LES simulation, we are particularly concerned with the eddies that are just larger than the mesh size. These eddies are too large to be broken down by molecular viscosity. Therefore, we need to find a way to model and remove these eddies from the grid. If we do not, the turbulent kinetic energy predicted in our large eddy simulation will be too high. These eddies are removed by applying an additional stress term to the Navier-Stokes equations, known as the sub-grid scale (SGS) stress. While various SGS models such as Smagorinsky model [113], dynamic Smagorinsky model [55,56], WALE model [137], and others are available, none are universally applicable to all types of flows, leading to potential inaccuracies in specific scenarios. Different SGS models all they do is provide different methods for calculating the subgrid kinematic viscosity, ${\mathsf{\nu }}_{sgs}$.

One of the shortcomings of the Smagorinsky Subgrid Scale model which is the most frequently adopted in CFD simulation of PCSI combustion is that it contains a model constant Cs called Smagorinsky coefficient that is not universal and depends on the local flow conditions and the fraction of the cell size that gives the sub-grid length scale. Cs for homogeneous isotropic turbulence is around 0.17. This value is very crucial for accurate simulation and has come from an analytical mathematical derivation based on a homogeneous isotropic turbulence. It is noted that modern CFD codes use different values of Cs. In STAR-CD [55], Fire [56], Fluent [131] Cs=0.1 and in PHOENICS [132], Cs=0.17. However, no matter what Cs values are, this is not true in case of rotation or near wall because of too much dissipation near the wall. Additionally, in Smagorinsky model, sub-grid stress is not damped close to wall. Hence, some kind of modification to the model should be needed. Therefore, various types of sub-grid scale models have been developed. The most frequently used near-wall treatment of Smagorinsky model is the Van Driest damping function which is a damping function for proper results in wall-bounded flows [55,56,131]. Figure 9 illustrates the necessary adjustments to the Smagorinsky model for accurate near-wall treatment. To address this, the sub-grid scale stress should ideally approach zero as the wall is approached. This can be achieved through various methods. One option is to implement an entirely different model, such as one based on sub-grid scale kinetic energy, replacing the traditional Smagorinsky model. Alternatively, we can modify the length scale near the wall, reducing it to zero to account for sub-grid scale eddies. Another approach is to adjust the velocity scale, making it less dependent on strain rate. The primary goal of these adjustments is to decrease the sub-grid kinematic viscosity near the wall, which will, in turn, reduce the sub-grid scale stress, effectively simulating the damping effects on the eddies.

Instead of a single user-defined constant Cs, modern CFD codes like Fire, STAR-CD implemented the Dynamic Smagorinsky Subgrid Scale model (134, 135) which computes a local time varying Cs value by test-filtering the flow field at a length scale greater than the grid length scale, which allows it to compute the correct result for wall-bounded flows without the use of damping functions.

The WALE (Wall-Adapting Local-Eddy Viscosity) Subgrid Scale model (136) is a more modern subgrid scale model that uses a novel form of the velocity gradient tensor in its formulation and widely adopted by CFD codes like STAR-CCM+ and OpenForm []. Similar to the Smagorinsky Subgrid Scale model, it suffers from the limitation that the model coefficient Cw is not universal. It is known that the WALE model is seemingly less sensitive to the value of this coefficient than the Smagorinsky model. Another advantage of the WALE model is that it does not require any form of near-wall damping—it automatically gives accurate scaling at walls [136].

Figure 5 schematically shows various options for correcting sub-grid scale kinematic viscosity near the wall. Here, S_ij is the strain rate of the resolved eddies on the CFD mesh.

Additionally, a Coherent structure model(CSM) as a subgrid-scale model is applied [64,67] was adopted in Fire code [56]. It is reported that the CSM gives good predictions and is almost the same performance the dynamic Smagorinsky model for various complex geometries [137].

Additionally, applying appropriate boundary conditions for LES is complex. Inflow boundary conditions, in particular, need to accurately represent turbulent fluctuations, which is challenging to achieve in practice. Due to the aforementioned factors, achieving dependable LES simulations necessitates a higher level of expertise and experience compared to the requirements for utilizing the RANS model [138].

Bolla et al. [100,102] executed numerical studies of RANS-LES comparison using an RCEM to analyze an automotive-sized scavenged pre-chamber, aiming to compare the two turbulence models’ ability to predict turbulence and fuel-air mixing. In this study, a Smagorinsky-type sub-grid scale model [139] was used to compute the unresolved turbulent scales for the LES turbulence model, while the RANS turbulence model employed the time scale Bounded k−ε Turbulence Model in VECTIS, which is an enhanced version of the standard k−ε model for high strain rates or strong adverse pressure gradients. [57]. The outcomes indicated that the RANS-based model could effectively reproduce the key ensemble-averaged flow patterns seen in LES for two pre-chamber setups. However, in RANS often falls short in accurately capturing the radial fuel-air mixing compared to LES.

3.2. Physical Phenomena and Combustion Models of PCSI Engines

The most commonly used combustion models for engine simulations are based on the flamelet combustion regime, particularly for spark ignition (S.I.) engines. In this regime, turbulence can distort and increase the surface area of the flame front while maintaining its inner laminar structure and flame speed. This condition is characterized by a Damköhler number (Da) greater than 1 and a Karlovitz number (Ka) less than 1.

Numerous previous studies [17,25,29,45,140,141,142,143,144,145,146] have investigated the operating mechanisms of TJI in PCSI engines, which are schematically illustrated in the above figure. The characteristics of TJI can be categorized into enrichment, thermal effects, and chemical effects. Key heat transfer phenomena during flame propagation include thermal quenching, which occurs due to rapid heat transfer to solid surfaces as the flame passes through the nozzle, and hydrodynamic quenching, which happens when the flame mixes with the cool, lean mixture as it enters the main chamber. For a detailed explanation of TJI, please refer to the referenced literature [45]. Recent studies [95,141,142,143,144,145,146] have validated that the two-stage combustion process in the main chamber comprises a jet-dominant phase and mixing controlled, which is influenced by the combustion intensity within the prechamber, and a flame propagation phase, which depends on the reactivity of the mixture in the main chamber driven by the in-cylinder bulk flow and the associated turbulence.

As described above, capturing the flame dynamics across the pre-chamber nozzles is a significant challenge due to the complexity of flame propagation through the orifices [147]. Most previous premixed combustion models were based on the corrugated flame zone; however, recent studies have demonstrated that PCSI engines operate at highly diluted condition, which approach the thickened flame regime. Consequently, flamelet-based combustion models, conventionally used for SI engine combustion analysis, are not suitable for application in PCSI engines operating under these ultra-lean conditions. These models fail to capture the combustion behavior of flame quenching and stretching through the orifices, potentially shortening ignition delays in the main chamber [116].

In the studies conducted thus far on 3D CFD combustion analysis of PCSI engines, the combustion models applied have been flamelet-based premixed combustion models originally used for premixed combustion in homogeneously operated spark ignition engines. However, the PCSI engines operate under very lean conditions (λ > 2) and in high Karlovitz number (Ka > 1) regimes due to the presence of turbulent jets and turbocharging, which result in extremely high turbulence intensities. Therefore, to accurately predict the combustion characteristics of PCSI engines, a combustion model must effectively capture multiple and distributed ignition points within the main chamber, covering both premixed and partially premixed combustion regimes. Furthermore, under the typical operating conditions of TJI systems, characterized by lean mixtures and highly turbulent flow fields, the foundational assumptions of flamelet-based models may no longer be applicable [93]. Therefore, recent studies have raised questions about the predictive performance of flamelet-based combustion models and have made significant efforts to find countermeasures. In this section, we discuss the application cases, comparative predictive performance, and limitations of combustion models that have been applied to the combustion analysis of PCSI engines. Additionally, we review various methods to enhance predictive performance and overcome the limitations of flamelet-based models.

3.2.1. Flamelet Assumption

As mentioned above, combustion models for PCSI engines are most often based on the so-called flamelet assumption [54,55,56,57,147,148], which has frequently been adopted for premixed combustion in homogeneously operated spark ignition engines [148]. The combustion models based on flamelet assumption simplify the complex interactions between turbulence and chemical reactions by assuming that the flame can be represented as an ensemble of thin, locally laminar flame structures, or “flamelets,” embedded within the turbulent flow. They assume a clear separation of scales between the turbulent eddies and the flame thickness. This allows the detailed chemical reactions to be precomputed and stored in a database, which can be accessed during the simulation. The chemical reactions are solved in a laminar flame configuration under varying conditions of temperature, pressure, and mixture composition. The results are stored in flamelet libraries, which provide information about species concentrations, temperature, and reaction rates as functions of mixture fraction and scalar dissipation rate.

In summary, flamelet models strike a balance between accuracy and computational cost, making them suitable for capturing the essential features of premixed combustion in S.I. engines. However, it’s essential to recognize their limitations, especially when dealing with non-premixed or partially premixed combustion regimes. Researchers continue to refine these models and explore more detailed approaches to improve engine combustion simulations.

3.2.2. G-Equation Model

A flamelet-based combustion model, the G-equation, is one of the widely adopted combustion models for simulating the combustion processes of PCSI engines in CFD simulations within the engine modeling community. This approach utilizes a level-set method, which represents moving interfaces or boundaries on a fixed computational mesh. It is particularly useful for problems where the computational domain is divided into two regions separated by an interface. The Level Set modeling technique allows the fluid-fluid interface to move within any given velocity field [1,148,149]. Detailed information on G-equation model can be found in the literature [61,71,81,86,90,91,94,100,101,102,103,104], and only the brief descriptions are provided here. In order to obtain a formulation that is consistent with the well-established use of Favre averages in premixed turbulent combustion, we split G and the velocity vector v into Favre means and fluctuations. Using a number of additional closure assumptions described in Peters [4], one finally obtains governing equations for $\tilde{G}$and its variance $\widetilde{G^{\text{'}\text{'}2}}$ are defined [1,2,3,4]. In order to solve $\tilde{G}$equation, a model for the turbulent flame velocity must be provided. The turbulent flame speed is a key parameter and is typically modeled as a function of the laminar flame speed and turbulence characteristics such as turbulent intensity and length scales. This allows the model to incorporate the effects of turbulence on flame propagation without directly solving the detailed turbulence-chemistry interactions.

Within the G-equation context, several correlations for turbulent burning velocity from literatures [4,56] are presented and evaluated in the previous studies [152]. Among the most popular correlation is Peters’s correlation [150] which is valid for both large-scale and small-scale turbulence as:

$$S_T=S_L+u^{\text{'}}\left\{-\frac{a_4{b}_3^2}{2b_1}Da+{\left[{\left(\frac{a_4{b}_3^2}{2b_1}Da\right)}^2+a_4{b}_3^{2}Da\right]}^{1/2}\right\}$$

(1)

where S_L is laminar flame speed, u^’ the fluctuating turbulent component, δ_L laminarflame thickness, $l_t$ the integral length scale, b₁ and b₃ are model constants corresponding to large and small-scale turbulence enhancement, respectively, and Da Damkohler number which is a ratio between the flow time scale($l_t/u^{\text{'}})$ over the chemical time scale(${\delta }_L/S_L).$The laminar flame speed, S_L depends upon the local pressure, the fresh gas temperature, the local unburned fuel-air equivalence ratio using the Metghalchi and Keck correlation [103,107] and the chemical time arising due to the flame stretching [56,151,153,154,155]. These common correlations for S_L are equations derived from fitting forms based on combustion experiments conducted over various temperature and pressure ranges. Therefore, outside the range of these correlations, the S_L is calculated using extrapolation methods, which inherently introduce relevant input errors into any combustion model.

Another way to get S_L is using tabulated values which was created based on the 30-species skeletal mechanism developed by Lu and Law [91].

The laminar flame thickness is calculated from the temperature profile along the normal direction of the flame front and also from the chemical time. The chemical time is calculated from the characteristic time of the laminar flame using Zeldovich Number which depends on the activation temperature of the fuel oxidation.

3.2.3. The Extended Coherent Flame Model(ECFM)

The Extended Coherent Flame Model (ECFM) [56,151] builds upon the basic principles of the Coherent Flame Model (CFM) by incorporating additional features to handle more complex combustion scenarios which focused on turbulence and flame interaction and includes detailed modeling of how turbulent eddies affect flame stretch and flamelet behavior. This involves correcting for the effects of turbulence on flame stretch, considering the turbulence intensity, and adjusting for the curvature and thermal expansion of flamelets [156]. Thus, the flame stretch is influenced by turbulence, as well as the ratios of turbulent to laminar flame velocities and lengths. It is adjusted for curvature and thermal expansion effects caused by laminar combustion in flamelets, based on the assumption of isotropic flame distribution [56].

The model calculates the rate of fuel consumption based on the flame surface density (FSD) and the reaction rate per unit flame surface area.

In the case of the coherent flame model, flame surface area per unit volume, defined as

$$\mathsf{\Sigma }=\frac{A_l}{V}$$

(2)

Using flamelet assumption, the mean turbulent reaction rate is computed as the product of the flame surface density Σ and the laminar burning velocity S_L via:

$$\overline{\rho \dot{r_{fu}}}=-\omega \mathrm{L}\Sigma ={\rho }_{\mathrm{fu},\mathrm{fr}}{\mathrm{S}}_{\mathrm{L}}$$

(3)

In this ω_L as the mean laminar fuel consumption rate per unit surface along the flame front, ρ_fu,fr the partial fuel density of the fresh gas, ρ the density of the fresh gas and y_fu,fr is the fuel mass fraction in the fresh gas.

When combustion starts, several new terms have to be computed. Amongst them are source terms and two quantities in order to use equation (Eq. 1): Σ and S_L.

$$\frac{\partial \Sigma }{\partial t}+\nabla \bullet \left(\Sigma u\right)=S_{production}-S_{destruction}$$

(4)

The first term of left hand side is the time dependent component and second term is the convective transport of the FSD. The first term of right hand side is source term which represents the production of flame surface density comes essentially from the turbulent net flame stretch, the second term is sink term which represents the quenching effect referring to the local extinction or reduction of the flame surface density due to unfavorable conditions, such as excessive strain, heat loss, or insufficient reactants. Hence, the FSD transport equation incorporates these effects into a combined source term , S_∑ which includes both the creation and destruction mechanisms:

S_∑ = S_production -S_destruction = S_stretch + S_quenching + S_other

(5)

Here, S_other represents additional source terms that might be relevant depending on the specific combustion scenario.

This approach allows for detailed tracking of how turbulence affects the flame surface and, consequently, the combustion process.

Recently, the ECFM-3Z model has been extensively adopted for 3D CFD analysis of PCSI engines using AVL’s Fire CFD code [56]. The ECFM-3Z model is is an extension of the ECFM combustion model based on FSD transport equation and mixing model that can describe inhomogeneous turbulence premixed and diffustion combustion and operates within three distinct zones or regions: fuel, air, and the air-fuel mixture. In this model, the fuel can be represented as a mixture of various components. Both the burnt and unburnt gases are categorized into these three zones.

The extent of mixing among these zones is determined using a characteristic time scale, which is computed based on the k-zeta-f turbulence model [61,65,68,69]. The ECFM-3Z model assumes that the composition of unburnt gases, including air and Exhaust Gas Recirculation (EGR), remains consistent across both mixed and unmixed zones. The properties of the burnt gases are calculated based on the reaction progress variable.

In conclusion, this combustion model focuses on flame propagation and the interaction between turbulent flow and flame, while utilizing simplified global kinetics for chemical kinetics.

3.2.4. The Multizone Well-Stirred Reactor (MZ-WSR) Model

A homogeneous reactor-type combustion model, MZ-WSR, operates on the premise that sharp gradients in temperature and density are unlikely to occur within a cell and models combustion as an ignition-based phenomenon. This model divides the reactor into several well-mixed zones, each of which is assumed to be perfectly mixed with uniform composition. Additionally, chemical reactions are assumed to occur instantaneously within each zone. The governing equation of the MZ-WSR (Multi-Zone Well-Stirred Reactor) combustion model can be represented as follows:

$$\frac{\mathrm{d}{Y}_i}{\mathrm{d}\mathrm{t}}=\gamma \dot{{\omega }_i}+\frac{\dot{m}}{V}(Y_{i,in}-Y_i)$$

(6)

where, γ is the multiplier, Y_i the mass fraction of species i, $\dot{{\omega }_i}$ the reaction rate of species i, $\dot{m}$ the mass flow rate, and $Y_{i,in}$ the mass fraction of species i in the inflow.

The MZ-WSR model, which is coupled with detailed kinetic calculations, is particularly suitable for modeling the jet ignition process. Therefore, most of the studies [83,84,87,88,89,90,91,93,94,106] investigating the PCSI system using WSR models have considered detailed chemical mechanisms derived from GRI Mech 3.0 by Lu and Law [158] by utilizing the integrated chemistry solver known as SAGE [78,88,157], which is included in the commercial code CONVERGE [54]. However, the use of the SAGE model for combustion involves a considerable simplification. In the simulations, turbulence models are not applied to the mean chemical production terms in the governing equations. This means the potential impacts of turbulent fluctuations on these terms are not accounted for in the simulations. It should be noted that the influence of turbulence, as modeled using the RANS k-ɛ turbulence approach discussed earlier, is applied exclusively to the transport equations of mass, energy, and momentum in the averaged equations. Prior numerical investigation into the jet ignition process [157] have demonstrated that the WSR assumption can lead to cooler temperatures in cells where thin flamelets are present. This issue can become more pronounced in ultra-lean mixtures. Recently, several papers [90,93] have evaluated the prediction performance of this model by comparing it with the G-equation model for the combustion processes in an both passive [93] and active PCSI engine [90] fueled by natural gas. The results of these literatures showed that MZ-WSR model predicts faster combustion rates in the main chamber than the G-equation model fails to match the pre-chamber combustion phase in the right place. However, the prediction of the combustion rate in the main chamber by this model matched well with the experimental results. The MZ-WSR model requires that each cell be treated as an individual well-stirred reactor. For this to be accurate, the characteristic time of the turbulence in the cell must be significantly smaller than the characteristic time of the combustion chemistry. In other words, the Damköhler number must be much less than 1 to assume a well-stirred reaction. Therefore, in the case of PCSI combustion, where lean and highly heterogeneous turbulent flows exist, this model can result in significant errors.

Despite its limitations, this model has been extensively used for analyzing the combustion processes of PCSI engines within the RANS framework up to the present [26,27]. This is because the MZ-WSR model excels in providing a detailed and chemically accurate representation of combustion processes, and allows for more flexibility in adjusting the model to account for different fuels such as natural gas [27,28] and gasoline [26] and dual-fuel combustion conditions [] due to its detailed chemical kinetics. This makes it particularly useful for research and development where precise emission predictions and understanding of combustion chemistry are crucial. The experimental validation of this model was conducted using in-cylinder pressure and heat release profiles, and it demonstrated good predictive accuracy for both pressure and the combustion process within range of 1.6 < λ < 2. [78,88,157].

In recent, When the G-equation is utilized, the MZ-WSR model is integrated before, during, and after the flame front to calculate the intermediate and post-reaction species instead of using simplified global kinetics [91,93,104]. However, its effects have not been quantitatively proven.

3.3. Turbulence-Chemistry Interaction

The interaction between turbulence and chemistry-related quantities plays a fundamental role in determining combustion characteristics of PCSI applications, including ignition, flame propagation speed in both the prechamber and main chamber, and the burn rate trend. This is because, compared to traditional SI engines, the turbulent flow field in PCSI systems is highly inhomogeneous, has large spatial gradients at the jet boundaries, and exhibits rapid temporal evolution. Additionally, a high level of turbulent fluctuation near the spark plug, due to jet-to-jet interaction produced by the throttling effect during the scavenging process, is one of the key factors distinguishing the flame evolution characteristics of the PCSI system from conventional SI engines. To simulate these complex turbulence-chemistry interactions, previous literatures utilized the source terms of combustion models that were developed for pre-mixed or non-premixed conventional gasoline S.I engines.

As previously explained, the ECFM-3Z and G-equation combustion models are the most widely and frequently used for combustion analysis of I.C engines including PCSI system. These two models have different approaches to representing turbulence-chemistry interactions. Therefore, they calculate different laminar flame speeds, turbulence intensities, and spatiotemporal scales, resulting in different turbulent flame speeds. Turbulent flame speed significantly affects engine output, fuel consumption rate, and pollutant formation. Thus, the predictive accuracy of these two combustion models ultimately depends on the accurate calculation of turbulent flame speed 159,160]. The G-equation combustion model calculates the turbulent flame speed using explicit correlations, as shown in Equation (1), whereas the ECFM-3Z model derives it from the FSD equation. By examining the changes in velocity scales ratio($u^{,}/S_L)$ and length scale ratio($l_t/{\delta }_L)$ during the combustion process for an optically accessible GDI engine under stoichiometric conditions using the two combustion models(G-equation with three different turbulent flame speed correlations), as depicted in the Borghi-Peters diagram in Figure 8 below, clear differences can be observed as combustion progresses [161]. The differing results from the two combustion models can be attributed to their distinct approaches in handling the flame brush, particularly in terms of the spatial distribution of the reaction volume predicted by each model [161]. Specifically, the differing flame brush volumes and notable deviations between two models are due to different turbulence intensities and scales, namely mean turbulence-flame interaction.

The flamelet-based combustion models undergo significant differences in the combustion process in terms of turbulence/flame interaction. However, validation results compared to experimental ensemble-averaged in-cylinder and burn rate traces do not show substantial differences. The figure below compares pressure traces on the left and heat release rate curves on the right from experiments and various combustion models under the same engine and operating conditions [161]. From an engineering perspective, the results show overall good agreement with the experimental ensemble average traces for all models, exhibiting only minor deviations. Moreover, a noteworthy point from these results is that the results of the ECFM-3Z and the G-equation with Bradley correlation are almost identical. Consequently, it is concluded that comparing the predictive performance of various combustion models using ensemble-averaged pressure traces and burning rate curves is not appropriate. Thus, the overall burn rate only partially reflects the accuracy of the simulation framework concerning turbulent flame speed (S_T). This is because combustion models must consider a broad range of interactions between flow and chemical reactions [159,160,161].

The FSD transport equation indeed includes source and sink terms that account for both the stretching and quenching effects. These terms are crucial for accurately modeling the dynamic behavior of the flame front in both laminar and turbulent combustion conditions. By incorporating these effects, the FSD equation can provide a comprehensive description of the flame surface evolution, which is essential for predicting combustion performance and stability. As shown in equation (7), the stretching effect in the FSD transport equation accounts for the deformation of the flame surface due to the flow field. This can be caused by both laminar and turbulent strain rates. In turbulent flows, this term can be more complex and might include contributions from both large-scale and small-scale turbulent eddies that stretch the flame front. To simulate stretching and quenching conditions in term S_∑, of equation(7), the intermittent turbulent net flame stretches (ITNFS) model [159] was used. The ITNFS model concept involves using an advanced model to characterize the interaction between a single vortex and the flame front, and then extrapolating this to represent the complete turbulent flow. This model postulates that each turbulence scale influences the flame independently, implying no interaction between different turbulence scales. It is based on the idea that the overall effects of turbulent fluctuations can be predicted from the behavior of individual scales in the unburned gas mixture as shown in Figure 9.

By applying this model to a complete turbulent flow, it is assumed that the cumulative effect of all turbulent fluctuations can be inferred from the behavior of each individual scale. The limitation of this model is that it cannot account for vortex interactions. It is clear that this model has limitations when applied to PCSI engines, where strong turbulence intensity is distributed near the flame front of the hot turbulent jet and highly inhomogeneous turbulence exists. The production of flame surface density primarily results from the net flame stretch due to turbulence. This flame stretch is expressed as the large-scale characteristic strain (ε/k), adjusted by a function C_t, which considers the size of turbulence scales, and viscous and transient effects [160]. C_t depends on turbulence parameters and the properties of the laminar flame. Hence, the right-hand side of Equation (3) can be expressed as follows.

$$S_{production}=\alpha K_{eff}\Sigma and S_{destruction}=\beta \frac{{\rho }_{fu,fr}{S}_L}{{\rho }_{fu}}{\Sigma }^2$$

(7)

$$\mathrm{H}\mathrm{e}\mathrm{r}\mathrm{e},\mathrm{}{K}_{eff}=K_t=\frac{\epsilon }{k}{C}_t$$

(8)

K_t is a very important property since it influences the source term for the flame surface and therefore the mean turbulent reaction rate. The coefficients of α known as stretching factor and β in equation(7) are arbitrary tuning constants used in ECFM.

However, it is important to better understand the contexts in which the ITNFS model can be effectively applied and where additional considerations or alternative models may be necessary, particularly due to the lack of consideration for nonlinear effects in the interaction of turbulence scales with the flame front, flame-generated turbulence, and reignition of fresh gases crossing a locally quenched flame front [*].

*159 C. Meneveau and T. Poinsot, “ Stretching and quenching of flamelets in Premixed turbulent combustion”, Combustion and Flame, Vol.86, pp.311-332, 1991.

[**] 160 Meneveau, C. and Sreenivasan, K.R. “The Multifractal Nature of Turbulent Energy Dissipation.” Journal of Fluid Mechanics 224 (1991): 429-484.

Several previous studies [11,12,150,162,163] have reported the average turbulence/flame interaction experienced by engine flames in PCSI engines and recently developed highly downsized engines on the Borghi-Peters diagram, confirming the occurrence of multiple regimes during S.I. combustion. The results of these studies emphasize the necessity for combustion models to be predictive across all potential combustion regimes [150,160,163]. Namely, the high levels of turbulence generated by the pre-chamber (PC), combined with the reduction in laminar flame speed due to dilution, significantly impact the combustion regime of PC-initiated combustion systems. These interactions are illustrated in the schematic Borghi diagram in Figure 10, where the Karlovitz (Ka) and Damköhler (Da) numbers are used to compare relevant timescales of turbulent combustion to determine the combustion regime. Consequently, the PC-initiated jet ignition combustion regime shifts into the thin reaction zone, bringing it much closer to the stability limits.

The two figures below are based on the same passive PCSI engine with a port fuel injection system in a gasoline engine. Figure 10(a) shows the progression of the combustion region at 2000 rpm on the Borghi-Peters diagram for λ = 1 and λ = 2 cases [12]. Figure 10 (b) illustrates the change in flame structure during the combustion period at 4500 rpm by increasing the air-fuel ratio using air and EGR [11]. From these results, it is clear that lean combustion shifts the combustion characteristics to the thickened flame regime. As is well known, this occurs because the eddies become smaller than the flame thickness, allowing some eddies to penetrate the pre-heat zone of the flame. Another important piece of information evident from Figure 10 (b) is that EGR dilution exhibits higher sensitivity compared to air dilution [11,12].

In very recent, gasoline-fueled passive pre-chamber engine was numerically modelled using RNG k-ε turbulence model and MZ-WSR with GAGE combustion model and investigated combustion charateristics for λ=1.0 and 1.2 at 4000rpm [83]. Figure 11 shows the evolution of the turbulent regimes in pre-chamber and main chamber for two operating conditions on the Borghi-Peters diagram. Noteworthy feature of this figure is that almost entire pre-chamber combustion, even in lean case, evolve in the thin reaction regime due to high level of turbulent turbulence intensity produced by strong jet to jet interaction during compression stroke. As shown in Figure 11(a), during the initial stages of combustion, the velocity scale ratio decreases due to the weakening of the initial turbulence intensity. Subsequently, it recovers as the residual gas decreases. The length scale ratio shows that the integral length scale(l_l) remains constant, and the laminar flame thickness(δ_L) gradually decreases and then stabilizes. Under lean conditions, the laminar flame thickness increases and the laminar flame speed decreases, causing the curve to shift upward. From these results, it can be deduced that the MZ-WSR model is suitable for pre-chamber combustion. Figure 11(b) shows that compared to pre-chamber combustion, main chamber combustion occupies a wider area on the Borghi–Peters diagram, starting from the border of the broken reaction zone, passing through the thin reaction zone, and moving into the wrinkled flame zone. This occurs because the high turbulence intensity is distributed across the flame front when the hot turbulent jet is ejected. As the intensity of the jet and turbulence decreases, the combustion quickly transitions through the thin reaction zone and linearly moves into the corrugated reaction zone. Therefore, combustion in the main chamber experiences three combustion regimes due to the rapid changes in velocity scale ratio and length scale ratio, indicating the presence of complex turbulence-flame interactions. This figure suggests that for PCSI engines, a multi-mode combustion model that can cover a wide range of the Borghi–Peters diagram needs to be developed. It is important to note that Figure 10 and Figure 11 are based on the RANS type k-epsilon turbulence model using Boussinesq’s isotropic eddy viscosity assumptions and flamelet-based combustion models, which may introduce some level of error.

3.2.5. The well-Tuned Versions of Combustion Models

As previously explained, to accurately predict the combustion in a PCSI engine, it is essential to develop a model capable of calculating multi-mode combustion. Currently, utilizing a flamelet-based combustion model is a practical alternative. Therefore, recent research has focused on tuning the model coefficients b₁, b₃ in equation (1), which represent the large and small-scale turbulence enhancements in the G-equation model described by equation (1) [91,93,104]. This tuning process is performed ad hoc.

In the case of combustion models using Flame Surface Density (FSD) or MZ-WSR, efforts have been made to adjust the flame stretch factor, α in equation (7) for FSD and multiplier, γ in equation (6) ad hoc to align with experimentally obtained in-cylinder pressure traces and heat release rate profiles [164]. However, these treatments could not be adequate to achieve a precise correlation between simulations and experimental results in both the pre-chamber and the main combustion chamber [93,104,164]. Such studies, as described above, attempt to address the limitations of describing the turbulence-flame interaction using RANS-type turbulence models and flamelet-based combustion models by calibrating the model coefficients included in the equations that represent combustion speed or burn rate.

Kim et al. [93] tuned the multiplier of the MZ-WSR model within the range of 1.0 to 1.4. Additionally, for the G-equation model, they adjusted the coefficients b₁ and b₃ between 1.0 and 2.5. for two air-fuel ratios at 1200rpm. The engine used in this study is a 1.86-liter passive PCSI PFI engine fueled by natural gas. They validated the improvement in prediction accuracy by comparing the in-cylinder pressure traces and heat release rate profiles with experimental data. Figure 5 presents a comparison between the cylinder pressure and heat release rate obtained from experimental data and the simulation results using the tuned models. The experimental data includes results from 300 cycles (represented by light gray lines) and the averaged pressure (depicted by a black bold line). Figure 12 and 13 compares the simulation results without tuning the model coefficients to the corresponding experimental results. From the results, it is evident that without tuning the model coefficients, it is impossible to accurately predict the combustion characteristics of the PCSI engine. Additionally, the discrepancies are more pronounced in the lean region, likely due to the inability to account for the transition to the thin reaction zone as the flame thickness increases and the laminar burning velocity decreases.

The tuning of model coefficients significantly improves the predictive accuracy of the G-equation model. However, the MZ-WSR model fails to improve predictive accuracy with coefficient tuning, as it does not account for the effects of small-scale turbulence on the reaction rates.

Silver M et al. [104] conducted a detailed analysis on a 2.1-liter active pre-chamber engine using natural gas as fuel, focusing on the impact of turbulent jets ejected from the pre-chamber on the burn rate in the main chamber. This study utilized the G-equation combustion model and RANS-type turbulence model for two different orifice diameters. In this work, the turbulent flame speed equation was tuned with coefficients b₁ and b₃ set to 0.78 and 2.0, respectively. The figure below compares the pressure traces in the pre-chamber and main chamber at 1200 rpm for the two orifice diameters between the calculated and experimental results. As observed, there is a excellent degree of agreement between the experimental and predicted data.

In this computation, a highly sophisticated 0D/1D wave model [**] was employed to obtain initial conditions and temperatures at the walls, as well as temperature and pressure values at the inflow and outflow boundaries for the 3D-CFD analysis. However, this model requires input conditions such as lift curves of intake and exhaust valves, pressure pulsations in the intake and exhaust pipes, and fuel lines. Additionally, experimentally obtained pressure traces are needed to analyze combustion phenomena in the pre-chamber and main chamber.

Therefore, without preceding experiments, it is challenging to perform precise 3D CFD analysis of the PCSI engine. Moreover, the tuning of model constants in the combustion model based on experimental values is necessary, which varies depending on the engine displacement, fuel used, pre-chamber geometry, and fuel supply system.

To the best of the authors’ knowledge, there have been only a few studies utilizing well-calibrated combustion models [46,62]. Although the physical phenomena of TJI combustion are not yet fully understood, model constants in the correlation or source term for various combustion models adjusted to correct the laminar flame speed and align with experimental pressure traces [41,48,49]. However, the ability to accurately predict all phenomena involved in TJI combustion remains uncertain [46].

Figure 15 shows the recent share of (a)combustion modelsused and fuels for PCSI simulations. The majority, 45.9% studies, relied on the MZ-WSR model, while 32.4% studies opted for G-equation model. Subsequently, the frequency of usage for ECFM(3Z) model accounted for 16.2%, while less frequently used combustion models, such as the Weller model, were grouped as ‘etc.’ Additionally, the fuels used in the CFD simulation for the PCSI engine were categorized and shown in Figure 15(b). As illustrated in the figure, 50% of PCSI engine studies utilized natural gas fuel, followed by gasoline fuel. Recently, research has emerged on using carbon-free fuels such as ammonia and hydrogen in PCSI engines. These studies employ a dual fueling strategy, utilizing natural gas and diesel fuel for initial stage of combustion, and have been grouped under dual fuel.

3.4. Numerical Grid of Pre-Chamber Engine

The process of subdividing a physical domain into smaller subdomains is known as grid generation or spatial discretization which divides the spatial domain into a mesh or grid. Each grid point, cell, or control volume represents a discrete location in the domain. In this context, meshing plays a crucial role in defining the discrete elements where physical laws are applied. These meshes establish the cells or elements over which flow calculations are conducted. The primary goal of the simulation is to generate numerical values for key variables, such as velocity, pressure, and temperature, at specific points within the mesh. As the grid size decreases, the solution to the discretized problem converges to the solution of the continuous problem. Thus, the grid significantly influences both the rate of convergence and the accuracy of the solution. Therefore, the grid generation step is a critical component of the simulation process [165] and generating a sufficiently fine and uniform grid for complex geometries is both challenging and time-intensive. It is difficult to ensure that the grid adequately captures all relevant scales of turbulence while keeping the computational requirements manageable.

3.3.1. Mesh Generation for LES Turbulence Model

In LES turbulent model, the mesh is limiting the size of eddy that we can physically represent on the mesh. LES resolves scales from the domain size L down to the filter size Δ, requiring substantial resolution of high wave number turbulent fluctuations. This necessitates either high-order numerical schemes or fine grid resolution if low-order numerical schemes are employed. Δ is the length scale or grid filter width and determined as follows: The integral length scale,$l_0$, size of eddy is directly related to the computational cell volume, V as follows:

$$\mathsf{\Delta }={\left(\mathrm{V}\right)}^{1/3}=l_0/5$$

(9)

where, the integral length scale,$l_0$, size of eddy is directly related to the computational cell volume, V

It is known that good LES should resolve at least 80% of the turbulent kinetic energy which is dependent on mesh size. As mentioned above, the remaining kinetic energy is modelled in the sub-grid model. Finally, the finer mesh, the more turbulent kinetic energy is resolved. As shown in Figure 5, the sub-grid scale viscosity is also a function of the mesh size, hence we are actually solving a set of different equations on each level of refined mesh. Consequently, a traditional mesh independence study, which aims to ensure that simulation results are not influenced by mesh density, becomes impractical. Instead, alternative validation techniques must be employed to ensure the reliability of LES results.

Accordingly, the LES turbulence model requires representing the integral length scale of each eddy with at least five grid cells, which demands an enormous number of grids. Additionally, there is the difficulty of meticulously ensuring that the conditions of equation (9) are satisfied.

3.3.2. Mesh Generation for URANS Turbulence Model

It is well known fact that using excessively large grid sizes in RANS simulations can lead to under-resolved mean fields. This discrepancy creates a significant difference between the actual RANS field and the computationally resolved field. Such under-resolution has a substantial impact on the accuracy of chemical simulations, often more so than the direct effects of turbulent fluctuations on the chemistry, which turbulent combustion interaction models are typically designed to account for. Consequently, in scenarios of poor resolution, typical turbulent combustion models must compensate for under-resolution to maintain accuracy, even though they were not originally intended for this purpose. Therefore, increasing the resolution in under-resolved areas of the simulation can minimize the errors that the unresolved field has on the chemistry. Once the resolution-related errors are mitigated, it can be observed that the remaining effects of turbulent interactions on combustion chemistry are of the same order of magnitude as the accuracy of the detailed kinetic mechanism being used. Additionally, it should be noted that in the case of an under-resolved flame front, the second derivative in the species conservation equation (representing diffusion) and the second derivative in the energy equation (representing conduction) would be under-estimated. This under-estimation reduces the mixing rates and consequently results in a lower calculated laminar flame speed.

All commercial CFD software provides grid generation tools in their pre-processors. When using the Converge CFD code, the Adaptive Mesh Refinement (AMR) algorithm is employed to increase the grid refinement level only in regions characterized by high velocities and large temperature gradients (e.g., flame fronts), without excessively slowing down the simulation with a uniformly refined grid. Namely, The mesh was re-generated at each iteration and dynamically refined. This approach resulted in the total number of cells in the computational domain ranging from 300,000~700,000 during the compression stroke to ranging from 1.23~3 [83] million to 9 million [78]at the end of the combustion phase. Studies using the Converge CFD code have set the base grid size to 1 [93]~4 [13,78]mm, and through a series of localized mesh refinements, achieved a minimum cell size of 0.125 [83]~0.5 [88,93]mm. By varying the grid size in the pre-chamber, the y+ value is monitored to ensure that the simulation results fall either within the recommended range (30 < y+ < 100) or within a range (20 < y+ < 30) that can be handled by a wall function [93]. Additional fixed mesh refinement is applied locally around the spark plug gap to accurately capture the ignition phase [48,83,88,90,93]. Additionally, both the pre-chamber area, including the nozzles, and the regions traversed by turbulent jets entering the main chamber are locally refined to a resolution of 0.125 [48,88]~0.25mm [78,93].

3.3.3. Differencing Process

A differencing scheme is a specific method used to approximate the derivatives in the discretized equations. These schemes define how the values at the discrete points in the grid are used to approximate the derivatives of the flow variables. There are several differencing schemes for convection and diffusion terms and frequently adopted one for PCSI engine CFD simulation is upwind differencing scheme [166] which uses values from upstream points to approximate the derivatives. It introduces numerical false diffusion [166] but is more stable for convection-dominated problems. There are several variations, such as first-order upwind and higher-order upwind schemes. Although higher-order differencing schemes such as the QUICK(Quadratic Upwind Interpolation for Convective Kinematics) scheme TVD(Total Variation Diminishing) scheme [167] can achieve high-accuracy solutions, they often compromise solution stability. As a result, for CFD analysis of PCSI (Pre-Chamber Spark Ignition) engines, which requires the simultaneous resolution of complex multi-mode combustion coupled with strong turbulence and conjugate heat transfer, most studies do not employ differencing schemes with accuracy higher than second order. Finally, it is always big challenge for CFD engineers to choose appropriate differencing scheme that provides a better balance between accuracy and stability. In the realm of automotive and mechanical engineering, the sensitivity of discretization schemes for convection terms in Large Eddy Simulation (LES) has been thoroughly explored, as documented [168].

The main function of a Sub-Grid Scale (SGS) model in Large Eddy Simulation (LES) is to dissipate resolved turbulent fluctuations effectively. The SGS model is specifically designed to provide the appropriate level of dissipation. Therefore, in LES, we should employ central differencing schemes, as they do not introduce any numerical dissipation due to their symmetric treatment of fluxes, ensuring that no artificial viscosity is introduced into the system. Central differencing schemes are second-order accurate, providing a good balance between accuracy and computational cost. However, they can be less stable than upwind schemes.

Central Differencing Scheme: Uses the average of values at surrounding points to approximate derivatives. It is second-order accurate but may introduce numerical dispersion and can be unstable in convection-dominated flows. Opting for the central difference scheme is generally favored to mitigate numerical inaccuracies; nonetheless, a common practice involves blending the outcomes of upwind schemes with those of the central scheme to enhance computational stability.

3.4. Time Discretization

To obtain a numerical solution, partial differential equations must be discretized in both spatial and temporal domains. Time-dependent variables, also known as transient terms, are mathematically represented as derivatives with respect to time. However, from a physical standpoint, these terms require special handling. Transient terms describe the variation of a specific variable over time within an infinitesimally small control volume, adhering to conservation principles while preserving generality. When addressing the discrete treatment of transient factors, it is generally preferable to seek a time-varying solution, as this approach directly influences the accuracy of the numerical results. For flows that are predominantly steady, it is advisable to first determine a time-dependent approximation and subsequently transition to a steady-state approximation. Selecting the appropriate time step is crucial in simulating turbulent jet ignition and spark ignition combustion in in PCSI combustion (IC) engines. Accurately capturing pressure fluctuations due to deflagration within the cylinder or other phenomena in two divided chamber IC engine, depends heavily on the time step resolution used in the simulation. A variable time-step approach, utilizing the Courant-Friedrichs-Lewy (CFL) number, has been implemented in almost all previous literature []. The maximum CFL number differed across various regions of the domain, spanning from 1 to 5, and the time step was automatically adjusted within a range of 0.001 to 1 CAD [49] or 2.5e-05 to 1e-08 sec [78].

3.5. Appropriate Initial and Boundary Conditions

The most straightforward method for obtaining initial conditions (such as intake and exhaust temperatures, pressures, injection quantity, and ignition timing) and boundary conditions are to utilize experimental data. However, this approach cannot be employed until a prototype engine is constructed, so it is only applied in a few specialized research studies [Novella]. In addition, one of the primary difficulties in 3D CFD simulating PCSI combustion systems is the scarce availability of boundary and initial conditions within not only main chamber but also the pre-chamber. Typically, installing measurement instruments in the confined space of a pre-chamber is extremely difficult. Therefore, in most previous studies related to 3D CFD simulation, initial values of all required thermodynamic parameters (including composition) were directly transferred from the 0-D/1-D model(s). This section describes the 0D/1D simulations that provide the necessary initial and boundary conditions for the 3D CFD simulation of the PCSI engine, based on the findings of previous studies.

3.5.1. Inflow and Outflow Boundary Conditions of Intake and Exhaust

In the CFD analysis of a PCSI engine using RANS models, the most influential boundary conditions are the inflow and outflow conditions. To achieve this, it is necessary to obtain time-dependent variations in velocity, pressure, temperature, and chemical composition at the cross-sections of the engine’s intake and exhaust pipes. The value of CFD simulations lies in their ability to provide precise foundational design data before the prototype is created. Therefore, an additional simulation is required to obtain the aforementioned boundary conditions before the prototype engine is built. Therefore, time-varying pressure, temperature and chemical compositions profiles, derived from a 0D/1D model built using commercially available 1-D gas dynamic codes such as GT-SUITE [169], GT-POWER [71,72,73,74,75,76,77,78,79,80], and WAVE [61], validated with experimental data, were implemented as the inflow and outflow boundary conditions. These same 1D-CFD models were utilized to specify the chemical composition.

These one-dimensional gas dynamic commercial programs thermodynamically model the components comprising the engine and consider the intake and exhaust pipes, which connect each component, as one-dimensional reacting compressible. They utilize the Method of Characteristics (MOC) [170] to calculate the pulsating flow characteristics that occur during the intake, compression, combustion, and exhaust processes. This approach allows for fast and accurate calculations of these characteristics at the system level of the entire engine. The gas dynamic simulation consists of intake, exhaust system and cylinder units. Geometrical data, flow losses, and heat transfer between the pipe wall and the gas are considered. The pre-chamber combustion engine model incorporates two different thermodynamic descriptions of the main chamber and pre-chamber operating cycle, integrated with the gas exchange processes through not only the intake and exhaust valves but also orifices in pre-chamber. The main chamber, pre-chamber and the intake and exhaust pipes are treated as a control volume, ensuring that energy and mass balances are maintained during the intake and exhaust phases, which utilize two distinct systems, each containing two thermodynamic zones, to independently model the combustion processes within the pre-chamber and main chamber. These systems are interconnected through nozzles, facilitating the exchange of mass and enthalpy between them. During the compression, combustion, and expansion phases of both chambers, the energy balance of the combustion, heat release rate is obtained from either experimental information [11] or two-zone model like approach [171,172] which requires the 3D-derived turbulent intensity profiles in both chambers because preliminary 3D CFD analyses or experiments should be carried out in motored conditions to tune constants of turbulence sub-model included in this model [8,9]. Another straightforward way of modeling the combustion in pre-chamber is using the Wiebe function [11,170] which were obtained from heat release rate from experimental pressure trace [11].

In very recent, very delicate 1D model based on GT-Powe is used to give further refined initial and boundary conditions of pressure and temperature to 3D CFD model [71,72,73,74,75,76,77,78,79,80,81,82,83]. However, this 1 D model requires intake and exhaust valve profile curves, intake, exhaust, fuel lines, pre-chamber and main chamber pressure traces. Combustions in pre- and main chamber are expressed by heat release based calibration(HRBC) where the heat release rates are computed from the measured pressure traces. Additionally, the discharge coefficient through the pre-chamber nozzle is calibrated to accurately replicate the pre-chamber pressure rise observed in the experimental data [104,173].

Figure 18 illustrates the single-cylinder pre-chamber engine modeled using GT-Power. As shown in the figure, the pre-chamber is modeled as a constant volume chamber, a check valve, and a nozzle. To obtain more refined 1D model, calibration must be performed using experimentally obtained pressure traces from both the pre-chamber and main chamber to get accurate flow coefficients of the nozzle and check valve [104,173].

Applying appropriate boundary conditions for Large Eddy Simulation (LES) is complex, especially since inflow conditions must accurately represent turbulent fluctuations. This is critical due to the essential requirement of supplying a dynamic velocity field that fluctuates across both temporal and spatial dimensions. However, studies using LES [64,97,105] do not provide detailed information regarding the boundary conditions.

3.5.2. Wall Heat Transfer Modeling

Given the challenges in experimentally measuring instantaneous gas-to-wall heat fluxes, 3D-CFD simulations of in-cylinder processes have become indispensable. These simulations are crucial for assessing not only the overall heat transfer to the combustion chamber walls and its spatial distribution but also the heat losses through the boundaries of the computational domain. For precise simulation of the uneven temperature profiles in both the pre-chamber and main chamber, and to effectively capture their effects on the turbulent jet quenching and stretching phenomena, a conjugate heat transfer model is indispensable. While numerous studies have emphasized the importance and necessity of developing models suitable for PCSI engines, to the best of my knowledge, no related research results have been reported to date.

In the 3D CFD analysis of IC engines, the heat flux to the walls significantly impacts the temperature inside the combustion chamber, thus greatly influencing the accuracy of predicting emissions such as NOx [179]. Therefore, numerous heat transfer models and wall heat flux correlations exist in the literature [174,175,176,177,178], many of which were developed and validated through experiments conducted in research laboratories over the past decades, typically under low-load and low-speed engine conditions. Recently, due to engine downsizing, operating conditions in the low-speed high-load range have become increasingly important, resulting in a significant increase in the thermal loads on engine components facing the combustion chamber. Recent study [179] have revealed that existing wall heat transfer models(Angelberger’s [177] and the Han and Reitz’s ones [178]) tend to overestimate wall heat transfer, as evidenced by experimental engine thermal surveys and temperature measurements conducted on four currently produced engines. Consequently, existing heat transfer models require new calibration. Recent studies have proposed alternative heat transfer models for the optimization of wall heat transfer in such highly downsized spark-ignition engines [179] and for engines using carbon-free fuels like hydrogen [180].

As can be inferred from previous research findings, all these commonly used wall heat transfer models have been found to be sufficient under certain conditions. However, their adequacy for PCSI engines has yet to be evaluated and remains uncertain. For PCSI engines, the validation and calibration of wall heat transfer models are crucial. This is because the conjugate heat transfer between the pre-chamber surface and the gas greatly affects flame quenching, thereby significantly impacting the prediction of flame speed. The importance of developing wall heat transfer models suitable for PCSI engine modeling has been highlighted in previous research. Chinnathambi et al. [78] underlined that a specialized wall heat transfer treatment will be necessary for the pre-chamber wall and the nozzles to accurately address the uneven temperature distribution and the quenching effects.

The wall heat transfer model by O’Rourke and Amsden [58] and Angelberger [80,177] are the most frequently adopted wall heat transfer model for PCSI engines. These model is one of the most popular isothermal and non-isothermal wall heat transfer models. The model by O’Rourke and Amsden assumes a constant near-wall flow temperature and density, and a constant Prandtl number within the boundary layer. On the other hand, Angelberger’s model allows the near-wall flow temperature and density to vary within the modelling range [181] Surface temperatures were based on simplified predictive FEM model from simulations using calibrated 0-D/1-D models [67].

Furthermore, utilizing the lumped model based on a thermal resistor network [170] with the experiments can yield the wall temperatures for the piston, liner and cylinder head within the 3D CFD model [29,182].

In thermal resistor network model, the film coefficient, necessary to calculate heat loss between the gas and the wall, is obtained empirical correlations such as Woschni’s or Annand’s correlations [170,183]. Recent research [184] suggested the use of a modified Woschni’s correlation for pre-chamber heat transfer analysis to precisely simulate the characteristics of PCSI combustion. It is important to note that the two heat transfer empirical formulas require calibration of constants included the correlations to accurately describe the conjugate heat transfer between the cylinder wall and in-cylinder gases, based on factors such as fuel type, combustion method, bore, stroke, and turbocharging type [183,185]. However, a specific heat transfer model for pre-chamber walls is currently unavailable [184], and calibrating these empirical formulas for PSCI engines is time-consuming. Consequently, the uncertainties related to this aspect cannot be entirely eliminated at present. Figure * illustrates a schematic diagram of the single-cylinder pre-chamber system modeled based on a 0D/1D approach [64].

4. Experimental Validation

Experimental validation plays a vital role in CFD modeling as it determines the accuracy of a model in representing real-world scenarios for specific engineering applications. Conducting validation involves comparing CFD simulation outcomes with empirical data to ensure precision in engineering analyses. For PCSI engines, most numerical studies have used the in-cylinder pressure traces and heat release rate curves in the main chamber as validation data. As mentioned in Figure 19, even though different combustion models and turbulent flame speed relations show different post-ignition combustion paths on the Borghi-Peters diagram, the pressure traces appear similar. This trend is also observed in the heat release rate curves. Because the HRR is derived from the pressure signal in the main chamber. It’s important to recognize that the experimental HRR provides a simplified estimation of the energy released during combustion. Therefore, for validation purposes, only the ignition onset, the peak combustion rate, and the duration of combustion should be rigorously considered.

Additionally, as seen in Figure 12 and Figure 13, in the case of the G-equation, it was found that the model constants related to large and small-scale turbulence enhancement, which represent flame propagation speed, can significantly impact the pressure and heat release rate behavior inside the cylinder of a PCSI engine. However, this does not necessarily guarantee an accurate representation of turbulence-chemistry interaction. Therefore, before fine-tuning the physical models used in predictive PCSI engine modeling, it is crucial to identify and understand the major sources of uncertainties. Many researchers have pointed out various uncertainties, including the lambda value inside the pre-chamber, excessively rich air-fuel ratios [186], and the correct calculation of laminar flame speed under ultra-lean conditions [91].

Syrovatka et al. [67] performed multi-cycle simulations on a single-cylinder active PCSI engine fueled by natural gas using the high-fidelity turbulence model LES and the ECFM combustion model. The results were compared to experimentally obtained heat release rate curves, as shown in the figure below. The operating conditions were set at a speed of 1800 rpm with an air excess ratio of 1.05. To analyze the cycle-by-cycle variation (CCV) effect, heat release rate curves for three cycles, excluding the first cycle calculated from the CFD model, were displayed for two engines with different pre-chamber volumes. For experimental validation, the maximum and minimum heat release rate curves obtained from the experimentally measured pressure traces over 120 cycles were shown. Figure 19 indicates that the LES modeling approach effectively captures the cycle-by-cycle variation of the PCSI engine, with all CFD traces falling within the limits of the experimental data. The figure also highlights the significant CCV in heat release rate curves due to lean combustion in the main chamber.

From these results, it is evident that analyzing model uncertainty by validating pressure and heat release rate curves averaged over 300 cycles using RANS-type turbulence models is very challenging. These findings suggest that the development of high-fidelity turbulence models capable of accounting for vortex-to-vortex interaction and non-uniform distribution of vortices, along with combustion models that can simulate multi-mode combustion, will be pivotal in establishing CFD analysis technology as a vital process in design engineering.

5. Conclusions and Future Research Directions

CFD holds the potential to serve as more than just a predictive tool; it can evolve into an advanced design instrument. To fully realize this potential, it is essential to enhance CFD with precisely tailored, physics-based numerical models. While CFD is extensively used in industrial and academic research, as well as in optimizing pre-chamber designs, current commercial CFD software still exhibits significant limitations in its predictive capabilities. These limitations are most pronounced in areas such as turbulence modeling, which affects mixing and combustion, wall heat transfer modeling, and the accurate depiction of combustion processes. In the specific context of pre-chamber spark ignition, there are numerous unresolved challenges related to flame and jet-wall interactions, ignition and combustion dynamics within the main chamber, and the ability of simplified combustion models in commercial CFD tools to accurately represent these complex processes.

This paper provides a comprehensive review of how CFD technology has been applied and developed over the past 20 years despite above mentioned challenges. It also thoroughly analyzes the models used for 3D CFD simulations of the emerging PCSI technology and their interactions. In particular, it offers an in-depth review of the limitations and errors associated with these complex physics-based numerical models, and explores how CFD analysis technology must evolve to transition from an academic research tool to an essential design tool. Additionally, this paper also discusses the difficulties in obtaining initial and boundary conditions, which is one reason why CFD currently fails to serve as a useful tool in the initial design phase, and explains technically why relying on 0D/1D analysis cannot fulfill the role of a predictive tool.

Ultimately, for high-fidelity CFD modeling of pre-chamber combustion engines, the primary sources of uncertainty must be identified and corrected before proposing fine-tuning and improvements related to various turbulence and combustion models.