Christelle Kouamé

Gérard Loiseau
email: gerard.loiseau@supagro.fr

Joël Grabulos

Renaud Boulanger

Christian Mestres

Development of a model for the alcoholic fermentation of cocoa beans by a Saccharomyces cerevisiae strain

Keywords: yeast, starter, chocolate, gamma concept, predictive microbiology

ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Introduction

Flavor is the essential attribute of the quality of chocolate. It is directly linked to that of the raw cocoa beans that depends essentially on the post-harvest processing of beans, namely fermentation and drying. Post-harvest processing is performed on-farm or in nearby small-scale processing units and is poorly controlled leading to variable and/or poor quality cocoa beans [START_REF] Aprotosoaie | Flavor chemistry of cocoa and cocoa products-an overview[END_REF]. Obtaining good quality marketable raw cocoa beans requires good agricultural practices (selection of plants, selection of cocoa pods/beans, maintenance of fermentation equipment, removal of the placenta) and good fermentation practices. Fermentation indeed plays a key role in developing chocolate flavor precursors [START_REF] Afoakwa | Flavor formation and character in cocoa and chocolate: a critical review[END_REF][START_REF] Beckett | Is the taste of British milk chocolate different?[END_REF][START_REF] Camu | Fermentation of cocoa beans: influence of microbial activities and polyphenol concentrations on the flavour of chocolate[END_REF][START_REF] Owusu | Influence of raw material and processing on aroma in chocolate[END_REF][START_REF] Pereira | Great intraspecies diversity of Pichia kudriavzevii in cocoa fermentation highlights the importance of yeast strain selection for flavor modulation of cocoa beans[END_REF][START_REF] Ravelomanana | The yeast flora of cocoa bean fermentation in the Ivory Coast[END_REF][START_REF] Schwan | The microbiology of cocoa fermentation and its role in chocolate quality[END_REF]. Before fermentation, cocoa beans are removed from the pod. They are surrounded by a white mucilaginous pulp; rich in pectin, citric acid and sugars [START_REF] Pettipher | Analysis of cocoa pulp and the formulation of a standardised artificial cocoa pulp medium[END_REF]. Generally, in Côte d'Ivoire, the beans are piled in heaps or placed in wooden boxes covered with banana leaves. The cocoa beans then begin to ferment. Fermentation occurs spontaneously and usually lasts from five to seven days. It has two main phases. The first is anaerobic and lasts an average of two days during which yeasts convert the sugars in the pulp into ethanol by alcoholic fermentation with an accompanying slight increase in the temperature that can reach 35 °C. At the same time, the pulp is hydrolyzed; this allows the entry of air into the mass, creating aerobic conditions maintained by periodic turning of the heaps of beans that favor the start of the second phase. During the second phase, oxygen is consumed by acetic bacteria to oxidize ethanol into acetic acid, which can be over-oxidized to carbon dioxide and water. As these reactions are exothermic, the mass temperature increases to reach or even exceed 50 °C. At the same time, ethanol and acetic acid diffuse in the beans leading to the death of the embryo but also to biochemical changes in the cotyledons such as the reduction of the sapid compounds responsible for bitterness and astringency, changing the color of the beans and the production of the precursors of the chocolate flavor [START_REF] Kadow | Identification of main fine flavour components in two genotypes of the cocoa tree (Theobroma cacao L.)[END_REF][START_REF] Trognitz | Diversity of cacao trees in Waslala, Nicaragua: Associations between genotype spectra, product quality and yield potential[END_REF]. The presence of a large lactic population during the fermentation of cocoa beans has been reported for many years. In addition, some authors identify a lactic phase, and lactic acid bacteria would play a relevant positive or negative role. Overall, if the yeasts S. cerevisiae and H. opuntiae, the lactic acid bacteria L. fermentum and the acetic acid bacteria A. pasteurianus lead the fermentation, good marketable beans are generally obtained [START_REF] Camu | Fermentation of cocoa beans: influence of microbial activities and polyphenol concentrations on the flavour of chocolate[END_REF][START_REF] Camu | Dynamics and biodiversity of populations of lactic acid bacteria and acetic acid bacteria involved in spontaneous heap fermentation of cocoa beans in Ghana[END_REF][START_REF] Daniel | Yeast diversity of Ghanaian cocoa bean heap fermentations[END_REF][START_REF] Nielsen | The microbiology of Ghanaian cocoa fermentations analysed using culture-dependent and cultureindependent methods[END_REF]Papalexandratou et al., 2011;Papalexandratou and De Vuyst, 2011). For these authors, lactic acid bacteria by consuming more fructose than glucose and citric acid to produce lactic acid, acetic acid, ethanol and mannitol would be essential for the proper progress of acetic fermentation and the aromatic quality of fermented beans. However, for [START_REF] Gálvez | Study on the microflora and biochemistry of cocoa fermentation in the Dominican Republic[END_REF] and [START_REF] Ho | The effect of lactic acid bacteria on cocoa bean fermentation[END_REF], lactic acid bacteria would compete for glucose with yeasts which would lead to poorly ethanolic fermentation of beans and to beans with rich lactic acid content; this excessive acidification can depreciate chocolate quality.

Inoculation with selected microorganisms has long been practiced to control the quality of fermented foods including wine, beer and bread [START_REF] Leroy | Lactic acid bacteria as functional starter cultures for the food fermentation industry[END_REF]. Recently several authors tested selected starters for the fermentation of cocoa beans. Inoculated fermentation at laboratory or pilot scale levels revealed the potential of this approach. The inoculation levels are generally 10 5 or 10 6 CFU/g for yeasts. These levels lead to regular fermentation and the development of desirable characteristics of the cocoa beans; higher inoculation levels would be prohibitive for large-scale commercial use. Inoculated fermentation could also reduce the fermentation period to between three and five days instead of five to seven days for natural fermentation [START_REF] Lefeber | On-farm implementation of a starter culture for improved cocoa bean fermentation and its influence on the flavour of chocolates produced thereof[END_REF][START_REF] Meersman | Breeding strategy to generate robust yeast starter cultures for cocoa pulp fermentations[END_REF][START_REF] Moreira | Impact of a microbial cocktail used as a starter culture on cocoa fermentation and chocolate flavor[END_REF][START_REF] Sandhya | Inoculum of the starter consortia and interactive metabolic process in enhancing quality of cocoa bean (Theobroma cacao) fermentation[END_REF][START_REF] Schwan | Cocoa fermentations conducted with a defined microbial cocktail inoculum[END_REF]. The quality of the inoculated beans was generally better than the quality obtained from spontaneous fermentation. [START_REF] Leal | Fermentation of cacao (Theobroma cacao L.) seeds with a hybrid Kluyveromyces marxianus strain improved product quality attributes[END_REF], during controlled fermentation of 45 kg cocoa beans inoculated with a strain of K. marxianus with high pectinolytic activity observed that this strain was present throughout the fermentations characterized by a faster and more intense pulp hydrolysis and a higher (+2°C) increase in temperature of beans mass. Inoculated cocoa beans showed a higher level of cotyledonary storage proteins degradation, a lower level of titratable acidity and a better cut-test score, indicating a higher proportion of well fermented beans. Sensorial analysis of the chocolate produced from cocoa beans fermented with K. marxianus gives better results in terms of the flavor attributes and global acceptability. [START_REF] Crafack | Influencing cocoa flavour using Pichia kluyveri and Kluyveromyces marxianus in a defined mixed starter culture for cocoa fermentation[END_REF] used a mixed starter culture composed of P. kluyveri or K. marxianus associated with Lb. fermentum and A.

pasteurianus. During fermentations of 20 kg of cocoa beans inoculated with a mixed starter a faster consumption of citric acid was observed compared to spontaneous fermentation, whereas no major differences in the pulp concentrations of acetic acid, lactic acid and mannitol were observed between the inoculated and spontaneous fermentations. After sensory analysis, chocolate produced from inoculated beans was found to be sweeter, fruity and with a full cocoa flavor. However, tests are still required at industrial scale to define the required conditions for inoculated fermentation as the raw material is subject to great variability (the origin of the plants, selection of the cocoa pods/beans, their maturity, the delay between pod opening and fermentation) and the environmental conditions. Proper inoculation with different strains (at different levels and at different times) and fermentation thus remain difficult to achieve without understanding the trophic competition and the succession of microbial strains, and the impact of the changes in environmental conditions, especially temperature, on the development of the strains. A modeling approach could help take all these parameters into account and help control the fermentation process.

The use of tools and methods of predictive microbiology [START_REF] Augustin | Modélisation de la croissance microbienne et gestion de la sécurité sanitaire des aliments[END_REF] has enabled the development of mathematical models of microbial growth under different environmental conditions, including the fermentation of Gowé, by [START_REF] Mestres | Modeling cereal starch hydrolysis during simultaneous saccharification and lactic acid fermentation; case of a sorghum-based fermented beverage, gowé[END_REF] and [START_REF] Munanga | Modeling lactic fermentation of Gowé using Lactobacillus starter culture[END_REF].

Modeling makes it possible to run simulations to evaluate the ability of the strains to control fermentation, to define the relative rates of the different microbial species that are used in the starter, the conditions of inoculation, and finally to define the parameters to control the process including temperature, pH, and substrate concentration [START_REF] Arellano-Plaza | Unstructured kinetic model for tequila batch fermentation[END_REF][START_REF] Mestres | Modeling cereal starch hydrolysis during simultaneous saccharification and lactic acid fermentation; case of a sorghum-based fermented beverage, gowé[END_REF][START_REF] Munanga | Modeling lactic fermentation of Gowé using Lactobacillus starter culture[END_REF].

This study was carried out as part of a project to control the production of premium chocolate characterized by floral and fruity notes. To achieve this goal, it was planned using a mixed starter consisting of selected strains of yeast and acetic bacteria and to develop a global model of the cocoa beans. The expected global model should be able to describe the growth, the use of substrates and the production of metabolites under various environmental conditions (temperature, pH, O2, organic acids content and ethanol) of strains constituting the starter. The global model will result of assembling two sub-models, that of alcoholic fermentation, and that of acetic fermentation. This global model will be built progressively, step by step, by progressively assembling sub-models recognized as being robust without their integration in the global model modifying the architecture of this one. The model of alcoholic fermentation developed here accounts for the most important intrinsic (glucose, ethanol, and nitrogen contents, pH) and extrinsic (temperature, oxygen level) fermentation parameters and can predict the growth of the yeast strain under the increasing temperature observed during cocoa bean fermentation in real conditions. This model provides insights into the parameters that drive the relay of the yeast population to that of acetic bacteria observed during cocoa beans fermentation.

Material and methods

Microbial strain

The Saccharomyces cerevisiae LM strain was provided by Lallemand S.A.S. (Blagnac, France) as active dry yeast at a viable cell concentration of 10 10 CFU/g. The yeast was stored at 4 °C.

Before experiments, yeast was reactivated in sterile water at 35 °C for 30 minutes.

Batch culture

