Anne Gaudin

Véronique Ansan

Jean-Pierre Lorand

Sylvain Pont

Schists At

Armorican Massif

Genesis of a florencite-bearing kaolin deposit on ordovician schists at Saint-Aubin-des-Châteaux, Armorican Massif, France

Keywords: Kaolinite, Al phosphate, REE, weathering, hydrothermal alteration, Armorican 15 Massif 16 17 18

archive ouverte pluridisciplinaire HAL, est destinée au dépôt

135 suggesting a genetic relationship between both events [START_REF] Chauris | Metallogeny of Armorican Massif[END_REF]Brown and 136 Dallmeyer, 1996;[START_REF] Bouchot | Late Variscan mineralizing 864 systems related to orogenic processes: the French Massif Central[END_REF][START_REF] Gloaguen | Palaeozoic oolitic ironstone of the French Armorican Massif: a chemical and structural trap for 1041 orogenic base metal-As-Sb-Au mineralisation during Hercynian strike-slip deformation[END_REF]. At the end of Late 137 Carboniferous, a network of N160°E trending fractures and faults crosscut oblique folds and 138 granite intrusions [START_REF] Gapais | Is the Hercynian belt of Brittany a major shear zone?[END_REF].

139

Since the end the Paleozoic era, the Armorican Massif underwent intense sub-tropical 140 wet climatic conditions, favoring weathering and the formation of lateritic soils (e. g. Estéoule-141 Choux, 1983;[START_REF] Le Corre | Le Massif Armoricain[END_REF]. 150 m 2 , with a mean depth of ~ 5 to 7 m (Fig. 2c). The pit face shows four superposed lithological 151 levels (Table 2 and Fig. 3): A lowermost level (TRS-0) of greyish schist preserving sub-vertical 152 (WNW-ESE/85°S) cleavage (S) with veinlets and nodules of iron oxides associated with 153 centimetric thick quartz veins (Fig. 3a). Above this level, a ~4 m thick darker claystone schist 154 level containing a few trilobite fossils [START_REF] Dubreuil | Saint-Mars-la-Jaille (421)[END_REF] 429 (2014). With respect to most predominant and trace elements, the mudstone composition 430 sampled in the Tertre Rouge schists is similar to that of quoted averages for post-Archean shales 431 [START_REF] Wedepohl | Environmental influences on the chemical composition of shales and clays[END_REF]. Transition and non-transition metal contents agree reasonably well with 432 worldwide average shale compositions [START_REF] Johnson | Secular distribution of Highly Metalliferous Black Shales corresponds 1072 with Peaks in Past Atmosphere Oxygenation[END_REF]. Several elements, however, do 458 Light Rare Earth Elements), P, Th, and Sr budget is clearly controlled by accessory phosphates 459 (namely florencite; predominant over monazite and xenotime-Y). Florencite is a phosphate 460 member of a phosphate subgroup of the alunite supergroup of minerals having the general 461 formula AB 3 (XO 4 ) 2 (OH) 6 , where the A-site is filled with Ce, La, and Nd; the B-site with Al;

462 and the X-site with P [START_REF] Bayliss | Alunite supergroup: recommended 845 nomenclature[END_REF]. The structure of florencite also accommodates 463 varying amounts of Ca 2+ , Pb 2+ , Hg 2+ , K + , Ba 2+ , Sr 2+ , Rb + (in A-site), Fe 3+ , Cu 2+ , Zn 2+ , Sn 2+ , V 3+ , 464 Cr 3+ , and Ga 3+ (in B-site) (Dill, 2001;[START_REF] Gaboreau | Aluminum phosphate-sulfate 1012 minerals associated with Proterozoic unconformity-type uranium deposits in the East Alligator 1013 River Uranium Field, Northern Territories, Australia[END_REF][START_REF] Repina | Florencite (Sm) -(Sm, 1154 Nd)Al 3 (PO 4 ) 2 (OH) 6 : A New Mineral Species of the Alunite -Jarosite Group from the Subpolar 1155[END_REF] ). They may be authigenic in shales and sandstone (e.g. early diagenetic 479 REE, Al phosphate minerals in marine sandstones are a major sink for oceanic biogenic 480 phosphorus [START_REF] Pouliot | Florencite: A first occurrence in Canada[END_REF]Rasmussen, 1996 and references therein)). They may 481 be a weathering product in argillized igneous rocks, resulting from pervasive chemical 482 alteration of primary phosphates [START_REF] Nriagu | Phosphate Minerals[END_REF][START_REF] Schwab | The formation of aluminous 1174 phosphates through lateritic weathering of rocks[END_REF] and coexisting

  pit located to the SE of Saint-Aubin-des-châteaux (Fig. 2b) has been 146 working since 1948. Today the B.N.Z. company still extracts the kaolinized schist for the 147 manufacture of firebricks. It is excavated into the Middle Ordovician slaty schists ("Calymene-148 bearing schists" in ancient geological maps) from the "Angers-Traveusot" formation located 149 into the south limb of WNW-ESE trending anticlinal fold. The pit covers a surface of 20,000

  408 the part played by rutile and florencite in the trace element budget. A more complex behavior 409 can be inferred from Fig.12for the transition elements (Sc, V, and Cr) Y and HREE that do 410 not show systematic concentration changes with respect to the schists. Their behavior likely 411 reflects a multimineral partitioning effect and a stronger preference for the crystalline structure 412 of clay minerals (cf Linnen et al., several kaolin deposits described thus far, the Tertre Rouge deposit shows 417 evidence of a polygenic formation process of clay minerals, such as (i) the pervasive 418 replacement of the schists by kaolinite and other secondary minerals (Fe oxyhydroxides); (ii) 419 the crystallization of nacrite-dickite-pyrophyllite inside veins/veinlets that were unrelated to the 420 pervasive pit scale formation process of kaolinite. Both processes involved a peculiar carbon-421 rich schist parent material characterized by aluminum phosphates (namely Sr-and Ca-bearing 422 florencite (Ce)). This protolith, ascribed to Middle Ordovician continental platforms, will be 423 discussed first, before addressing the timing and petrogenesis of clay formation. 424 425 7.1 Florencite-bearing Ordovician sedimentary parent rock composition 426 The bulk-rock analyses of the Tertre Rouge schists are indicative of an upper crust 427 mudstone/shale composition characterized by a two-fold metal enrichment (As, Be, Bi, Cr, Ga, 428 Cu, In, Pb, Sc, Sn, Ta, Th, V, P, and Ti) compared with the estimate made by Rudnick and Gao

Table 4 ,

 4 TableS1). If close to 1, their coefficients correlations with the major index elements 388 (TiO 2 , K 2 O, P 2 O 5 ) indicate the major mineral concentrator for each element, in addition to SEM 389 EDS analyses of Figs 9 and 10. It can be deduced that LILE (Ba, Rb, and Cs) and Pb are 390 preferentially hosted in mica/illite. Sr-bearing florencite is the main concentrator for Sr, Th, and Fig.12). Compared to the schists, all the veins are 405 strongly depleted in LILE (Cs, Rb, Ba, Sr, Th, U, and LREE) and HFSE (Zr, Hf, Nb, and Ta); 406 this agrees with the near 0 modal abundance of mica, florencite, and rutile in those veins. Their 407 HFSE contents positively correlate with TiO 2 , while LREE correlate with P 2 O 5 , thus outlining

	(Figs. 3b-e) exists; this level was

391 LREE, with LREE+MREE and Th being also hosted in monazite (SEM EDS data). Xenotime 392 controls Y and HREE. The homogeneous distribution of florencite in all of the three analyzed 393 schists perfectly accounts for the very flat and highly reproducible LREE segment in the crust-394 normalized pattern of Fig.

12

. It also explains the coherent geochemical behavior of P, V, and 395 Cr; all three elements may replace PO 4 in the anionic group of phosphates (see

Chang et al., 396 1997)

. By contrast, the variability of MREE and HREE likely reflects more complex 397 distribution patterns between two trace phases (xenotime + zircon). Zircon is the main 398 concentrator of Sc, Zr, Hf, and U while rutile carries the budget of Nb, Ta, and W. Note that 399 multiple positive correlations between chalcophile elements Cu-Co-Zn and Sn, yet no Cu-Zn 400 sulfides were detected with the SEM. The good positive correlation between Sb and Mo is 401 characteristic of euxinic shale sediments (residual altered pyrite)

[START_REF] Samalens | The role of black shales as a source of sulfur and 1170 semimetals in magmatic nickel-copper deposits: Example from the Partridge River Intrusion, 1171 Duluth Complex, Minnesota, USA[END_REF] 402 Coveney and Pasava, 2004;[START_REF] Large | Multistage Sedimentary and Metamorphic Origin of 1084 Pyrite and Gold in the Giant Sukhoi Log Deposit, Lena Gold Province[END_REF]

