Roméo Barnabé Bohounton

Oswald Yédjinnavênan Djihinto

Sedjro-Ludolphe Oronce

Dedome

Mahudro Réné

Laurette Yovo

Koffi Djossou

Aristide Koba

Pierre Adomou

Luc Salako Villeneuve

Fidèle Djogbénou

Paul

Rome ´o Barnabe ´bohounton

Oswald Ye ´djinnavênan Djihinto

Oronce Sedjro- Ludolphe

Re ´ne ´mahudro Yovo

Laurette Djossou

Koffi Koba

Aristide Adomou

Pierre Villeneuve

Salako Luc

Djogbe ´nou

Fidèle Paul Tchobo
email: fideletchobo@gmail.com

Euclasta condylotricha flowers essential oils: A new source of juvenile hormones and its larvicidal activity against Anopheles gambiae s.s. (Diptera: Culicidae)

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

Introduction

Essential oils (EOs), also called secondary metabolites, are volatile liquids with a strong odour obtained from various aromatic plant parts. EOs contain diverse bioactive compounds generated by plants to protect themselves against pathogens and insects [START_REF] Lee | Effects of Essential Oil from Hinoki Cypress, Chamaecyparis obtusa, on Physiology and Behavior of Flies[END_REF][START_REF] Isman | Essential Oil-Based Pesticides: New Insights from Old Chemistry[END_REF]. The bioactive compounds possess a wide range of biological activities such as toxicity and repellence to insects, antimicrobial, antioxidant, anticancer and antimalarial [START_REF] Muthiah | Botanical essential oils and uses as mosquitocides and repellents against dengue[END_REF][START_REF] Valdivieso-Ugarte | Antimicrobial, Antioxidant, and Immunomodulatory Properties of Essential Oils: A Systematic Review[END_REF][START_REF] Lesgards | Anticancer Activities of Essential Oils Constituents and Synergy with Conventional Therapies: A Review[END_REF][START_REF] Milhau | In Vitro Antimalarial Activity of Eight Essential Oils[END_REF]. Since EOs are mixture of different bioactive compounds, those with proven larvicidal and adulticidal properties were assumed to offer less chance of resistance development in mosquito vectors and could be an alternative to synthetic insecticides [START_REF] Sarma | Combinations of Plant Essential Oil Based Terpene Compounds as Larvicidal and Adulticidal Agent against Aedes aegypti (Diptera: Culicidae)[END_REF].

In Benin, vector-borne infectious diseases, such as malaria, acquired through the bite of an infected female arthropod, are still a major public health concern. Mass distribution of insecticide-treated nets (ITNs) and indoor residual spraying (IRS) with chemical insecticides are the main strategies implemented in malaria vector control programs [START_REF] Carnevale | Les anophèles: Biologie, transmission du Plasmodium et lutte antivectorielle[END_REF]9]. Unfortunately, vectors resistance to the different classes of insecticides currently used in public health (pyrethroids, carbamates, and organophosphates) could jeopardize the vector control interventions [START_REF] Soonwera | Adulticidal, larvicidal, pupicidal and oviposition deterrent activities of essential oil from Zanthoxylum limonella Alston (Rutaceae) against Aedes aegypti (L.) and Culex quinquefasciatus (Say)[END_REF]. Indeed, several studies have reported that natural mosquito populations in Benin showed resistance to pyrethroids, DDT, carbamates and organophosphates insecticides [START_REF] Benoı ˆt | Phenotypic effects of concomitant insensitive acetylcholinesterase (ace-1R) and knockdown resistance (kdrR) in Anopheles gambiae: a hindrance for insecticide resistance management for malaria vector control[END_REF][START_REF] Fr | Impact of Insecticide Resistance on the Effectiveness of Pyrethroid-Based Malaria Vectors Control Tools in Benin: Decreased Toxicity and Repellent Effect[END_REF]. In addition, N'guessan et al. [START_REF] N'guessan | Reduced efficacy of insecticide-treated nets and indoor residual spraying for malaria control in pyrethroid resistance area, Benin[END_REF] have demonstrated that, pyrethroid resistance in the primary malaria vector An. gambiae was threatening the effectiveness of ITN and IRS in the areas of high resistance in Benin. The widespread of the insecticide resistance phenomenon in addition to the environmental pollution from chemical insecticides and its high operational cost have led to the need to develop alternative malaria control approaches. Developing effective and eco-friendly tools for reducing the burden or eliminating malaria is highly recommended in this context.

The plant extracts properties as insecticide, repellent or fumigant against mosquitoes and other pests and insect vectors have been well reported. Especially, EOs of several plants have been shown to exhibit significant insecticidal and repellent properties against adult mosquitoes [START_REF] Bassole | Ovicidal and larvicidal activity against Aedes aegypti and Anopheles gambiae complex mosquitoes of essential oils extracted from three spontaneous plants of Burkina Faso[END_REF][START_REF] Nonviho | Insecticidal activities of essential oils extracted from three species of poaceae on Anopheles gambiaespp, major vector of malaria[END_REF][START_REF] Seyoum | Repellency of live potted plants against Anopheles gambiae from human baits in semi-field experimental huts[END_REF], as well as against mosquito larvae [START_REF] Tchoumbougnang | Activite ´larvicide sur Anopheles gambiae Giles et composition chimique des huiles essentielles extraites de quatre plantes cultive ´es au Cameroun[END_REF][START_REF] Madhu | Efficacy of bioactive compounds from Curcuma aromatica against mosquito larvae[END_REF][START_REF] Kweka | Insecticidal activity of the essential oil from fruits and seeds of Schinus terebinthifolia Raddi against African malaria vectors[END_REF][START_REF] Kalaivani | Biological activity of selected Lamiaceae and Zingiberaceae plant essential oils against the dengue vector Aedes aegypti L. (Diptera: Culicidae)[END_REF]. Therefore, it is worth investigating the Beninese flora to find out the plant species whose essential oils could display insecticidal activity against malaria vectors. The grass genus Euclasta Franch (Poaceae) comprises two species, Euclasta clarkei (Hack.) Cope which is distributed in South-western Asia, and Euclasta condylitricha (Steud) (Eu. condolytricha) which is found in Africa, southern Asia and America [START_REF] Piątek | A molecular phylogenetic framework for Anthracocystis (Ustilaginales), including five new combinations (inter alia for the asexual Pseudozyma flocculosa), and description of Anthracocystis grodzinskae sp[END_REF]. Eu. condolytricha, also known as Andropogon condylotrichus, is an annual species, and the height of mature plants can reach up to 2 meters [START_REF] Akoègninou | Flore analytique du Be ´nin[END_REF]. To our knowledge, there is no report on the chemical constituents of Eu. condyloticha essential oil nor its insecticidal activity. Considering the concerted efforts to develop plant products-based insecticides as an excellent alternative to synthetic insecticides, the current study aims to: (i) investigate the chemical profile of the essential oils obtained from Eu. condyloticha flowers harvested in Benin; (ii) evaluate the larvicidal activity of the EOs on Anopheles gambiae s.s. mosquitoes. ´. The collected flowers were dried at 25˚C ±2˚C for 72 hours. Three batches of 100 g of dried flowers were submitted to hydro-distillation in the Clevenger apparatus at 100˚C for 150 min. The distilled oil was dried using anhydrous sodium sulphate, transferred into an airtight amber-colored vial, and stored at 4˚C until further use. The yields were averaged over the three extraction experiments per plant material.

Experimental

Plant material and extraction

Analysis of volatile compounds

The chemical composition of the EOs was evaluated by gas chromatography coupled with mass spectrometry (GC-MS) according to a methodology previously reported [START_REF] Bohounton | Chemical composition and the insecticidal activity of Aeollanthus pubescens leaf essential oil against Anopheles gambiae sensu stricto[END_REF]. Compounds identifications were performed using published data [START_REF] Adams | Identification of essential oil components by gas chromatography/quadrupole mass spectroscopy[END_REF] and a comparison with the NIST mass spectral library.

Mosquito strains

Three An. gambiae s.s. laboratory strains (Kisumu, Acerkis and Kiskdr) from the insectary of the laboratory of Vector-Borne Infectious Diseases at the Institut Re ´gional de Sante ´Publique Alfred Quenum (IRSP-AQ) of the University of Abomey-Calavi in Ouidah (Benin) were used in this study. Kisumu strain originating from Kenya is a reference strain susceptible to all insecticides [START_REF] Shute | A Method of Maintaining Colonies of East African Strains of Anopheles Gambiae[END_REF]. Acerkis, resistant to organophosphate and carbamate-based insecticides, is homozygous for (G119S) mutation [START_REF] Djogbe ´nou | Characterization of insensitive acetylcholinesterase (ace-1R) in Anopheles gambiae (Diptera: Culicidae): resistance levels and dominance[END_REF]. Kiskdr is homozygous for kdr R (L1014F) allele and confers resistance to pyrethroids and DDT [START_REF] Alout | Insecticide resistance alleles affect vector competence of Anopheles gambiae s.s. for Plasmodium falciparum field isolates[END_REF]. Both AcerKis and Kiskdr were supposed to share the same genetic background as the Kisumu strain but differ by the presence of resistance alleles. The colonies of the three strains were maintained in the insectary under optimum conditions (25-27˚C temperature, 70-80% relative humidity and 12:12 light and dark period).

Larval bioassay

Third instar larvae of each mosquito strain were used for the bioassays. Bioassays were performed according to the standard method recommended by the World Health Organization as previously described in [START_REF] Bohounton | Chemical composition and the insecticidal activity of Aeollanthus pubescens leaf essential oil against Anopheles gambiae sensu stricto[END_REF], using seven concentrations of each EO sample as follow: 10, 20, 30, 40, 50, 60 and 70 ppm.

Data analysis