One liter Lambda Minifor fermenters equipped with a cold water condenser on air exit pipes.

(LAMBDA Instruments GmbH, Baar-Switzerland), were used. Agitation, temperature, pO2, airflow, and pH were monitored and controlled. 1 M hydroxide sodium and 0.5 M sulfuric acid was used to adjust the pH. Yeast growth was monitored using an infrared in-line turbidity sensor (Optek FC20-ASD19-, Elscolab, Arcueil, France).

Microplate assays

Culture tests were carried out with a Bioscreen C MBR (Labsystems, Helsinki, Finland), in sterile microplates under controlled stirring and temperature (30 °C). Turbidity was measured at wavelength 600 nm. The growth rate was determined according to [START_REF] Membré | Experimental protocols and strain variability of cardinal values (pH and aw) of bacteria using Bioscreen C: microbial and statistical aspects[END_REF].

Plate counting methods

Tenfold serial dilution of samples (0.1 mL) was carried out in 9 g/L NaCl water and plated on Sabouraud agar (Biokar-diagnostics, Beauvais, France). Plates were incubated at 30 °C for 48 h.

Counts are based on triplicate plating.

Chemical analyses

Glucose, ethanol and glycerol contents were measured by HPLC (Dionex Ultimate 3000, Thermo Electron SAS, Courtaboeuf, France) with separation on an Aminex HPX87H column (Biorad, Hemel Hempstead, UK) and detection with both refractive-index and ultraviolet detectors (210 nm). Elution was performed at 30 °C with 5 mM sulfuric acid at a flow rate of 0.6 mL/min.

Modeling

Primary model

The logistic model with delay and rupture [START_REF] Rosso | Modélisation et microbiologie prévisionnelle: Élaboration d'un nouvel outil pour l'agroalimentaire[END_REF] was used to describe microbial growth, according to Equation (1).

= 0 , ≤ μ 1 - , > (1)
where λ is the lag phase (h), μmax the maximum specific growth rate (h -1 ), Nt and Nmax (CFU/mL) the microbial count at time t and at the stationary growth phase.

Modeling the lag phase (λ)

In our conditions, i.e. direct inoculation with the Saccharomyces cerevisiae LM strain in the form of active yeast, the duration of the latency phase (λ) in hours, increases linearly with temperature (T, °C) according to empirical Equation ( 2).

λ = a * T + b
(2)

Secondary model