.

  Nd 0.13 Ce 0.28 Sr 0.24 Fe 0.05 Ca 0.17 Al 3.32(P 1.87 O 8 ) (OH) 6 (calculated on the basis of 11 oxygen 468 atoms). Their EDS spectra systematically detected Ca and Sr, indicating solid solutions toward 469 goyazite and crandallite, which are two other types of aluminum phosphate (Fig.10). Thorium 470 was detected in half of the EDS spectra as well, suggesting the substitution of 2REE 3+ by Ca 2+

									. The Tertre
	465 Rouge florencites are generally too fine-grained to provide reliable in-situ compositions; only
	466 one	particle	allowed	us	to	obtain	its	structural	formulae:	-
	467 La 0.15 471 + Th 4+ . Sulfur was also detected, suggesting a significant solution toward sulfates (a property
	472 of Al phosphates in the APS association, see Dill, 2001 and reference therein).
	473	APS minerals of the goyazite (SrA1 3 (PO 4 ) 2 (OH) 5 -H 2 0) group, including florencite-Ce,
	474 have been identified in a number of Australian coals, as well as in kaolinite claystones
	475 (tonsteins) associated with coal seams in New South Wales (Ward et al., 1996). These
	476 phosphates typically represent a range of solid solutions containing goyazite, gorceixite
	477 (BaAI 3 (PO 4 ) 2 (OH) 6 • H 2 0), crandallite (CaA1 3 (PO 4 ) 2 (OH) 5 • H 2 0), and florencite ((La, Ce,
	478 Nd)Al 3 (PO 4 ) 2 (OH) 6

Table 3 .

 3 Major element contents of Tertre Rouge samples (wt%).

	Element	TRS-0b	TRS-0c	TRS-1d	TRV-1	TRV-2a	TRV-2b
	SiO 2	54.99	58.72	53.23	55.19	45.9	46.36
	TiO 2	1.06	1.21	1.14	0.01	0.26	0.08
	Al 2 O 3	26.52	25.28	29.5	34.05	37.45	38.14
	Fe 2 O 3	4.48	1.98	2.99	0.21	2.94	2.77
	MgO	0.79	0.52	0.27	0.13	0.05	0.02
	MnO	0.006	n.d.	n.d.	n.d.	n.d.	n.d.
	CaO	0.13	0.13	0.16	0.07	0.06	0.11
	Na 2 O	0.82	0.92	0.85	0.21	0.09	0.02
	K 2 O	3.93	3.32	3.24	1.07	0.34	0.1
	P 2 O 5	0.25	0.25	0.4	0.03	0.13	0.07
	L.O.I.	7.46	7.86	8.97	9.62	13.77	13.51
	Total	100.44	100.18	100.75	100.58	100.98	101.18
	Organic C	0.23	0.89	0.58	n.a.	n.a.	n.a.

L.

O.I. : loss on ignition; n.d. : not detected; n.a.: not analyzed

L'

Highlights

 The Tertre Rouge kaolin deposit replaces ordovician schists of the Armorican Massif  Al-phosphates are major concentrators of LREE and Sr  Kaolinization was mainly produced by Paleogene weathering  A minor, older hydrothermal event is suggested by high-T clay mineral veinlets

Declaration of interests

+ The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

☐The authors declare the following financial interests/personal relationships which may be considered as potential competing interests:


Florian Naudet
email: florian.naudet@chu-rennes.fr

Clara Locher

Alain Braillon

André Gillibert

Letter to the editor: Covid-19 And Vit-D: Disease Mortality Negatively Correlates With Sunlight Exposure

The manuscript "Covid-19 And Vit-D: Disease Mortality Negatively Correlates with Sunlight Exposure" [START_REF] Lansiaux | Covid-19 And Vit-D: Disease Mortality Negatively Correlates With Sunlight Exposure[END_REF] held our attention as we found fatal shortcomings that invalidates the analyses and conclusions.

General considerations

First, the title itself is misleading as the manuscript does not report any information on Vitamin D since it was not measured in the context of this study.

Second it relies on unfounded hypotheses as the authors refer to a website [2] without any scientific validity rather than to peer-reviewed articles.

Third, the reporting is poor. The manuscript does not conform to any reporting guidelines such as STROBE [START_REF] Von Elm | The Strengthening the Reporting of Observational Studies in Epidemiology (STROBE) statement: guidelines for reporting observational studies[END_REF]. The public data sources citations (references 13 to 15) are hyperlinks to main pages of websites, insufficient to find the actual data, hampering any reproducibility effort. From authors' words, it is not clear whether "mortality rate" refer to the number of confirmed COVID-19 deaths divided by the number of inhabitants of each region (i.e. incidence of lethal COVID-19) or to the number of confirmed COVID-19 deaths divided by the number of COVID-19 confirmed cases (i.e. lethality of confirmed COVID-19). The former would be related to the incidence of COVID-19 while the latter is related to the prognosis of COVID-19. Trying to reproduce the results from Table 1, we eventually found that the "mortality rate" refer to the latter and so, is related to the prognosis of COVID-19.

Major statistical flaws

The authors claimed that they have analyzed data of 64,553,275 French citizens to explore the correlation between Sunlight exposure and Mortality. The main finding is a correlation between mortality rate and Sunlight exposure (that is considered as a surrogate marker of Vit-D) with an impressive "p-value of 1.532x10-32) correlated to the COVID-19 mortality rate, with a Pearson coefficient of -0.636".

To explore this correlation, authors used a "Pearson correlation" (that explores association between two quantitative variables), and refer to aggregate data reported by various health agencies and institutions, i.e. Santé Publique France, INSEE, Meteo France). Upon request, the primary author shared on social media [4] the scatterplot (not published in the manuscript) displaying the correlation they studied, suggesting an ecological analysis on data aggregated at regional level, confirmed by our further analysis (see below). Indeed, there are only 12 points in the figure, corresponding to the 12 regional areas reported in Table 1. Exploring such a correlation leads to weak evidence as it is prone to ecological fallacy, a bias that refers to inappropriate conclusions at an individual scale based on the analysis of aggregate data. Moreover, a careful look at the manuscript, shows that the analysis done by the authors is totally flawed. From aggregate regional data provided in Table 1, we found that all correlation coefficients have been computed with 12 data points -one for each region -but that the statistical test has been performed as if there were as many data points as individuals (n=64,553,275). Indeed, the usual Student's t statistic associated to a Pearson's correlation coefficient is t=r×√[(n-2)/(1-r 2 )]. For n=64,553,275 and r=0.663, one can easily compute that t=7115.601 as is displayed on the first line of table 3 for correlation between male sex and COVID-19 mortality rate. Performing such an analysis violates the logic of Pearson's correlation. Indeed, Pearson's correlation is based on the empirical variance observed between the twelve data points that were given. The reduction of variance due to large regional sample sizes is already taken in account in this empirical variance and the Student's t statistic must be based on the sample size that has been used for the empirical variance calculation (i.e. n=12). This first error makes the estimation of statistical error highly biased. The actual error is at least √((64553275-2)/(12-2)) higher, i.e. 2540.7 times higher. After correction of this error, from Table 1 data, the significance level is P=0.03. Moreover, this statistic is based on two other dubious or erroneous assumptions. The first error is assumption of binormal distribution; to avoid this assumption, Spearman's correlation would be better than Pearson's correlation (P=0.04). The second assumption is that the 12 data points are independent. This second assumption is erroneous at least for sunlight exposure, since two neighboring regions have more similar climates than two distant regions. Consequently, even the P-value of Spearman's correlation is erroneous and the actual unbiased P-value would be higher if spatial correlation is taken in account! Furthermore, we are surprised that all regions of metropolitan France were included but Corsica (Corse in French). Authors state that they excluded Corsica "because of poorer access there than on the continent", a statement made without any reference. In contrast, the Insee statistics indicate that the density of healthcare providers per inhabitants is in the same range as on the continent [START_REF]Professionnels de santé au 1ᵉʳ janvier 2018[END_REF]. As the authors point in the manuscript, they worked from a protocol. A time stamped protocol (e.g. registered on the Open Science Framework before any data extraction and analysis) is warranted to demonstrate that exclusion of Corsica was an a priori choice. Indeed, an analysis based on Spearman correlation updated with data from Santé Publique France (Data from April 25 th , 2020, extracted on July 24 th , 2020 [START_REF] Santé_Publique_France | Données hospitalières relatives à l'épidémie de Covid-19 2020[END_REF]) and Méteo France extracted on July the 24 th 2020 [START_REF] Météo_France | Climat France / Ensoleillement[END_REF], show that the correlation observed between intrahospital mortality (cumulative hospital deaths divided by cumulative hospital deaths plus cumulative discharged patients plus number of currently hospitalized patients) and Sunlight without Corsica is not significant anymore when Corsica is added. With our new data extraction, Spearman's R=-0.47 (p=0.11) with Corsica versus R=-0.71 (p=0.01) without Corsica (Figure 1 and Table 1 that includes data gathered from [START_REF] Santé_Publique_France | Données hospitalières relatives à l'épidémie de Covid-19 2020[END_REF][START_REF] Météo_France | Climat France / Ensoleillement[END_REF][START_REF] Ined | Estimation de population au 1er janvier, par région, sexe et grande classe d'âge 2020[END_REF]). Note that our data extraction provided slightly different figures than author's ones, as we could not found the exact same data sources. Therefore, exclusion of Corsica may influence the results. If author cannot prove that this choice was made a priori, this is suspect of p-Hacking. Authors' specify that they had assessed confounding factors. Indeed, with their erroneous analyses, they found that many potential confounding factors have significant effects (see Table 3) on the mortality rate, but no statistical adjustment is performed to cancel the effects of these factors. Moreover, climate confounding, such as temperature have not even been analyzed. Let's now apply reductio ad absurdum. Applying the same methods as authors of this article, we can calculate the association between sunlight exposure and number of senior assisted living beds for 100 000 inhabitants of more than 75 years. With authors' raw data (table 1 of their article) we find a Pearson's R=-0.693 a Student's t statistic equal to -7720.827 and a Pvalue inferior to 10 -30 . After fixing the Student's t calculation, we find a P-value equal to -3.039 and a P-value equal to 0.01. Using updated data, including Corsica, and Spearman's correlation gives a R=-0.69 and a P-value equal to 0.009. Therefore, this method provides better evidence for the fact that sunlight exposure makes people build nursing homes (maybe via increased serological vitamin D in stakeholders building nursing homes, influencing their minds) than for the fact that vitamin D improves the prognosis of COVID-19; indeed, it is robust to the inclusion of Corsica. Another explanation would be the fact that the climate is influenced by nursing homes.

Conclusion

For all these reasons the manuscript has no informative value at all concerning any association between "Covid-19 And Vit-D". Therefore, we think that the article methods and conclusions are too flawed to have any value. Table1: Updated data gathered [START_REF] Santé_Publique_France | Données hospitalières relatives à l'épidémie de Covid-19 2020[END_REF] for COVID-19 data, [START_REF] Météo_France | Climat France / Ensoleillement[END_REF] (Yearly regional climate → Reporting of each station in each region), and [START_REF] Ined | Estimation de population au 1er janvier, par région, sexe et grande classe d'âge 2020[END_REF] for number of male and female inhabitants in 2019.

Figure 1 :

 1 Figure 1: Scatterplot and Spearman's correlations of intra-hospital mortality (number of COVID-19 deaths in hospital divided by number of COVID-19 hospitalizations) and Sunlight exposure (with and without Corsica)


Margaux Lhuissier

Annabelle Couvert

Abdoulaye Kane

Abdeltif Amrane

Jean-Luc Audic

Pierre-Francois Biard

Keywords: Absorption, Hydrodynamics, waste oils, PDMS, packed column, modeling

come

Experimental evaluation and modeling of the hydrodynamics in structured packing operated with viscous waste oils

Introduction

Volatile Organic Compounds (VOCs) have harmful effects on the environment through their global warming contribution. Disturbing the Chapman cycle in the atmosphere by reacting with radicals, VOCs cause tropospheric ozone accumulation and thus intensify global warming [START_REF] Cloirec | COV (composés organiques volatils)[END_REF].

Besides, some VOCs could be toxic towards human health. Their impact on liver, blood and nervous systems have been demonstrated [START_REF] Kampa | Human health effects of air pollution[END_REF].

Several mature technologies can be implemented to remove VOCs from air such as adsorption on activated carbon filters, thermal and catalytic oxidations, absorption (gas-liquid scrubbing) or biofiltration [START_REF] Ruddy | Select the best VOC control strategy[END_REF]. To target hydrophobic VOCs, absorption in a non-aqueous phase (NAP) using silicone oils, phtalates, adipates or even ionic liquids, would be a promising treatment [START_REF] Hadjoudj | Absorption of Chlorinated VOCs in High-Boiling Solvents: Determination of Henry's Law Constants and Infinite Dilution Activity Coefficients[END_REF][START_REF] Heymes | Hydrodynamics and mass transfer in a packed column: Case of toluene absorption with a viscous absorbent[END_REF][START_REF] Bourgois | Solubilities at High Dilution of Toluene, Ethylbenzene, 1,2,4-Trimethylbenzene, and Hexane in Di-2-ethylhexyl, Diisoheptyl, and Diisononyl Phthalates[END_REF][START_REF] Bourgois | Study of mass transfer of VOCs into viscous solvents in a pilot-scale cables-bundle scrubber[END_REF]Darracq et al., 2010a;Darracq et al., 2010b;[START_REF] Guihéneuf | Absorption of Hydrophobic Volatile Organic Compounds in Ionic Liquids and Their Biodegradation in Multiphase Systems[END_REF][START_REF] Biard | A simple and timesaving method for the mass transfer assessment of solvents used in physical absorption[END_REF][START_REF] Biard | Volatile organic compounds absorption in packed column: theoretical assessment of water, DEHA and PDMS 50 as absorbents[END_REF][START_REF] J O U R N A L P R E -P R O O F Rodriguez Castillo | Physicochemical properties of some hydrophobic room-temperature ionic liquids applied to volatile organic compounds biodegradation processes: Physicochemical properties of some hydrophobic room-temperature ionic liquids[END_REF][START_REF] Rodriguez Castillo | Assessment of VOC absorption in hydrophobic ionic liquids: Measurement of partition and diffusion coefficients and simulation of a packed column[END_REF]. Nonetheless, all these commercial or laboratory made solvents are very expensive, which would increase the CAPEX of an industrial process. Consequently, several studies have J o u r n a l P r e -p r o o f recently investigated the potential of cheaper NAP, such as waste oils [START_REF] Bay | Absorption of Volatile Organic Compounds in Biodiesel[END_REF][START_REF] Lalanne | Absorption of a mixture of volatile organic compounds (VOCs) in aqueous solutions of soluble cutting oil[END_REF][START_REF] Ozturk | Absorptive Removal of Volatile Organic Compounds from Flue Gas Streams[END_REF].

In a previous paper investigating the physico-chemical properties (viscosity, volatility, partition coefficient for toluene and dichloromethane) of several industrial waste oils [START_REF] Lhuissier | Characterization and selection of waste oils for the absorption and biodegradation of VOC of different hydrophobicities[END_REF], lubricant and transformer oil were proved to be potential NAP candidates for VOC absorption.

Nevertheless, these waste oils are from 10 to 70 times more viscous than water and have a low surface tension, which can affect significantly their hydrodynamics, especially in irrigated packed columns which are the most commonly used contactors at industrial scale for gas scrubbing (Le [START_REF] Cloirec | COV (composés organiques volatils)[END_REF][START_REF] Minne | Effect of liquid and gas physical properties on the hydrodynamics of packed columns[END_REF].

In a counter-current packed column, it is recommended to operate in the loading zone located between the loading and flooding points [START_REF] Billet | Prediction of mass transfer columns with dumped and arranged packings: updated summary of the calculation method of Billet and Schultes[END_REF]. Indeed, for gas velocities below the loading points (in the pre-loading zone), the downward flow of liquid, and consequently the liquid holdup, depend only on the liquid rate and are not influenced by the gas velocity, leading to low interfacial area and gas-liquid interactions [START_REF] Billet | Prediction of mass transfer columns with dumped and arranged packings: updated summary of the calculation method of Billet and Schultes[END_REF][START_REF] Billet | Fluid dynamics and mass transfer in the total capacity range of packed columns up to the flood point[END_REF], 1993;[START_REF] Heymes | Hydrodynamics and mass transfer in a packed column: Case of toluene absorption with a viscous absorbent[END_REF][START_REF] Mackowiak | Fluid Dynamics of Packed Columns -Principles of the Fluid Dynamic Design of Columns for Gas/Liquid and Liquid/Liquid Systems[END_REF][START_REF] Stichlmair | General model for prediction of pressure drop and capacity of countercurrent gas/liquid packed columns[END_REF]. In this pre-loading zone, the pressure drop increases linearly with the gas superficial velocity (Piché et al., 2001a). In the loading zone, the shear forces in the gas flow support some of the descending liquid allowing to increase advantageously the liquid holdup and the interfacial area [START_REF] Billet | Prediction of mass transfer columns with dumped and arranged packings: updated summary of the calculation method of Billet and Schultes[END_REF]Piché et al., 2001a). In the loading zone, the pressure drop increases with the gas superficial velocity to a power of around 1.8 (Piché et al., 2001a).

Then, when the flooding point is reached, the shear stress of the gas flow is high enough to maintain the liquid at the top of the column, i.e. the liquid starts to overflow and the pressure drop increases sharply [START_REF] Billet | Prediction of mass transfer columns with dumped and arranged packings: updated summary of the calculation method of Billet and Schultes[END_REF][START_REF] Billet | Fluid dynamics and mass transfer in the total capacity range of packed columns up to the flood point[END_REF], 1993;Piché et al., 2001b). Consequently, for a given gas flowrate, the choice of a working point in the loading zone allows subsequently to determine the diameter of the contactor [START_REF] Roustan | Transferts gaz-liquide dans les procédés de traitement des eaux et des effluents gazeux[END_REF]. More precisely, a working gas superficial velocity equal to 70-80% of the flooding velocity is recommended by several authors [START_REF] Billet | Prediction of mass transfer columns with dumped and arranged packings: updated summary of the calculation method of Billet and Schultes[END_REF][START_REF] Maćkowiak | Pressure drop in irrigated packed columns[END_REF], even if a lower boundary of 60% was also proposed [START_REF] Roustan | Transferts gaz-liquide dans les procédés de traitement des eaux et des effluents gazeux[END_REF].