The analysis of dose-mortality responses in larval bioassays was performed using the BioRssay script version 6.2 [START_REF] Roussel | Activite ´compare ´e des hormones juve ´niles en C-18 (JH-I) et en C-16 (JH-III) chez Locusta migratoria[END_REF] in R software Version 3.0 [START_REF] Toong | Isolation of insect juvenile hormone III from a plant[END_REF]. This script calculates the mortality-dose regression using a generalised linear model (GLM). To assess the adequacy of the model, a chisquare test between the observed dead numbers and the dead numbers predicted by the regression is used. Using a likelihood ratio test (LRT), it also tests whether the mortality-dose regressions are similar for the different strains, using a likelihood ratio test (LRT). If there are more than two strains test, it also computes the pairwise test and corrects it using sequential Bonferroni correction [START_REF] Milesi | BioRssay: A R Script for Bioassay Analyses[END_REF]. Finally, it computes the lethal concentrations inducing 50% (LC 50 ) and 95% (LC 95 ) mortality recorded in each mosquito strain and the associated confidence intervals; the resistance ratios, i.e. RR 50 or RR 95 (LC 50 or LC 95 in each strain, divided respectively by the LC 50 or LC 95 of the reference strain) and their 95% confidence intervals. Susceptible or resistant status was defined according to Mazzarri & Georghiou [START_REF] Mazzarri | Characterization of resistance to organophosphate, carbamate, and pyrethroid insecticides in field populations of Aedes aegypti from Venezuela[END_REF] and Bisset et al. [START_REF] Bisset | Temephos resistance and esterase activity in the mosquito Aedes aegypti in Havana, Cuba increased dramatically between 2006 and 2008[END_REF] criteria: RR 50 �1 and RR 50 >1 indicate respectively susceptibility and phenotypic resistance response against the essential oil.

Results and discussion

Essential oils chemical composition

This study characterised the EOs chemical composition of Eu. condylotricha flowers harvested in three localities (Ouessè, Parakou and Sinende ´) of the Benin Republic. The dried flowers of Eu. condylotricha were hydrodistilled and produced an EOs with a characteristic odour. The extractions yields were 0.097%, 0.070% and 0.085%, respectively for EOs from Ouessè, Parakou and Sinende ´plant samples.

Gas Chromatography-Mass Spectrometry (GC-MS) applied to the EOs revealed fifty-five (55) chemical constituents, representing 94.95% to 97.78% of the total essential oils (Table 1). However, four compounds have not been identified, including one compound common to the three EO samples (1.27 to 2.12%) and three compounds in EO from Sinende ´(0.31% to 1.69%). The EO sample from Sinende ´contained forty-four (44) compounds representing 95.86% of the oil. Forty (40) chemical constituents were identified in the EO from Parakou, accounting for 94.95% of the essential oil's mass. The lowest number of phytochemicals, thirty-five [START_REF] Bede | Insect Juvenile Hormones in Plants[END_REF] components, representing 95.86% of all the compounds, was recorded in the EO from Ouessè. There was a wide range of variations in EOs constituents. Each of the three EOs was dominated by juvenile hormones. Indeed, the primary class of components identified was juvenile hormones with 74.87%, 61.29% and 54.92%, respectively in EOs samples from Ouessè, Parakou and Sinende ´. Among the juvenile hormones, the most abundant chemical constituents were the juvenile hormone C16 (48.72% in EO from Ouessè; 35.97% in EO from Sinende ´; and 32.08% in EO from Parakou), followed by methyl farnesoate 10,11-diol (28.73%, 25.67%, and 18.56% respectively in EO from Parakou, Ouessè, and Sinende ´) (Table 1). The sesquiterpenes were the second representative class of components identified in three EOs samples. The major group of sesquiterpenes was the hydrocarbons sesquiterpenes (15.17% in Ouessè, 27.77% in Parakou and 15.36% in Sinende ´), while the oxygenated sesquiterpenes accounted for 1.66% in Ouessè, 2.23% in Parakou and 21.44% in Sinende ´. The hydrocarbons sesquiterpene β-eudesmene was found in EO from Parakou (24.12%), Ouessè (13.36%), and Sinende (12.75%). Among the oxygenated sesquiterpenes, the major compound was Tau-cadinol (18.54%), found mainly in the EO from Sinende ´(Table 1).

PLOS ONE

The EOs analysed in this study were enriched in juvenile hormones. These compounds were known as insect hormones involved in vital physiological processes in larvae and adult insects [START_REF] Goodman | 8-The Juvenile Hormones[END_REF][START_REF] Roussel | Activite ´compare ´e des hormones juve ´niles en C-18 (JH-I) et en C-16 (JH-III) chez Locusta migratoria[END_REF]. However, to date, only two studies reported the presence of juvenile hormones in plant extracts. The first by Toong et al., [START_REF] Toong | Isolation of insect juvenile hormone III from a plant[END_REF] in the whole plants of Cyperus iria L and C. aromaticus Ridl. These authors showed that the excess of the juvenile hormone could distort the wings, change the colour, and induce infertility in adult migratory grasshopper Melanoplus sanguinipes females. The second study has shown the existence of juvenile hormone in aqueous methanolic extract of Cananga latifolia stem bark [START_REF] Yang | Plant-derived juvenile hormone III analogues and other sesquiterpenes from the stem bark of Cananga latifolia[END_REF]. Moreover, juvenile hormones were reported to affect the development of nematodes and represent valuable phytocompounds for insect pest control [START_REF] Goodman | 8-The Juvenile Hormones[END_REF][START_REF] Bede | Insect Juvenile Hormones in Plants[END_REF].

We observed a significant variation in the yield and composition of the essential oils of Eu. condylotricha flowers collected in three different localities in the Benin Republic. This variation can be attributed to the different environmental conditions, such as altitude, solar exposure, and soil composition [START_REF] Figueiredo | Factors affecting secondary metabolite production in plants: volatile components and essential oils[END_REF]. Such variations in EOs yield and phytochemical constituents dependent on the plant's geographical area have been reported for several other plant species, demonstrating that distinct plants from different locations display different chemotypes [START_REF] Saei-Dehkordi | Chemical composition of essential oils in Zataria multiflora Boiss. from different parts of Iran and their radical scavenging and antimicrobial activity[END_REF][START_REF] Yayi-Ladekan | Variation diurne de la composition chimique et influence sur les proprie ´te ´s antimicrobiennes de l'huile essentielle d'Ocimum canum Sims cultive ´e au Be ´nin[END_REF].

The current study is the first work on Eu. condylotricha (Poaceae) flowers' EOs composition. However, several studies have been carried out on the chemical composition of essential oils from other plants of the Poaceae family. Indeed, it was shown that the essential oils of Bothriochloa spp.; Cymbopogon spp. and Chrysopogon zizanioides are dominated by the presence of sesquiterpenes [START_REF] Scrivanti | Essential oil composition of Bothriochloa Kuntze (Poaceae) from South America and their chemotaxonomy[END_REF][START_REF] Verma | Essential oil composition of the sub-aerial parts of eight species of Cymbopogon (Poaceae)[END_REF][START_REF] Campos | Acaricidal properties of vetiver essential oil from Chrysopogon zizanioides (Poaceae) against the tick species Amblyomma cajennense and Rhipicephalus (Boophilus) microplus (Acari: Ixodidae)[END_REF][START_REF] Assaeed | Sesquiterpenes-Rich Essential Oil from Above Ground Parts of Pulicaria somalensis Exhibited Antioxidant Activity and Allelopathic Effect on Weeds[END_REF]. As sesquiterpenes have been found as major constituents of the EOs analysed in this study, EOs from the Poaceae family is enriched in sesquiterpenes.

Larvicidal activity

For several decades, plant derivatives like EOs from aromatic plants have been screened for effective larvicidal activities against An. gambiae larvae [START_REF] Gnankine | Essential Oils as an Alternative to Pyrethroids' Resistance against Anopheles Species Complex Giles (Diptera: Culicidae)[END_REF][START_REF] Osanloo | Plant-Derived Essential Oils; Their Larvicidal Properties and Potential Application for Control of Mosquito-Borne Diseases[END_REF][START_REF] Dantanko | Bioactivity of Essential Oils of Laggera pterodonta and Laggera aurita against Larvae of Anopheles gambiae, Malaria Vector[END_REF]. The EOs of Eu. condylotricha The RR 50 < 1 suggests that the pyrethroid-resistant (Kiskdr) and the carbamate/organophosphate resistant (Acerkis) mosquito strains were susceptible to the Eu. condylotricha flowers essential oils, except Kiskdr strain to the sample from Ouessè. This finding indicates that the presence of insecticide resistance mutations (L1014F in Kiskdr and G119S in Acerkis) did not induce a cross-resistance to the EO extracts. Although, there were no previous studies on the larvicidal activity of Eu. condylotricha flowers EO, some EOs from the plants of the Poaceae family with proven bioactivity against Anopheles species larvae have been reported [START_REF] Govindarajan | Larvicidal and repellent properties of some essential oils against Culex tritaeniorhynchus Giles and Anopheles subpictus Grassi (Diptera: Culicidae)[END_REF][START_REF] Tchoumbougnang | Activite ´larvicide sur Anopheles gambiae Giles et composition chimique des huiles essentielles extraites de quatre plantes cultive ´es au Cameroun[END_REF][START_REF] Manimaran | Larvicidal and knockdown effects of some essential oils against Culex quinquefasciatus Say, Aedes aegypti (L.) and Anopheles stephensi (Liston)[END_REF]. The larvicidal effect of Eu. condylotricha EOs might be caused by the secondary metabolites contained in the OEs. The high amount of juvenile hormones and other sesquiterpenes could be the basis of the observed larvicidal activity. Indeed, several ethnobotanical studies have reported that many sesquiterpene-rich EOs demonstrated excellent larvicidal activities. This is the case of Murraya exotica EO with 61.5% of sesquiterpene hydrocarbons and 6.01% of oxygenated sesquiterpenes that exhibited strong larvicidal activity against Anopheles stephensi larvae (LC 50 = 31.3 ppm) [START_REF] Krishnamoorthy | Identification of chemical constituents and larvicidal activity of essential oil from Murraya exotica L. (Rutaceae) against Aedes aegypti, Anopheles stephensi and Culex quinquefasciatus (Diptera: Culicidae)[END_REF]; Zingiber nimmoni EO that contains 51.9% of sesquiterpene hydrocarbons and 16.2% of oxygenated sesquiterpenes which showed high larvicidal activity against Anopheles stephensi larvae (LC 50 = 41.2 ppm) [START_REF] Govindarajan | Larvicidal and repellent potential of Zingiber nimmonii (J. Graham) Dalzell (Zingiberaceae) essential oil: an eco-friendly tool against malaria, dengue, and lymphatic filariasis mosquito vectors?[END_REF]; and Chloroxylon swietenia EO that has 53.9% of sesquiterpene hydrocarbons and 3.0% of oxygenated sesquiterpenes which was very active on Anopheles stephensi larvae (LC 50 = 14.9 ppm) [START_REF] Ravi Kiran | Composition and larvicidal activity of leaves and stem essential oils of Chloroxylon swietenia DC against Aedes aegypti and Anopheles stephensi[END_REF]. Moreover, some sesquiterpenes exhibited larvicidal activity against mosquito larvae. Sesquiterpenes such as humulene (LC 50 = 6.19 ppm); caryophylene (LC 50 = 41.66 ppm); elemene (LC 50 = 10.26 ppm); germacrene D-4-ol (LC 50 = 6.12 ppm) and α-cadinol (LC 50 = 10.27 ppm) were demonstrated to be toxic to Anopheles subpictus larvae [START_REF] Govindarajan | Larvicidal and repellent potential of Zingiber nimmonii (J. Graham) Dalzell (Zingiberaceae) essential oil: an eco-friendly tool against malaria, dengue, and lymphatic filariasis mosquito vectors?[END_REF][START_REF] Pavela | Essential oils for the development of eco-friendly mosquito larvicides: A review[END_REF]. Caryophyllene oxide (LC 50 = 49.81 ppm) and germacrene D (LC 50 = 49.46 ppm) were also shown toxic to Anopheles anthropophagus larvae [START_REF] Zhu | Chemical composition and larvicidal activity of essential oil of Artemisia gilvescens against Anopheles anthropophagus[END_REF].