The gamma concept secondary model (Zwietering et al., 1993) was used to describe the impact of temperature, pH, and glucose, ethanol, nitrogen content on the growth rate (μmax), according to Equation (3). μ = !"T°C% * !"pH% * !"(Glucose0% * !"(Ethanol0% * !"(Nitrogen0%

with ! values ranging from 0 to 1.

Impact of temperature and pH on yeast growth

The values of gamma factor for temperature and pH, namely the cardinal temperature and pH model (CTPM), proposed by [START_REF] Rosso | Convenient model to describe the combined effects of temperature and pH on microbial growth[END_REF], were calculated according to equation (4)

9:;< = ">% = ? @ @ A @ @ B 0 , > ≤ > C= "DED FGH %"DED IJ % J KD LMN ED IJ % JOP QKD LMN ED IJ RKDED LMN REKD LMN ED RS"=ET%D LMN UD IJ E=DVW 0 , > ≥ > , > C= < > < > (4)
where 9:;< = values range between 0 and 1. X corresponds to temperature or pH and the indices min, opt, max, to the minimum, optimum and maximum values for temperature or pH. n is a shape parameter that equals 1 for pH and 2 for temperature.

Temperature test cultures were performed in a one-liter glass fermenter using yeast nitrogen base (YNB) without amino acids (Biokar-Diagnostics, Beauvais, France) supplemented with 20 g/L glucose at a controlled pH of 4. Temperatures from 10 °C to 45 °C were tested.

pH test cultures were carried out Bioscreen C microplates according to [START_REF] Membré | Experimental protocols and strain variability of cardinal values (pH and aw) of bacteria using Bioscreen C: microbial and statistical aspects[END_REF], in YNB without amino acids supplemented with 20 g/L glucose. The pH was adjusted from 3 to 8 (in progressive steps of 0.5) by adding sterile 1 M sodium hydroxide and 0.5 M sulfuric acid. The experimental temperature was set at 30 °C and all the tests were conducted in duplicate.

Effect of glucose

In agreement with the substrate inhibition model of [START_REF] Andrews | A mathematical model for the continuous culture of microorganisms utilizing inhibitory substrates[END_REF], the value of the gamma factor for glucose and was calculated according to equation (5).

!"(S0% = ([0 / "]^+ ([0 + ([0 _ /Ki% (5)

where γ([S]) is the gamma factor of glucose ranging between 0 and 1, [S] is the initial glucose concentration (g/L), Ks (g/L) is the saturation constant for glucose and Ki, the inhibition constant (g/L).

Glucose culture assays were performed in sterile Bioscreen C microplates at pH 4 and 30 °C, in YNB without amino acids, supplemented with concentrations of glucose ranging from 0 to 50 g/L.

Two tests were performed per glucose concentration.

Nitrogen content on yeast growth rate

In agreement with the Monod model [START_REF] Liu | A simple thermodynamic approach for derivation of a general Monod equation for microbial growth[END_REF], the value of the gamma factor for nitrogen was calculated according to equation ( 6).

!"(S0% = ([0 /"]a + ([0 % (6)

where γ([S]) is the gamma factor of nitrogen, ranging between 0 and 1, [S] is the initial nitrogen concentration (mg/L) and Kn (mg/L) is the saturation constant of the Monod model for nitrogen.

Nitrogen culture tests were carried out on sterile Bioscreen C microplates at pH 4 and 30 ° C, in YNB broth without amino acids and without ammonium sulfate added with 20 g / L of glucose and nitrogen at concentrations ranging from 0 to 20 mg / L by adding di ammonium sulfate salt.

Effect of ethanol content on yeast growth

The effect of ethanol on yeast growth was determined using equation ( 7) proposed by [START_REF] Coroller | Modelling the influence of single acid and mixture on bacterial growth[END_REF].

!"(SI0% = 1 - (cd0 efg (7)
where [SI] is the ethanol concentration in g/L and MIC is the minimum inhibitory concentration of ethanol.

Ethanol culture assays were performed in sterile Bioscreen C microplates at pH 4 and 30 °C, using YNB broth without amino acids supplemented with 20 g/L of glucose and with ethanol concentrations ranging from 0 to 75 g/L, in duplicate.

Kinetics of ethanol production

The kinetics of ethanol production was directly proportional to yeast growth by equation ( 8), using the model of Luedeking and Piret in which maintenance energy is disregarded [START_REF] Van Impe | Towards a novel class of predictive microbial growth models[END_REF].

h(ijkl0 hj = m "n op/ % * q (8)
where [EtOH] is the concentration of ethanol (g/L) and m "n op/ % is the yield of ethanol in relation to the yeast population (CFU/mL) and μ is the specific growth rate (h -1 ). After integration, this equation gives:

(r st0 = m "n op/ % * "u %

Where N0 (CFU/mL) is the initial the population

Glucose consumption

Similarly, glucose consumption is directly proportional to the yeast population and glucose content can be modeled by equation ( 10).

(Glucose0 = (Glucose0 u -m "vwxyoz{ ⁄ % * " -u % (10)
where (Glucose0 u is the intial glucose content and m "vwxyoz{/ % , is glucose consumption in relation to the yeast population.

Effect of oxygen

The effect of the oxygen concentration in the culture medium on the Saccharomyces cerevisiae LM strain growth was studied in one-liter Lambda glass fermenter containing YNB broth without amino acid supplemented with 20 g/L glucose. The pH was set at 4 and the temperature at 30 °C.

Three levels of oxygenation (0, 4 and 8 mg/L) were tested.

Validation of the yeast growth model

The model was validated using the one-liter Lambda glass fermenter with YNB broth without amino acid supplemented with 60 g/L glucose, at pH 4. A temperature profile (from 32 °C to 45 °C) was programmed to follow as faithfully as possible the same increases in temperature as those recorded during spontaneous semi-industrial fermentation of cocoa beans in Côte d'Ivoire (Fig. 1).

Two inoculation rates (10 4 or 10 6 CFU/mL) were tested in duplicate. Infrared absorbance was monitored throughout fermentation and samples were collected at hourly intervals to assess the metabolites.

Table 1: Parameters of the growth model in YNB medium taken into account and environmental conditions.

Parameters

• µopt • tlag (λ) • Tmin, Topt, Tmax • pHmin, pHopt, pHmax • MIC EtOH • KS glucose • Ygluc/CFU/mL • Yethanol/CFU/mL
Environmental conditions

• [Glucose] initial = 60 g/L • pH = 4 • PO2 = 0 mg

Data processing

The parameters of the models and their standard errors were adjusted with Origin® (OriginLab, Northampton, MA, United States) software using linear or non-linear curve fitting with the Levenberg and Marquardt iteration algorithm. Simulink® (MATLAB and Statistics Toolbox Release 2012b, The MathWorks, Inc., Natick, Massachusetts, United States.) software was used to build the global model of yeast growth. As shown in Figure 2, the CTPM secondary model fit the experimental data well (RMSE 0.04).

Results and discussion

Modeling growth kinetics on a synthetic medium

The optimal growth temperature of the Saccharomyces cerevisiae LM strain was 33.3 °C, with µopt of 0.46 h -1 (Fig. 2a). These values are in agreement with those in the literature. [START_REF] Walsh | Growth of Saccharomyces cerevisiae and Saccharomyces uvarum in a temperature gradient incubator[END_REF] The Saccharomyces cerevisiae LM strain appeared to be acid-tolerant with a pHopt value of 4.5 and an almost constant maximum growth rate of 0.32 to 0.38 h -1 with an increase in pH from 2.9 to 5.7. Absence of variation in the growth rate was also reported by for S. cerevisiae. [START_REF] Charoenchai | Effects of temperature, pH, and sugar concentration on the growth rates and cell biomass of wine yeasts[END_REF] showed that a pH between 3 and 4 did not significantly affect the growth rate and cell biomass of wine yeasts. At pH 6.5 and above, there was a marked decrease in the growth rate of the Saccharomyces cerevisiae LM strain and at pH 8, growth stopped (Fig. 2b). During the alcoholic fermentation of cocoa beans, the pH value is fairly stable, i.e. between 3.5 and 4, at the end of the fermentation but sometimes reaches 5.5 [START_REF] Gálvez | Study on the microflora and biochemistry of cocoa fermentation in the Dominican Republic[END_REF]. This suggests that pH has little effect on the growth of the Saccharomyces cerevisiae LM strain during fermentation.

The Ks and Ki values for glucose are respectively of 0.27 g/L and 209 g/L were determined by adjusting the Andrew model (Eq. 5) to the experimental data (Fig. 3a). These results indicate that the growth of the Saccharomyces cerevisiae LM strain is very slightly impacted by the glucose concentration, which has been reported for glucose concentrations above 200 g/L by many authors including [START_REF] Charoenchai | Effects of temperature, pH, and sugar concentration on the growth rates and cell biomass of wine yeasts[END_REF]. Considering the low level of inhibition of the growth of the Saccharomyces cerevisiae LM strain for initial glucose concentrations close to 60 g / L, this effect will not be taken into account in the model.

Considering the impact of the concentration of mineral nitrogen on the growth rate of Saccharomyces cerevisiae LM strain, Kn value for the Monod model (Eq. 6) was very low (fitted at 0.04 mg/L of free nitrogen, with a quite large uncertainty) (Fig. 3b). There are more free amino acids in pulp and their content increase from 1 to 5 mg/g upon fermentation [START_REF] Camu | Dynamics and biodiversity of populations of lactic acid bacteria and acetic acid bacteria involved in spontaneous heap fermentation of cocoa beans in Ghana[END_REF].

This shows that the amount of free nitrogen is not a limiting factor in the fermentation of cocoa beans and this parameter will consequently not be taken into account in the model.

Figure 3: Effect of the initial glucose concentration (a) and the nitrogen concentration (b) on the growth rate of the Saccharomyces cerevisiae LM strain. Experimental data (black dots) and fitted (dashed line).

Equation ( 7) satisfactorily describes the decrease in μmax with an increase in the ethanol concentration (Fig. 4).

Figure 3: Inhibitory effect of ethanol on growth of the Saccharomyces cerevisiae LM strain. Experimental data (dots), fitted data (dashed line).

A MICEthanol of 92.5 g/L was calculated for the Saccharomyces cerevisiae LM strain at 30 °C and pH 4. The Saccharomyces cerevisiae LM strain appears to be fairly resistant to ethanol. This result is consistent with the values reported in the literature by [START_REF] Rose | Recent research on industrially important strains of Saccharomyces cerevisiae, in: Biology and Activities of Yeast[END_REF] ranging from 55 g/L to 103 g/L, by [START_REF] Arroyo-López | Susceptibility and resistance to ethanol in Saccharomyces strains isolated from wild and fermentative environments[END_REF] values ranging from 95.6 to 141.4 g/L and by [START_REF] Antoce | Quantitative study of yeast growth in the presence of added ethanol and methanol using a calorimetric approach[END_REF] values ranging from 88.8 to 107.7 g/L. [START_REF] Samagaci | Growth capacity of yeasts potential starter strains under cocoa fermentation stress conditions in Ivory Coast[END_REF] isolated 18 out of 36 pectinolytic strains isolated from cocoa able to grow at 30 °C at an ethanol concentration of 79 g/L and 7 strains able to grow up to 94.8 g/L. The maximum concentration of ethanol is less than 20 g/L [START_REF] Camu | Fermentation of cocoa beans: influence of microbial activities and polyphenol concentrations on the flavour of chocolate[END_REF][START_REF] Camu | Dynamics and biodiversity of populations of lactic acid bacteria and acetic acid bacteria involved in spontaneous heap fermentation of cocoa beans in Ghana[END_REF][START_REF] Lefeber | On-farm implementation of a starter culture for improved cocoa bean fermentation and its influence on the flavour of chocolates produced thereof[END_REF][START_REF] Schwan | Cocoa fermentations conducted with a defined microbial cocktail inoculum[END_REF] during cocoa bean fermentation and is in any case much lower than the MICEthanol of the Saccharomyces cerevisiae LM strain; the ethanol concentration consequently has very little impact on the growth of this strain. Ethanol can also have an antagonistic effect of temperature as high ethanol concentrations decrease the optimum growth temperature of yeast strains and increase the lag phase [START_REF] Medawar | Yeast growth: lag phase modelling in alcoholic media[END_REF]. [START_REF] Sá-Correia | Temperature profiles of ethanol tolerance: Effects of ethanol on the minimum and the maximum temperatures for growth of the yeasts Saccharomyces cerevisiae and Kluyveromyces fragilis[END_REF] reported that with an ethanol concentration close to 11% (v/v) the maximum growth temperature of a S. cerevisiae strain drops from 42 °C to 27 °C. [START_REF] Samagaci | Growth capacity of yeasts potential starter strains under cocoa fermentation stress conditions in Ivory Coast[END_REF] reported that 36 yeast pectinolytic strains isolated from fermented cocoa beans showing both tolerance to alcohol (10-12% alcohol) and to temperature (up to 40 °C) failed to grow at only 35 °C with less than 8-10% alcohol. The maximum level of ethanol recorded during cacao fermentation being 5 to 10 times lower (less than 2%), the interaction between ethanol level and temperature can be neglected in the model.

The value of the maximum growth rate was between 0.41 and 0.43 h -1 (results not shown) for dissolved oxygen concentrations 0; 4 and 8 mg/L. This factor was thus not taken into account in the model.

Figure 5 shows the consumption of glucose and the production of ethanol by the Saccharomyces cerevisiae LM strain. The mass yield of ethanol from glucose (Yp/s) is 0.46 g/g, which is close to the theoretical value 0.51 g of ethanol per gram of glucose consumed. This theoretical yield is never reached as part of the glucose is used for the synthesis of biomass and secondary metabolites such as glycerol. In S. cerevisiae under oenological conditions, the observed yields are between 0.47 g/g and 0.48 g/g [START_REF] Lonvaud | Microbiologie du vin : bases fondamentales et applications[END_REF]. Similarly, [START_REF] Laopaiboon | Ethanol production from sweet sorghum juice in batch and fed-batch fermentations by Saccharomyces cerevisiae[END_REF] found values ranging to 0.42 g/g to 0.48 g/g by fermenting sweet sorghum juice.

3.1.Validation of the model on synthetic medium

The global alcoholic fermentation model was developed using the values of the parameters and environmental conditions listed in Table 1. During the beans fermentation, some amount of glucose is also consumed by lactic acid bacteria, this is not included in this model, only built to predict the behavior of the yeast strain. The Saccharomyces cerevisiae LM strain was cultured in yeast nitrogen base without amino acids, supplemented with 60 g/L glucose at pH 4, in conditions approaching those observed at the beginning of cocoa bean fermentations characterized by absence of sucrose, but an equal amount of glucose and fructose [START_REF] Ho | The effect of lactic acid bacteria on cocoa bean fermentation[END_REF]Papalexandratou et al., 2011). We only added glucose as the Saccharomyces cerevisiae

LM strain can grow indifferently on glucose or fructose (results not shown). The temperature was regulated to mimic the measured profile observed in a wooden box during spontaneous fermentation in Côte d'Ivoire (Fig. 1) and close to the temperature profile observed in Ghana and the Dominican Republic [START_REF] Camu | Dynamics and biodiversity of populations of lactic acid bacteria and acetic acid bacteria involved in spontaneous heap fermentation of cocoa beans in Ghana[END_REF][START_REF] Gálvez | Study on the microflora and biochemistry of cocoa fermentation in the Dominican Republic[END_REF][START_REF] Pettipher | Analysis of cocoa pulp and the formulation of a standardised artificial cocoa pulp medium[END_REF]. Two levels of inoculation, 10 4 and 10 6 CFU/mL, were tested to assess the influence of the initial yeast population on the course of the fermentation, given that this range is often used for the inoculation of cocoa beans [START_REF] Meersman | Breeding strategy to generate robust yeast starter cultures for cocoa pulp fermentations[END_REF][START_REF] Ramos | Impact of different cocoa hybrids (Theobroma cacao L.) and S. cerevisiae UFLA CA11 inoculation on microbial communities and volatile compounds of cocoa fermentation[END_REF]. In comparison, the initial populations in natural cocoa fermentations are generally between 10 3 and 10 7 CFU/mL [START_REF] Ardhana | The microbial ecology of cocoa bean fermentations in Indonesia[END_REF][START_REF] Batista | Dynamic behavior of Saccharomyces cerevisiae, Pichia kluyveri and Hanseniaspora uvarum during spontaneous and inoculated cocoa fermentations and their effect on sensory characteristics of chocolate[END_REF][START_REF] Camu | Dynamics and biodiversity of populations of lactic acid bacteria and acetic acid bacteria involved in spontaneous heap fermentation of cocoa beans in Ghana[END_REF][START_REF] Nielsen | The microbiology of Ghanaian cocoa fermentations analysed using culture-dependent and cultureindependent methods[END_REF]. It is indeed likely that pod ripeness and postharvest pod age (due to transport and temporary storage) determine the initial amounts of yeasts and LAB [START_REF] Camu | Dynamics and biodiversity of populations of lactic acid bacteria and acetic acid bacteria involved in spontaneous heap fermentation of cocoa beans in Ghana[END_REF].

Figure 6 compares the experimental population measured during duplicated fermentation with the population predicted by the model under the same conditions.

The model predicted a lag phase of nearly 7 h for both fermentations whereas it was longer, close to 10 h with an inoculation level of 10 4 CFU/mL and shorter, close to 3 h with an inoculation level of 10 6 CFU/mL. The lag phase is indeed always difficult to predict [START_REF] Augustin | Mathematical modelling of the growth rate and lag time for Listeria monocytogenes[END_REF][START_REF] Medawar | Yeast growth: lag phase modelling in alcoholic media[END_REF][START_REF] Munanga | Modeling lactic fermentation of Gowé using Lactobacillus starter culture[END_REF] and the model presented a RSE of 2.1 h for the lag phase. This difficulty in measuring and predicting the lag phase has several possible explanations, the heterogeneity of the population, the balance between the individual latency of each cell and the latency of the global population, or the phenomena of intermediate latency explained by [START_REF] Swinnen | Predictive modelling of the microbial lag phase: a review[END_REF]. Nevertheless, the differences between the predicted and experimental populations are less than one log unit, and one can thus consider that the model predicts the growth of the Saccharomyces cerevisiae LM strain fairly well. However, the model could be improved by a better determination of the lag phase.

With an inoculation level of 10 4 CFU / mL, the population increased for approximately 23 h and reached a stationary phase of 33 h, while the model predicts growth arrest at 33.6 h when the temperature reached 41 °C, i.e. the Tmax of the Saccharomyces cerevisiae LM strain and when the γtemperature became null (Fig. 6a). These results are consistent with those of [START_REF] Mensonides | The metabolic response of Saccharomyces cerevisiae to continuous heat stress[END_REF], who showed that increasing the culture temperature of a strain of S. cerevisiae from 28 °C to 43 °C resulted in a sharp decrease in the growth rate at 40 °C, a complete end to growth at 42 °C, and total loss of viability above 43 °C. On the other hand, with an inoculation level of 10 6 CFU/mL, the growth of the population stopped after 21 h according to the model and at 26 h in the experiment when the temperature was close to 37.5 °C, which corresponds to a γ-temperature of 0.7 (Fig. 6b). Consequently, another parameter was responsible for stopping growth.

Figure 7 shows the changes in glucose concentration predicted by the model and those observed during the growth of the yeast strain. With an inoculation level of 10 4 CFU/mL, the model predicted changes in the glucose concentration fairly well up to 33 h when the temperature reached 41 °C. While the model predicted glucose consumption would stop at that time, in the experiment, the consumption of glucose continued after 33 h and was only complete after 48 h (Fig. 7a).

With an inoculation level of 10 6 CFU/mL, the consumption of glucose by the Saccharomyces cerevisiae LM strain was well predicted by the model up to 18 h when the temperature approaches 36 °C. The model predicted total exhaustion of glucose at 20 h, while in the experiment; glucose was exhausted at 25 h. Despite this slight discrepancy, in these fermentation conditions, the exhaustion of glucose clearly stopped the growth of the Saccharomyces cerevisiae LM strain, both in the model (γ-glucose of 0) and in real conditions.

Concerning ethanol production, with an inoculation level of 10 4 CFU/mL the experimental and predicted data matched very well (Fig. 8a). With an inoculation level of 10 6 CFU/mL (Fig. 8b), the model predicted over-production of ethanol after 18 h (by 10 g/L), as it predicted over consumption of glucose. These discrepancies appear to be linked to slight overestimation of the growth and final population by the model (5.10 7 CFU/mL); the slight difference in the populations indeed represents an extra ethanol production of 9 g/L. Anyway, the final concentrations of ethanol in the validation tests were in the range of the peak ethanol production (between 10 and 24 g/L at 24-48 h) observed in natural or inoculated fermentations of cocoa [START_REF] Crafack | Influencing cocoa flavour using Pichia kluyveri and Kluyveromyces marxianus in a defined mixed starter culture for cocoa fermentation[END_REF].

Gamma glucose

The model predicts yeast growth, glucose consumption and ethanol production fairly well, except at the end of the exponential phase. This is particularly true with the lower inoculation level when the temperature reaches Tmax of the strain. This can be attributed to the fact that we did not take decoupling between the production of ethanol and the consumption of glucose into account in the model, which reduces the yield of biomass and ethanol. This phenomenon is widely reported in the literature as being the result of the adaptation of yeast to heat stress; as a consequence, the energy required by the cell for maintenance increases due to the production of heat shock proteins and metabolites such as trehalose and glycerol [START_REF] Aldiguier | Synergistic temperature and ethanol effect on Saccharomyces cerevisiae dynamic behaviour in ethanol bio-fuel production[END_REF][START_REF] Kukec | Influence of temperature shock on the glycerol production in cv[END_REF][START_REF] Mensonides | A kinetic model of catabolic adaptation and protein reprofiling in Saccharomyces cerevisiae during temperature shifts[END_REF][START_REF] Omori | Enhancement of glycerol production by brewing yeast (Saccharomyces cerevisiae) with heat shock treatment[END_REF]. The concentrations of glycerol measured after total glucose consumption were 4.15 g/L and 4.76 g/L respectively for inoculation levels of 10 4 and 10 6 CFU/mL, i.e. nearly one quarter of the total amount ethanol produced. The ratio of glycerol produced to glucose consumed (60 g/L) was between 0.069 g/g and 0.079 g/g. This ratio is higher than that calculated for non-botrytized wines in which between 4 and 9 g/L of glycerol is produced from an initial total sugar contents in grape juice of 200 g/L giving a ratio of between 0.02 g/g and 0.045 g/g [START_REF] Remize | Re-assessment of the influence of yeast strain and environmental factors on glycerol production in wine[END_REF]. To improve the accuracy of the model, it thus appears necessary to take into account the maintenance energy and the production of glycerol, particularly when the model temperature increases. This is all the more important as, according to [START_REF] Illeghems | Assessment of the contribution of cocoa-derived strains of Acetobacter ghanensis and Acetobacter senegalensis to the cocoa bean fermentation process through a genomic approach[END_REF], glycerol could allow the growth of A. pasteurianus strains often found during the acetic fermentation of the beans.

The model and the validation sets, performed in conditions that mimic cocoa fermentation, clearly demonstrate that the inoculation level and the speed of increase of temperature drive the yeast growth while other factors such as pH and ethanol, free nitrogen and oxygen levels do not have a significant impact on yeast growth. If the inoculation level is high and/or the increase in temperature is low, glucose is completely exhausted after 1 day, while with a low yeast inoculation level and/or a rapid increase in temperature, yeast stops growing before the glucose is exhausted.

Similar observations have been made for cocoa fermentation, i.e. complete exhaustion of glucose after high level of yeast inoculation and 1-2 day fermentation has been reported [START_REF] Crafack | Influencing cocoa flavour using Pichia kluyveri and Kluyveromyces marxianus in a defined mixed starter culture for cocoa fermentation[END_REF][START_REF] Lefeber | On-farm implementation of a starter culture for improved cocoa bean fermentation and its influence on the flavour of chocolates produced thereof[END_REF][START_REF] Sandhya | Inoculum of the starter consortia and interactive metabolic process in enhancing quality of cocoa bean (Theobroma cacao) fermentation[END_REF] while residual glucose has been reported after natural fermentation or with a low level of yeast inoculation. The former leads to high transient ethanol content during cocoa fermentation, but to acidic, astringent chocolate with a quite low hedonic rating [START_REF] Sandhya | Inoculum of the starter consortia and interactive metabolic process in enhancing quality of cocoa bean (Theobroma cacao) fermentation[END_REF] while the latter can lead to chocolate with a high hedonic rate.

In any case, the death of yeasts at the end of the anaerobic phase is largely documented (De Vuyst and Weckx, 2016) and attributed to depletion of fermentable sugar, an increase in temperature and in pH and/or ethanol inhibition [START_REF] Crafack | Influencing cocoa flavour using Pichia kluyveri and Kluyveromyces marxianus in a defined mixed starter culture for cocoa fermentation[END_REF][START_REF] De Vuyst | The cocoa bean fermentation process: from ecosystem analysis to starter culture development[END_REF]. The model developed here and the validation tests confirmed that the increase in the concentration of ethanol and the change in pH do not have a significant impact on yeast growth or death during cocoa fermentation, since most yeasts showed the same behavior as the Saccharomyces cerevisiae LM strain. The increase in temperature and glucose exhaustion thus clearly drive yeast growth and death during the anaerobic phase of cocoa fermentation. It should be noted that we did not observe yeast death as the yeast population was calculated from measured absorbance and this increases with biomass whether the biomass is alive or dead. As already mentioned above, maintenance metabolism and energy consumption during the stationary phase were not included in the model, The next step is developing, with a similar approach, a model acetic fermentation using a strain of acetic bacteria isolated from an industrial fermentation in Côte d'Ivoire and to combine both models to design a global model of the fermentation of the pulp of cocoa for predicting the evolution of substrates (glucose, fructose, ethanol) and metabolites (ethanol, acetic acid) that can have a direct impact of cocoa quality.

Figure 1 :

 1 Figure 1: Temperature programmed for validation tests during model validation (dashed curve). Temperatures recorded during spontaneous fermentation in an ivory wooden box Coast, (solid curve).

Figure 2 :

 2 Figure 2: Effect of the temperature at pH =4 (a) and of the pH (b) at 30 °C on the maximum growth rate of Saccharomyces cerevisiae LM strain.

Figure 4 :

 4 Figure 4: Relation between the concentration of glucose (•), ethanol (▲) and the yeast population size.

Figure 5 :

 5 Figure 5: Saccharomyces cerevisiae LM strain growth ( rep1 and • rep2); simulated growth of yeasts (---); temperature profile (-.-.-.-) and γ-temperature ( ) with inoculation level at 10 4 CFU/mL (a) and 10 6 CFU/mL (b).

Figure 6 :

 6 Figure 6: Changes in glucose concentration ( rep1 and ■ rep 2); simulated changes in glucose (---) and γ-glucose ( ) during fermentation with inoculation levels of 10 4 CFU/mL (a) and 10 6 CFU/mL (b).

Figure 7 :

 7 Figure 7: Changes in ethanol concentrations ( rep1 and ▲rep 2); simulated changes in ethanol (---) and γ-ethanol ( ) during fermentation with inoculation levels of 10 4 CFU/mL (a) and 10 6 CFU/mL (b).

  neither was yeast adaptation or the kinetics of death depending of the speed of the temperature increase. Taking these phenomena into account could improve the predictability of the model, particularly the consumption of glucose at the end of the growth phase when temperature is the limiting growth factor, for example, with low inoculation levels. This work is a first step for developing a global model of cocoa pulp fermentation that can help for controlling cocoa fermentation. The development of the model of alcoholic fermentation of the pulp surrounding cocoa beans by Saccharomyces cerevisiae LM strain was based on the use of predictive microbiology models and particularly by using the gamma concept. The next step is under development; a similar approach is being used to model acetic fermentation using a strain of acetic bacteria isolated from an industrial fermentation in Côte d'Ivoire. Once validated, it will be combined with the alcoholic fermentation model to obtain the overall model of cocoa bean fermentation. It will thus be validated at a laboratory scale, with a pilot scale fermenter, and on an industrial scale in a fermentation center in Côte d'Ivoire. In addition, to take into account the consumption of the hexoses by LAB strains, an existing model of the lactic fermentation[START_REF] Mestres | Modeling cereal starch hydrolysis during simultaneous saccharification and lactic acid fermentation; case of a sorghum-based fermented beverage, gowé[END_REF][START_REF] Munanga | Modeling lactic fermentation of Gowé using Lactobacillus starter culture[END_REF]) can be integrated into the global model of the fermentation of cacao.ConclusionThis study is the first step of development of a global model of the cocoa beans fermentation of a project whose goal is to control the aromatic quality of chocolate by using a mixed microbial starter under controlled fermentation conditions. A model of alcoholic fermentation of the pulp surrounding cocoa beans by Saccharomyces cerevisiae LM strain, based on the use of the gamma concept, has been developed. Simulations and validations tests on a laboratory scale carried out with the model during fermentations under conditions mimicking those recorded during industrial fermentations in Côte d'Ivoire reveal the potential of modeling for the control of the anaerobic phase of cocoa fermentation. The model clearly show that the level of inoculation and the speed of the increase in temperature clearly drive yeast growth, while other factors including pH and ethanol, free nitrogen and oxygen levels have no significant impact on yeast growth.

  Table 2 lists the values of growth kinetics parameters for the Saccharomyces cerevisiae LM strain. The cardinal temperature values and pH were estimated by adjusting the cardinal model (Eq.

4) to the experimental data sets by minimizing the residual standard root mean square error (RMSE).

Table 2 : Secondary model parameter values ± standard errors for Saccharomyces cerevisiae LM strain Modeled Variable Equation Parameter Value

 2

			pHmin	2.20 ± 0.14
			pHmax	7.92 ± 0.19
			pHopt	4.54 ± 0.14
	µmax (h -1 )	(4)	µopt	0.46 ± 0.01
			Tmin (•C)	8.2 ± 2.5
			Tmax (•C)	41.8 ± 0.3
			T opt (•C)	33.3 ± 0.5
	λ(h)	(2)	a (h.°C -1 )	-0.61 ± 0.06
			b (h)	27.5 ± 1.8

  reported a Topt range of 30 °C to 35 °C for S. cerevisiae and[START_REF] Arroyo-López | Effects of temperature, pH and sugar concentration on the growth parameters of Saccharomyces cerevisiae, S. kudriavzevii and their interspecific hybrid[END_REF], 34.1 °C. According to[START_REF] Torija | Effects of fermentation temperature on the strain population of Saccharomyces cerevisiae[END_REF], the optimal growth temperatures for oenological yeast strains are between 30 °C and 35 °C. The Saccharomyces cerevisiae LM strain appeared slightly

thermotolerant, with a Tmax of 41.6 °C. In a comparative study of 27 strains of the genus Saccharomyces and non-Saccharomyces,

[START_REF] Salvadó | Quantifying the individual effects of ethanol and temperature on the fitness advantage of Saccharomyces cerevisiae[END_REF]

reported that S. cerevisiae was the

prediction of microbial spoilage in foods. J. Ind.Microbiol. Biotechnol. 12,[324][325][326][327][328][329]

Acknowledgments: This work was supported by the Frenchoc Premium project (20 th Fond Unique Interministériel) and its three actors: the French Center for Agricultural Research for International Development (CIRAD), the companies CEMOI and Lallemand S.A.S (Lallemand S.A.S -B.P. 59 -31702 Blagnac Cedex -France). We thank Lallemand for providing the Saccharomyces cerevisiae LM strain.

Author's contributions: Gérard Loiseau and Christian Mestres designed the protocols, reagents and facilities; Christelle Kouame and Joël Grabulos carried out the experiments; Christelle Kouamé, Christian Mestres and Gérard Loiseau analyzed the data. They wrote the article before submission and revised it according to the comments of the editor and reviewers.

Conflicts of interest:

Authors declare they have no conflicts of interest.


Amaury Gendron

Philippe Renaud
email: philippe.renaud@freescale.com

Marise Bafleur

Nicolas Nolhier

N Nohier

Analytical description of the injection ratio of self-biased bipolar transistors under the very high injection conditions of ESD events

come

Analytical description of the injection ratio of self-biased bipolar transistors under the very high injection conditions of ESD events

Introduction

Protecting high voltage smart power technologies against electrostatic discharges (ESD) is a great challenge [START_REF] Kirk | A theory of transistor cutoff frequency (fT) falloff at high current densities[END_REF] C. T. Kirk,"A theory of transistor cutoff frequency (fT) falloff at high current densities", IRE Trans. on electron devices, pp [164][165][166][167][168][169][170][171][172][173]1961. [START_REF] Markus | ESD Protection Considerations in Advanced High-Voltage Technologies for Automotive[END_REF]. Typical protection structures are based on self-biased NPN bipolar transistors given their good ESD robustness. The main drawback of the NPN transistor is its strong snapback behavior [START_REF] Mergens | Analysis and Compact Modeling of Lateral DMOS Power Devices under ESD Stress Conditions[END_REF] [4] [START_REF] Bertrand | Analysis and Compact Modeling of a Vertical Grounded-Base n-p-n Bipolar Transistor Used as ESD Protection in a Smart Power Technology[END_REF], that requires stacking several structures to increase the clamping voltage above the power supply value and then fulfill the requirement of latch-up free operation [START_REF] Ker | Double Snapback Characteristics in High-Voltage nMOSFETs and the Impact to On-Chip ESD Protection Design[END_REF].

To avoid stacking stacking structures, which is detrimental to both silicon area and on-state resistance, a solution consists in controlling this strong snapback effect (mostly design approach [START_REF] Reynders | Electrostatic Discharge Protection Device[END_REF]) or using reduced gain PNP bipolar transistors that do not exhibit this effect [START_REF] Gendron | Area-Efficient Reduced and No-Snapback PNP-based ESD Protection in Advanced Smart Power Technology[END_REF] [9] . In both cases, a deep understanding of the involved physical mechanisms under the very high injection conditions that occurs during an ESD event is required to improve the efficiency of such protection devices.

The analytical study of the physical mechanisms involved in a bipolar transistor is a topic that has been extensively studied; with the Boltzmann approximation as a corner stone of the analytical explanation of bipolar effect. However, under very high injection condition, this approximation does not remain valid and some corrective terms need to be added. Other derivations of semiconductor devices fundamental equations may be better suited to treat the case where the current densities create a space charge approaching the doping level of a considered region.

Historically, very high injection effects were treated for current densities up to 10 4 A.cm -2 , because, in most of applications, higher current densities values lead to silicon melting, hence device destruction. However, during the very short duration of an ESD event, protection devices can briefly sustain current densities as high as 10 6 A.cm -2 without destruction. Since the current peak lasts only for a few hundreds of nanoseconds, the heat generation in the protection device may not be sufficient to induce a melting of the silicon, and to damage the structure.

This paper describes the effect of such level of current densities and the impact of the resulting induced charges on the behavior of a self-biased bipolar transistor. In a self-biased configuration, base and emitter are connected or even the base can be left floating. The bias current is supplied by the avalanche current generated at the reverse biased base-collector junction (Figure 1). Starting from the same basic equations the cases of NPN and PNP will be treated separately.

Theoretical consideration for the determination of current gain at very high current densities

For the sake of simplicity, and in the aim of developing analytical models, we choose here to firstly consider a one-dimensional model of an NPN transistor in open-base configuration.

This latter choice is motivated by the fact that a 1D approach does not allow a realistic base biasing. We also consider constant doping levels and abrupt junctions. The recombinations will be neglected since, in most of the advanced technologies (in particular, CMOS and smart power technologies), the emitter and base widths are much smaller than the diffusion lengths of their respective minority carriers. Finally, we will consider isothermal conditions, that means we will not take into account the self-heating of ESD protections. The implications of these simplifying hypothesis will be discussed at the end of this paper, where we will qualitatively study their effects on the results' accuracy.

The numerical applications will be based on physical parameters' values at 300K. In particular for the electrons and holes mobilities, respectively µ n and µ p , and saturation velocities, respectively v sat,n and v sat,p , we will use the values:

1 6 , 1 , 1 1 2 1 1 2 . 10 6 . 8 . 7 . 10 . . 480 . . 1417 - - - - - - = = = = s cm v s cm v V s cm V s cm p sat n sat pB nB µ µ
Previous studies established that the current gain β of a bipolar transistor decreases for high current densities in the following way [START_REF] Webster | On the Variation of Junction-Transistor Current-Amplification Factor with Emitter Current[END_REF] [11] [START_REF] Rey | Some Aspects of Current Gain Variations in Bipolar Transistors[END_REF] :

H C B C J J J J + = = 1 0 β β (1)
Where J C and J E are respectively the collector and base current densities, β 0 is the well known maximal current gain expressed as:

B B pE iE E E nB iB W N D n W N D n 2 2 0 = β (2)
Where n i is the intrinsic carrier density, D n and D p the minority carrier diffusion lengths, N the doping concentration, W E and W B the emitter and base widths, respectively. Indexes E and B indicate whether the emitter or the base is considered. J H is the critical current density that can be expressed [START_REF] Rey | Some Aspects of Current Gain Variations in Bipolar Transistors[END_REF] with the base parameters only:

B B nB H W N qD J 2 = (3)
Where q is the electron charge.

The previous analytical expressions are only valid under the assumption of a high transistor gain. That is to say the emitter injection current J e-should be much greater than the base injection current J h+ . However, the current gain drastically decreases with the increase in current density and the previous assumption is no more fulfilled. This is why, we propose to re-derive the basic analytical expressions by taking into account both emitter and base injection currents. This way, we should be able to describe the mechanisms involved in bipolar transistors at very low current gain values, as it is experimentally encountered during ESD events.

On the Injection ratio J e-/ J h+

Starting from the formalism and expressions set by J. P. Bailbé [START_REF] Bailbé | Contribution à l'Etude Physique des Transistors Bipolaires[END_REF] it seems appropriate to treat the problem in terms of the total stored charges within the emitter and base. The charges are provided both by the extrinsic doping and by the current flow. This formalism leads to the following expression of the injection ratio:

2 2 2 2 2 2 2 2 iE nE SE iB nB SB iB nB B iE pE SE iB pB SB iE pE E BC E i n BC E i p h e n D Q n D Q n D Q n D Q n D Q n D Q dx n D p dx n D n J J + + + + = = = ∫ ∫ + - γ (4)
The integral boundaries are set by the emitter contact and the transition from the base quasi-neutral zone to the base-collector space charge region. Hence, taking into account all different contributions, Q E is the total ionized doping charge in the emitter, Q B the total ionized doping charge in the base, Q SB the total charge induced by the electron current within the base, Q SE the total charge induced by the hole current within the emitter.

The very general expression (4) is valid for any current value. In particular, its validity domain is not limited by any current dependent assumption, such as the Boltzmann approximation, for instance.

We remind that it is important to differentiate the values of both diffusion coefficients and intrinsic carrier densities depending on whether the emitter or base is considered, as these two parameters are greatly dependant on the doping concentration.

Evolution of the injection ratio under very high injection a. Preliminary considerations

Under very high injection, the charge Q B is negligible in comparison with the charge Q SB .

Thus the term depending on Q B can be systematically suppressed, and expression (4) becomes:

2 2 2 2 2 iE nE SE iB nB SB iE pE SE iB pB SB iE pE E n D Q n D Q n D Q n D Q n D Q + + + = γ (5)
None of the other terms can be neglected without any careful evaluation. First of all, for low injection ratio values, electron and hole currents are in the same range. Hence Q SB and Q SE are also in the same range, unless the difference between W B and W E favors the charge storage into either the emitter or the base. Besides, they should not be neglected in front of Q E , to preserve the expression validity even at very high current.

b. Generic transit time expressions

Let us now express the different types of charges as a function of the standard transistor parameters. The determination of the charges resulting from the ionized impurities within the emitter is trivial:

E E E W qN Q = (6)
The charges induced by the current into the base and the emitter are:

- = e B SB J Q τ (7) + = h E SE J Q τ (8)
Where τ B and τ E are the transit times, respectively of the holes within the base and of the electrons within the emitter. There is not general expression of the transit times that can be analytically derived, in function of the transistor parameters. However, we can use approximated expressions under some particular operating conditions [START_REF] Rittner | Extension of the Theory of the Junction Transistor[END_REF]:

For low injection condition:

€ τ B = W B 2 2D nB (9)
For high injection condition:

€ τ B = W B 2 4D nB ( 10
)
Equivalent expressions can be derived for the transit time of the electrons within the emitter,

τ E .
Expression [START_REF] Webster | On the Variation of Junction-Transistor Current-Amplification Factor with Emitter Current[END_REF] for the high injection is only valid under the high gain assumption, which is not the case during an ESD event. Thus, we need to derive a specific expression for the transit times that applies in our particular case.

c. Hole transit time into the base

Under very high injection conditions, the drift current into the base is no more negligible compared to the diffusion current. To take into account both the diffusion and drift currents, a method consists in defining an effective diffusion coefficient D nB,eff [START_REF] Webster | On the Variation of Junction-Transistor Current-Amplification Factor with Emitter Current[END_REF]. When the low gain assumption is fulfilled, diffusion and drift currents become equal, and the resulting effective diffusion coefficient is then two times the carrier diffusion coefficient (D nB,eff =2D nB ); as exemplified in expressions ( 9) and [START_REF] Webster | On the Variation of Junction-Transistor Current-Amplification Factor with Emitter Current[END_REF].

Using the same methodology as developed in [START_REF] Webster | On the Variation of Junction-Transistor Current-Amplification Factor with Emitter Current[END_REF], let us recalculate the expression of this effective diffusion coefficient without neglecting the base injection current J h+ .

Starting from the continuity equations:

      + = - dx dn D E n q J nB nB e µ (11)       - = + dx dp D E p q J pB pB h µ (12)
Where n and p are the electron and hole concentrations, gradients, E is the electric field within the base quasi-neutral zone. For high injection, the base neutrality implies n=p, thus [START_REF] Rittner | Extension of the Theory of the Junction Transistor[END_REF] and ( 12) can be combined as:

dx dn D q J J nB h pB nB e 2 + = + - µ µ (13)
Using the definition of the injection ratio (

€ γ = J e- J h +
), an expression of J e-proportional to the electron gradient can be obtained:

dx dn D q J pB nB nB e γ µ µ - = - 1 2 (14)
This expression allows defining the effective diffusion coefficient:

γ µ µ pB nB nB effectif nB D D - = 1 2 , (15)
Compared to the previous effective diffusion coefficient (equals to 2D nB ), this new expression evidences a further increase by a factor of

€ 1 1- µ nB µ pB γ
. This anticipates that drift current becomes the dominant conduction type with respect to diffusion current.

Using [START_REF] Reynders | Electrostatic Discharge Protection Device[END_REF], the electron transit time within the base becomes:

nB B pB nB B D W 4 1 2         - = γ µ µ τ (16)
And the base current induced charge:

nB e B pB nB SB D J W Q 4 1 2 -         - = γ µ µ (17)

d. Electron transit time into the emitter

The determination of the electron transit time within the emitter is not as straightforward.

Since the emitter doping level can be rather high, the high injection condition may not be satisfied over the full duration of the ESD pulse. In particular, the approximation n=p cannot be applied, in the same way as previously done for the base. Regarding the emitter, both low and high injection conditions can be encountered during an ESD event. As previously mentioned, a unique transit time analytical expression cannot be determined. Then, there is no way to have a precise analytical description of the effects resulting from the current induced charge into the emitter Q SE . However, most bipolar transistors used as ESD protections, exhibit a base that is significantly wider than the emitter, making possible to neglect Q SE compared to Q SB . Using this new condition on the emitter and base dimensions, an injection ratio expression without the Q SE terms can be derived from expression (5):

2 2 2 iB nB SB iB pB SB iE pE E n D Q n D Q n D Q + = γ (18)
This expression will be preferentially considered in the following.

Besides, we can take advantage of expression [START_REF] Bertrand | Analysis and Compact Modeling of a Vertical Grounded-Base n-p-n Bipolar Transistor Used as ESD Protection in a Smart Power Technology[END_REF] to give a rough evaluation of the impact of Q SE on the injection ratio. Two values of Q SE should be considered, one under low injection and high gain conditions with transit time expression [START_REF] Renaud | Efficient High Voltage No-Snapback and Low Ron ESD Protection Device for Smart Power Technologies[END_REF]:

+ = h pE E SE J D W Q 2 2 (19)
and one under high injection and low gain conditions with transit time expression (10) :

2 ' SE SE Q Q = (20)
As the purpose is only to determine the magnitude of the modification induced on the injection ratio, it is not worthwhile to make a more complex computation under high injection and low gain conditions.

e. Injection ratio expression

First, the current induced charge into the emitter Q SE is not taken into account. The computation is based on expression [START_REF] Amerasekera | Self-Heating Effects in Basic Semiconductor Structures[END_REF]. By inserting expressions ( 6) and ( 17) of Q E and Q SB , respectively, we obtain:

E E n D J W n D D J W n D W qN - -         -         - + = γ µ µ γ µ µ γ (21)
This expression can be re-written as:

pB nB e pB nB Cr J J µ µ γ µ µ γ +         - = - 1 2 (22)
Where the parameter J Cr , homogenous to a current density, is expressed as:

2 2 2 2 2 B iE pE E E iB nB Cr W n D W N n qD J = (23)
The factorization of γ in [START_REF] Ker | The Impact of Low-Holding-Voltage and the Design of Latch-up Free Power-Rail ESD Clamp Circuit for LCD Drivers Ics[END_REF] leads to the following second degree equation:

0 2 2 2 =         +         + - - pB nB e Cr pB nB J J µ µ γ µ µ γ (24)
The reduced determinant of the corresponding binomial is: The reduced determinant is positive; equation [START_REF] Gendron | Deep Trench NPN transistor for Low-RON ESD Protection of High-Voltage I/Os in Advanced Smart Power Technology[END_REF] has then two mathematical solutions. One solution increases with the electron current and has no physical meaning in our study, thus only the decreasing solution is valid:

0 2 2 2 ' >         + =         -         + = Δ - - - e
2 1 2 1 1 r r J J r J J e Cr e Cr -       + + + = - - γ (26) With: pB nB r µ µ = 1 (27)
Expression (26) shows that a modification of the behavior occurs when € J Cr J e- becomes close to the mobility ratio r 1 . Then, J Cr can be considered as a critical current density to which J e-has to be compared to characterize the injection ratio variations. To evaluate J Cr , expression [START_REF] Li | PMOSFET-based ESD Protection in 65 nm Bulk CMOS Technology for Improved External Latchup Robustness[END_REF] can be re-written as:

H Cr J J 0 γ = (28)
Where J H is the critical current density describing the beginning of the gain fall, given by (3), and γ 0 the same expression as the maximum injection ratio at low injection:

B B iE pE E E iB nB W N n D W N n D 2 2 0 = γ (29)
The trend of the injection ratio versus the normalized electron current (J e-/J Cr ) is plotted in Figure 2. At 300K, r 1 , the ratio of the electrons mobility on the holes one, is equal to 3.0.

Mobilities Ratio

f. Asymptotic behavior at extreme high injection

For J e->>J Cr , the injection ratio expression (26) has an asymptotic minimum limit equal to the mobility ratio:

1 r l = γ (30)
A simple manipulation of expression [START_REF] Bertrand | Analysis and Compact Modeling of a Vertical Grounded-Base n-p-n Bipolar Transistor Used as ESD Protection in a Smart Power Technology[END_REF], in which Q SE is neglected, allows analyzing the physical mechanisms associated with this asymptotic behavior:

1 2 2 2 2 2 r n D Q n D Q n D Q n D Q n D Q iB nB SB iE pE E iB nB SB iB pB SB iE pE E + = + = γ (31)
With this expression, it is obvious that the injection ratio reaches its asymptotic limit when the current induced charge within the base becomes greater than the doping charge within the emitter (Q SB >>Q E ). Considering [START_REF] Rittner | Extension of the Theory of the Junction Transistor[END_REF] and ( 12) continuity equations when n=p, an injection ratio limit equal to the mobility ratio could only be reached if, for both electron and hole currents within the base, the drift current due to the electric field becomes dominant compared to the diffusion current. This result has an important implication since, until today, it has been considered that under high injection conditions, the drift current reaches only the value of the diffusion current. To our knowledge, it has never been ever reported that for extreme high injection conditions, the drift current could be the main component of a bipolar transistor current.

We should outline that we do not know at the present stage of the study if this limit can be achieved. As J Cr depends on the current density, via the modification of the internal base width W E under high injection, we are not sure that the condition J e->>J Cr can be fulfilled.

Further investigations on that aspect will be performed in the following of the paper.

g. Evaluation of Q SE contribution

The purpose is to determine the corrective terms in the injection ratio expression (26), when Q SE is taken into account. As stated before, there is no analytical expression of Q SE valid over the full duration of the ESD stress and we will only make a rough evaluation of its contribution.

At first, the low injection within the emitter is considered. By inserting expressions ( 6), [START_REF] Jensen | Coupled Bipolar Transistors as Very Robust ESD Protection Devices for Automotive Applications[END_REF] and [START_REF] Arora | Electron and Hole Mobilities in Silicon as Function of Concentration and Temperature[END_REF] into the injection ratio expression (5), we obtain:

2 2 2 2 2 2 2 2 2 2 2 2 4 1 2 4 1 iE pE nE h E iB nB e B pB nB iE pE h E iB pB nB e B pB nB iE pE E E n D D J W n D J W n D J W n D D J W n D W qN + - + - +         - +         - + = γ µ µ γ µ µ γ (32)
As previously done when Q SE was neglected, similar mathematical manipulations lead to a new expression of the injection ratio:

( ) ( ) 2 2 1 2 1 1 1 1 E E e Cr E e Cr r r r r J J r r J J + -       - + + - + = - - γ (33)
Where r E1 and r E2 are two parameters without dimensions given by:

2 2 2 2 2 1 B iE pE nE E iB nB E W n D D W n D r = (34) 2 2 2 2 2 2 2 2 B iE pE E iB nB E W n D W n D r = (35)
This expression of the injection ratio still exhibits an asymptotic limit when J e->>J Cr :

( ) ( ) 2 2 1 2 1 1 1 1 E E E l r r r r r r + - - + - = γ (36)
To evaluate the effect of Q SE , we should compare r E1 to r 1 and r E2 to r 1 2 . This effect can be neglected if:

   << << 2 1 2 1 1 r r r r E E (37)
Besides, the coefficients r E1 and r E2 are proportional to

€ W E W B       2
, which confirms that Q SE is negligible for a base significantly wider than the emitter.

Regarding the case of high injection with low gain, the study shows that the terms r E1 and r E2 should be divided by two. The injection ratio is then given by:

2 2 2 2 2 1 2 1 1 1 1 E E e Cr E e Cr r r r r J J r r J J + -             - + +       - + = - - γ (38)
This reduction of the corrective terms is due to the decrease of the transit time at high injection, which leads to lower charge Q SE .

In this case, the asymptotic limit for J e->>J Cr is given by:

2 2 2 2 2 1 2 1 1 1 1 E E E r r r r r r + -       - +       - = γ (39)

Application to a self-biased bipolar transistor during an ESD stress. a. On-state resistance (R ON ) of an ESD protection

The on-state resistance, R ON , is one of the most critical parameters of an ESD protection; as it must be sufficiently low to prevent over voltage reaching oxide breakdown. For a protection based on a self-biased bipolar transistor, the R ON results mainly from high current effects at the base-collector junction [START_REF] Gendron | Area-Efficient Reduced and No-Snapback PNP-based ESD Protection in Advanced Smart Power Technology[END_REF]. The charge induced by the current adds up to the base doping, which is virtually higher, and subtracts to the collector doping, which becomes virtually lower [START_REF] Kirk | A theory of transistor cutoff frequency (fT) falloff at high current densities[END_REF]. This modulation of the space charge region (SCR) leads to a field modification resulting in an equivalent resistance. During an ESD stress, the current induced charge in the SCR could reach 10 18 cm -3 . In smart power technologies, the doping levels encountered are lower or in the range of this value, meaning that these high injection effects have a significant impact.

Calculation of the current induced charge in the SCR indirectly leads to an estimate of the SCR width variation and of the voltage drop in this region, if we consider:

-The maximum value of the electric field at the base-collector junction approximately constant with current.

-An abrupt collector-base junction with a collector doping level much higher than the base one to approximate a static junction position.

-The current induced charge within the base is much higher than the doping charge.

These conditions are not too restrictive, in particular for high voltage ESD protection. By this way, we can determine a first approximation of the R ON variations without developing an analytical expression for the SCR width.

b. Current induced charge in the base-collector SCR

Within the SCR, the carriers reach their saturation velocity. Thus, the electrons n and holes p concentrations in the region can be expressed as: The resulting total charge within the SCR becomes:

n sat e p sat h csr v J v J , , - + - = ρ (42)
Using J=J e-+J h+ , equation (42) can be expressed as a function of the total current J:

-         + - = e p sat n sat p sat scr J v v J v , , , 1 1 1 ρ (43)
An analytical expression of € J e-as a function of J can be extracted by using J=J e-+J h+ and

€ γ = J e- J h + : - - - = e e J J J γ (44)
By combining it with the expression (26) of the injection ratio, we obtain a second order equation in J e-:

( ) ( )( ) 0 2 2 1 2 2 1 = + + - - - - - - - e Cr e e e J J J r J J J J r (45)
whose solution is:

( ) ( ) ( ) 2 1 2 1 1 1 1 1 2 _ 1 r J J J r J J r r J Cr Cr Cr e + + + + + = - (46)
Inserting expression (46) into (43), we obtain an analytical expression of the SCR charge versus the total current density J:

( ) ( ) ( ) 2 1 2 1 1 1 , , , 1 1 2 _ 1 1 1 1 r J J J r J J r r v v J v Cr Cr Cr p sat n sat p sat scr + + + + +         + - = ρ (47)

c. Variations of the current induced charge in the base-collector SCR

To plot the evolution of ρ SCR , it is convenient to consider a value normalized on a critical charge ρ Cr defined as:

        + = p sat n sat Cr Cr v v J , , 1 1 ρ (48)
From (47), it is straightforward to obtain the expression of this normalized value of ρ SCR :

( ) ( ) ( ) 2 1 1 1 1 2 1 1 1 2 1 1 r J J r J J r r J J r Cr Cr Cr Cr scr + + + + - + - = ρ ρ (49) With p sat n sat n sat v v v r , , , 2 + = (50)
The trend of the ρ SCR normalized on ρ Cr versus the current density normalized on J Cr is plotted in Figure 3, both for a NPN and a PNP. The SCR charge in an NPN and a PNP is equivalent as far as the current density is smaller than the critical current J Cr . For higher values, the charge in an NPN monotonously increases, whereas the charge in a PNP decreases when a threshold current density J I is exceeded.

For the NPN transistor, the charge variation results in a continuous SCR narrowing, when the current increases. This SCR narrowing tends to reduce the R ON , which is beneficial for the performance of ESD protections. Concerning the injection ratio, the related internal base widening leads to a decrease of J Cr with current. If J Cr is lower or similar to the current density during an ESD stress, the condition J e->>J Cr is fulfilled and the value of the injection ratio tends to the mobility ratio.

In the case of a PNP, the charge decrease results in a SCR widening, which has a detrimental impact on the R ON . Besides, the normalized current density value keeps lower than seventeen (Figure 3), so that the condition J e->>J Cr may not be satisfied. This means the injection ratio value here, may not tend to the mobilities ratio. For further investigations, numerical applications are now undertaken., both for NPN and PNP transistors.

• NPN

For the numerical application, the parameters r 1 and r 2 are calculated for the values of mobilities and the saturation velocities at 300K defined in the section 2 :

€ r 1 = µ nB µ pB = 3.0 r 2 = v sat, n v sat, n + v sat, p = 0.56
The variations of ρ SCR can be estimated through the derivative of expression (47):

€ dρ scr dJ = 1 v sat, p - 1 1+ r 1 1 v sat, n + 1 v sat, p         r 1 + J Cr 2 1+ r 1 ( )⋅ J Cr J + J Cr 2         (51)
To evaluate the evolution with the current density, the limits of expression (51) are calculated for J<<J Cr and J>>J Cr :

n n Sat J J l v dJ d Cr , 1 , 1 =   →  << ρ (52) € dρ dJ J >>J Cr  →   1 1+ r 1 ( ) ⋅ r 1 1 v Sat, n - 1 v Sat, p         = l 2;n (53)
The ratio of these two limits is:

€ l 2, n l 1, n = 1 1+ r 1 ( ) ⋅ r 1 - v Sat, p v Sat, p         = 0.44 (54)
The slope of the charge increase is approximately reduced by two when the current density increases. However, this should not result in a strong degradation of the ESD protection robustness.

• PNP Analytic developments for the PNP are equivalent to those of the NPN by inverting the indexes of holes (p) and electrons (n).

In the case of a PNP transistor, the parameters r 1 and r 2 are given by :

€ 2 1+ r 1 ( )⋅ J J Cr + 1 = 1 1+ r 1 ( )⋅ r 2 -r 1 [ ] (55)
Equation ( 55) has only a solution if:

€ 1+ r 1 ( ) ⋅ r 2 -r 1 ≥ 0 (56)
Using the chosen parameters for the numerical application, it results that this condition is fulfilled:

€ 1+ r 1 ( ) ⋅ r 2 -r 1 = 0.26 ≥ 0
and the solution of equation ( 55) is given by:

€ J I = J Cr 1-1+ r 1 ( )⋅ r 2 -r 1 ( ) 2 2 1+ r 1 ( )⋅ 1+ r 1 ( )⋅ r 2 -r 1 [ ] 2 = 5.40 J Cr (57)
To characterize the PNP behavior, it is also important to determine the injection ratio at extremely high current. As we have to consider that the total charge present in the SCR remains positive:

0 , , > - = - + n sat e p sat h csr v J v J ρ (58)
a new condition on the injection ratio can be defined:

n sat p sat e h v v J J , , ≥ = - + γ (59)
For a PNP, the ratio of the saturation velocities (≈0.8) is greater than the ratio of mobilities (≈0.33), and hence is the limit of the injection ratio. This physical behavior results from the effect of the base narrowing that induces an increase of the carrier gradient into the base. As a consequence, the diffusion current remains significant compared to the drift current.

Therefore, the current induced charge within the base Q SB never becomes dominant compared to the doping charge within the emitter Q E . The maximum ratio of the drift current to the diffusion current can be determined from the following generic expression:

€ J Drift J Diffusion = 1 V T ⋅ pE dp dx (60)
Where V T is the thermodynamic potential. From the continuity equations ( 11) and ( 12) where, under high injection conditions n=p, we obtain:

€ J h + J e- = µ pB µ nB ( pE -V T dp dx ) ( pE + V T dp dx ) (61)
Given this injection ratio expression, the condition (59) becomes:

€ J Drift J Diffusion ≤ v sat, p µ nB v sat,n µ pB + 1 v sat, p µ nB v sat,n µ pB -1 = 2.4 (62)
This condition means that, in a PNP bipolar transistor, the drift current cannot exceed 70 % of the total current. By contrast, in an NPN one, the drift current can become negligible under very high injection.

Comparison with numerical simulation

To validate the proposed analytical description, we have compared it to 1D finite elements simulation results obtained with ISE ® simulation tool. The simulated devices are an NPN and a PNP with the same dimensions and doping types simply inverted to switch from NPN to PNP. As for the analytical approach, the doping levels are constant and the junctions abrupt.

Besides, carriers' recombination is neglected. The dimensions and doping level are chosen in the same range as the ones used for a high voltage ESD protection in smart power technology.

An emitter width of 0.25µm, a base width of 3µm and a collector width of 1µm are chosen.

The emitter, base and collector doping levels are 10 18 cm -3 , 10 16 cm -3 and 10 20 cm -3 , respectively. With these parameters, the maximum gain will be high, resulting in important dynamic fluctuations. To avoid any interference of very high injection effects with thermal effects, the simulations are carried out under isothermal conditions that properly reflect the thermal conditions under very fast ESD pulses.

a. NPN

The simulations are performed for current densities from 10 4 A.cm -2 , at which the high injection assumption is valid, up to 10 6 A.cm -2 , which is the maximum reached in an ESD protection. The effective base widths W B are extracted and the critical current densities J Cr [START_REF] Li | PMOSFET-based ESD Protection in 65 nm Bulk CMOS Technology for Improved External Latchup Robustness[END_REF] are calculated, for several current densities (Figure 4). As expected, they increase monotonously. Since the order of magnitude of the critical current density is 10 4 A.cm -2 , the asymptotic behavior should be reached before the end of the simulation at 10 6 A.cm -2 . The electron and hole current densities within the base are extracted to determine the injection ratio (Figure 5). At high current densities, its value tends to the mobility ratio, as predicted by (30). Then, the injection ratio is calculated with the expressions (26), in which Q SE is neglected, and (38), in which Q SE is taken into account considering high injection in the emitter (Figure 5). The agreement with the simulation is acceptable, in particular at high current densities where the assumptions made are completely fulfilled. Besides, the effect of Q SE results in a slight variation that can be neglected. For the proposed analytical description, the mobilities are considered constant in the quasineutral base region, even if they rapidly decrease in a very narrow region closed to the basecollector SCR. This approximation may explain the difference between the simulated and calculated injection ratios.

Mobilities Ratio Mobilities Ratio

b. PNP

The analytical expressions obtained for the NPN will be used as reference, keeping in mind that electron and hole indexes have to be inverted to obtain the expression for the PNP.

The simulations are performed from 10 3 A.cm -2 , from which the high injection assumptions become valid, and up to 2 10 5 A.cm -2 , above which the approximation on the electric field gradient does not remain valid. As for the NPN, the effective base widths W B and the critical current densities J Cr [START_REF] Li | PMOSFET-based ESD Protection in 65 nm Bulk CMOS Technology for Improved External Latchup Robustness[END_REF] are calculated for several current densities (Figure 6). As theoretically predicted, a base narrowing is observed at high current. The comparison of the electric field profiles before and after the occurrence of this base narrowing effect (Figure 7)

shows a gradient reduction, hence a reduction of the charge, resulting in a SCR extension increase. Using expression (57), we predict that this behavior occurs at a current density of As for the NPN, the injection ratio is directly calculated from the electron and hole current densities and theoretical expressions ( 26) and (38) (Figure 8). In accordance with (59), its value does not tend to the mobilities ratio, but to the saturation velocity ratio. This confirms the difference between the behaviors of NPN and PNP at extremely high current densities.

Above 4 10 3 A.cm -2 , the accuracy of the analytical expression is acceptable compared to the direct calculation of ratio between hole and electron current densities.

Saturation Velocities Ratio Saturation Velocities Ratio

Hypothesis Discussion

This study of the injection ratio under very high injection was possible by the virtue of some hypothesis. The most limiting ones are the 1D approach with floating base and the isothermal boundaries. For both of them, the discrepancy with the actual behavior of an ESD protection will be qualitatively evaluate. This complementary study is mandatory to determine how accurate are the results of the analytical derivation.

a. 1D description

As a constraint of the 1D description, the base has to be floating, although it is usually connected to the emitter in self-biased bipolar transistors. Nevertheless, at very high current, this hypothesis is not restrictive as the current flowing to the base contact become negligible compared to the current flowing through the emitter contact. More precisely, the ratio of the avalanche-generated carriers collected at the base contact on those assuring the reverse injection at emitter-base junction tends to zero [START_REF] Esmark | Device Simulation of ESD Protection Elements[END_REF]. In an ESD protection, the only purpose of the base contact is to provide a control on the trigger voltage V t1 [START_REF] Amerasekera | An Investigation of BiCMOS ESD Protection Circuits Elements and Applications in Submicron Technologies[END_REF] and the standby leakage current.

The 1D description assumes that the paths followed by the bias current and by the emitter injection current are the same. This hypothesis is valid for most of the self-biased ESD protections, in which the bias current is generated by the avalanche multiplication of the carriers injected by the emitter (this ensures the maximum of the impact ionization is located along the injection current path) [16] [17]. In particular, this hypothesis is valid for the widely-used vertical NPN transistors. In this type of component, the emitter, the base and the collector are respectively formed by a shallow high-doped N diffusion, a P-well and a Nburied layer (NBL) biased via a deep N-well (Figure 9). If the distance between the emitter and the NBL is shorter than the distance between the emitter end the deep N (i.e. the vertical base is shorter than the lateral base), the electrons injected by the emitter flow vertically under the emitter. Besides, the high avalanche region is located at the P-well-NBL junction so that, the bias current flow vertically from this junction to the emitter. It is important to notice the 1D hypothesis restricts the study to the case of the self-biased transistor. The results cannot be extended to the more general case where the bias current is directly provided at the base contact. For this latter case, the base contact is required, and the bias current does not follow the same path as the emitter injection current. The bias current comes from the base contact whereas the emitter injection current flows toward the basecollector junction.

b. Isothermal condition

Besides the very high injection conditions, one of the major specificities of an ESD protection is its strong self-heating. The main limitation in the reduction of the protection's size is the temperature, which could be high enough to lead to a second breakdown (between 700 K and 1250 K depending on the structure) [START_REF] Amerasekera | Self-Heating Effects in Basic Semiconductor Structures[END_REF]. This second breakdown induces current filaments where the temperature sharply increases up to the melting point of the silicon and the degradation of the protection [START_REF] Esmark | Device Simulation of ESD Protection Elements[END_REF]. The protection's behavior is influenced by these temperature variations through the modifications of the silicon physical parameters. For the injection ratio of a self-biased bipolar transistor, the most important parameters are the mobility, the saturation velocity (both involved at very high currents), and the impact ionization coefficient (as the bias current is provided by an avalanching junction).

• Mobility and saturation velocity

Both the mobility and the saturation velocity decrease when the temperature increases. At 600K, their values are typically equal to:

€ µ nB = 310 cm 2 .s -1 .V -1 µ pB = 120 cm 2 .s -1 .V -1
v sat, n = 5.85 10 6 cm.s -1 v sat, p = 5.85 10 6 cm.s -1

These variations does not modify the analytical expression of the injection ratio limits, (30) for an NPN and (59) for a PNP, as the ratio of the electrons mobility to the holes one (based on the model proposed by Arora and al [START_REF] Arora | Electron and Hole Mobilities in Silicon as Function of Concentration and Temperature[END_REF]) remains higher than the ratio of the electrons saturation velocity to the holes one (based on the model proposed by Canali and al. [START_REF] Canali | Electron and Hole Drift Velocity Measurement in Silicon and their Empirical Relation to Electric Field and Temperature[END_REF]) for the validity range of these models up to 600K (Figure 10). As a result, the physical effects described in this paper remain valid, in particular the base narrowing for a PNP. Only the values of the injection ratio limits are modified. At 600 K, these limits are equal to 2.6 for an NPN and 1.0 for a PNP. These values remain close to those at 300K, corresponding to variations of respectively 13% and 25%.

• Impact ionization coefficient

In an ESD protection, the avalanche phenomenon is strongly influenced by the temperature, as it is located in a region where both the current density and the electric field are high and as a result, the heat generation is also high. The amount of carriers generated by avalanche depends on the impact ionization coefficients, which decrease when the temperature increases, for both electrons and holes [START_REF] Reggiani | Experimental extraction of the electron impactionization coefficient at large operating temperatures[END_REF]. To maintain the avalanche current required to bias the bipolar transistor, the electric field maximum should increase. This effect tends to widen the base-collector SCR, and hence to reduce the effective base width, whatever the transistor type, NPN or PNP. The consequence on the I-V characteristic of the protection is a R ON degradation. To counteract this degradation and to obtain a low-R ON protection, it is mandatory to optimize the design of the structure so as to limit or delay the effects leading to the base narrowing.

In the case of an NPN, this optimization can be achieved by strongly decoupling the avalanche injection from the bipolar current gain mechanism that will tend to induce current filamentation and local heating [START_REF] Tremouilles | Design Guidelines to Achieve a Very High ESD Robustness in Self-Biased NPN[END_REF]. Such decoupling can be obtained by designing a very wide base bipolar transistor. Using such design guidelines, it was shown that introducing a deep trench into the base to widen it, allows meeting a very high robustness (>8kV HBM) together with a very low-R ON (>>1Ω) [START_REF] Gendron | Deep Trench NPN transistor for Low-RON ESD Protection of High-Voltage I/Os in Advanced Smart Power Technology[END_REF]. In the case of a PNP, since this base narrowing effect is much more detrimental and current gain cannot be degraded since already very low, the efficient design approach consists in tuning the base and collector doping profiles so as to enhance the space charge region narrowing under very high injection. The proposed analytical model was used to support the design of high voltage PNP-based ESD protection structures before running TCAD simulations. Similarly to NPN-based protections, these structures exhibited a very high ESD robustness (>10kV HBM) and a greatly improved R ON (4 times) that allowed significantly reducing the ESD protection area of high voltage I/Os in an advanced smart power technology [START_REF] Gendron | Area-Efficient Reduced and No-Snapback PNP-based ESD Protection in Advanced Smart Power Technology[END_REF].

Conclusion

An analytical description of the current injection ratio of a self-biased bipolar for very high current densities such as encountered during an ESD event have been successfully derived. It demonstrates that the decrease of the injection ratio under very high injection conditions has an asymptotic limit that is almost independent on the transistor doping level and dimensions.

It is also shown that this limit is equal to the ratio of the mobilities for the case of an NPN transistor and to the ratio of the saturation velocities for the PNP. This asymptotic behavior is reached when the induced current charge within the base becomes in the range of the emitter doping level (for instance a current density of 10 5 A.cm -2 in the quasi-neutral base region leads to a carrier concentration of 10 18 cm -3 ). At such high current density levels, the drift current within the base becomes greatly significant, which is not the case at standard low current levels. Moreover, the expression of the injection ratio allows determining the current induced charge within the base-collector charge space region and then anticipating the influence of such injection levels on the effective SCR and base widths. By this way, we have demonstrated that under very high injection conditions, in contrast to the NPN transistor, the PNP transistor is subjected to a base narrowing effect (inversion of the standard "base pushout"). This effect is reported for the first time. The results of this analytical approach were fully validated by a comparison with numerical simulation. Finally, the discussion of the most restrictive hypothesis, the 1D analysis and the isothermal conditions, has shown that our conclusions should remain valid for an actual ESD protection.

Besides, this analytical study allows explaining several specificities in the development of ESD protection. Despite the fall of the bipolar gain under high injection, the ESD designers do not optimize their protections to reduce this effect. Sometimes, they even use very low gain transistors. Such approaches are efficient because ESD protections are mainly based on self-biased bipolar transistors. In this configuration, the limit of the injection ratio at high current allows preserving the bipolar effect. This is a major difference compared to others domains of the electronics, such as applications for power management or high frequency. In these applications, the transistors are biased in a regular configuration (bias current provided

Figure 1 :

 1 Figure 1: Configuration of a self-biased NPN bipolar (a) and its equivalent circuit (b).

Figure 2 :

 2 Figure 2: Injection ratio versus normalized electron current drawn from (26).

Figure 3 :

 3 Figure 3: Normalized charge (ρ SCR /ρ Cr ) versus normalized current density (J/J Cr ), plotted for (49).

µ

  To characterize the PNP behavior, the most important parameter is the current density at which the charge begins to decrease. This optimum is reached for € dρ scr dJ = 0. Given expression (51), the solution of this equation should verify:

Figure 4 :

 4 Figure 4: Effective base width (W B ) and critical current density (J Cr ) versus the total current density, for an NPN.

Figure 5 :

 5 Figure 5: Injection ratio versus current density, obtained from simulation and analytical expressions (26) and (38).

4. 1 2 Figure 6 :Figure 7 :

 1267 Figure 6: Effective base width (W B ) and critical current density (J Cr ) versus the total current density, for a PNP. The calculated current density J I , where the base variation changes, is pointed out.

Figure 8 :

 8 Figure 8: Injection ratio versus current density, obtained directly for the simulation and from the analytical expressions (26) and (38).

Figure 9 :

 9 Figure 9: Cross-section of a vertical NPN transistor.

Figure 10 :

 10 Figure 10: Electrons to holes mobilities and saturation velocities ratio, and vice versa.

through the base contact), and the designers have often to work very hard to keep a high gain value under high injection. Another specificity in the field of ESD is the lack of interest in PNP transistors. The few publications about PNP transistors usually conclude on poor ESD capabilities due to a high R ON or a low robustness [22] [23]. These observations may be in relation with the modification of the base-collector SCR in a PNP, which is detrimental for the R ON , in particular if it leads to a base narrowing. The low robustness can be a direct consequence of a high R ON , as a higher voltage drop across the protection results in higher energy dissipation. To develop PNP-based ESD protection, it is mandatory to apply specific design rules in order to reduce the effects of the high injection on the base-collector SCR.

Based on this thorough analysis of very high injection effects, specific ESD design guidelines were withdrawn and successfully applied to the design of very efficient and robust PNP-based ESD protections.