Several semi-empirical correlations have been developed in the literature to determine the loading and flooding points as-well-as liquid holdup and pressure drop in the pre-loading and loading zones of either random and structured packing [START_REF] Maćkowiak | Determination of flooding gas velocity and liquid hold-up at flooding in packed columns for gas/liquid systems -Maćkowiak -1990[END_REF][START_REF] Maćkowiak | Pressure drop in irrigated packed columns[END_REF][START_REF] Miyahara | Flow dynamics in low height packed column having large fractional void space[END_REF][START_REF] Rocha | Distillation columns containing structured packings: a comprehensive model for their performance. 1. Hydraulic models[END_REF] J o u r n a l P r e -p r o o f 1993; [START_REF] Billet | Prediction of mass transfer columns with dumped and arranged packings: updated summary of the calculation method of Billet and Schultes[END_REF]S. Piché et al., 2001aS. Piché et al., , 2001bS. Piché et al., , 2001c;;[START_REF] Heymes | Hydrodynamics and mass transfer in a packed column: Case of toluene absorption with a viscous absorbent[END_REF]. These correlations take into account the column diameter, packing characteristics (porosity, particle diameter, specific surface area, etc.), liquid and gas properties (viscosity, density, superficial tension) and liquidto-gas mass flowrate ratio L/G. Besides, the pressure drop in the loading zone depends also on the selected liquid and gas velocities. Most of these correlations were established using data gathered with water-like solvents having low or moderate viscosities. The Billet-Schultes correlations were developed in the 90 ties for both random and structured packing with solvents having dynamic velocities up to around 100 mPa s (Billet andSchultes, 1999, 1993). In 2006, Heymes et al. showed that these correlations were the most accurate to predict the pressure drop of dumped packing using DEHA. In 2017, Minne also proved the good predictive capacity of these correlations using solvents with different viscosities, densities and surface tension (ethylene glycol, silicone oil, water and Isopar G).

The procedure developed by Billet-Schultes involves the calculation of liquid holdup, which is the ratio of the liquid within the column volume to the packing volume.

Only a few studies have been performed to assess the hydrodynamics of viscous solvents in packed columns, and up to now, no one has been dedicated to waste oils [START_REF] Brunazzi | Hydrodynamics of a gas-liquid column equipped with MellapakPlus packing[END_REF][START_REF] Darracq | Couplage de l'absorption dans une phase organique et de la biodégradation dans un réacteur multiphasique[END_REF][START_REF] Guillerm | Absorption of toluene in silicone oil: Effect of the solvent viscosity on hydrodynamics and mass transfer[END_REF][START_REF] Heymes | Hydrodynamics and mass transfer in a packed column: Case of toluene absorption with a viscous absorbent[END_REF][START_REF] Minne | Effect of liquid and gas physical properties on the hydrodynamics of packed columns[END_REF][START_REF] Minne | The effect of fluid properties and packing size on the hydrodynamics of packed columns[END_REF][START_REF] Tsai | Influence of viscosity and surface tension on the effective mass transfer area of structured packing[END_REF].

Thus, the purpose of this work is to study the hydrodynamics of two waste oils, a lubricant and a transformer oil, in a laboratory scale packed column, by comparison to another reference solvent: a commercial silicone oil (PDMS 20). A modern 4 th generation structured packing made of metallic corrugated sheets was used (Flexipac® 500Z HC). Both the loading and flooding points for several liquid-to-gas mass flowrate ratios (L/G) have been determined, as well as the pressure drop evolution with the gas and liquid velocities. Besides, the reliability of the Billet-Schultes correlations was evaluated allowing to subsequently simulate different conditions and to assess the scale-up of the process at industrial scale.

Material and method

Chemical products

J o u r n a l P r e -p r o o f

The silicone oil used was composed of polydimethylsiloxane (PDMS), a synthetic oil named Rhodorsil 47V20 provided by Bluestar Silicones (Table 1). The waste oils were selected among those collected by Chimirec (Javené, France).

Pilot-scale packed column

The irrigated packed column (120 mmID (Dcol) and packing height Z of 1 m) was operated at counter current (Fig. 1). The air flow (G = 1.17 kg m -3 at 25°C and 1 atm) was introduced at the bottom of the column and the liquid solvent was introduced at the top of the column. The column was filled with structured packing (Flexipac®, provided by Koch Glitsch-USA) made out of corrugated sheets arranged in a crisscrossing configuration to create flow channels for the gas phase (Table 2).

The air flowrate was regulated by means of a membrane valve located after a fan and measured by a rotameter (GF Type SK 20 CH-8201, Switzerland). The column was fed with the liquid absorbent by means of a centrifugal pump (Iwaki MD100, Japan). The liquid flowrates were regulated with a valve and measured by a previously calibrated rotameter (GF type SK 11 CH-8201, Switzerland). The experimental conditions and the range of operating parameters are shown in Table 3. For each liquid flowrate selected, the gas flow was increased incrementally up to the flooding point (determined by visual observation), and pressure drops (ΔP) were measured three times for each point using a vertical U-shaped tube filled with water. Only the pressure drop corresponding to the packing was measured. The temperature of the liquid phase ranged from 294 to 298 K, and the gas temperature was kept constant at 298 K by means of a heat exchanger. Considering all the experiments, the difference between the inlet temperatures of the gas and liquid phases was always lower than 3 K.

The accuracy and reliability of this set-up was previously checked by Darracq and Guillerm (Darracq, 2011;[START_REF] J O U R N A L P R E -P R O O F Guillerm | Optimisation du couplage de l'absorption et de la biodégradation pour l'élimination de Composés Organiques Volatils hydrophobes[END_REF].

Billet-Schultes correlations Billet and Schultes (B-S) developed a set of correlations allowing to determine (i) both the loading and flooding gas superficial velocities (UG,Lo and UG,Fl) for a given L/G ratio and (ii) the pressure drop, liquid holdup and the interfacial area for the selected working gas superficial velocity (UG) in both the preloading and loading zones [START_REF] Billet | Prediction of mass transfer columns with dumped and arranged packings: updated summary of the calculation method of Billet and Schultes[END_REF][START_REF] Billet | Fluid dynamics and mass transfer in the total capacity range of packed columns up to the flood point[END_REF], 1993[START_REF] Maćkowiak | Pressure drop in irrigated packed columns[END_REF]. These correlations, presented in Table 4, involve several constants which depend specifically on the packing (CP, CLo, CFl and Ch). The structured packing used in this study was not previously tested by Billet-Schultes.

J o u r n a l P r e -p r o o f

Therefore, the constants used in our calculations were deduced from the experimental results (part 3.3). For a given packing, the loading and flooding points depend on the selected L/G ratio, the liquid and gas properties and the packing properties. They can be deduced from the equations 1-7 and 8-11, respectively. These sets of equations need to be solved by iteration (for example using the Excel® Solver). The determination of the loading points involved the determination of the ratio of the hydraulic to the geometric surface area (ah/a) according to Eqs. 6 or 7. This ratio must not be confused with the ratio of the specific interface between the phases to the geometric area (a°/a), determined with Eqs.

18-20. Besides, these equations involve several resistance coefficients () whose calculation depends on the nature of the dispersed phase (Eqs. 3-4 and Eqs. 10-11).

Results and Discussion

Experimental pressure drop of the irrigated column

The experimental pressure drop of the irrigated packed column increases logically with both the liquid and gas superficial velocities according to Fig. 2. The gas velocity can be easily converted to the vapor capacity factor FV, often used to report data dealing with the hydrodynamics of an irrigated packed column [START_REF] Minne | Effect of liquid and gas physical properties on the hydrodynamics of packed columns[END_REF], by multiplication with the square root of the air density (which was almost constant in this study and equal to 1.17 kg m -3 ). Lower liquid velocities were applied for the lubricant since the flooding was reached at a lower L/G ratio for a given gas velocity than for silicon and transformer oils. This observation is consistent with the fact that both loading and flooding points appear at lower gas velocities for viscous solvents at similar liquid densities [START_REF] Minne | Effect of liquid and gas physical properties on the hydrodynamics of packed columns[END_REF].

Nonetheless, even with a high viscosity, it would be feasible to use this lubricant in a packed column.

In Fig. 2, experimental points located in both pre-loading and loading zones are represented. Some values of the pressure drop were also measured in the flooding zone but are not represented in Fig. 2 for clarity. In the pre-loading zone, the pressure drop increases linearly with the gas velocity and is minimally influenced by the nature and the flowrate of the liquid. [START_REF] Minne | Effect of liquid and gas physical properties on the hydrodynamics of packed columns[END_REF] and [START_REF] Brunazzi | Hydrodynamics of a gas-liquid column equipped with MellapakPlus packing[END_REF] also observed this behavior for various liquids with similar densities but very different viscosities for random and structured packing, respectively. It is consistent with the fact that gas and liquid have almost no interaction in the pre-loading zone. However, this behavior completely changes J o u r n a l P r e -p r o o f when the loading point is reached. Indeed, the pressure drop at given liquid and gas velocities is significantly higher for the lubricant in the loading zone than for the transformer oil and the silicon oil.

Furthermore, the pressure drop increases with the gas velocity to the power of 1.8. All these observations in the loading zone are in agreement with the literature [START_REF] Minne | Effect of liquid and gas physical properties on the hydrodynamics of packed columns[END_REF][START_REF] Minne | The effect of fluid properties and packing size on the hydrodynamics of packed columns[END_REF]Piché et al., 2001a). The results emphasize that the pressure drop values stay advantageously under 450 Pa.m -1 in both the pre-loading and loading zones for each of the viscous liquid phases, showing that the implementation of viscous solvents would not affect significantly the operating costs (OPEX).