Conclusion

Local plant derivatives with insecticidal properties constitute a promising alternative for malaria vector control. The EOs of Eu. condylotricha samples collected in different regions of Benin Republic were analyzed and found to contain high concentrations of sesquiterpenoid compounds. The juvenile hormone C16 was the major constituent of the phytochemical compounds of all EOs samples. The first report of this hormone in Eu. condylotricha provides solid background for a wide range of applications for malaria control. The OEs showed significant larvicidal activity against resistant An. gambiae strains. Before translating these research findings into operational interventions, further investigations on mechanisms by which EOs mediate their bioinsecticidal activity may require evaluations of the major EOs' constituents separately. This study findings contribute to the dissemination of knowledge regarding the chemical composition and larvicidal activity of this Beninese plant species, which is almost incipient in the literature.

Fig 1 .

 1 Fig 1. Euclasta condylotricha Steud harvested in Benin. A: The whole plant with stem, leaves and flowers. B: Areal part of the plant showing the flowers. https://doi.org/10.1371/journal.pone.0278834.g001

  flowers were tested against third instar larvae of two resistant (Kiskdr, Acerkis) and one susceptible (Kisumu) laboratory strains of An. gambiae. The larvicidal activity of the EOs under a series of concentrations ranging from 10 to 70 ppm were evaluated. As shown in Fig 2, toxicity of the three EOs is concentration-dependent. The EOs exhibited significant larvicidal activity, with LC 50 values ranging from 35.21 to 52.34 ppm (Table 2) after 24 hours of exposure. The EO sample from Sinende ´was the most active, with LC 50 of 35.21 ppm on Kiskdr larvae, followed by the EO from Ouessè with LC 50 of 38.10 ppm and 38.46 ppm on Acerkis and Kisumu larvae, respectively. Acerkis strain larvae (LC 50 = 38.10 ppm) were significantly more susceptible to EO from Ouessè compared to Kiskdr (LC 50 = 46.30 ppm, p < 0.001), while no significant difference was observed compared to Kisumu larvae (LC 50 = 38.46 ppm, p = 0.33). However, for the Sinende śample, Kisumu and Kiskdr larvae susceptibility was significantly different (p < 0.001). The LC 50 value recorded in Kiskdr larvae (35.21 ppm) was significantly lower than that in Acerkis (44.45 ppm, p < 0.001) and Kisumu larvae (52.34 ppm, p < 0.001). This is the first report on the larvicidal activity of Eu. condylotricha flowers OEs against An. gambiae larvae. In agreement with the criteria established by Cheng et al., [46] plant EOs showing LC 50 values < 50 ppm within 24 hours were very active, and LC 50 (24h) values < 100 ppm were declared active. Based on those classifications, Eu. condylotricha flowers EOs could represent an inexpensive and environmentally benign agent for controlling malaria vectors.