However, the results clearly show that lower L/G ratios must be selected when the solvent viscosity increases to avoid loading and flooding points that would be too low, highlighting that the determination of the loading and flooding points would be determinant when designing a column. To reach this goal, robust models must be implemented.

Determination of the experimental loading and flooding points

The log-log plots corresponding to Fig. 2 are provided as supplementary material (Fig S .1). For the three solvents, the experimental loading and flooding superficial velocities, summarized in Table 5 andTable 6, were determined from the break-up of the slopes between the pre-loading zone, the loading zone and the flooding zone on these log-log plots.

The loading and flooding points for a given liquid-to-gas mass flowrate ratio (L/G) were close to the results obtained previously by [START_REF] Guillerm | Absorption of toluene in silicone oil: Effect of the solvent viscosity on hydrodynamics and mass transfer[END_REF] with the same packing using PDMS 50 (viscosity of 50 mPa.s). The observed loading and flooding velocities are relatively low, in the range 0.39-0.65 and 0.56-1.06 m s -1 , respectively. This behavior is justified by both the nature of the packing, which is a structured packing with a high interfacial area leading to a high liquid capacity, and the high viscosity of the solvents investigated. Indeed, in agreement with the experimental observations and theoretical considerations of other authors [START_REF] Billet | Fluid dynamics and mass transfer in the total capacity range of packed columns up to the flood point[END_REF][START_REF] Brunazzi | Hydrodynamics of a gas-liquid column equipped with MellapakPlus packing[END_REF][START_REF] Minne | The effect of fluid properties and packing size on the hydrodynamics of packed columns[END_REF], both loading and flooding velocities decreased with the solvent viscosity because of higher frictional forces. Thus, lower liquid loads were operated for the most viscous solvent, i.e. the lubricant. Nonetheless, having similar J o u r n a l P r e -p r o o f viscosities, densities and surface tensions, the loading and flooding velocities for the transformer oil and the silicon oil are not significantly different at a given liquid load. On the one hand, the influence of the density on the hydrodynamics is rather well described in the literature. Indeed, the gravitational forces are counterbalanced by the frictional forces at lower velocities when the density decreases [START_REF] Minne | Effect of liquid and gas physical properties on the hydrodynamics of packed columns[END_REF]. On the other hand, the possible influence of the surface tension is still controversial in the literature and seems to depend on the nature of the packing and of the liquid load [START_REF] Minne | Effect of liquid and gas physical properties on the hydrodynamics of packed columns[END_REF].

Anyway, the Billet-Schultes model, applied in the next part, does not take the influence of the surface tension into account.

Application of the Billet-Schultes: determination of the specific constants for Flexipac® 500Z HC packing.

Eqs 8-11 show that the determination of the flooding point at a given L/G ratio depends on the properties of the liquid phase (viscosity  and density ) and on the properties and the nature of the packing through its void fraction , specific area a and through a specific constant CFl. By trying to minimize the sum of the relative errors (RE) between the experimental gas flooding velocities and those deduced from the B-S correlations, the value of CFl was determined (Table 7). The value found is in agreement with the range previously determined by Billet and Schultes [START_REF] Billet | Prediction of mass transfer columns with dumped and arranged packings: updated summary of the calculation method of Billet and Schultes[END_REF]. The low average relative error (ARE) of 7.5% (calculated for 11 values) between the experimental and theoretical gas flooding velocities (summarized on Table 6) shows that B-S correlations are efficient for the prediction of the flooding points even for solvents with very different properties. [START_REF] Minne | Effect of liquid and gas physical properties on the hydrodynamics of packed columns[END_REF] found similar RE for a silicon oils using two dumped packing. The parity plot is provided as supplementary material and does not put in evidence any systematic error (Fig. S-2(b)).

From Eqs 8-11, the liquid holdup at the flooding points was also determined (Table 6).

Eqs. 1-7 show that the determination of the loading point at a given L/G ratio depends on the properties of the liquid phase (viscosity  and density ) and on the properties and the nature of the packing through its void fraction , specific area a and through two specific constant CLo and Ch.

Several couples (CLo ;Ch) allowed to minimize the ARE between the experimental gas loading velocities and those deduced from the B-S correlations. Thus, the best couple summarized on Table 7 was determined through the simultaneous minimization of the ARE of both the experimental and theoretical pressure drops and loading points. Only the pressure drops measured in the loading zone J o u r n a l P r e -p r o o f were considered. Indeed, the pressure drop in the pre-loading zone is low and suffers from higher experimental uncertainties. Besides, Minne (2017) pointed out the relatively low performance of the Billet-Schultes model to predict the pressure drop on the pre-loading zone.

The values of Ch and CLo found are in agreement with the ranges previously determined by Billet and Schultes [START_REF] Billet | Prediction of mass transfer columns with dumped and arranged packings: updated summary of the calculation method of Billet and Schultes[END_REF], even if the value of CLo is slightly lower than the range of Table 7.

The low ARE of 8.0% (calculated for 10 values) between the experimental and theoretical gas loading velocities (summarized on Table 5) shows that the B-S correlations are efficient for the prediction of the loading points, even for solvents with very different properties. The parity plot is provided as supplementary material and shows that the model slightly overestimates the loading points for the lubricant and slightly underestimates them for the transformer oil (Fig. S-2(a)), but these discrepancies remain acceptable.

The results summarized in Table 5 and Table 6 clearly confirm that the liquid holdup would be higher at a given L/G ratio for the lubricant, i.e. for a higher viscosity, for both the loading and flooding points.

Thus, even if viscous solvents must be operated at a lower superficial gas velocity and/or a lower L/G ratio (i.e. lower absorption rate [START_REF] Biard | Volatile organic compounds absorption in packed column: theoretical assessment of water, DEHA and PDMS 50 as absorbents[END_REF]), a higher interfacial area would be expected.

Finally, the value of Cp (Table 7) was determined through the minimization of the ARE between the experimental and theoretical (Eqs 12-17) pressure drops of the points located in the loading zones.

The value of Cp found is in agreement with the range previously determined by Billet and Schultes [START_REF] Billet | Prediction of mass transfer columns with dumped and arranged packings: updated summary of the calculation method of Billet and Schultes[END_REF]. The low average relative error of 15.1% (calculated for 134 values: 42 for PDMS 20, 32 for transformer oil, 60 for lubricant) between the experimental and theoretical pressure drops in the loading zone shows that B-S correlations estimate with a good confidence level the pressure drop even with viscous solvents. On the pre-loading zone, the agreement is lower, with an ARE of 34.1%, which is consistent with the observations of [START_REF] Minne | The effect of fluid properties and packing size on the hydrodynamics of packed columns[END_REF]. Nonetheless, it is not recommended to operate on the pre-loading zone [START_REF] Billet | Prediction of mass transfer columns with dumped and arranged packings: updated summary of the calculation method of Billet and Schultes[END_REF]. The agreement between the experimental and theoretical pressure drops on the pre-loading and loading zones can be assessed Fig. 2. The detailed ARE for the three solvents in each zone are summarized in Table 8. It shows that the best agreement is found for PDMS 20 and transformer oil. For the lubricant, the ARE is higher, around 27%, but remains fully acceptable. Furthermore, the parity plot of the pressure drop is provided as supplementary material and does not put in evidence any systematic error (Fig. S-3).

J o u r n a l P r e -p r o o f

Discussion and simulations

From the constants determined previously (Table 7) and the B-S correlations, it was possible to estimate the loading and flooding points for the two waste oils for an increasing L/G ratio in the 0.5-15 range (Fig. 3).

The loading and flooding velocities obviously decrease with the L/G ratio. For the lubricant (i.e. the most viscous solvent), the loading velocity is from 6% (at L/G = 0.5) to 20% higher (at L/G = 15) and the flooding velocity is from 10% (at L/G = 0.5) to 30% higher (at L/G = 15) than for the transformer oil.

Because of almost identical properties, the data for PDMS 20 are similar to those of the transformer oil and are not represented. For the sake of comparison, a simulation for water was also carried out (Fig. S-4). The results show that for this solvent, the loading and flooding points would be respectively from 17 to 52% higher and from 25 to 48% higher than for the transformer oil. Fig. 3 shows that the choice of a working velocity at 80% of the flooding velocity allows operating in the loading zone in any case. However, the choice of a working velocity at 60% of the flooding velocity, proposed as a lower boundary by some authors [START_REF] Roustan | Transferts gaz-liquide dans les procédés de traitement des eaux et des effluents gazeux[END_REF], would be inappropriate, especially at high L/G ratio and for the lubricant. The breakthrough of the curves observed at a L/G ratio around eleven is due to a shift of the dispersed phase (from the liquid to the gas phase) corresponding to a value of 𝐿 𝐺 × √ 𝜌 𝐺 𝜌 𝐿 higher than 0.4 (Table 4). Fig. 3 also highlights that the selection of an excessive L/G ratio would be inappropriate for several reasons. Indeed, for the viscous solvents investigated, it would impose to select a low superficial gas velocity, i.e. to increase the column diameter for a given gas flowrate to treat. Moreover, a low superficial gas velocity would affect the turbulences in the gas phase, decreasing the gas-phase mass transfer rate. Since a significant part of the mass transfer resistance can be located in the gas phase during hydrophobic VOC absorption, even using viscous solvents [START_REF] Biard | Volatile organic compounds absorption in packed column: theoretical assessment of water, DEHA and PDMS 50 as absorbents[END_REF][START_REF] Rodriguez Castillo | Assessment of VOC absorption in hydrophobic ionic liquids: Measurement of partition and diffusion coefficients and simulation of a packed column[END_REF], would decrease the overall mass transfer rate. Finally, for a high L/G ratio, the loading zone is narrower. Thus, low variations of the gas velocity may induce significant variations of the hydrodynamics within the column, which would complicate the monitoring of a packed column at industrial scale.

J o u r n a l P r e -p r o o f

The evolution of the pressure drop DP/Z (orange curves), the liquid holdup hL (yellow curves) and the ratio a°/a was simulated for several liquid loads (Fig. 4). hL and the ratio a°/a remain constant in the pre-loading zone in agreement with the experimental observations of several authors with structured packing [START_REF] Billet | A physical model for the prediction of liquid hold-up in twophase countercurrent columns[END_REF][START_REF] Brunazzi | Hydrodynamics of a gas-liquid column equipped with MellapakPlus packing[END_REF][START_REF] Zakeri | Experimental investigation of liquid holdup in structured packings[END_REF]. In the pre-loading zone, hL and the ratio a°/a increase by 100% for the transformer oil when the liquid load increases by a factor of 4 (from 5 to 20 m h -1 ). Moreover, hL and the ratio a°/a are around 30% higher for the lubricant in the pre-loading zone. [START_REF] Brunazzi | Hydrodynamics of a gas-liquid column equipped with MellapakPlus packing[END_REF] also showed that the liquid holdup increases with the solvent viscosity even in the pre-loading zone. Thus, a viscous solvent exhibits advantageously a higher interfacial area although the pressure drop is almost insensitive to the nature of the solvent, at least in the pre-loading zone.

When the loading zone is reached (around 0.60 m s -1 for UL = 5 m h -1 and around 0.45 m s -1 for UL = 20 m h -1 for the transformer oil; around 0.55 m s -1 for UL = 5 m h -1 for the lubricant), the gas and the liquid start to interact. Thus, both hL and the ratio a°/a increase until the flooding point is reached (around 1.20 m s -1 for UL = 5 m h -1 and around 0.78 m s -1 for UL = 20 m h -1 for transformer oil; around 1.02 m s -1 for UL = 5 m h -1 for the lubricant). A rather high liquid holdup, from 25 to more than 45%, would be reached at the flooding point. High liquid hold ups are concomitant with a good wetting of the column and with high interfacial area a°. The interfacial area can be even higher that the geometrical area a of the packing according to the experimental observations of [START_REF] Tsai | Influence of viscosity and surface tension on the effective mass transfer area of structured packing[END_REF] on similar structured packing fed with solvents having viscosities up to 14 mPa s. [START_REF] Cherif | Comparison of models for the prediction of hydrodynamic parameters in structured packing columns for biogas purification[END_REF] showed the good predictive capacity for a° of the Billet-Schultes correlations for a structured packing. In the loading zone, contrarily to the pre-loading zone, the pressure drop is particularly sensitive to the liquid load and the viscosity of the solvent. Thus, for example, at a liquid load of 5 m h -1 and a gas velocity of 0.8 m s -1 , the interfacial area and the pressure drop for the lubricant would be around 40% and 15% higher than for the transformer oil, respectively.

Packing scale-up

J o u r n a l P r e -p r o o f

According to part 3.3, the Billet-Schultes correlations provide fair estimations of the pressure drop and loading and flooding points. These parameters are crucial for the evaluation of the operating cost through the determination of the packed column diameter and the selection of the L/G ratio (part 3.4).

In this part, in order to study the scale-up of a column with this specific packing, a realistic case involving a flowrate of 4000 Nm 3 .h -1 of air to be treated was considered (where N stands for the standard temperature and pressure (STP) conditions, meaning 1 bar and 0°C according to the IUPAC). The transformer oil was selected among the three solvents considering the observations of the parts 3.3 and 3.4. [START_REF] Biard | Volatile organic compounds absorption in packed column: theoretical assessment of water, DEHA and PDMS 50 as absorbents[END_REF][START_REF] Rodriguez Castillo | Assessment of VOC absorption in hydrophobic ionic liquids: Measurement of partition and diffusion coefficients and simulation of a packed column[END_REF] already studied the hydrodynamics of DEHA (12.5 mPa s), PDMS (50 mPa s) and two ionic liquids (50-60 mPa s) in a packed column (Dcol = 1 m) filled with bulk Pall® rings for with the same gas flowrate using the Billet-Schultes correlations. Apart from a different packing nature (corrugated sheets vs. random rings), Pall® rings offer a geometric surface area of only 139.4 m 2 m -3 , around 3 times lower than Flexipac®500 Z HC. Thus, both the computed loading and flooding superficial velocities for similar L/G ratio (in the range 1-5) are significantly lower with the Flexipac® packing. Thus, it was impossible to select the same working velocity of 1.52 m s -1 as Biard et al. [START_REF] Biard | Volatile organic compounds absorption in packed column: theoretical assessment of water, DEHA and PDMS 50 as absorbents[END_REF][START_REF] Rodriguez Castillo | Assessment of VOC absorption in hydrophobic ionic liquids: Measurement of partition and diffusion coefficients and simulation of a packed column[END_REF] since this one is higher than the computed flooding velocity for the Flexipac® packing whatever the L/G ratio considered.

Thus, for this scale-up assessment, a L/G ratio in the range 1-5 to avoid flooding velocities that were too low (Fig. 3 (a)) was selected. Then, a working gas velocity taken at 80% of the flooding velocity as recommended in the literature was considered [START_REF] Billet | Prediction of mass transfer columns with dumped and arranged packings: updated summary of the calculation method of Billet and Schultes[END_REF]. A working gas velocity (which decreases with the L/G ratio), in the range 0.54-0.86 m s -1 , and a column diameter (which increases with the L/G ratio), in the range 1.33-1.68 m, were calculated. Fig. 5 presents the evolution of the liquid holdup hL, interfacial area a° and pressure drop DP/Z calculated with the B-S correlations for increasing L/G ratios.

The calculated pressure drops would be in a reasonable range (95-175 Pa m -1 ) for an industrial application. The pressure drop would be lower at industrial scale at given liquid and gas velocities than at the lab-scale, because of negligible wall effects. These wall effects are taken into account by the Billet-Schultes correlations through the column diameter . The pressure drop decreases with the L/G ratio and with the working gas velocity UG, even if the liquid holdup and the liquid load J o u r n a l P r e -p r o o f increase. The low values of the pressure drop must be balanced with the low gas velocity tolerated by this packing. It confirms that viscous solvents do not necessarily lead to high pressure drop (Biard et al., 2018). Advantageously, the computed interfacial area is significant, higher than 200 m 2 m -3 . These results emphasize the good hydrodynamic performance offered by this packing. Nonetheless, the use of a low gas velocity would be detrimental for the mass transfer rate in the gas phase, in which the main part of the resistance can be located, especially during the absorption of hydrophobic VOCs for which the selection of viscous absorbents is justified [START_REF] Rodriguez Castillo | Assessment of VOC absorption in hydrophobic ionic liquids: Measurement of partition and diffusion coefficients and simulation of a packed column[END_REF].

Conclusion

The purpose of this work was to study the hydrodynamic behavior of two viscous waste oils (a transformer oil and a lubricant characterized by viscosities of 19 mPa s and 79 mPa s, respectively) and a silicone oil (20 mPa s) in a laboratory-scale packed column (Dcol = 0.12 m) filled with structured packing (Flexipac® 500Z HC). First, the results showed that it was possible to successfully apply such viscous solvents in structured packing made of corrugated sheets. The loading and flooding points were determined for L/G ratios between 0.3 and 14.0. Pressure drop in the loading zone values were in the range from 50 to 450 Pa.m -1 . The pressure drop values were significantly higher in the loading zone at the same gas and liquid velocities for the most viscous solvent (lubricant) than for the other solvents. Besides, lower liquid loads (i.e. liquid velocities) must be applied for the lubricant, since the flooding point was reached at lower gas velocities for a given L/G ratio. Considering that the high lubricant viscosity will also lower the mass transfer rate (directly, through lower turbulences, and indirectly, through lower diffusion coefficients [START_REF] Rodriguez Castillo | Assessment of VOC absorption in hydrophobic ionic liquids: Measurement of partition and diffusion coefficients and simulation of a packed column[END_REF]), the selection of such a viscous solvent would be justified only for pollutant/solvent systems characterized by low partition coefficient (i.e. high affinities).

In order to apply waste oils in a scrubbing process, the accurate prediction of the hydrodynamic performance is a prerequisite. The Billet-Schultes correlations [START_REF] Billet | Prediction of mass transfer columns with dumped and arranged packings: updated summary of the calculation method of Billet and Schultes[END_REF] were used to predict the loading points and flooding points, as well as the pressure drop, the liquid holdup and the interfacial area at the working velocity. The specific constants of the packing used were determined by numerical resolution. The ARE between the experimental and theoretical pressure drops in the loading J o u r n a l P r e -p r o o f zone was around 15%. Furthermore, the Billet-Schultes correlations successfully determined the loading and flooding points with ARE around 7-8%. Several simulations were performed with these correlations to assess the influence of the operating conditions. It highlighted the high hydrodynamic performance of this packing. It also demonstrated that the use of viscous solvents does not necessarily lead to excessive pressure drop in agreement with the experimental results. A fortiori, high liquid holdup and interfacial area would be advantageously obtained. Nonetheless, low loading and flooding velocities were computed, even for a low L/G ratio, which leads to select low working gas velocities, leading to high diameter columns and lowered mass transfer rate in the gas phase. Finally, the scale-up of a realistic packed column for gas treatment was assessed and showed that that it would be possible to apply these kinds of solvents with the Flexipac®500Z HC packing even at industrial scale. For a complete techno-economic assessment, the determination of the mass transfer coefficients in both the liquid and gas phases will be necessary, allowing to evaluate concomitantly the pressure drop, the overall mass transfer coefficient and thus, the resulting removal efficiency in a given packed column [START_REF] Biard | Volatile organic compounds absorption in packed column: theoretical assessment of water, DEHA and PDMS 50 as absorbents[END_REF][START_REF] Rodriguez Castillo | Assessment of VOC absorption in hydrophobic ionic liquids: Measurement of partition and diffusion coefficients and simulation of a packed column[END_REF].

U L (m h -1 ) U L (m h -1 ) U L (m h -1 )
J o u r n a l P r e -p r o o f 0.021 0.027 0.031 * Lubricant and transformer oil have complex chemical structures. Several samples were investigated by 1 H-NMR spectroscopy and showed quite similar chemical natures between lubricant and transformer oil and allowed to roughly estimate an average chain length of 15 carbons which is coherent with literature data for mineral oils [START_REF] Aluyor | Biodegradation of mineral oils-A review[END_REF].

P T D Z FG FL UG UL bar °C m m Nm 3 h -1 L h -1 m s -1 m s -1 1
25 0.12 1 9 -50 15 -750 0.2 -1.2 3.10 -4 -2.10 -2

J o u r n a l P r e -p r o o f Eq. 8

Eq. 9

Eq. 10

Eq. 11 Step 1: hL,Lo from Eqs. 5-7 and calculated with UL

Step 2: ℎ 𝐿,𝐹𝑙 = 2.2ℎ Eq. 15

Eq. 16

Eq. 17 Eq. 18

Eq. 19

Eq. 20

  Margaux Lhuissier a , Annabelle Couvert a , Abdoulaye Kane b , Abdeltif Amrane a , Jean-Luc Audic a , Pierre-François Biard a a Univ Rennes, Ecole Nationale Supérieure de Chimie de Rennes, CNRS, ISCR -UMR 6226, F-35000 Rennes, France b UniLaSalle-Ecole des Métiers de l'Environnement, Campus de Ker Lann, 35170 Rennes, France r n a l P r e -p r o o f Nomenclature a: specific surface area of packing (m 2 m -3 ) ah : hydraulic surface area (Eq. 6 and Eq. 7, m 2 m -3 ) a°: interfacial area relative to the packing volume (m 2 m -3 ) ARE: Average Relative Error CP, CLo, CFl, Ch : constants relative to each commercial packing according to Billet-Schultes dh: hydraulic diameter = 4/a (m) Dcol : column diameter (m) F: volume flowrate (m 3 h -1 or m 3 s -1 ) g: specific gravity constant (9.81 m s -2 ) G: gas mass flowrate (kg s -1 ) hL: liquid holdup (-) L: liquid mass flowrate (kg s -1 ) Re: Reynolds number (-) U: superficial velocity (m s -1 ) Z: packing height (m) RE: Relative error Greek letters ΔP/Z : Linear pressure drop (Pa m -1 ) : Packing void fraction (-) J o u r n a l P r e -p r o o f : resistance coefficient (-)  dynamic viscosity (Pa s)  density (kg.m -3 ) Subscripts L: Relative to the liquid phase G: Relative to the gas phase Lo: At the loading point Fl: At the flooding point

Figure 1 :

 1 Figure 1: Experimental set-up.

Figure 2 :

 2 Figure 2 : Evolution of the pressure drop with the gas velocity for PDMS (a), transformer oil (b) and lubricant (c) for increasing values of the liquid load UL (in m 3 m -2 h -1 ). The straight lines correspond to the calculations according to the Billet-Schultes correlations (part 3.3).