Fig 2 .

 2 Fig 2. Percentage of dead An. gambiae larvae after 24 h exposition at different concentrations of EOs. 0 ppm indicates control solution with water and ethanol without EO. https://doi.org/10.1371/journal.pone.0278834.g002

  LaBEV) of the University of Abomey-Calavi, Benin, and the voucher specimens are kept at the National Herbarium of the University. The voucher specimen numbers are AAC 190/HNB, AAC 191/HNB and AAC 192/HNB respectively for samples from Ouessè, Parakou and Sinende

	Euclasta condyloticha Steud flowers (Fig 1) were collected in November 2018 from the locali-
	ties of Ouessè (middle) (8˚38'19.3"N, 2˚38'04.7"E), Parakou (northeast) (9˚22'55.2"N, 23
	6'48.1"E) and Sinende ´(north-west) (10˚06'52.7"N, 2˚22'58.7"E) of Benin Republic. The taxo-
	nomic identification of the plants was carried out at the 'Laboratoire de Botanique et E ´cologie
	Ve ´ge ´tale (

Table 1 . Chemical composition of the Euclasta condylotricha Steud essential oils.

 1

	No	Rt	Compounds		Percentage of total composition
				Ouessè	Parakou	Sinende 1
		13.87	α-Pinene	-	tr	-
	2	15.08	β-Myrcene	tr	tr	tr
	3	17.08	o-Cymene	-	tr	-
	4	17.36	β-Terpinyl acetate	-	tr	-
	5	22.5	Linalool	-	0.15	tr
	6	26.76	trans-Borneol	tr	-	tr
	7	26.8	3-Pinanylamine	-	-	tr
	8	27.6	Terpinen-4-ol	tr	tr	-
	9	28.59	L-α-Terpineol	-	tr	-
	10	28.93	n-Hexyl butanoate	tr	tr	tr
	11	29.2	trans-2-Hexenyl Butyrate	tr	tr	tr
	12	32.97	2-Methylbutyl caproate	tr	tr	tr
	13	41.67	(-)-β-Elemene	0.72	1.15	1.02
	14	41.93	(2Z)-2-Hexenyl butyrate	tr	0.15	-
	15	42.73	β-Gurjunene	tr	tr	-
	16	43.19	β-Caryophyllene	0.36	0.44	0.33
	17	44.5	L-Alloaromadendrene	tr	tr	tr
	18	45.28	α-Humulene	0.17	0.17	0.23
	19	45.43	E.E-8.10-Dodecadien-1-ol	-	0.22	tr
	20	45.52	α-Selinene	0.16	0.15	0.15
	21	46.08	1.6-Dimethylhepta-1.3.5-triene	1.03	0.57	0.84
	22	46.73	(+)-Aromadendrene	0.11	0.20	0.15
	23	47.32	β-Eudesmene	13.46	24.12	12.75
	24	47.83	γ-Gurjunene	0.94	2.45	1.18
	25	48.4	Aciphyllene	0.13	0.25	0.16
	26	49	γ-Cadinene	-	-	0.56
	27	49.12	(-)-α-Panasinsene	tr	0.14	tr
	28	49.63	NI	-	-	0.31
	29	50.38	α-Cadinene	-	-	tr
	30	51.15	Elemol	0.13	tr	0.10
	31	52.02	Geranyl isobutyrate	-	0.10	-
	32	52.12	D-Nerolidol	0.11	-	0.26
	33	52.53	β-copaene	-	-	tr
	34	52.63	(+)-Spathulenol	tr	tr	-
	35	52.87	β-Caryophyllene oxide	tr	tr	tr
	36	53.41	Epiglobulol	-	-	0.42
	37	53.88	Ledol	tr	tr	tr
	38	54.66	(-)-Neointermedeol	0.34	0.61	0.28
	39	54.79	Epicubenol	-	-	0.23
	40	55.10	7-epi-cis-sesquisabinene hydrate	-	tr	-
	41	55.12	Germacrene D-4-ol	-	-	tr
	42	55.54	Di-epi-1.10-cubenol	-	-	tr
	43	55.74	(-)-β-Longipinene	tr	tr	tr
	44	56.36	tau.-Cadinol	-	-	18.54
	45	56.67	β-Eudesmol	tr	0.10	-
	46	56.98	Globulol	0.52	0.96	1.08
						(Continued )

Table 1 .

 1 (Continued)

	No	Rt	Compounds		Percentage of total composition
				Ouessè	Parakou	Sinende 47
		58.89	Shyobunol	-	-	tr
	48	62.59	NI	-	-	0.49
	49	64.27	Methyl farnesoate 10.11-diol	25.67	28.73	18.56
	50	66.58	2.3-Dimethylpentanal	tr	tr	tr
	51	66.99	Juvenile hormone C18	0.48	0.48	0.39
	52	68.36	Juvenile hormone C16	48.72	32.08	35.97
	53	68.55	NI	2.12	1.27	1.56
	54	69.23	NI	-	-	1.69
	55	69.93	trans-Z-α-Bisabolene epoxide	0.69	0.46	0.53
			Oxygenated monoterpenes	0.54	1.3	0.9
			Sesquiterpene hydrocarbons	15.17	27.77	15.36
			Oxygenated sesquiterpenes	1.66	2.23	21.44
			Juvenile hormones	74.87	61.29	54.92
			Others	3.62	2.49	5.29
			Total identified (%)	95.86	94.95	97.78

Compounds are listed in order of elution from a HP-5MS column. Rt: retention times; tr: trace amount (<0.1%); (-): not detected; NI: not identified. https://doi.org/10.1371/journal.pone.0278834.t001

Table 2 . Toxicity against An. gambiae larvae after 24 h exposure. EOs samples Strains LC 50 (ppm) 95% C.I RR 50 95% CI [LCL-UCL] LC 95 (ppm) 95% CI [LCL-UCL] Slope ±S.E Intercept ±S.E Chi(p) value

 2

		Kisumu	38.46	36.94-39.96	-	-	74.45	69.65-80.67	5.73 ± 0.45	-9.09 ± 0.47	0.98
	Ouessè	Acerkis	38.10	36.43-39.75 0.99	0.86-1.14	79.33	73.42-87.20	5.16 ± 0.28 -8.16 ± 0.45	0.95
		Kiskdr	46.30	44.15-48.59 1.20	1.04-1.38	92.05	83.43-104.65	5.51 ± 0.38	-9.18 ± 0.63	0.71
		Kisumu	48.06	45.99-50.31	-	-	108.81	98.07-124.02	4.63 ± 0.27	-7.79 ± 0.45	0.92
	Parakou	Acerkis	40.25	38.09-42.55 0.84	0.73-0.95	92.20	82.94-105.51	4.57 ± 0.29	-7.33 ± 0.48	0.64
		Kiskdr	47.41	45.36-49.60 0.99	0.87-1.21	102.88	93.09-116.72	4.88 ± 0.29	-8.19 ± 0.49	0.59
		Kisumu	52.34	50.36-54.50	-	-	100.72	92.10-112.85	5.78 ± 0.35	-9.94 ± 0.6	0.69
	Sinende ´Acerkis	44.45	42.23-46.79 0.85	0.74-0.97	93.49	84.18-107.20	5.09 ± 0.35	-8.4 ± 0.58	0.45
		Kiskdr	35.21	31.61-38.63 0.67	0.59-0.77	103.31	85.28-138.55	3.51 ± 0.38	-5.44 ± 0.61	0.44

No mortality was observed in the control LC 50/95 : lethal concentrations; S.E: standard error; C.I: Confidence interval; RR 50 is resistance ratio at LC 50 : LC 50 (resistant strain)/ LC 50 (Kisumu). LCL: Lower confidence limit; UCL: Upper confidence limit Chi(p) is indicated to judge whether the data are well fitted to the regression or not. The fits are acceptable when the p-value is over 0.05. https://doi.org/10.1371/journal.pone.0278834.t002

PLOS ONE | https://doi.org/10.1371/journal.pone.0278834January 23, 2023

This study received partial support from the Wellcome Trust in the form of an International Intermediate Fellowship in Public Health and Tropical Medicine grant to LSD [N˚109917/Z/15/ Z]. The funders had no role in study design, data collection and analysis, publication decision, or manuscript preparation.

All relevant data are within the paper and its Supporting Information files.

Supporting information


Dmytro Kandaskalov
email: dmytro.kandaskalov@im2np.fr

Liangzhao Huang

Johnathan Emo

Philippe Maugis
email: philippe.maugis@im2np.fr

begindocument/before

published or not. The documents may come

Carbon diffusion in bcc-and bct-Fe: Influence of short-range C-C pair interactions studied from first-principles calculations

is usually not affected by the carbon-carbon interactions. However, carbon-rich regions are commonly observed in out-of-equilibrium structures such as supersaturated bainite or as-quenched martensite, with solute carbon fractions up to a few at% [START_REF] Garcia-Mateo | Low temperature bainitic ferrite: Evidence of carbon super-saturation and tetragonality[END_REF]. Carbon segregation at structural defects (grain boundaries and dislocations) can also raise the solute fraction up to ∼8 at% [START_REF] Wilde | Three-dimensional atomic-scale mapping of a Cottrell atmosphere around a dislocation in iron[END_REF]. Lastly, low-temperature spinodal decomposition of martensite leads to carbon-rich body-centered regions ranging from 10 to 14 at% C [START_REF] Zhu | Carbide Characterisation in Low Temperature Tempered Steels[END_REF][START_REF] Haasen | The early stages of the decomposition of alloys[END_REF][START_REF] Rochegude | Interstitial Atom Ordering in Binary Fe-N Solid Solutions Studied by Mössbauer Spectrometry[END_REF][9][START_REF] Taylor | Spinodal decomposition during aging of Fe-Ni-C martensites[END_REF].

By itself the phenomenon of carbon diffusion is important in many ways. Technological processes such as surface hardening by carburizing or nitriding are realized by interstitial diffusion [START_REF] Rajesh | Evaluation and Diffusion Assessment for Surface Hardening Processes[END_REF]. Molecule-surface interactions, such as those in corrosion and heterogeneous catalysis, can also involve diffusion of adsorbed atoms into subsurface layers with further diffusion into the bulk of the material. The local stress field has a strong effect on the interstitial diffusion in iron, due to the induced deformation of the interstitial sites [START_REF] Maugis | Stress-controlled carbon diffusion channeling in bct-iron: A mean-field theory[END_REF][START_REF] Svoboda | An innovative concept for interstitial diffusion in stressed crystals[END_REF]. Furthermore, interstitial diffusion causes phenomena such as internal friction [START_REF] Bagramov | Internal friction in a martensitic high-carbon steel[END_REF], strain ageing [15], embrittlement [START_REF] Rellick | Intergranular embrittlement of iron-carbon alloys by impurities[END_REF] and steel erosion [START_REF] Green | Influence of microstructure on the erosion of plain carbon steels[END_REF]. Finally, a better understanding of the interstitial diffusion mechanisms in iron will aid controlling phase transformations in steel and modelling microstructure formation in the material.

The classical measurement of carbon diffusion in ferrite was done in 1950 by Wert [START_REF] Wert | Diffusion Coefficient of C in -Iron[END_REF]. He showed that in the temperature range of 238-473 K, the data agree very well with a linear Arrhenius plot, which yields a diffusion barrier of 0.87 eV. After Wert, carbon diffusion in ferrite was measured with various methods and in different temperature ranges [19]. For the low-temperature linear Arrhenius regime, the assumed atomic mechanism involves carbon migration from an octahedral O-site to a nearest-neighbour O-site via a saddle point (S) (i.e., the tetrahedral T-site for an isolated carbon atom) [19][20][START_REF] Beshers | Diffusion in Body-Centered Cubic Metals American Society for Metals[END_REF][START_REF] Condit | Interstitial Diffusion in the B.C.C. Lattice[END_REF][START_REF] Mclellan | The thermodynamics of dilute interstitial solid solutions with dual-site occupancy and its application to the diffusion of carbon in alpha iron[END_REF][START_REF] Da Silva | Diffusion of carbon and nitrogen in BCC iron[END_REF]. Those measurements are however limited to very dilute alloys. Diffusion measurements in carbon-supersaturated alloys are difficult, owing to the instability of the solid solution, this is why results in the literature are only indirect [9,[START_REF] Decristofaro | The kinetics of carbon clustering in martensite[END_REF]. They nevertheless evidence a significant decrease in carbon diffusivity in martensite when the carbon content is increased. Theoretical analyses showed that C-C long-range elastic interactions are responsible for this decrease, via an increase of the migration barriers [START_REF] Maugis | Stress-controlled carbon diffusion channeling in bct-iron: A mean-field theory[END_REF][START_REF] Lawrence | Carbon diffusion in supersaturated ferrite: a comparison of mean-field and atomistic predictions[END_REF].

In more carbon-rich regions, short-range C-C interactions cannot be neglected. There, migration of a C atom is expected to be affected by the presence of neighbouring C atoms. This short-range interaction is difficult to study by experimental techniques and thus theoretical studies at the atomic level are necessary. In 2003, Jiang et al. studied by density functional theory (DFT) the diffusion of C in ferrite [START_REF] Jiang | Carbon dissolution and diffusion in ferrite and austenite from first principles[END_REF]. Their results on single carbon atom migration correspond well to the experimental data in dilute alloys [START_REF] Wert | Diffusion Coefficient of C in -Iron[END_REF][START_REF] Da Silva | Diffusion of carbon and nitrogen in BCC iron[END_REF]. This is also true for recent DFT [START_REF] Maugis | Revisiting the pressure effect on carbon migration in iron[END_REF] and embedded atom method (EAM) approaches [START_REF] Becquart | Atomistic modeling of an Fe system with a small concentration of C[END_REF]. However, those studies do not take into account the effects of short-range C-C interactions on diffusion. There are several theoretical works on C-C interactions when both atoms sit on octahedral interstitial sites (OO-configurations) [START_REF] Jiang | Carbon dissolution and diffusion in ferrite and austenite from first principles[END_REF][START_REF] Chentouf | Effect of interstitial carbon distribution and nickel substitution on the tetragonality of martensite: A first-principles study[END_REF][START_REF] Ruban | Self-trapping of carbon atoms in -Fe during the martensitic transformation: A qualitative picture from ab initio calculations[END_REF][START_REF] Domain | Ab initio study of foreign interstitial atom (C, N) interactions with intrinsic point defects in -Fe[END_REF]. There is no similar study of transition-state OS-configurations, although they are essential to computing diffusion migration energies.

This article provides a multi-scale approach to study the influence of short-range C-C interactions on C migration in body-centered iron by combining ab initio calculations and kinetic Monte Carlo (KMC) simulations. We discuss the energetic properties of pair configurations with C-C distances up to 6 Å. One of the C atoms (called the static atom) occupies an O-site. The other C atom (the migrating atom) is positioned in an O-site (stable state) or a S-site (saddle-point state). We propose a shell-nomenclature to classify all the C-C configurations. The migration paths of a carbon atom in the presence of a neighbouring static atom are discussed and compared to the case of an isolated migrating atom. The obtained migration data are used to compute carbon order parameter and diffusion coefficients in concentrated alloys by KMC simulations.

II. METHODOLOGY A. DFT computational techniques

In the present work, all calculations have been carried out using the Vienna Ab-initio Simulation Package (VASP) [START_REF] Kresse | Ab initio molecular dynamics for liquid metals[END_REF] based on the density functional theory (DFT). To resolve the self-consistent Kohn-Sham equations, the projectedaugmented wave (PAW) method was employed [START_REF] Kresse | From ultrasoft pseudopotentials to the projector augmented-wave method[END_REF]. The PAW method is a generalization of the pseudo-potential and linear augmented-plane-wave methods (within the frozen-core approximation) with the computational efficiency of pseudo-potential DFT algorithms. The atomic basis for PAW were 2s 2 2p 2 and 4s 1 3d 7 for C and Fe atoms respectively.

The exchange-correlation functional PBE [START_REF] Perdew | Generalized Gradient Approximation Made Simple[END_REF] was used to calculate the electronic structure within the generalized gradient approximation (GGA). All calculations are presented with spin polarization. The plane-wave energy cut-off was established at 400 eV for all studied systems. The supercells were Fe 128 C n , 0 ≤ n ≤ 2. A 4 × 4 × 4 Monkhorst and Pack grid [START_REF] Monkhorst | Special points for Brillouin-zone integrations[END_REF] of k-points has been used for integration in the irreducible part of the Brillouin zone. These parameters were sufficient to converge the total energies down to 50 meV per supercell for studied systems.

All systems have been fully relaxed with respect to the atomic positions with constant parameters of the supercell.

The atomic relaxations were performed with a conjugate-gradient algorithm. The computed cell parameter of iron (a 0 =2.825 Å) corresponds well to the other theoretical works [START_REF] Domain | Ab initio study of foreign interstitial atom (C, N) interactions with intrinsic point defects in -Fe[END_REF][START_REF] Fang | Stability and crystal structures of iron carbides: A comparison between the semi-empirical modified embedded atom method and quantum-mechanical DFT calculations[END_REF][START_REF] Connetable | First principle calculations of the -Fe3AlC perovskite and iron-aluminium intermetallics[END_REF][START_REF] Kandaskalov | Study of multivacancies in alpha Fe[END_REF][START_REF] Kandaskalov | Maugis A first-principle study of the structural, elastic, lattice dynamical and thermodynamic properties of a"-Fe16C2 and a"-Fe16N2 phases[END_REF]. The computed magnetic moment of the Fe bulk is 2.16 µ B /at which in a good agreement with experimental data 2.12µ B /at [START_REF] Kittel | Introduction to Solid State Physics[END_REF].

To determine the migration minimum energy paths (MEPs), the nudged elastic band (NEB) method was used [START_REF] Henkelman | Nudged-elastic band method with two climbing images: Finding transition states in complex energy landscapes[END_REF][START_REF] Henkelman | A climbing image nudged elastic band method for finding saddle points and minimum energy paths[END_REF]. The transition-state searches were carried out with the convergence criteria of 10 -5 eV for total energies of ionic step iterations and 0.01 eV/ Å for forces acting on the atoms. We used 3 images in these NEB analyses (1/4, 1/2 and 3/4 points on reaction coordinate). Then, all energy landscapes were fitted with polynomials of degree 6.

B. Stable and saddle-point configurations

To study the atomic jump of a migrating carbon atom located in the vicinity of another carbon atom, we defined two types of C-C pair configurations: an OO-configuration consists of two carbon atoms in octahedral interstitial sites, while an OT-configuration consists of one carbon atom in the octahedral site and the other in the tetrahedral site. The OO-configurations are stable and correspond to initial and final states of an elementary carbon jump. In OT-configurations, the migrating atom is located between the initial and final sites (OO-configurations). In most cases, this T site corresponds to a saddle point (S) as for single C atom migration in iron. We ranked the OO-and OT-configurations according to the distance between the carbon atoms, represented by the vector 1 2 a 0 [u, v, w]. Notice that for a given C-C distance, several non-equivalent configurations exist. The set of all configurations sharing the same distance is referred to as a shell (first, second, third...). The equivalent configurations in the same shell form a subshell. The subshells are labelled with small letters (a, b, c...) for OT-configurations and capital letters (A, B, C...)

for OO-configurations. For example, the 7th shell of OT-configuration has three subshells: 7a, 7b and 7c (see Figure 1); the 7th shell of OO-configuration has only two subshells: 7A and 7B.

In our model we studied the configurations up to the 14th shell, which represents a reasonable limit for our cubic supercell of 128 Fe atoms. The shell number N is directly related to the distance vector 1 2 a 0 [u, v, w]: For OT-configurations the shell number is equal to

N OT = u 2 + v 2 + w 2 + 3 4 . (1)
In the case of OO-configuration the shell number is

N OO = u 2 + v 2 + w 2 -A, (2)
where A = 0 for

u 2 + v 2 + w 2 ≤ 6, A = 1 for 8 ≤ u 2 + v 2 + w 2 ≤ 15 and A = 2 for 16 ≤ u 2 + v 2 + w 2 ≤ 23.
The ranking of subshells is based also on the vector 1 2 a 0 [u, v, w]. Let us take the same 7th OT-shell as an example (N OT = 7). This shell has 3 subshells defined by the vectors 1 2 a 0 [2, 1.5, 0], 1 2 a 0 [0, 1.5, 2] and 1 2 a 0 [2.5, 0, 0]. We attribute the labels a, b and c in order of increase of |w|. Thus, the subshell with the vector 1 2 a 0 [0, 1.5, 2] is 7c. If |w| is the same for several subshells, as we have in this example, we rank them by increasing the sum |u| + |v|. In the present case, the labels are as follows: 1 2 a 0 [2.5, 0, 0] is 7a and 1 2 a 0 [2, 1.5, 0] is 7b (see Table I).

In our work we label OT-configurations as OS-configurations taking into account that our interest consists specifically in the finding of saddle-points, however for some topological discussions when we will discuss precisely the mid position between initial and final state we will return to using OT-configurations.

We also need to introduce sublattice information for the interstitial atoms of the same shell. A carbon atom on an octahedral site can be of three types (O x , O y or O z ) depending in which sublattice (x, y or z) is located this atom.

Mechanically, the presence of a carbon atom creates a prolate tetragonal distortion along the direction corresponding to its sublattice.

A carbon atom on a tetrahedral site creates an oblate tetragonal distortion and can be of three types (T x , T y or T z ) depending on direction of the oblate effect [START_REF] Maugis | Stress-controlled carbon diffusion channeling in bct-iron: A mean-field theory[END_REF]. Figure 2a illustrates a carbon jump starting from an octahedral site of type z, through a tetrahedral site of type x, and ending at an octahedral site of type y, i.e. a "z -x -y" jump. During this path, the bcc lattice is successively tetragonally distorted along directions z, x and y. For any kind of jump, the lattice is distorted in three successive directions [START_REF] Maugis | Revisiting the pressure effect on carbon migration in iron[END_REF]. Each octahedral interstitial site has 4 nearest-neighbour octahedral sites. Thus any carbon atom can migrate towards 4 possible directions (Figure 2b). The 5 sites involved (central site and its 4 neighbours) are located in the same plane, which is perpendicular to the axis of the central site. Without loss of generality, the static atom was positioned on the z sublattice (O z type of site). The migrating atom, for its part, can occupy any of the three sublattices x, y, z, which causes tetragonal distortions in one of three possible directions. Notice that these distortions do not generate additional shells because each shell is made of the same type of interstitial carbon atoms: z-type, or x-with y-type (x and y are equivalent in this respect).

C. Carbon migration barrier

For an isolated carbon atom, the migration path is symmetric since initial and final states are geometrically identical, and the saddle point is at mid-path. For this case, there is only one migration energy E m0 , since the direct and reverse jumps are equivalent (Figure 3a ).

In the case of a pair of carbon atoms with 128 iron atoms per supercell, we used as reference energy the state of two isolated carbon atoms positioned on octahedral interstitial sites. The reference energy is therefore twice the insertion energy of carbon in O-site. The interaction energy E C-C of a carbon-carbon pair is calculated as

E C-C = E ins (C -C) -2E ins (C) (3)
where E ins (C) is the insertion energy of one carbon atom in octahedral site and E ins (C -C) is the insertion energy of two neighbouring carbon atoms. E C-C > 0 means repulsion. The insertion energies are calculated as

E ins (C) = E(Fe 128 C 1 ) -E(Fe 128 ) -E(C) E ins (C -C) = E(Fe 128 C 2 ) -E(Fe 128 ) -2E(C) (4)
where E(Fe 128 C 2 ), E(Fe 128 C 1 ) and E(Fe 128 ) are the energies of the Fe 128 C n systems with n =2, 1, 0 in correspondence.

In this work we are interested in three energetic parameters: C-C interaction energy E C-C mentioned above and the direct and reverse migration energies E m1 and E m2 (Figure 3b). The direct jump is defined as the case of C-C pair dissociation (i.e. the C-C distance increases), while the reverse jump corresponds to the pair association (the C-C distance decreases).

Let us define E C-C (OS) as the interaction energy in the transition state, i.e. when the migrating atom is at the saddle point. E C-C (OO init ) and E C-C (OO fin ) are the interaction energies at the initial and the final state respectively.

Then the migration energies are simply expressed as

E m1 = E C-C (OS) -E C-C (OO init ) E m2 = E C-C (OS) -E C-C (OO fin ) (5)
At variance to the case of an isolated carbon atom, presence of the static carbon atom breaks the symmetry of the direct and reverse jumps. Thus, as long as the interaction between the two carbon atoms exists, we expect that:

(a) The energies of the initial, transition and final states will be modified;

(b) The interaction energies E C-C (OO init ) and E C-C (OO fin ) will be different from zero.

In the following, we will use the energetic quantities E C-C , E m1 and E m2 to discuss the influence of the geometric parameters of various C-C pair configurations on the interaction and migration energies. Em1 and Em2 are different also, the saddle point is slightly shifted from mid-path.

D. Mean-field elastic model

Zener ordering of carbon leads to a stress-strain field that affects carbon migration. Using a mean-field elastic model [9,[START_REF] Maugis | Stress-controlled carbon diffusion channeling in bct-iron: A mean-field theory[END_REF], we can express the strain components as a function of the carbon content. First, we define the atom fractions c i = n i /n, where n i is the number of carbon atoms on the octahedral sites of type i = x, y, z, and n is the total number of iron atoms. The total fraction c = i c i is the number of carbon atoms per iron atom in the crystal. In the case of Zener ordering along the z-direction, the order parameter is defined as η = (c z -c x or y )/c. In the mean-field approximation, the volume density of dipole tensor p is given by

p = 1 V 0 i c i P Oi , (6)
where P Oi is the dipole tensor of a carbon atom on an octahedral site of type i and V 0 is the atomic volume of iron.

The carbon-induced strain tensor at equilibrium writes

= S p, ( 7
)
where S is the compliance tensor of iron. They are allowed to jump during the simulation. The C-C pairwise interactions that are given by the DFT calculations are included in the carbon jump frequency:

Γ = ν 0 exp - E m0 + (E S pair -E O pair ) -(P S -P O ) : k B T , (8)
with

E S pair = k∈S z k E k C-C (OS), E O pair = k∈O z k E k C-C (OO init ), (9)
where P S (resp. P O ) denotes dipole tensor associated to carbon atom in the saddle-point state (resp. octahedral site) and ν 0 is the attempt frequency that is assumed to be strain-independent. The values of these parameters can be found in [START_REF] Maugis | Stress-controlled carbon diffusion channeling in bct-iron: A mean-field theory[END_REF]. S and O are respectively the collections of all the investigated tetrahedral and octahedral shells, and z k is the number of carbon atoms at the k-th shell of the considered carbon atom. In order to estimate the order parameters, the occupancies of the carbon atoms are measured after each atomic jump. After the system reaches equilibrium, 100 measurements of the mean squared displacement (MSD) of the carbon atoms after 3 × 10 5 atomic jumps (i.e., a total of 3 × 10 7 carbon jumps for each simulation) are performed to compute the tracer diffusion coefficients. The relative convergence of the obtained carbon diffusivities is within ±1% with these settings. The tracer diffusivities D i (i = x, y, z) are related to the MSD R 2 i along the x, y, z directions during a period of time τ via the Einstein-Smoluchowski equation: R 2 i = 2D i τ . The MSD is given by

R 2 i = N n=1 R n,i 2 = f i N n=1 R 2 n,i , ( 10
)
where N is the total number of jumps in the period of time τ and R n,i is the displacement of a carbon atom along the i-direction due to the n-th jumps.

N n=1 R 2 n,i
is the uncorrelated MSD corresponding to the random diffusion path, and f i is the correlation factor characterizing the deviation between the tracer diffusion and a random walk.

Note that the displacements of about 1000 C atoms are simultaneously measured during the simulation. Therefore, the final MSD results from the average of all carbon atoms.

III. RESULTS AND DISCUSSION

A. Isolated carbon migration

To test our methodology, we calculated the migration energy of carbon in pure iron and compared it with experimental and theoretical literature. The experimental migration energy values are 0.86 eV [START_REF] Wert | Diffusion Coefficient of C in -Iron[END_REF] by Wert and 0.874 eV by Da Silva et al. [START_REF] Da Silva | Diffusion of carbon and nitrogen in BCC iron[END_REF]. These results were obtained by study of lnD evolution as 1/T function.

Theoretically, by EAM (Embedded Atom Method) calculations the migration energy of 0.81 eV was found by ab initio NEB approach which is more precise than EAM calculations. In the case of Domain [START_REF] Domain | Ab initio study of foreign interstitial atom (C, N) interactions with intrinsic point defects in -Fe[END_REF] and Bialon [START_REF] Bialon | Ab initio study of boron in -iron: Migration barriers and interaction with point defects[END_REF] the overestimation of migration energy to the small underestimation of T-site energy stability because of symmetry imposed in their methodology.

We obtained a migration energy E m0 = 0.869 eV, which agrees precisely with experimental literature and most theoretical works. Our value, obtained from volume-fixed calculations, is close to 0.875 eV calculated with cell relaxation [START_REF] Maugis | Revisiting the pressure effect on carbon migration in iron[END_REF]. These results allow us to conclude once more that, for 128 iron atoms per supercell, cell relaxation has a negligible effect on the migration energy.

B. OO-interactions

DFT studies of some OO-configurations have been published by Ruban [START_REF] Ruban | Self-trapping of carbon atoms in -Fe during the martensitic transformation: A qualitative picture from ab initio calculations[END_REF] and by Chentouf [START_REF] Chentouf | Effect of interstitial carbon distribution and nickel substitution on the tetragonality of martensite: A first-principles study[END_REF] under constant supercell volume. In the present study, we recalculated a full set of OO-configurations, as they constitute input data for our NEB calculations. We also re-adapted the shell nomenclature for OO-configurations and complemented published energetic results with some extra data. In his work, Ruban computed the same configurations as ours, but the exact values for energies were given only for a few structures. We extracted approximate values out of his paper, when possible. In Chentouf's work, the precise values of energies are given, but some configurations are missing. Our energetic results correspond well to these two works [START_REF] Chentouf | Effect of interstitial carbon distribution and nickel substitution on the tetragonality of martensite: A first-principles study[END_REF][START_REF] Ruban | Self-trapping of carbon atoms in -Fe during the martensitic transformation: A qualitative picture from ab initio calculations[END_REF].

Our configuration nomenclature introduces shells and subshells: each shell corresponds to a specific distance between the carbon atoms while a subshell represents equivalent configurations for the same distance. In the case of OOconfigurations the maximum number of subshells is two (labeled A and B), except for the 15th shell, which has 3 subshells (A, B and C).

We studied 24 configurations, which represent 14 shells plus configuration 15C (Table I andFigure The repulsive character of configurations 3A, 4B, 7B, 8B, 10B and 14B can be interpreted in terms of elastically interacting pairs of Fe-C-Fe dipoles. In those configurations, the two carbon atoms sit in the same O z sublattice, such that both dipoles are parallel to the z-axis. This creates cumulative tetragonal distortions along the z-axis.

Furthermore, this may shift neighbouring iron atoms in either identical or opposite directions (see Figure 6). Several cases can be distinguished. The least stable configuration 4B shows an iron atom "compressed" between 2 carbon atoms. In configuration 14B, the strong interaction corresponds to two dipoles in line. Interestingly, repulsion is not present in configuration 4A despite the proximity of the dipoles, because the shifts of neighbouring Fe atoms are synchronized and the system comports similar to two isolated carbon atoms. Conversely, configuration 7B shows antisynchronized iron atom shifts. Compared to configuration 4A, the antisynchronized shift in configuration 7B destabilizes the structure and creates a repulsion. This is also true for configurations 3A and 10B. Notice that configuration 14B is similar to that existing in the α -Fe 16 C 2 phase (labeled 2NN' in Refs. [START_REF] Kandaskalov | Maugis A first-principle study of the structural, elastic, lattice dynamical and thermodynamic properties of a"-Fe16C2 and a"-Fe16N2 phases[END_REF][START_REF] Kandaskalov | Thermodynamic stabilities in the Fe -Fe 16C2 system: Influence of carbon-carbon interactions studied by DFT[END_REF]).

C. OS-interactions

Using the same method, we studied the OT-configurations, which correspond to the transition states for carbon migration. The first 14 shells, corresponding to 29 configurations, are gathered in Figure 7 and Table II. The 1st shell was not considered, since its very small C-C distance ( 1 4 a 0 0.7 Å) is more than twice smaller than the C-C distance For OS-configurations, the shells can contain up to 5 subshells (labelled a, b, c, d and e). Since OS-configurations are not stable because one of the carbon atoms occupies a saddle point, it is not possible to obtain the C-C interaction energies for this type of configurations by a standard total energy minimization procedure. Therefore, we used the NEB method to compute the C-C interaction energy evolution during migration paths, to obtain the energies for all OS-configurations. We remind that for two isolated OS carbon atoms the E C-C energy is 0.87 eV and thus if OS-interaction energy is higher than 0.87 eV it is considered as repulsive interaction (Figure 3a). The results (Table II) show us that mostly all interactions are repulsive in this way. We could see a group of interactions highly repulsive with energies from 1.38 eV up to 1.74 eV. These configurations correspond to small distances between 2 carbon atoms (2a, 3a, 3b, 5c) and also repulsive configuration 8b which will be discussed further. Other configurations have energies in range 0.86 -1.04 eV showing small repulsion effect. And only one configuration (14b) is clearly attractive with 0.76 eV energy.

The results on OS E C-C energies would be most efficient from the discussion of migration paths.

As it was shown in the section IIC, the migration energy of dissociation (E m1 ) is the difference between the energy of OS-configuration and OO-configuration of initial state and the migration energy of association (E m2 ) is the difference between the energy of OS-configuration and OO-configuration of final state. Thus, the migration energies (E m1 and E m2 ) depend on the relative stability of corresponding OO-and OS-configurations.

The static carbon atom (which is located in z-site) influences the energy path of the migrating carbon atom via the initial OO-site, the transition OS-site and the final OO-site together. Since these three atomic positions occupy three different sublattices, the effect is contrasted. As mentioned before, a C atom in O-site creates a prolate distortion whereas in S-site it creates an oblate distortion. The prolate and oblate distortions may compensate each other when they act in the same direction. In this case, the OS-configuration is "stabilized" and the interaction energy is smaller than 0.87 eV. Conversely, if the oblate and prolate distortions act in opposite directions, the OS-configuration is destabilized and the C-C interaction energy is larger than 0.87 eV.

D. Analysis of migration paths

We present our results of the migration path in the dissociation way, i.e. the C-C distance is closer in the initial state than in the final state. During an elementary migration path, the initial and final states are always OO-configurations, even if a given OO-configuration is by itself very repulsive (such as configurations 1A or 4B). All elementary migration paths were arranged by axes. Figures 8 and9 present 4 migrations paths along z axis (we recall that the migrating atom initially sits on a z-site). Figure 10 presents 7 migration paths in the x -y plane. These 11 paths cover all possible migration scenarios within a 2a 0 radius around the static carbon atom.

It is to be noted that the transition-state OS-configuration is in most cases located at the mid transition path, i.e.

at a crystallographic T-site. However, for the most repulsive configurations (2a, 3a, 4a and 5c) the saddle point is shifted towards the initial state. The migration paths 1A → 2a → 2A (Figure 8), 1A → 3a → 4A and 4B → 5c → 5B (Figure 10) have their maximum at 1/4 migration path, i.e. in the middle between initial OO-configuration and OT-configuration. For the migration path 2A → 4a → 5B the maximum is at 1/3 migration path. All 4 migration energies are very small: 0.14-0.27 eV which is related first of all to a large repulsion of two carbon atoms in the initial OO-configuration, and second, to the displacement of Fe atoms around which perturbs the topology of O-and S-sites in corresponding OO-and OS-configurations. Note that even very repulsive OO-configurations such as 1A, 2A, 3A have an activation barrier towards dissociation. In details, the migration paths will be discussed in the next subsection.

We might expect that during a migration path, when the migrating C atom is moving away from the static carbon atom along a linear path O -S -O -S -O, the OO-configuration energy should gradually decrease down to 0 eV as the C-C repulsion decreases. For its part, the OS-configuration energy should gradually decrease down to 0.87 eV.

Our results reveal a more complex behaviour.

Let's analyse the migration in [1 0 k] direction (Figure 8b). If we look at the evolution of OO-configuration energy along the migration path, we have the expected picture: the energy decreases gradually as 1.96 → 0.83 → 0.04 → 0.07 → -0.04 eV. However, during this path the energy of OT-configurations doesn't decrease, but it oscillates: 1.60 → 1.06 → 1.33 → 0.76 eV. In this case we need to discuss the OT-configurations to have the moving atom topologically in the T site to discuss energy trend variation. For the last two migrations in [1 0 k] direction the OT-configuration coincides with OS-configuration, but for the first two configurations they are different because of high repulsion between the carbon atoms. This particular difference between OS and OT-sites exists only for 2a, 3a, 4a and 5c configurations.

In the initial stage, this energy (1.60 eV) is almost twice larger than 0.87 eV, due to the large C-C repulsion.

The energy of the second OS-configuration is 1.06 eV (i.e. the energy is decreased 1.5 times), but for the 3rd OSconfiguration the energy increases up to 1.33 eV. It finally reaches the value of 0.76 eV, smaller than the expected final value (0.87 eV). For [k 2 2] axis we see that the energies of configurations 7B, 9b, 8B and 11d are shifted by +0.15 eV with respect to isolated systems, thus the migration energy remains close to 0.87 eV. Only the association of 11A to 8B has larger activation energy (0.95 eV) as the configuration 11A is completely isolated. For [k 2 0] which is more close to the static carbon atom than the previous axis [k 2 2]. We observe for the first migration smaller migration energies 0.79/0.76 for dissociation/association and the next two migrations have the behaviour of isolated C atom. The smaller migration energies for the first migration corresponds to the fact that OS-configuration 5a is more stable than isolated S-with O-configurations as synchronized "prolate-oblate" effect is present. Along axis [k 0 2] we observe that 4B is very repulsive but then stabilization the system is similar to isolated. Migration along the remaining axes [k 3 1] and [k 1 3] proceeds as if C atoms were isolated except for the repulsive 10B configuration, which produces slight decrease in migration energies 0.73/0.78 eV for dissociation of 10B to 13B and association 10B to 9B.

E. Order parameter at equilibrium

The KMC simulation results of the Zener order parameter at different carbon contents are plotted in Fig. 11. Note that if the C-C pair interactions are neglected, the order parameter can be given by the mean-field elastic model presented in [START_REF] Maugis | Stress-controlled carbon diffusion channeling in bct-iron: A mean-field theory[END_REF]. The good agreement between the analytical results and the KMC simulation results neglecting the C-C pair interactions validates the KMC approach. The mean-field approach predicts an order/disorder transition at c = 1.27 at% and T = 300 K. This transition carbon content is slightly modified due to the introduction of C-C pair interactions. According to the KMC simulations, the transition content should be between 1.27 at% and 1.30 at% at

We investigated as well the carbon order parameter at different temperatures (cf. Fig. 12). We set the total fraction of carbon atoms to 2.50 at%. Here again, KMC simulation results neglecting the C-C pair interactions correspond well to the analytical results from the mean-field approach. Following the KMC simulation results, the order/disorder transition temperature is affected by the C-C pair interactions: the predicted transition temperature is about 590 K if the interactions between carbon atoms are neglected and about 570 K if not.

F. Tracer diffusion in ferrite and martensite

We present in Fig. 13 the KMC simulation results of the carbon tracer diffusion coefficients and the correlation factors at different carbon contents. Here again, if the C-C pair interactions are neglected, the diffusion coefficients can be given by an analytical approach presented in [START_REF] Maugis | Stress-controlled carbon diffusion channeling in bct-iron: A mean-field theory[END_REF]. We have a good agreement between the analytical results and the KMC simulation results neglecting the C-C pair interactions. As predicted in previous studies [START_REF] Maugis | Stress-controlled carbon diffusion channeling in bct-iron: A mean-field theory[END_REF], carbon diffusion is highly anisotropic due to carbon ordering. The carbon diffusivity along the z-axis (i.e., the Zener ordering direction) is much lower than that along x and y directions (perpendicular to the ordering direction). Comparison between the KMC simulation results with and without considering the C-C pair interactions highlights the effect of these interactions on the kinetic properties. The carbon diffusivity along the z-axis is significantly increased by the C-C repulsive interactions (up to 10 times at high carbon contents). The increase in carbon diffusivity along x and y directions are less important (e.g., about 20% when c = 4.0 at% at 300 K). Moreover, the interactions between carbon atoms lead to important kinetic correlations: the correlation factors are smaller than 0.2 when c > 3.0 at% at 300 K.

Also, we studied the effect of C-C pair interactions on the kinetic properties at different temperatures (cf. Fig. 14).

According to the simulation results, the increase of the carbon diffusivity is significant at low temperatures. For example, the acceleration of the carbon diffusivity along the z direction is over 10 times at T < 300 K. A relevant kinetic correlation occurs at low temperatures: the correlation factors are smaller than 0.4 at T < 400 K.

As presented in Table II, most of the migration barriers between the investigated configurations are smaller than the isolated carbon migration barrier (0.87 eV). This decrease of migration barriers explains the increase of the carbon diffusivity. We consider the migration paths along z-direction (cf. Fig. 9 and Fig. 10). As mentioned in Section III D, the carbon migration path along [1 2 k], [3 0 k] and [3 2 k] is similar to the isolated one because the activation energy is close to 0.87 eV. However, a "prolate-oblate" effect occurs along the [1 0 k]-migration path. The carbon atom can be "trapped" in 5A configuration by performing forward and backward jumps connecting to 2A configurations. This effect should be responsible for the important kinetic correlation along the z-direction. As for the migration paths in the x-y plane, there is no important oscillation of the migration barriers as seen along the z-direction. Hence, it is unlikely for the carbon atoms to be trapped in the x-y plane. However, the migration paths along [k 0 0], [k 1 1] and [k 0 2] are particularly different from the isolated one (i.e., the random diffusion path) due to the important C-C repulsive interactions in 1A, 2A and 4B OO-configurations-very small carbon migration barriers for an "avoid" from these configurations (0.17 eV, 0.55 eV and 0.89 eV, respectively). This difference leads to the kinetic correlations along

x-and y-directions because the carbon atoms tend to "escape" each other.

The tracer diffusion coefficients of carbon atoms in Fe-C martensite were previously calculated by Lawrence et al. [START_REF] Lawrence | Carbon diffusion in supersaturated ferrite: a comparison of mean-field and atomistic predictions[END_REF] by applying the molecular dynamic (MD) simulation based on an EAM potential. This simulation was limited to high temperature (> 800 K) so that the statistic to measure the MSD was sufficient in the time scale of MD.

As predicted in their studies, our results also showed slow carbon diffusion in martensite. However, we cannot quantitatively compare the diffusion coefficients presented in the present study with those obtained in [START_REF] Lawrence | Carbon diffusion in supersaturated ferrite: a comparison of mean-field and atomistic predictions[END_REF] because, on one hand, the investigated carbon concentrations and temperatures are not the same; and on the other hand, our simulation is based on a more recent energy database of the elastic dipoles and the migration barriers given by the DFT calculation. Therefore, the diffusion properties obtained in this paper must be different from those of Lawrence et al. [START_REF] Lawrence | Carbon diffusion in supersaturated ferrite: a comparison of mean-field and atomistic predictions[END_REF]. Compared with their studies, our approach based on KMC simulation allows to compute the tracer diffusion coefficient at a wider range of temperature (> 200 K in the present study). Moreover, we provide a more complete investigation into the anisotropy of the carbon diffusion, as well as the kinetic correlations in martensite.

IV. CONCLUSION

We studied by a DFT method short range C-C interactions in ferrite and their influence on the migration energies of carbon. Stable octa-octa (OO) and unstable octa-saddle state (OS) configurations were investigated.

We found that the configuration stability is not always a simple function of C-C distance and there is no direct correlation in energy between neighbours OO-and OS-configurations. The absence of correlation is related to the fact that during a simple migration path the migrating C atom goes through 3 sublattices (x, y, z) while a static carbon atom remains always in the same sublattice. Concerning OO-configurations, in most cases the interaction energy is close to 0.0 eV. However, when 2 carbon atoms occupy the same sublattice the interaction energy rises to 0.10-0.30 eV.

In the case of small C-C distances the interaction energy increases up to 2.0 eV, indicating strong repulsion.

For the first time the OS-configurations were studied in detail by the NEB technique. These configurations are unstable, being in the vicinity of tetrahedral sites. For OS-configurations we observed an interesting effect of "oblateprolate" interactions: when two carbon atoms are located in the same sublattice some configurations reduce the interaction energy (0.60-0.80 eV) while others increase it, compared to an isolated carbon atom (0.95-1.10 eV).

Finally the migration energy was studied for all configurations where C-C distance is in the limit of 2a 0 . We observed a range of small migration energies (0.14-0.27 eV) and a larger range of medium migration energies (0.60-0.80 eV).

Most configurations have migration energies (0.87±0.05 eV) close to the isolated system. Some configurations have migration energies larger than 0.92 eV.

FIG. 1 :

 1 FIG.1:The 7th shell of OT-configurations contain 20 sites distributed in 3 subshells (7a, 7b and 7c) (a) and the 7th shell of OO-configurations contain 12 sites distributed in 2 subshells (7A and 7B). In these figures, the static carbon atom is in the center of the polyhedron, in a Oz site. For OT-configurations each subshell includes the geometrically equivalent T -sites for carbon atoms: 8 sites for 7a, 4 sites for 7b and 8 sites for 7c and for OO-configurations each subshell includes the geometrically equivalent O-sites for carbon atoms: 4 sites for 7A and 8 sites for 7B

FIG. 2 :

 2 FIG. 2: (a) Tetragonal distortion effect during carbon migration from an O-site to a nearest-neighbour O-site. In the case shown, the migration path is Oz -Tx -Oy. (b) A carbon atom in a Oz-site has four possible jumps: two towards Ox-sites through Ty-sites, and two towards Oy-sites through Tx-sites.

FIG. 3 :

 3 FIG. 3: (a) Migration of C atom in ferrite. The two values over the saddle point (0.87, 0.87) are right and left migration energies Em1 and Em2 (0.87 eV each). The interstitial energy of C in O-site is the energy reference. (b) An example of migration path of C atom in presence of static C atom. Here the initial and final state are different and as a consequence migration energies

E

  . Kinetic Monte Carlo simulations Kinetic Monte Carlo (KMC) simulations based on the residence-time algorithm are performed on the body-centered cubic (bcc) lattice with the octahedral sites. Periodic boundary conditions are applied to the simulation box. About 1000 carbon atoms are inserted. The size of the simulation box is adjusted by the atomic fraction of carbon atoms.

FIG. 4 :

 4 FIG.4: The first 14 shells for OO-configurations, drawn at the same scale and the same space orientation. Different subshells inside a shell are represented by different colors and labelled in case if there is more than one subshell per shell.

FIG. 5 :FIG. 6 :

 56 FIG.5: Key OO-configurations illustrating the elastic dipole-dipole interaction. Static carbon atom (green), carbon atoms in the same plane (brown), carbon atoms out-of-plane in configurations 3A and 10B (violet). Each OO-configuration is labelled (in blue) and the C-C interaction energy in eV is given (in pink). The iron atoms are represented by large golden spheres.

FIG. 7 :

 7 FIG.7: The first 14 shells for OS-configurations, drawn at the same scale and the same space orientation. Different subshells inside a shell are represented by different colors and labelled in case if there is more than one subshell per shell.

FIG. 8 :

 8 FIG. 8: Prolate-oblate interaction and its influence on the C-C interaction energy of OS-configurations. Complete migration path along [1 0 k] vector with indication of OO-and OS-configurations along the path (b). The couple of migration energies (Em1, Em2) are indicated under the saddle points. Top: the atomic configurations corresponding to the tetragonal positions k = 0.5, 1.5, 2.5 and 3.5 (a). Blue arrows corresponds to Fe atom displacement by static C atom and red arrows -to displacement by dynamic C atom.

FIG. 9 :

 9 FIG. 9: Three migration paths parallel to the z-axis: [1 2 k], [3 0 k] and [3 2 k].

FIG. 10 :FIG. 12 :FIG. 14 :

 101214 FIG. 10: All 7 migration paths for C atom developing along x axis: [k 0 0], [k 1 1], [k 2 2], [k 2 0], [k 0 2], [k 1 3] and [k 3 1].

  Garruchet et al.[START_REF] Garruchet | Modelling the carbon Snoek peak in ferrite: Coupling moleculardynamics and kinetic Monte-Carlo simulations[END_REF]. The DFT calculations give results in the range of 0.85 to 0.94 eV: Fors et al. 0.85 eV [45], Jiang et al. 0.87 eV [27], Souissi et al. 0.865 eV [46], Domain et al. 0.90 eV [32] and Bialon et al. 0.94 eV [47]. These results are based on NEB approaches as our or on a simple comparison of energy difference between interstitial T and O sites. Domain et al. demonstrated that smaller supercells give higher activation energies (0.92 eV for 54 atoms vs 0.90 eV for 128 atoms). We observe that generally theoretical results could precisely simulate experimental data by

TABLE I :

 I DFT results for octahedral-octahedral carbon-carbon configurations (OO-configurations): configuration label, number of equivalent configurations (# config.), interaction energy EC-C compared to literature (in eV), vector components [u v w] (values in parentheses can be permuted) and bond length dC-C (in a0 units).

	Configuration # config. EC-C ref. [30] ref. [31]	Vectors	dC-C
	1A	4	1.96	1.95	2.0	[(±1 0) 0]	√ 4/4
	2A	8	0.83	0.83	0.8	[(±1 0) ±1]	√ 8/4
	3A	8	0.26	0.27	0.25	[±1 ±1 ±1]	√ 12/4
	4A	4	0.01	0.02	0.0	[(±2 0) 0]	√ 16/4
	4B	2	2.00	1.96	2.0	[0 0 ±2]	√ 16/4
	5A	8	0.04	0.05	0.0	[(±2 ±1) 0]	√ 20/4
	5B	8	0.04	-	0.0	[(±1 0) ±2]	√ 20/4
	6A	16	0.07	0.07	0.0 [(±2 ±1) ±1]	√ 24/4
	7A	4	0.03	0.03	0.0	[±2 ±2 0]	√ 32/4
	7B	8	0.16	0.15	0.1	[(±2 0) ±2]	√ 32/4
	8A	4	0.06	0.06	0.0	[(±3 0) 0]	√ 36/4
	8B	16	0.17	0.17	0.1 [(±2 ±1) ±2]	√ 36/4
	9A	8	0.07	-	0.0	[(±3 0) ±1]	√ 40/4
	9B	8	0.07	0.07	0.0	[(±1 0) ±3]	√ 40/4
	10A	16	0.02	0.02	0.0 [(±3 ±1) ±1]	√ 44/4
	10B	8	0.30	0.30	0.3	[±1 ±1 ±3]	√ 44/4
	11A	8	0.04	0.03	0.0	[±2 ±2 ±2]	√ 48/4
	12A	8	0.00	0.00	0.0	[(±3 ±2) 0]	√ 52/4
	12B	8	-0.01	-	0.0	[(±3 0) ±2]	√ 52/4
	13A	16	0.00	-	0.0 [(±3 ±2) ±1]	√ 56/4
	13B	16	0.01	0.00	0.0 [(±2 ±1) ±3]	√ 56/4
	14A	4	-0.02 -0.02	0.0	[(±4 0) 0]	√ 64/4
	14B	2	0.83	0.83	0.8	[0 0 ±4]	√ 64/4
	15C	16	-0.04	-	-	[(±1 0) ±4]	√ 68/4

in graphite (1.42 Å) or diamond (1.54 Å). It is thus highly improbable in the bcc structure.

TABLE II :

 II DFT results for octahedral-saddle carbon-carbon configurations (OS-configurations): configuration label, number of equivalent configurations (# config.), interaction energy EC-C (in eV), migration path, dissociative and associative migration energies Em1 and Em2, vector components [u v w] and bond length dC-C (in a0 units).

	Configuration # config. EC-C	Migration path	Em1	Em2	Vectors	dC-C
	1a	4	-	-	-	-	[(±0.5 0) 0]	√ 1/4
	2a	8	2.12	1A→ 2a → 2A	0.27	1.30	[(±1 0) ±0.5]	√ 5/4
	3a	4	2.09	1A → 3a → 4A	0.17	2.12	[(±1.5 0) 0]	√ 9/4
	3b	16	1.38	2A → 3b → 3A	0.55	1.12 [(±0.5 ±1) ±1]	√ 9/4
	4a	8	1.03	2A→ 4a → 5B	0.23	1.02	[(±1 0) ±1.5]	√ 13/4
	5a	8	0.81	4A → 5a → 5A	0.79	0.76	[(±0.5 ±2) 0]	√ 17/4
	5b	16	0.88	3A → 5b → 6A	0.62	0.80 [(±1.5 ±1) ±1]	√ 17/4
	5c	8	2.14	4B → 5c → 5B	0.14	2.10	[(±0.5 0) ±2]	√ 17/4
	6a	16	0.90	5A → 6a → 6A	0.86	0.83 [ (±2 ±1) ±0.5]	√ 21/4
	7a	4	0.81	4A → 7a → 8A	0.81	0.75	[(±2.5 0) 0]	√ 25/4
	7b	8	0.96	5A → 7b → 7A	0.80	0.80	[(±1.5 ±2) 0]	√ 25/4
	7c	8	0.84	5B → 7c → 7B	0.90	0.83	[(±1.5 0) ±2]	√ 25/4
	8a	16	1.00	6A → 8a → 8B	0.93	0.83 [(±2 ±1) ±1.5]	√ 29/4
	8b	8	1.37	5B → 8b → 9B	1.33	1.30	[(±1 0) ±2.5]	√ 29/4
	9a	16	0.92	6A → 9a → 10A	0.85	0.90 [(±2.5 ±1) ±1]	√ 33/4
	9b	16	1.04	7B → 9b → 8B	0.89	0.88 [(±0.5 ±2) ±2]	√ 33/4
	10a	8	0.95	8A → 10a → 9A	0.89	0.88	[(±3 0) ±0.5]	√ 37/4
	11a	8	0.89 7A → 11a → 12A 0.86	0.90	[(±2.5 ±2) 0]	√ 41/4
	11b	16	0.91 9A → 11b → 10A 0.84	0.89 [(±0.5 ±3) ±1]	√ 41/4
	11c	8	0.94	7B → 11c → 12B	0.77	0.93	[(±2.5 0) ±2]	√ 41/4
	11d	16	0.99 8B → 11d → 11A 0.81	0.95 [(±1.5 ±2) ±2]	√ 41/4
	11e	16	1.09	9B → 11e → 10B	1.02	0.78 [(±0.5 ±1) ±3]	√ 41/4
	12a	8	0.87 9A → 12a → 12B 0.81	0.88	[(±3 0) ±1.5]	√ 45/4
	12b	16	0.95 8B → 12b → 13B 0.78	0.94 [(±2 ±1) ±2.5]	√ 45/4
	13a	4	0.90 8A → 13a → 14A 0.83	0.91	[(±3.5 0) 0]	√ 49/4
	13b	16	0.87 10A → 13b → 13A 0.85	0.86 [(±1.5 ±3) ±1]	√ 49/4
	13c	16	1.03 10B → 13c → 13B 0.73	1.03 [(±1.5 ±1) ±3]	√ 49/4
	14a	16	0.86 12A → 14a → 13A 0.86	0.85 [(±3 ±2) ±0.5]	√ 53/4
	14b	8	0.76 9B → 14b → 15C 0.69	0.79	[(±1 0) ±3.5]	√ 53/4

Acknowledgments

The Agence Nationale de la Recherche (contract C-TRAM ANR-18-CE92-0021) is thanked for its financial support. Centre de Calcul Intensif d'Aix-Marseille is acknowledged for granting access to its high performance computing resources. This work was also supported by the HPC resources of IDRIS under allocation number 100524 made by GENCI.

Thus, this work showed that the role of static carbon atom is not negligible for the diffusion phenomena and migration energy varies from 0.15 eV up to 2.12 eV while for the single C atom migration the energy is 0.87 eV. The computed energy database was used to calculate the carbon order parameter and diffusion coefficients by kinetic Monte Carlo simulations. According to our simulations, C-C pair repulsive interactions lead to (1) the modification of the order/disorder transition conditions depending on the carbon content and temperature; (2) the increase of the carbon diffusivity, especially along the Zener ordering direction (up to ten times); (3) the important kinetic correlations. These results clearly show that while calculated by mean-field model diffusion coefficient D x is not highly depend on C-C interactions (at 4 % atomic content of C, D x is about 10% higher when C-C interactions are included), the diffusion coefficient D z is highly depended of C-C interactions and the diffusion coefficient D z is one order of magnitude higher