Figure 3 :

 3 Figure 3: Evolution of the loading and flooding gas velocity for transformer oil (a) and lubricant (b) for increasing values of the L/G ratio (Flexipac® 500Z HC). The green frame corresponds to the loading zone (located between the loading and the flooding velocities). The orange and blue lines correspond to 80%and 60% of the flooding velocity, respectively. These data were simulated using the B-S correlations.

Figure 4 :Figure 5 :

 45 Figure 4: Evolution of the pressure drop (orange curves), the liquid holdup (yellow curves) and the ratio a°/a (green curves) for increasing values of the gas velocity for transformer oil (a) and lubricant (b) (Flexipac® 500Z HC). Two liquid loads (5 and 20 m h -1 ) were simulated for the transformer oil and one liquid load (5 m h -1 ) was simulated for the lubricant. These data were simulated using the B-S correlations.

  Determination of the interfacial area (a°) at the working point J o u r n a l P r e -p r o o f

Table 1 : Properties of the studied liquid absorbents. PDMS Transformer oil Lubricant Chemical composition

 1

		Siloxanes Linear alkanes
	Dynamic viscosity L (mPa s) at 25°C 20	19	79
	Molar mass (g mol -1 )	3000	212*
	Density L (kg m -3 )	900	865	875
	Superficial tensionL (N m -1 )

Table 2 : Packing characteristics

 2

	Name	dh (m)	a (m -1 ) ε (-)
	Flexipac® 500Z HC 7.60.10 -3 500	0.95

Table 3 : Operating conditions

 3

Table 4 : Billet-Schultes correlations developed for the determination of the parameters relative to the hydrodynamics of an irrigated counter-current packed column (Billet and Schultes, 1999, 1995, 1993, 1991). Determination of the loading point (UG,Lo) for a given L/G ratio

 4 Determination

	𝑈 𝐺,𝐿𝑜 = √ 𝜓 𝐿𝑜 𝑔		× (𝜀 -ℎ 𝐿,𝐿𝑜 ) × √	ℎ 𝐿,𝐿𝑜 𝑎	× √	𝜌 𝐺 𝜌 𝐿	Eq. 1
	𝑈 𝐿,𝐿𝑜 =	𝜌 𝐺 𝜌 𝐿	×	𝐿 𝐺	× 𝑈 𝐺,𝐿𝑜	Eq. 2
	With:
	𝜓 𝐿𝑜 =	𝑔 𝐶 𝐿𝑜	2 × ( 𝐺. √ 𝐿	𝜌 𝐺 𝜌 𝐿	(	𝜇 𝐿 𝜇 𝐺	)	0.4	) 0,752	if	𝐿 𝐺	× √	𝜌 𝐺 𝜌 𝐿	≤ 0.4 (liquid dispersed in the gas)	Eq. 3
	𝜓 𝐿𝑜 =	𝑔 (0.695𝐶 𝐿𝑜 ( 𝜇 𝐿 𝜇 𝐺	)	0.1588	) 2 × ( 𝐺 √ 𝐿	𝜌 𝐺 𝜌 𝐿	(	𝜇 𝐿 𝜇 𝐺	)	0.4	) 1,46	if	𝐿 𝐺	× √ 𝜌 𝐺 𝜌 𝐿	> 0.4 (gas dispersed in the	Eq. 4
	liquid)
	And: ℎ 𝐿,𝐿𝑜 = ( 12 𝑔	𝜇 𝐿 𝜌 𝐿	𝑎 2 𝑈 𝐿,𝐿𝑜 )	1 3 ⁄	× (	𝑎 𝑎 ℎ	)	2 3 ⁄	Eq. 5
	With: ( 𝑎 ℎ 𝑎 ) = 𝐶 ℎ × ( ( 𝑎 ℎ 𝑎 ) = 0.85 × 𝐶 ℎ × ( 𝑈 𝐿,𝐿𝑜 .𝜌 𝐿 𝑎.µ 𝐿 ) 𝑈 𝐿,𝐿𝑜 .𝜌 𝐿 0.15 × ( 𝑎.µ 𝐿 ) 0.25 𝑈 𝐿,𝐿𝑜 𝑔 × ( 2 .𝑎 ) 𝑈 𝐿,𝐿𝑜 0.1 𝑔	if ( 2 .𝑎 ) 𝑈 𝐿,𝐿𝑜 .𝜌 𝐿 𝑎.µ 𝐿 0.1 if ( ) < 5 𝑈 𝐿,𝐿𝑜 .𝜌 𝐿 𝑎.µ 𝐿	) ≥ 5	Eq. 6 Eq. 7
	ℎ 𝐿,𝐹𝑙 3 . (3. ℎ 𝐿,𝐹𝑙 -𝜀) =	6 𝑔	. 𝑎 2 . 𝜀.	𝜇 𝐿 𝜌 𝐿	.	𝐿 𝐺	.	𝜌 𝐺 𝜌 𝐿	. 𝑈 𝐺,𝐹𝑙
	With:
	𝑈 𝐺,𝐹𝑙 = √	2. 𝑔 𝜓 𝐹𝑙	×	(𝜀 -ℎ 𝐿,𝐹𝑙 ) 1.5 𝜀 0.5	× √ ℎ 𝐿,𝐹𝑙 𝑎	× √	𝜌 𝐿 𝜌 𝐺
	And:
	𝜓 𝐹𝑙 =	𝑔 𝐶 𝐹𝑙	2 . ( 𝐺 √ 𝐿 𝜌 𝐺 𝜌 𝐿		(	𝜇 𝐿 𝜇 𝐺	)	0.2	) 0,388	if	𝐿 𝐺	× √ 𝜌 𝐺 𝜌 𝐿	≤ 0.4
	𝜓 𝐹𝑙 =	𝑔 (0.6244𝐶 𝐹𝑙 ( 𝜇 𝐿 𝜇 𝐺	)	0.1028	) 2 × ( 𝐺 √ 𝐿 𝜌 𝐺 𝜌 𝐿	(	𝜇 𝐿 𝜇 𝐺	)	0.2	) 1.416	if	𝐿 𝐺	× √	𝜌 𝐺 𝜌 𝐿	> 0.4

of the flooding point (UG,Fl) for a given L/G ratio

  𝐿,𝐿𝑜 (Step 3: ℎ 𝐿 = ℎ 𝐿,𝐿𝑜 + (ℎ 𝐿,𝐹𝑙 -ℎ 𝐿,𝐿𝑜 ) × (

														Eq. 12
											𝜇 𝐿 𝜌 𝑊 𝜇 𝑊 𝜌 𝐿	)	0.05	𝑈 𝐺 𝑈 𝐺,𝐹𝑙	) 13	(with hL,Lo and hL,Fl from steps 1 and 2)	Eq. 13 Eq. 14
	Determination of the pressure drop (DP/Z) at the working point
	𝛥𝑃 𝑍	= 𝜓 𝐿 ×	𝑎 (𝜀 -ℎ 𝐿 ) 3 ×	𝑈 𝐺 2 2 𝜌 𝐺	× (1 +	4 𝑎𝐷 𝑐𝑜𝑙	)
	With:
	𝜓 𝐿 = 𝐶 𝑃 × (	64 𝑅𝑒 𝐺	+	1.8 𝑅𝑒 𝐺 0.08 ) × 𝑒𝑥𝑝 (	13300 𝑎	. √ 𝑈 𝐿 2 𝑔	) × (	𝜀-ℎ 𝐿 𝜀	) 1.5	× ( ℎ 𝐿,𝐿𝑜 ℎ 𝐿	)	0.3	with hL and hL,Lo
	calculated from Eqs. 12 and 14
	𝑅𝑒 𝐺 =	6 𝑎	𝑈 𝐺 𝜌 𝐺 𝜇 𝐺	(1 +	4 𝑎𝐷 𝑐𝑜𝑙	)

-1

Table 5 . Comparison of the experimental loading gas velocities (UG,Lo) at a given liquid load (UL) and L/G ratio to the ones predicted using the Billet-Schultes correlations (Eqs 1-7). RE corresponds to the relative error between the experimental and the theoretical loading gas superficial velocities. The liquid holdup (hL,Lo) determined with the Billet-Schultes correlations is also provided.

 5

		UL (m h -1 ) L/G	UG,Lo experimental (m s -1 )	UG,Lo model (m s -1 )	RE	hL,Lo model
		3.06	1.02	0.64	0.64	0.0%	0.092
	PDMS	7.38 12.36	2.81 5.50	0.56 0.48	0.55 0.48	2.4% 0.5%	0.150 0.203
		23.65	12.95 0.39	0.36	8.2%	0.280
	Transformer oil	9.92 16.14	3.18 5.92	0.64 0.56	0.53 0.46	17.8% 0.160 17.6% 0.211
		1.40	0.45	0.65	0.66	1.2%	0.083
	Lubricant	2.69 3.97	1.00 1.72	0.56 0.48	0.59 0.54	5.0% 11.6% 0.158 0.122
		5.23	2.53	0.43	0.50	15.3% 0.188

J o u r n a l P r e -p r o o f

Acknowledgments

We are very grateful to the French Association Nationale de la Recherche et de la Technologie (ANRT) for the CIFRE PhD grant N° 2016/0238 attributed to Margaux Lhuissier, and also to the French governmental agency ADEME for the CORTEA funding n°1881C0001. We would like to thank our industrial partner Chimirec for providing the waste oils used in this study. We would like to give a warm thanks to Thomas HULL (UniLaSalle-Ecole des Métiers de l'Environnement) for proofreading this article.

Declaration of interests

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

J o u r n a l P r e -p r o o f J o u r n a l P r e -p r o o f


Member, IEEE Mahmoud Wafik Eltokhey

Senior Member, IEEE Mohammad Ali Khalighi

Senior Member, IEEE Zabih Ghassemlooy

Dimming-Aware Interference Mitigation for NOMA-Based Multi-Cell VLC Networks

Keywords: Visible light communications, multi-cell networks, non-orthogonal multiple access, inter-cell interference

This paper proposes inter-cell interference mitigation in a dimming-aware way for multi-cell visible-light communication networks, through efficient time-scheduling, scaling, and coordination of non-orthogonal multiple access (NOMA) transmissions at the access points. By this method, users are grouped and served in different time slots depending on whether they are at the center or at the edge of the cell. This reduces the number of NOMA users per time slot, which decreases the network computational complexity, by reducing the average number of successive interference cancellation steps. Comparison of the proposed scheme with the classical NOMA over different dimming constraints shows an improvement of up to 39% and 37% in the average sum-rate and fairness, respectively, for a 30% duty cycle transmission in a 4-cell scenario with 8 users.

I. INTRODUCTION

Visible-light communications (VLCs) have been receiving increasing attention for more than two decades, in particular, for providing wireless connectivity in indoor environments. The main factors driving this growing interest include the huge available bandwidth, immunity against radio-frequency (RF) induced interference, and inherent security [START_REF] Ghassemlooy | Visible Light Communications: Theory and Applications[END_REF]. For relatively large indoor spaces, multiple light emitting diode (LED) luminaires serve as access points (APs), which handle users in their coverage area, forming hence a multi-cell VLC network. There, users can be classified as cell-center or cell-edge users (CCUs or CEUs, respectively), depending on coverage by one or more APs. In such networks, the performances of CCUs are affected by inter-user interference (IUI), resulting from signals of users being in the same cell, whereas CEUs should mitigate both IUI and inter-cell interference (ICI), the latter resulting from signals of users in the neighboring cells. This reflects the importance of employing efficient multiple-access (MA) techniques in minimizing IUI and ICI effects.

Focusing on downlink transmission within this context, one popular scheme is the power-domain non-orthogonal MA (NOMA). By this scheme that we will call simply NOMA, the APs (i.e., the transmitters, Txs) multiplex users' signals in the power domain using superposition coding. At the user side (i.e., the receiver, Rx), to minimize multi-user interference (MUI), successive interference cancellation (SIC) is performed M. W. Eltokhey and M. A. Khalighi are with the Aix-Marseille University, CNRS, Centrale Marseille, Institut Fresnel, Marseille, France (e-mail: Mahmoud.Eltokhey@centrale-marseille.fr, Ali.Khalighi@fresnel.fr) Z. Ghassemlooy is with the Optical Communications Research Group, Faculty of Engineering and Environment, Northumbria University, Newcastle upon Tyne, NE1 8ST, UK (e-mail: z.ghassemlooy@northumbria.ac.uk)

This work has received funding from the European Union's Horizon 2020 research and innovation programme under the Marie Skłodowska-Curie grant agreement No. 764461 (VisIoN). [START_REF] Ding | Cooperative non-orthogonal multiple access in 5G systems[END_REF]. The particular interest of NOMA signaling in VLC networks arises from several facts: (i) APs typically need to handle a relatively small number of users, hence, requiring a small number of SIC detection steps; (ii) the limited mobility of users in indoor VLC scenarios and the overall slow channel time variations makes acquiring the users' channel state information (CSI) for SIC detection rather easy; (iii) the limited link distance enables signal detection at relatively high average signal-to-noise ratios (SNRs) [START_REF] Marshoud | Non-orthogonal multiple access for visible light communications[END_REF]. In addition, the advantage of NOMA over orthogonal MA techniques has been demonstrated in several previous works, e.g., [START_REF] Yin | Performance evaluation of non-orthogonal multiple access in visible light communication[END_REF], [START_REF] Eltokhey | Multiple access techniques for VLC in large space indoor scenarios: A comparative study[END_REF].

There is an apparent need for reducing the number of users in SIC detection. The reason is two fold: (i) the users with higher decoding orders (i.e., lower detection priority) require CSI of the preceding users; (ii) the users with lower decoding orders suffer from MUI arising from the higher decoding-order users. In this paper, we propose a NOMAbased scheme, called time-sliced (TS) NOMA, which adds a temporal dimension to the conventional NOMA. By TS-NOMA, CCUs and CEUs in a given cell are grouped and handled by NOMA in separate time slots. The duration and the position of these time slots are determined according to the requirements of light dimming and ICI mitigation. This solution offers reduced detection complexity and decreased IUI, compared with the conventional NOMA, while enabling efficient ICI mitigation and light dimming. It also facilitates network real-time adaptation in the case of users mobility.

The concept of adding a temporal dimension to NOMA has also been considered for RF networks. For instance, in [START_REF] Zhai | Energy-efficient user scheduling and power allocation for NOMA-based wireless networks with massive iot devices[END_REF], optimized power allocation (PA) and user scheduling for NOMA were considered for massive Internet-of-Things (IoT) networks. In the context of MA mobile edge computing, [START_REF] Wu | NOMAassisted multi-access mobile edge computing: A joint optimization of computation offloading and time allocation[END_REF] proposed optimized time allocation for mobile users for offloading their computational workload to edge servers using NOMA. Temporal fair user scheduling was also considered for NOMA signaling in [START_REF] Shahsavari | A general framework for temporal fair user scheduling in NOMA systems[END_REF]. Compared with these works, our proposed technique differs in considering dimming-compatible duty cycling of NOMA signals for multi-cell VLC networks.

II. MAIN ASSUMPTIONS

Consider a multi-cell VLC network where a central control unit coordinates N t APs and classifies N r users as CCUs or CEUs. Without loss of generality, we consider intensity modulation based on DC-biased optical orthogonal frequencydivision multiplexing (DCO-OFDM). Let the i th AP, AP i serves N i users, and denote by h ij the channel gain between AP i and the j th user U ij (handled by AP i ). Assume that (i) each AP (an LED luminaire) has a Lambertian pattern;

(ii) each Rx uses a PIN photo-detector (PD) with an optical concentrator; and (iii) any CEU is only associated with the AP corresponding to the strongest channel gain [START_REF] Feng | Multiple access design for ultra-dense VLC networks: Orthogonal vs non-orthogonal[END_REF]. Accounting for line-of-sight (LOS) propagation while neglecting non-LOS contribution, the channel DC gain is given by [START_REF] Ghassemlooy | Visible Light Communications: Theory and Applications[END_REF]:

h ij = R S (m + 1) A ij 2π l 2 ij cos m (φ ij ) cos(θ ij ), (1)
where φ ij is the Tx angle of emission, θ ij denotes the incident angle, and l ij is the distance between AP i and U ij . Also, R is the PD responsivity, S and m are the LED conversion efficiency and Lambertian order, respectively, and A eff is effective area of the Rx, which is given by

A eff = q 2 A PD / sin 2 (θ FOV ) [1]
, where A PD is the PD active area, q is the optical concentrator refractive index, and θ FOV is the Rx field-of-view (FOV), which are assumed to be the same for all Rxs for notation simplicity. To comply with the light dimming requirement, we set the duty cycle δ of the LED emission to δ = τ /T , where τ is the signal transmission time over a certain period T [START_REF] Lee | Modulations for visible light communications with dimming control[END_REF]. Denoting by P o (t) > 0 the LED output optical power, we define the dimming ratio γ as:

γ = 1 P lim T →∞ 1 T T 0 P o (t) dt , 0 ≤ γ ≤ 1 ( 2
)
where P is the average LED optical power in the case of no dimming, i.e., for γ = 1. The dimming ratio γ (in power domain) is related to the duty cycle δ (in time domain). Note that such a dimming control by duty cycling the transmitted signal can also be applied to DCO-OFDM signaling [START_REF] Wang | Performance of dimming control scheme in visible light communication system[END_REF], [START_REF] You | Efficient data transmission using MPPM dimming control in indoor visible light communication[END_REF]. In this case, the signal can be modeled as a zero-mean Gaussian random process based on the central limit theorem [START_REF] Mesleh | On the performance of different OFDM based optical wireless communication systems[END_REF], hence, dimming is mainly controlled by the DC bias.

III. CONVENTIONAL NOMA SIGNALING

At the Tx (i.e., the AP), NOMA users signals are multiplexed in the power domain using superposition coding. The assigned power to each user is set according to a PA scheme and its channel gain; users with higher channel gains are allocated less power. At the Rx, SIC detection is performed, where users with larger channel gains decode their signal in higher orders (IUI of lower order users being ideally eliminated by SIC [START_REF] Marshoud | Non-orthogonal multiple access for visible light communications[END_REF]). Assuming perfect CSI knowledge for all users, after excluding the DC component, the received signal by U ij from AP i is given by:

r ij = a ij P e h ij d ij + j-1 k=1 a ik P e h ij d ik + Ni k=j+1 a ik P e h ij d ik + z j , (3)
where P e and a ij represent the total transmit electrical power (excluding the DC offset) and the corresponding PA weight, respectively, where P o = S √ P e . Also, d ij is the desired data and z j denotes the Gaussian noise with variance σ 2 n . The first, second, and third terms in (3) represent the desired signal, the suppressed interference by SIC, and the residual interference, respectively. For the sake of simplicity and without loss of generality, we consider static PA such that a 2 ij = α a 2 i j-1 , where 0 < α < 1 represents the PA factor, i.e., the ratio between the power level of U ij to that of U ij-1 (the preceding user in the decoding order) [START_REF] Marshoud | Non-orthogonal multiple access for visible light communications[END_REF]. To ensure normalized P e , we set Ni j=1 a 2 ij = 1. The electrical signal-to-interference-plusnoise ratio (SINR) for U ij is then calculated as:

SINR NOMA,Uij = h 2 ij P e a 2 ij I ICI + h 2 ij P e k>j a 2 ik + σ 2 n , (4)
where I ICI denotes the ICI power. The upper bound on the achievable rate of U ij in terms of the duty cycle δ, corresponding to a dimming ratio of γ, is given as:

R NOMA,Uij = δ B 2 log 2 1 + SINR NOMA,Uij (bps), ( 5
)
where B is the system bandwidth and the division by the factor 2 is due to Hermitian symmetry constraint of DCO-OFDM.

IV. TS-NOMA SIGNALING A. Concept

The idea behind TS-NOMA is to adjust the duty cycle of the NOMA signal at the APs to eliminate ICI and to satisfy the light dimming requirement. Firstly, at the central control unit, users are classified as CCU or CEU depending on their channel gains. Afterwards, at each AP, NOMA is applied separately to signals of CEUs and CCUs, which are then transmitted in two different time slots, such that the duty cycle of the total transmission time satisfies the dimming requirement. In addition, CEU signals are transmitted in different nonoverlapping time slots to guarantee no ICI. To clarify better the idea, consider the example of Fig. 1(a) where a 4-cell VLC network is shown with the AP positioned at the center of each cell, and two CCUs within each cell and two CEUs in the intersecting areas of adjacent cells. Note that we define the cell boundaries based on the users CSI, according to a predefined threshold for channel gain [START_REF] Feng | Multiple access design for ultra-dense VLC networks: Orthogonal vs non-orthogonal[END_REF], which is set here to 9.74×10 -7 , corresponding to a cell radius of 2 m.

Here, we consider two schemes: (i) Scheme A, see Fig. 1(b), where in addition to separating CEU signals in time, all CCU signals (from different APs) are transmitted in non-overlapping time slots; and (ii) Scheme B, see Fig. 1(c), where the same time slot is attributed to the CCUs in different cells. The choice between these two schemes depends on the dimming requirement (see Subsection IV-B). In Fig. 1, d ij denotes the signal of U ij , and t c,i and t e,i refer to the time slots attributed to CCU and CEU groups served by AP i , respectively.

We denote by c int the number of coordinating APs, having intersecting coverage areas; c int = 4 in the example of Fig. 1. To mitigate ICI, the transmission cycle T is divided among the c int APs: each AP disposing of an interval of T /c int for signal transmission, where the maximum duty cycle is δ max = 1/c int . Note that, as CCUs are not affected by ICI, their signals could be transmitted in parallel by all the APs, as considered in Fig. 1(c). However, in this case, they must use the same t c,i , hence having less flexibility, as compared with Fig. 1(b).

Lets N c,i and N e,i denote the numbers of CCUs and CEUs served by AP i , respectively. Assuming perfect synchronization and excluding the DC signal component, the received signal for CCU j handled by AP i over time slot duration t c,i is:

r ij (t) = a ij √ P e h ij d ij + j-1 k=1 a ik P e h ij d ik + Nc,i k=j+1
a ik P e h ij d ik + z j . (6)
Similarly, for CEU j served by AP i in time slot t e,i we have:

r ij (t) = a ij √ P e h ij d ij + j-1 k=1 a ik P e h ij d ik + Ne,i k=j+1
a ik P e h ij d ik + z j . (7)
For both cases, the corresponding SINR is given by:

SINR TS,Uij = h 2 ij P e a 2 ij h 2 ij P e k>j a 2 ik + σ 2 n , (8)
and the approximate maximum achievable throughput is then:

R TS,Uij = δ i B 2 log 2 1 + SINR TS,Uij (bps), (9)
where δ i is the transmission duty cycle for U ij corresponding to AP i . Note that, δ i = t c,i /T if U ij is a CCU, and δ i = t e,i /T if it is a CEU. Also, (t c,i + t e,i )/T = δ.

B. TS-NOMA schemes

We consider three cases for the transmission scheme, depending on the required duty cycle δ. Note that we reasonably assume the same dimming level for all APs.

1) No dimming, δ = 1: In this case, we can use either Scheme A or B. However, to maximize network capacity, parallel transmission of CCUs signals is preferred, i.e., by Scheme B, Fig. 1(c). Note that, (i) in the intervals where no "signal" is transmitted, the corresponding APs still transmit the DC bias to ensure lighting; (ii) in case where there is no CEU in a cell, the corresponding CEU time slot is devoted to the transmission of CCU signals in all cells; and (iii) the absence of CCUs in a cell does not change the timing in order to avoid ICI with CCU signals from the other cells.

2) δ > δ max : In this case, δ cannot be satisfied using Scheme A because this latter can guarantee dimming with no ICI only for δ ≤ δ max . Using Scheme B, we transmit the signals of CCUs (which do not suffer from ICI) in parallel from all APs. This way, the time intervals for nonoverlapping transmission of CEU signals can be increased. We have, t c,i + cint i=1 t e,i = T . For simplicity, we consider the same t e,i for all APs. Since for dimming constraint we have (t c,i + t e,i )/T = δ, therefore:

δ i = (1 -δ)/(c int -1)
for CEUs, δ i = (δ c int -1)/(c int -1) for CCUs. [START_REF] Lee | Modulations for visible light communications with dimming control[END_REF] Note, in the absence of CEUs in a cell, the corresponding time slot is used by the CCUs in the same cell to comply with the dimming condition for all APs. If there is no CCU in a cell, the timing is still not changed to avoid ICI with CCU signals from the other cells.

3) δ ≤ δ max : Here, the required dimming condition can be satisfied by both Schemes A and B. However, as there is no room to increase the network capacity by improving time resources utilization, sequential transmission of CCU signals is preferred, i.e., Scheme A. This allows allocating different durations for t c,i , thus more flexibility in timing design. Similar to the previous case, if there is no CEU (CCU) in a cell, the corresponding time slot is used by the CCUs (CEUs) in the same cell.

C. Time-slot fixing strategies

For setting the different time intervals t c,i and t e,i for CCUs and CEUs, we consider six specific strategies that we will refer to as S1 to S6, as described in the following.

1) TS-NOMA-S1: As the simplest way, for all i, we use equal time slots, i.e., t c,i = t e,i . 2) TS-NOMA-S2: We set slot durations according to the associated number of users. For Scheme A, Fig. 1(b), this results in t c,i = N c,i /N i and t e,i = N e,i /N i .

For Scheme B, Fig. 1(c), t e,i is determined by N e,i , whereas t c is set considering the maximum number of CCUs N c,max per group among all CCU groups, i.e., N c,max = max(N c,1 , ..., N c,Nt ). This ensures setting the largest possible t c . We have then t c = N c,max /(N c,max + Nt i=1 N e,i ) and t e,i = N e,i /(N c,max + Nt i=1 N e,i ). 3) TS-NOMA-S3: We set t c,i and t e,i according to the average channel gain for each group of CCUs and CEUs (similar to the idea of static PA in NOMA [START_REF] Marshoud | Non-orthogonal multiple access for visible light communications[END_REF]). This way, for a larger average channel gain, we allocate a smaller duration to the corresponding time slot. The ratio between each two consecutive time slots ζ is a constant, called time allocation (TA) coefficient. This allows time slot allocation to be better tailored to the channel gains in every group. For Scheme B, t c is set by considering the minimum average channel gain per group (among all CCU groups), to ensure the largest possible t c . 4) TS-NOMA-S4: Setting slot duration is based on the same approach as in TS-NOMA-S3, while considering as criterion the average channel gain divided by the number of users in the group. This allows groups with a large number of users to be allocated a long time slot, hence, resulting in a still better fairness in the network. 5) TS-NOMA-S5: This strategy applies to the case of δ > δ max , where we use (10) to calculate t c,i and t e,i , to satisfy dimming and mitigate ICI. 6) TS-NOMA-S6: For the case of δ ≤ δ max , users are served as in conventional NOMA, while avoiding parallel transmission of APs (i.e., to avoid ICI). Note that the described strategies should be selected depending on the required dimming, or in other words, δ. For δ = 1, TS-NOMA-S1, -S2, -S3, and -S4 can be used. For δ ≤ δ max , these four strategies, as well as TS-NOMA-S6 can be used where t c,i and t e,i are calculated separately in each cell. For δ > δ max , only TS-NOMA-S5 can be used. In fact, TS-NOMA-S1 to -S4 do not satisfy the constraints of dimming and ICI mitigation in this case, and therefore cannot be used.

V. TS-NOMA PERFORMANCE

A. Network performance metrics

We consider two performance metrics, i.e., the maximum achievable throughput and the Jain's fairness index (FI) [START_REF] Eltokhey | Hybrid NOMA and ZF pre-coding transmission for multi-cell VLC networks[END_REF]. This latter is a measure of the throughput homogeneity among users, which for a total number of N r users is defined as:

FI = Nt i=1 Ni j=1 R ij 2 N r Nt i=1 Ni j=1 R 2 ij , (11)
where R ij denotes the maximum achievable throughput of U ij .

We have FI = 1 if all R ij are equal; as the difference between R ij increases, FI becomes smaller.

B. Main assumptions and considered scenarios

We consider a 4-cell network, as in Fig. 1(a), where APs are positioned at the center of the cells and the emitted optical power from each AP is set to P o = 1.584 W as in [START_REF] Eltokhey | Hybrid NOMA and ZF pre-coding transmission for multi-cell VLC networks[END_REF]. The considered heights of the APs and the Rxs are 2.5 and 0.85 m, 1 Number of LED chips per luminaire [START_REF] Eltokhey | Hybrid NOMA and ZF pre-coding transmission for multi-cell VLC networks[END_REF] 36 Maximum signal current per LED chip [START_REF] Eltokhey | Hybrid NOMA and ZF pre-coding transmission for multi-cell VLC networks[END_REF] 100 mA LED conversion efficiency S [START_REF] Eltokhey | Hybrid NOMA and ZF pre-coding transmission for multi-cell VLC networks[END_REF] 0.44 W/A PD responsivity R [START_REF] Yin | Performance evaluation of non-orthogonal multiple access in visible light communication[END_REF] 0.4 A/W PD area A PD [START_REF] Eltokhey | Hybrid NOMA and ZF pre-coding transmission for multi-cell VLC networks[END_REF] 1 cm 2 Rx's FOV θ FOV [START_REF] Eltokhey | Hybrid NOMA and ZF pre-coding transmission for multi-cell VLC networks[END_REF] 62 deg. Refractive index of optical concentrator q [START_REF] Eltokhey | Hybrid NOMA and ZF pre-coding transmission for multi-cell VLC networks[END_REF] 1.5 BW B 10 MHz Equivalent Rx noise power spectral density [START_REF] Yin | Performance evaluation of non-orthogonal multiple access in visible light communication[END_REF] 10 -21 A 2 /Hz respectively. Three different scenarios are considered, called Scenarios 1, 2, and 3, for which we set N r to 8, 12, and 16, and the corresponding number of CEUs to 4, 4, and 8, respectively. These represent an increased level of IUI and ICI in the case of conventional NOMA. For each scenario, 200 randomly generated users positions are used over which the performance is averaged. 1 The PA coefficient α is set to 0.3, providing the best compromise between sum-rate and FI. The same value is adopted for the TA coefficient Table I summarizes the other simulation parameters.

C. Numerical results

Here we compare the sum-rate and FI of TS-NOMA with conventional NOMA as shown in Fig. 2 for the case of no dimming, i.e., δ = 1, and the three considered scenarios described in Subsection V-B. We notice that, despite the ICI mitigation merits of TS-NOMA, NOMA achieves a better sum-rate, except for TS-NOMA-S4 case, which has a close performance to NOMA. The reason is that TS-NOMA-S1, -S2, and -S3 schemes have a poor sum-rate performance when a relatively large number of users are in a group, where they are allocated a relatively short time slot. The TS-NOMA-S4 strategy avoids such conditions, as described above, which results in a higher average sum-rate. As concerns FI, we notice from Fig. 2(b) that generally TS-NOMA schemes outperform conventional NOMA, except for TS-NOMA-S3, and for TS-NOMA-S4 in Scenario 3. The best FI are obtained for TS-NOMA-S1 and -S2, as expected: the former equally divides the time resources among the users, whereas the latter allocates them based on the number of users in each group. For TS-NOMA-S3 and -S4, FI is penalized by large variations in time slot durations due to considering a constant TA coefficient. Overall, compared with conventional NOMA, the presented results suggest TS-NOMA-S4 as the most appropriate scheme that makes a good compromise between sum-rate and FI. Finally, Fig. 3 compares TS-NOMA and NOMA over sumrate and FI, for δ < 1, namely 70%, 50%, 30%, and 25% for the three considered scenarios. Given that δ max = 1/4, for the three former cases where δ > δ max , we use the TS-NOMA-S5 strategy, whereas for the last case where δ = δ max , we use TS-NOMA-S4 (which achieved the best sum-rate and acceptable FI in Fig. 2) and TS-NOMA-S6. We notice from Figs. 3 (a), (b), and (c) that TS-NOMA-S5 outperforms NOMA in both sum-rate and FI in all scenarios. This superiority of TS-NOMA-S5 can be explained by the decreased IUI due to handling a smaller number of NOMA users per time interval, and the elimination of ICI. Similarly, from Fig. 3 (d) we notice a better sum-rate and FI for TS-NOMA-S4 and -S6, compared to conventional NOMA, due to minimized ICI. Although TS-NOMA-S6 achieves a higher sum-rate than TS-NOMA-S4 because of allocating more time resources per user, it needs the same number of SIC detection steps as for conventional NOMA. Overall, compared with conventional NOMA, we can conclude TS-NOMA-S4 and -S5 as the most appropriate schemes for δ ≤ δ max and δ > δ max , respectively, providing improved sum-rate and FI, as well as reduced Rx complexity. Note also that for all TS-NOMA schemes in Fig. 3, we notice a degradation of FI from Scenario 1 to 3 (from the smallest to the largest number of users), which is because of increased IUI. In addition, we note an increase in the average sum-rate for scenarios with a larger number of users. Yet, the increase in sum-rate from Scenario 2 to 3 is less significant than that from Scenario 1 to 2. This is due to the decrease in the available time resources per user.

VI. CONCLUDING REMARKS

We proposed the TS-NOMA scheme for multi-cell VLC networks and showed its efficiency in handling MUI, compared with the conventional NOMA. By separating in time domain the CCUs and CEUs in different NOMA groups, the proposed TS-NOMA offers ICI mitigation, dimming compatibility, as well as reduction in the number of users handled by NOMA per time interval, thus a reduced Rx computational complexity. We showed the advantage of TS-NOMA in terms of sum-rate and network fairness, especially for the case of light dimming.

As concerns the computational complexity of TS-NOMA, the main factors are: (i) power allocation (which is negligible due to using static PA); (ii) the number of SIC operations at the Rxs (with an advantage over conventional NOMA due to handling less users per time slot); and (iii) CSI acquisition (which is the same as for NOMA). Nevertheless, using TS-NOMA increases the network complexity due to synchronization management (which is rather low due to handling only two time slots per AP) and time-slot allocation. For this latter, the highest complexity in the case of dimming corresponds to TS-NOMA-S4 for δ ≤ δ max , which requires in the worst case N CCUs + N CEUs + c int -1 additions and 4 c int + 1 divisions. Overall, the increased complexity is quite justified, given the improvement achieved in the overall network performance.

Fig. 1 :

 1 Fig. 1: Illustration of the proposed TS-NOMA schemes. (a): top view of a 4-cell network example. (b) and (c): signal transmission timing for different users according to Scheme A and B, respectively.

Fig. 2 :

 2 Fig. 2: Comparison of sum-rate and FI for TS-NOMA and NOMA for the case of no dimming (δ = 1).

Fig. 3 :

 3 Fig. 3: Comparison of sum-rate and FI between TS-NOMA and NOMA for (a) δ = 70%, (b) 50%, (c) 30%, and (d) 25%.

TABLE I :

 I Simulation parameters

	Parameter	Value
	Room dimension	(7 m × 7 m × 3 m)
	LED luminaire Lamertian order m [14]

Note that we exclude in the simulations those scenarios for which CEUs are in coverage areas of more than two APs. Indeed, for such cases, the NOMA performance will be penalized considerably because of ICI, in contrary to TS-NOMA where CEUs of each AP are served in non-overlapping time slots