domain: chim spi

text:
Anna Voronova

Cristina Prieto

Maria Pardo-Figuerez

Jose Maria Lagaron

Amitav Sanyal

Bilal Demir

Thomas Hubert

Valérie Plaisance

Valérie Pawlowski

Séverine Vignoud-Despond

Valerie Plaisance

Valerie Pawlowski

Severine Vignoud

Alexandre Barras

Amar Abderrahmani

Rabah Boukherroub

Sabine Szunerits
email: sabine.szunerits@univ-lille.fr

Photothermal Activatable Mucoadhesive Fiber Mats for On-Demand Delivery of Insulin via Buccal and Corneal Mucosa

Keywords: electrospinning, on-demand release, reduced graphene oxide, mucoadhesive

Introduction

Electrospinning of fiber mats incorporating therapeutic agents have become widely explored for the development of fast-dissolving drug delivery systems. [START_REF] Huang | Stimuli-Responsive Electrospun Fibers And Their Applications[END_REF][START_REF] Celebioglu | Hydrocortisone/Cyclodextrin Complex Electrospun Nanofibers for a Fast-Dissolving Oral Drug Delivery System[END_REF][START_REF] Topuz | Electrospinning of Cylodextrin Functinal Nanobibers For Drug Delivery Applications[END_REF][START_REF] Xue | Electrospinning and Electropsun Nanofibers: Methods, Materials, And Applications[END_REF][START_REF] Li | Ultra-Fine Polyelectrolyte Fibers From Electrospinning Of Poly(Acrylic Acid)[END_REF] In contrast to other forms of drug carriers, such as liposomes and nanoparticles, electrospinning technology provides great flexibility in the therapeutic to be loaded into the fiber matrix and the material itself. [START_REF] Ning | Application Of Blocking And Immobilization Of Electrospun Fiber In The Biomedical Field[END_REF] While electrospun drug reservoirs have high drug entrapment efficiency, the delivery of the drug is guided by drug diffusion out of the matrix, as well as by the degradation/dissolution of the carrier polymer. A different concept is the controlled release of active principles at a certain time point, putting the patient into the center of action. pH-sensitive electrospun fiber are one of the first-generation stimuli-responsive fiber mats. Polymers with acid-labile functional groups as acetal and hydrazone groups induce drug release at low pH. 7 While experimentally validated, the use of pH sensitive drug release from fiber mats for medical applications is currently limited to examples of treating inflamed issue. [START_REF] Wallace | Challenges and Opportunities of pHIn Chronic Wounds Therap[END_REF] On the other hand, temperature responsive fiber mats such as poly(N-isopropyl acrylamide) (PNIPAM), take advantage of the polymer's lower critical solution temperature (LCST) at 31-35°C, where at temperatures above LCST it becomes insoluble. While such polymers are efficient drug loading matrixes at low temperature due the increased porosity, they are less adapted for temperature-controlled drug release in the form of external bandages on the human body. Electrospun fiber mats where an increase in temperature results in swelling of its 3D scaffold and drug release are more adapted for such applications. Some of us have demonstrated lately the utilization of reduced graphene oxide (rGO) loaded poly(acrylic acid) (PAA) based hydrophilic nanofiber matrix for photothermal release of antibiotics. PAA was electrospun in the presence of beta-cyclodextrin (20 wt.% of PAA) as a cross-linker. The formed PAA@rGO matrix features good photothermal heating properties reaching 51±2 °C upon laser (980 nm) irradiation at a power density of 1 Wcm -1 .

The specific features of electrospun polymer fiber mats such as high surface-volume ratio and interconnected nonwoven architecture provides them with mucoadhesive properties and should be considered for transmucosal drug carrier systems. [START_REF] Sofi | Electrospun Nanofibers For The Delivery Of Active Drugs Through Nasal, Oral And Vaginal Mucosa: Current Status and Future Perspectives[END_REF] In addition, transport of molecules across the buccal membrane is far simpler than the skin, alveoli or intestinal membranes. Therefore, mucoadhesive drug delivery is expected to improve drug bioavailability, to increase the patients' s therapy compliance, to sustain drug delivery and a relatively rapid onset of action with the possibility to remove the formulation if therapy is required to be discontinued. [START_REF] Zhang | Oral Mucosal Drug Deliveryclinical Pharmacokinetics And Therapeutic[END_REF][START_REF] Shaikh | Mucoadhesive Drug Delivery Systems[END_REF] Through less permeable that the sublingual area, the buccal oral cavity provides an ideal targeting zone. The buccal mucosa is well-vascularized, allowing drugs to be rapidly absorbed into the venous system underneath the oral mucosa. Being different from transdermal drug delivery route, mucosal surfaces are lacking the stratum corneum barrier, resulting in faster drug delivery. The inherent barriers for efficient buccal mucosal drug delivery are the keratinized tissue and the elimination of the drugs to the flashing action of saliva in the absence of uni-directional mucoadhesive mouth patches. [START_REF] Patel | Advances in Oral Transmucosal Drug Delivery[END_REF] We investigate, in this study, the usefulness of the PAA@rGO matrix as mucoadhesive films for buccal linings and insulin release (Figure 1). We speculated that the presence of ionizable groups such as carboxyl groups of PAA and rGO together with hydroxyl groups of the cyclodextrin cross-linker ensures the interaction with the mucosal membranes. Indeed, results of mucoadhesive studies confirmed excellent muco-adhesion. In vitro insulin release studies and ex vivo permeation studies using pig cheek lining and cow eyes demonstrated the possibility to deliver insulin via mucosal membranes and might open new ways for insulin therapy, which is so far based on multiple daily subcutaneous needle injections through pen or catheters connected to pumps. These methods are both painful and inconvenient, and significantly deteriorate the life quality of patients and lead patients to neglect or even give up the therapy.

Figure 1. Illustration of the formation of β-cyclodextrin cross-linked PAA@rGO nanofiber mats using electrospinning method and its use for on-demand insulin release via mucoadhesive interfaces such as cornea and buccal mucosa. Ex vivo examples: cow cornea and pig buccal lining.

Results and Discussion

2.1. Characteristics of mucoadhesive PAA@rGO fiber mats Mucoadhesive patches were fabricated by electrospinning a mixture of PAA, β-cyclodextrin (20 wt.% of PAA) and rGO (2.4 wt.%). The bandages consist of a dense fiber mat (Figure 2a) with fibers of 400±150 nm in diameter. Most importantly, PAA@rGO were water-insoluble and the fibrous nanostructures remained preserved even upon immersion for 1 week in water at 37°C (Figure 2b). This contrasts with rGO-free PAA fibers which dissolve immediate after immersion in water. The PAA@rGO films proved to be also stable upon immersion into saliva (Figure 2c). This makes it different to others works for which the nanofiber film disintegration was initiated upon contact with stimulated salivary fluid in a time frame of about 150 s. [START_REF] Illangakoon | Fast Dissolving Paracetamol/Caffeine Nanofibers Prepared By Electrosphinning[END_REF] Fiber mats assessed from four different batch fabrications were observe to have an average weight of 3.2±1.0 and thickness of 120±18 µm (Figure 2d). The degree of swelling was fast (Figure 2e), reaching up to 2000 % of its weight in the first 30 s and a steady state after 1-2 min. Swelling was stronger in alkaline pH. The results of swelling behavior correlate with other works, [START_REF] Li | Ultra-Fine Polyelectrolyte Fibers From Electrospinning Of Poly(Acrylic Acid)[END_REF] where the swelling ratio increased with increasing solution pH.

The presence of rGO was validated by the characteristic Raman bands at 1350 cm -1 (D-band) and 1580 cm -1 (G-band) (Figure 3a). The presence of rGO endows the mucoadhesive patch with photothermal properties (Figure 3b). Irradiation of PAA@rGO with a near infrared laser (980 nm, 500 mWcm -2 ) resulted in a surface temperature of 51±2 °C within 5 min. The mechanical properties of PAA@rGO were assessed in addition (Figure 3c). The maximal displacement at which the applied force dropped abruptly is identified as the breaking point [START_REF] Fonseca | Pullulan Microneedle Pathces For The Efficient[END_REF] between the patch and the mucoadhesive surface, buccal tissue in our case. PAA@rGO adheres well to the buccal tissue and shows a mucoadhesive force of 0.740.3 N cm -2 , in the order of reported mucoadhesive forces. [START_REF] Vaidyaa | Ionic Liquid-Mediated Delivery Of Insulin To Buccal Mucosa[END_REF] The biocompatibililty of PAA@rGO mats was assessed on the HeLa and ARPE-19 cell lines (Figure 3d). The toxicity was evaluated using the resazurin assay, based on the conversion of non-fluorescent dye to a fluorescent molecule by mitochondrial and cytosolic enzymes. No loss in cell viability and loss of metabolic activities was observed upon incubation of both cells for 24 h respectively. Lack of cytotoxicity suggests that no residual monomers, which could lead to potential cytotoxicity, are. Furthermore, eventual detrimental effect of thermal activation of PAA@rGO gel on cell viability was determined. Thermal activation of the of PAA@rGO mats for 10 min at 500 mW cm -2 did not induce any change in cell viability.

Insulin-loaded mucoadhesive PAA@rGO patches

We opted for using insulin as a model drug given that diabetes is one of the leading causes of early mortality worldwide. [START_REF] Wright | Life Expectancy And Cause-Specific Mortality In Type 2 Diabetes: A Population-Based Cohort Study Quantifying Relationships In Ethnic Subgroups[END_REF] As all current antidiabetics are unable to achieve a long-term glycemic control, insulin administration remains currently the only option to control blood glucose levels efficiently. To decrease blood glucose concentrations, 250-330 nM (145-192-µg mL -1 ) of insulin is needed. Human insulin (200 µg mL -1 ) was consequently integrated into PAA@rGO fiber mats through their simple immersion for 3 h under continuous shaking at 4 °C. The loading capacity for human insulin was evaluated through the determination of the concentration of insulin in solution before and after loading using HPLC, as reported by us recently. [START_REF] Pagneux | Electrothermal Patches Driving the Transdermal Delivery of Insulin[END_REF] Insulin could be efficiently loaded (89% loading efficiency) onto the mucoadhesive patches over a time span of 60 min (Figure 4a). These loading efficiencies are comparable to rGO loaded hydrogels (80%), as reported by some of us. [START_REF] Teodorescu | Photothermally Triggered On-Demand Insulin Release from Reduced Graphene Oxide Modified Hydrogels[END_REF] The loading capacity was found to be pH-dependent, with a higher insulin loading at pH 4 (Figure 4a). At pH 5, insulin is positively charged and interacts more strongly with the negatively charged backbones of PAA@rGO fibers. To determine the amount of insulin passively released from the mucoadhesive patch, the patch was immersed for 24 h in a solution of pH 7.4 (Figure 4b). About 5% of insulin was released, mostly due to some insulin being present outside the patch while most remained entrapped.

Most important, after 24 h, no further insulin was released. This contrasts with the amount of insulin released upon photothermal irradiation for 10 min (Figure 4c). The good heating capacity of the mucoadhesive patch results in an efficient heat triggered insulin release, reaching 150 µg (75%). Besides the sufficient photothermal release, PAA@rGO fiber mats featured appropriate reusability. Reloading of the mats with insulin followed by photothermal activation for 10 min could be performed with the same efficacy for 5 successive cycles.

Ex-vivo insulin permeation analysis

The distinct keratinization of the oral mucosa negatively impacts on the passage of drugs through the mucosa and limits therapeutic efficiency. [START_REF] Thirion-Delalande | Comparative Analysis of the Oral Mucosae From Rodents and Non-Rodents: Application To The Nonclinical Evaluation Of Sublingual Immunotherapy Products[END_REF] The absorption of insulin from the patches through buccal mucosa (Figure 5a) has been assessed in the following by ex vivo permeability studies, using pig-based models with close resemblance to the human buccal mucosa. [START_REF] Vaidyaa | Ionic Liquid-Mediated Delivery Of Insulin To Buccal Mucosa[END_REF][START_REF] Pinto | In Vivo,Ex Vivo And In Vitro Assessment Of Buccal Permeation Of Drugs From Delivery Systems[END_REF] For this, porcine buccal mucosa was separated from the underlining tissue using a scalpel and rinsed with PBS and fitted into Franz diffusion cells within 2 h after harvesting to preserve the membrane permeability and to ensure tissue integrity. To determine the permeability of insulin, electrospun fiber mats patch (1.5×1.5 cm 2 , 2 mg) was loaded with 500 µg mL -1 of insulin by immersion into 1 mL solution at 4 °C for 1 h. With a loading efficiency of 82% (Figure 4a), this accounts for 410 µg mL -1 insulin integrated to the patch. This patch was placed on the top of the porcine buccal mucosa (Figure 5b) and both passive and heatinitiated permeations of insulin were assessed (Figure 5c). As expected, without heat activation, the insulin-loaded fiber mat showed limited insulin permeation through the epithelial cells of the porcine mucosa. In contrast, photothermal laser activation for 10 min resulted in substantial insulin penetration across the buccal tissue. After 6 h, 122% of insulin (which accounts for 41.1 µg mL -1 insulin) has permeated the porcine mucosa model, with about 161 % of insulin having been detected in the porcine lining. In total, the heat activated patches were able to transport 282 % (82 µg mL -1 ) of the loaded insulin into the buccal lining. These findings are in line with recently reported buccal patches by Vaidya and Mitragotri using insulin loaded biodegradable polymeric patches using chitosan as mucoadhesive matrix and ionic liquids (ILs)/deep eutectic solvent (DES) as the transport facilitator. [START_REF] Vaidyaa | Ionic Liquid-Mediated Delivery Of Insulin To Buccal Mucosa[END_REF] An insulin flux of 16.6±2.9 μg cm -2 h -1 was recorded, being considerably elevated compared to that previously reported for photothermal activated insulin flux via porcine skin (as J=0.8±0.2 μg cm -2 h -1 ). [START_REF] Teodorescu | Photothermally Triggered On-Demand Insulin Release from Reduced Graphene Oxide Modified Hydrogels[END_REF] It is in comparable with results by some of us using an electrothermal activation approach with a flux of 15.6 ± 1.3 μg cm -2 h -1 17 and in line with insulin delivered from poloxamer-407 gels using transdermal iontophoresis in the presence of menthone as a skin enhancer. [START_REF] Pillai | Transdermal Delivery Of Insulin From Poloxamer Gel: Ex Vivo And In Vivo Skin Permeation Studies In Rat Using Iontophoresis And Chemical Enhancers[END_REF] As the buccal mucosa is deficient of a stratum corneum layer, insulin can pass unhindered into the entire buccal epithelium. The efficient mucosal insulin delivery under heat might further arises from the fluidization of the upper epithelium of the buccal tissue upon temperature increase, thereby indicating a paracellular mode of insulin transport. Fluorescence analysis of harvested buccal tissue before and after heat treatment and delivery of FITC-label insulin (Figure 5d) showed no structural damage of the tissue at the application site with the fluorescence signal remaining localized to the outside the tissue without heat application. After 10 min heat application and observing the fluorescence stained tissue after 6h, fluorescence of 500 µm in depth is observed, equal to the thickness of buccal nonkeratinized structures. Eventual, the increased temperature induces mucus thinning effect and the protective role of PAA@GO against saliva proteases could be reasons for the efficient insulin delivery. images of buccal tissue treated with a FTIC-insulin loaded PAA@rGO patch with (left) and with (10 min, 980 nm, 500 mW cm -2 ) activation.

In addition, comparable experiments were performed on pig cornea (Figures 6a,b), coated with a mixture of mucus, epithelial cells, proteinaceous and lipoidal material that adhere firmly to the corneal surface. Permeation of insulin via the pig cornea was comparable to that of pig buccal mucosa, with limited passive diffusion; 131% of insulin permeated the cornea and 251% remained on the cornea tissue after 6 h, corresponding to a total insulin of 371% i.e. 151 3 µg mL -1 and an insulin flux of 24.33.1 µg mL -1 h -1 (Figure 6b). This flux is severaltimes larger than that achieved using photothermal activation via porcine skin and 1.5 times using buccal mucosa.

2.4

In-vivo mucoadhesive properties and acceptability of PAA@rGO nanofibers mats by volunteers Pancreatectomized Gottingen minipigs are valuable models of insulin-dependent diabetes as they share many similarities with humans such as the morphology of the pancreas and the overall metabolic status of the two species. To investigate the efficiency of the nanofiber mat to lowering blood glucose level (BGL), insulin (100 μg mL-1 corresponding to 2.88 Insulin units) was loaded into the PAA@rGO fibers and applied inside the buccal cavity of the pig (Figure 7a). After thermal activation of the fiber mats for 10 min an immediate reduction of the BGL was observed (Figure 7b). In the control experiment without thermal activation, the BGL remained constant over the same time period. The decrease of the BGL agrees with an increase in the plasma insulin concentration, confirming the hypoglycemic effect of insulin released from the patch (Figure 7c). It furthermore coincides with the insulin plasma profile of subcutaneous injected insulin (2.88 IU). Figure 7d shows furthermore the buccal lining after removal of the insulin-loaded PAA@rGO fiber mat indicating no irritation and visual degradation of the mucus lining; Figure 7. In Vivo Permeation Study of insulin from PAA@rGO fiber mats via the buccal lines of mini-Pigs: (a) Image of pig buccal cavity in the presence of the PAA@rGO fiber mat.

(b) Effect of insulin (2.88 IU) loaded PAA@rGO fibers mats on blood glucose level (BGL) over time after 10 min thermal activation (black) and with no heat activation (grey). (c) Time dependent blood insulin level in minipig. In a control experiment (bright blue), insulin was injected subcutaneously (2.88 IU) and blood insulin was dosed over time. (d) Image of buccal cavity after PAA@rGO fiber mat removal.

2.5.

In-vivo mucoadhesive properties and acceptability of PAA@rGO nanofiber mats Finally, we applied PAA@rGO nanofibers mats to 6 healthy adult volunteers (3 male, 3 female) aged between 25 and 49 years to obtain information about in vivo residence time and patch acceptability (Figure 8a).The adhesion evaluated for a time span of 120 min was excellent, with all bandages residing for 120 min in the 6 cases f (Figure 8b). The volunteers rated the patch as excellent (note of 9.3) with most of volunteers feeling comfortable (note 8.3) while wearing the patches in the buccal cavity. 2 male volunteers considered it "mentally" difficult to have a nanomaterial in the mouth and the discomfort was related to this rather than the bandage itself. All volunteers agreed that the patches have no taste (note 10). The rating for removing the bandage is evaluated as medium (note of 6.1) indicating that it remains difficult to remove the bandage once installed removal easy. None of them reported interference with speak nor effects on saliva production and swallowing.

Conclusion

In summary, we have demonstrated that incorporation of reduced graphene oxide (rGO) into a water insoluble fiber mats, formed from poly(acrylic acid) (PAA) and β-cyclodextrin crosslinker by electrospinning, leads to the formation of water insoluble hydrogel with good photothermal characteristics. The resulting hydrogel fiber mat allowed efficient loading of a therapeutic peptide like insulin, while maintaining its good fibrous morphology. Using porcine buccal mucosa as well as porcine cornea as ex vivo models, the flux of insulin across the mucosal lining recorded upon passive or photothermal activation was investigated. As in the case of previously reported photothermal transdermal insulin delivery, thermal activation resulted in strongly increased insulin permeation rates reaching J=16.6±2.9 μg cm -2 h -1 for buccal mucosa and J= 24.33.1 µg mL -1 h -1 for porcine cornea, being several times more elevated to comparable transdermal delivery. This result confirms that buccal mucosa is suitable for the non-invasive delivery of insulin via electrospun fiber mats-based mucoadhesive patches. Efficient release of insulin by photothermal activation supports the idea of using mucoadhesive patches containing insulin for treating patients with diabetes. In type 2 diabetes, the most prominent form of diabetes, basal insulin therapy is often added to oral diabetes medication regimens that do not provide optimal glucose management. In general; this biotherapy does not require high dose of insulin. The quantity of insulin released from our mucoadhesive patches could be in line with their use in basal therapy. Therefore, compared to multiple insulin injection, this mode of delivery offers a huge advantage by increasing the compliance of patients to treatment, as it is resilient and can be easily self-administrated. In this perspective, our study allows future preclinical studies aiming to test the efficiency of this original platform to on-demand deliver insulin and to treat diabetes in animal models of diabetes. Although these dates cannot be directly related to human buccal mucosa since porcine mucosa epithelium is about 3-times thicker than in humans, insulin released from PAA@rGO fiber mats clearly show to penetrate buccal mucosa lining deeply.

EXPERIMENTAL SECTION

Materials

Polyacrylic acid (PAA, Mn 450 000 g mol -1 ), β-cyclodextrin, recombinant insulin, insulin FITC, artificial saliva, phosphate-buffered saline (PBS, pH 7.4), sodium hydroxide and hydrochloric acid were purchased from Sigma-Aldrich. Reduced graphene oxide (rGO) in powder form was purchased from Graphenea (Spain).

Electrospinning of nanofiber mats

The electrospinning of nanofiber mats was performed according to the procedure described by us previsouly. [START_REF] Altinbasak | Reduced Graphene Oxide Embedded Polymeric Nanofiber Mats: An 'On-Demand' Photothermally-Triggered Antibiotic Release Platform[END_REF] Briefly, dispersed aqueous rGO solution (1 mg mL -1 , 2.4 wt.%) was mixed with PAA (81.4 wt.%), and β-cyclodextrin (16.2 wt.%) and kept under magnetic stirring overnight.

The resulting clear solution was electrospun using a 1 mL syringe fitted with a 14-gauge blunt needle at the rate of 0.35 mL min -1 at 15 kV. The distance was kept at 15 cm during the electrospinning process. Resulting electrospun fibers were cross-linked under vacuum at 140 °C for 30 min.

Determination of swelling degree

The swelling degree of PAA@rGO nanofibers was determined following the method reported by Mittal et al. The nanofiber mat was cut into smaller pieces (1×1 cm) and dried at 70 °C for 24 h before weighting. The dried samples were then immersed in Milli-Q water at room temperature for 1 h and taken off at given time intervals (2s-10 min). The nanofiber pieces were dried with a cotton tissue in order to remove any remaining water and weighted using a microbalance. The swelling ratio of nanofiber mats was calculated as following:

Swelling degree (%) = × 100
where Wt is the weight of wet nanofibers at time t and Wd is the weight of dry nanofibers.

Determination of weight, thickness and pH of nanofiber mats

The measurements of weight, thickness and pH were performed in triplicates on randomly selected patches from independent batches. The thickness was determined using a thickness gage ID-C112XBS from Mitutoyo Corporation (Japan). For weight determination, the mats of 1×1cm 2 in size were weighted on a digital balance. The pH value was recorded using a pH meter (Hanna Instruments, Rhode Island, US) by immersion of the nanofiber mats into water for 2 min.

Drug loading

Insulin-FITC and insulin loaded patches were produced by immersion of 3.0±1.0 mg PAA@rGO electrospun mats into 1 mL (200 µg mL -1 ) of insulin-FITC or insulin aqueous solutions under continuous shaking at 150 rpm for 4 h at 4°C.

In vitro photothermal drug release studies

In vitro drug release experiments were performed in 1 mL 1× PBS solution (pH 7.4). The nanofiber mats were irradiated with a continuous mode laser (Gbox model, Fournier Medical Solution) with an output light at 980 nm at power densities 0.5, 0.7 and 1 W cm -2 for 10 min.

The temperature change was captured by an Infrared Camera (Thermovision A40) and treated using ThermaCam Researcher Pro 2.9 software. The amounts of insulin and rhodamine B released were evaluated by high-performance liquid chromatography (HPLC) of the solution collected after irradiation using calibration curves generated for each drug.

Cytotoxicity tests of PAA@rGO bandage

The ARPE-19 spontaneously immortalized cell line of human retinal pigment epithelium and the e HeLa cells, derived from cervical carcinoma from a 31-year old female [ATCC® CCL-2™, ECACC, Sigma Aldrich, Saint-Quentin Fallavier, France], were cultured and maintained in Dulbecco's Modified Eagle's medium (DMEM, Gibco®) high glucose with 2 mM Glutamine, DMEM/F12 media with 2 mM Glutamine and DMEM supplemented with 10% fetal bovine serum (FBS, Gibco®) and 1% penicillin-streptomycin (Gibco®), respectively, in a humidified incubator at 37 °C and 5% CO2. Cells were seeded at a density of 10 4 cells/well in a 96-well plate and grown for 24 h before assay.

The cell viability of PAA@rGO was evaluated using resazurin cell viability method. Briefly, 1 mL of the resazurin solution (11 µg mL -1 ) in complete medium was added to each well containing the gel and the plate was incubated for 4 h in the humidified incubator. The fluorescence emission of each well was measured at 593 nm (20-nm bandwidth) with an excitation at 554 nm (18-nm bandwidth) using a Cytation™ 5 Cell Imaging Multi-Mode Reader (BioTek Instruments SAS, France). Each condition was replicated three times

Ex-vivo permeation experiments

All animal studies were conducted at the Plateforme de Recherche Expérimentale (University of Lille, France) in accordance with the ethnic protocol entitled "Délivrance transdermique, de substances pharmacologiques, contrôlée par l'utilisation d'une énergie thermique chez la souris et le mini-porc" (Nr APAFIS#21080-2019111511521244). 3 female mini-pigs weighting 10-16 kg were used . Ex-vivo permeation experiments were performed using cornea and buccal mucosal lining from pigs. Fresh porcine buccal mucosa was separated from the underlining tissue using a scalpel and stored in PBS 1× at 4 °C not more than 12 h before the permeation experiment. Tissue with any visual damage was discarded. Intact mucosa was immersed in deionized water at 65 °C for 60 s in order to separate the epithelium from the connective tissue as reported previsouly. [START_REF] Franz-Montan | Evaluation Of Different Pig Oral Mucosa Sites As Permeability Barrier Models For Drug Permeation Studies[END_REF] Mucosa samples (thickness 310±80 µm) were gently cut into pieces (1cm 2 ) and placed between the chambers of Franz diffusion cell and ex-vivo experiment was performed. In the case of pig cornea, they were gently cut into 1.8 cm 2 circular pieces prior to their use in the Franz diffusion cell.

Insulin skin diffusion experiments were carried out using static Franz diffusion cells (Proviskin, France), exhibiting an effective area of 0.785 cm 2 . The receptor chamber was filled with 1× PBS (pH 7.4) solution and maintained at 37 °C using a circulating bath (Julabo) with magnetic stirring at around 500 rpm. The cornea pieces were carefully clamped between the donor and the receptor chambers (8 mL) and pre-incubated for 20 min. The drug loaded electrospun mats were irradiated with a continuous wave laser at 980 nm for 10 min (laser power=1 Wcm -2 ). At determined time intervals, 300 μL aliquots of diffused solution were removed from the receptor compartment and analyzed by HPLC. After each sampling, an equal volume of fresh diffusion medium was added to the receptor chamber.

Evaluation of insulin trapped in the ex vivo sample was performed by adding the species into water/ice mixture for 10 min and sonicated in the presence of ZnO2 beads (4 mm in diameter), before being centrifuged for 30 min at 13500 rpm using an ultracentrifuge (Mini Scan Fuge ORIGIO). The liquid phase was collected and filtrated through a 0.1 μm Nylon filter (Whatman Puradisc 13 mm) and the amount of insulin was determined by HPLC/UV.

In-vivo residence times

The residence time was evaluated on 6 volunteers (3 male and 3 females). Following the international standard of Good Clinical Practices, placebo patches (15×12 mm) were applied to the buccal mucosa for 5 s with applied pressure. Patch adhesion was monitored for 2h and residence time determined. The residence time was considered as the time until the patch dropped from the buccal mucosa. An acceptable questioner was completed by each volunteer to obtain information about parameters such as taste, comfort, dry mouth and salivation.

Data analysis

Results are presented as mean standard deviations unless overwise stated. Statistical analysis all date was carried out using GraphPad prism version 7.0 (GraphPad software Inc., San Diego, USA). All experiments were conducted at least in triplicate.

Characterization

Fourier transform infrared spectroscopy (FT-IR) analyses were performed on Thermo Fisher Scientific Inc, Nicolet 380. Raman spectroscopy measurements were recorded on a LabRam HR Micro-Raman system (Horiba Jobin Yvon, France) using a 473-nm laser diode as excitation source. Scanning electron microscopy (SEM) micrographs of prepared mats were performed with an electron microscope Hitachi S-4800 FE-SEM (Hitachi Technologies Corp., Tokyo, Japan). The samples were sputtered with a gold/palladium layer before scanning. The accelerating voltage was 10 kV, the probe current was 9800 nA. Mechanical properties of elastomeric PAA-@rGO fibers were investigated via a TAXT plus C Texture Analyzer. The fiber samples were soaked in water and squeezed between tensile grips; a tension test was applied using the following parameters with a pre-test speed of 1.0 mm/s, test speed of 0.5 mm/s and a post-test speed of 10 mm/s. Tensile test was performed until seeing a macroscopic damage of the fibers. Relative analysis was applied to calculate stress, strain and Young modulus of elasticity. The photothermal effect of the fibers was assessed using a 980 nm-continuous wave laser (Gbox model, Fournier Medical Solution) with an output light at 980 nm at power densities 0.5, 0.7 and 1 W cm -2 . This laser was injected into a 400 µm-core fiber and placed at 6 cm from the bottom of the wells. The output was not collimated and the resulting beam divergence allowed to illuminate uniformly up to 4 wells. The temperature changes on the surface of the nanofiber mats were captured by an infrared camera (Thermovision A40) and treated using ThermaCam Researcher Pro 2.9 software. The tissue cuts were done using The Leica CM3050 S cryostat (thickness of 10 µm) and visualized with Cytation™ 5 Imaging multi-mode reader 4.12. High-performance Liquid Chromatography (HPLC) for Quantification of insulin

The concentration of insulin loaded onto PAA@rGO was determined by an HPLC system (Shimadzu, Tokyo, Japan) equipped with a 5 μm C18 AQ Uptisphere X ® , 250 mm × 4.

Figure 2 .

 2 Figure 2. Characteristics of mucoadhesive PAA@rGO fiber mats: (a) SEM images at different magnifications. (b) SEM of PAA@rGO after immersion for 1 week at 37°C in water.

Figure 3 .

 3 Figure 3. Characterization of PAA@rGO: (a) Raman spectra of PAA and PAA@rGO mats. (b) Photothermal heating curves of the patch at 980 nm illumination as a function of power density (0.5-1 W cm -2 ). (c) Tensile tests via Force-displacement curve of PAA@rGO on porcine buccal mucosa. (d) Cell viability of ARPE-19 and HeLa cell in the presence of PAA@rGO without and with heat activation at 0.5 W cm -2 for 10 min.

Figure 4 .

 4 Figure 4. Insulin loaded PAA@rGO. (a) Human insulin (200 µg mL -1 ) loading capacity as a function of solution pH. (b) Passive release over time. (c) Photothermal release with time upon illumination at 980 nm using a laser power of 0.7 W cm -2 . (d) Loading and release of insulin into the bandage.

Figure 5 .

 5 Figure 5. Ex-vivo insulin permeation studies of insulin-loaded PAA@rGO nanofiber mats via porcine buccal mucosa: (a) Photographic image of porcine buccal tissue harvested and used for permeation studies. (b) Photograph of porcine buccal tissue fixed into Franzen diffusion cell and immobilization of insulin-patch on top. (c) Cumulative permeation profile of insulin through the buccal tissue after 6 h upon passive or photothermal activation (10 min action with laser at 980 nm) (n=3). Error bars represent means of ± SEM. (d) Fluorescence

Figure 6 .

 6 Figure 6. Ex-vivo insulin permeation studies of insulin-loaded PAA@rGO nanofibers mats via porcine cornea: (a) Photographic image of porcine buccal tissue harvested and used for permeation studies. (b) Photograph of porcine buccal tissue fixed onto Franzen diffusion cell and immobilization of insulin-patch on top. (c) Cumulative permeation profile of insulin through the buccal tissue after 24 h upon passive or photothermal activation (10 min irradiation at 980 nm, 500 mW cm -2 ) (n=3). Error bars represent means of ± SEM.

Figure 8 .

 8 Figure 8. In-vivo mucoadhesive performance of PAA@rGO nanofibers mats placebo: (a) Photograph of placement of mucoadhesive patches in humans attached by 5 s pressure application. (b) Acceptability evaluation of 6 adult volunteers on the adhesion time, overall bandage excellent, comfort, taste and removability in a scale of 0 (bad)-5 medium-10 (excellent)

4. 9 .

 9 In-vivo experiments using Mini-PigsInsulin delivery in vivo was tested on female domestic pigs (three) weighing 50-62 kg. The animals were subjected to total pancreatectomy a day before each experiment. A tunneled central venous catheter into the external jugular vein was placed in each animal facilitating a blood sampling. The insulin loaded patches (2.8 IU) were placed inside the buccal lining of the pig moth. No fixing is needed due to the excellent mucoadhesive character of the patch. The PDD@rGO patch was activated for 10 min giving the temperature ~50 °C. Glucose levels were measured directly from fresh blood via central catheter using a commercial glucometer. Blood samples were obtained at and different time intervals. Blood samples were centrifuged at 5.000 rpm for 10 min at 4 °C and plasma was immediately separated and stored at -80 ° C until analyzed. Radioimmunoassay kits were used for the measurements of insulin (Bi-Insulin RIA ® ; Bio-Rad, Elexience, Verrières-le-Buisson, France). All the manipulations in pigs, including total pancreatectomy, placing of the catheter, microneedles and patch activation as well as blood samplings, were performed under general anesthesia with a 2-3% concentration of isoflurane (Aerrane; Baxter, France) and under analgesia with transcutaneous fentanyl (Recuvyra ; Elanco; France). Animal studies were approved by the Institutional Ethics Committees for the Care and Use of Experimental Animals of the University of Lille (protocol n°21080). All experiments were performed in accordance with the guidelines for animal use specified by the European Union Council Directive of September 22, 2010 (2010/63/EU).

  6 mm column (Interchim, Montluçon, France) heated to 40 °C. The mobile phase consisted of a mixture of eluent A (formic acid 0.1% in water) and eluent B (formic acid 0.1 % in acetonitrile) at a flow rate of 1 mL/min. The samples were injected at a volume of 40 L and the detection was recorded at a wavelength of 215 nm.

ACKNOWLEDGMENTS

The Centre National de la Recherche Scientifique (CNRS), the University of Lille, the Hautsde-France region, and the CPER "Photonics for Society" are acknowledged for financial support. AV thanks the i-SITE foundation of the University of Lille for a PhD fellowship. This

We acknowledge the help and involvement of the Dhure (Département Hospitalo-Universitaire de Recherche Expérimentale) platform in work related to the mini-pig.


Yannick Hallez

Jacques Magnaudet

A numerical investigation of horizontal viscous gravity currents

We study numerically the viscous phase of horizontal gravity currents of immiscible fluids in the lock-exchange configuration. A numerical technique capable of dealing with stiff density gradients is used, allowing us to mimic high-Schmidt-number situations similar to those encountered in most laboratory experiments. Plane two-dimensional computations with no-slip boundary conditions are run so as to compare numerical predictions with the 'short reservoir' solution of Huppert (J. Fluid Mech., vol. 121, 1982, pp. 43-58), which predicts the front position l f to evolve as t 1/5 , and the 'long reservoir' solution of Gratton & Minotti (J. Fluid Mech., vol. 210, 1990, pp. 155-182) which predicts a diffusive evolution of the distance travelled by the front x f ∼ t 1/2 . In line with dimensional arguments, computations indicate that the self-similar power law governing the front position is selected by the flow Reynolds number and the initial volume of the released heavy fluid. We derive and validate a criterion predicting which type of viscous regime immediately succeeds the slumping phase. The computations also reveal that, under certain conditions, two different viscous regimes may appear successively during the life of a given current. Effects of sidewalls are examined through three-dimensional computations and are found to affect the transition time between the slumping phase and the viscous regime. In the various situations we consider, we make use of a force balance to estimate the transition time at which the viscous regime sets in and show that the corresponding prediction compares well with the computational results.

Introduction

Gravity currents are generated when a dense fluid (of density ρ 2 ) is released in a lighter one (of density ρ 1 ) or vice versa. Such flows are widely encountered in environmental and geophysical situations (storm outflows, avalanches and the like) [START_REF] Hoult | Oil spreading on the sea[END_REF][START_REF] Simpson | Gravity currents in the laboratory, atmosphere, and ocean[END_REF] as well as in engineering problems (e.g. warm water discharges in cold water, spreading of fire gases and ventilation; Linden 1999). While theoretical models and extensive laboratory experiments have been devoted to this class of flows for a long time (see [START_REF] Simpson | Gravity Currents[END_REF] for an excellent review of the field up to the last decade), use of direct numerical simulation (DNS) to investigate them in detail is relatively recent. Most of the computational studies to date have focused on the early slumping phase of inertia-dominated currents during which the front velocity is roughly constant (e.g. [START_REF] Härtel | Analysis and direct numerical simulation of the flow at a gravity-current head. Part 1. Flow topology and front speed for slip and non-slip boundaries[END_REF][START_REF] Necker | Mixing and dissipation in particle-driven gravity currents[END_REF][START_REF] Cantero | On the front velocity of density currents[END_REF]. This is of course partly because DNS can contribute to clarify several important open questions related to the internal structure of the flow during this phase. In contrast, the case of currents dominated by viscous effects or that of the late viscous phase of gravity currents initially dominated by inertia has not yet been widely addressed in DNS. As will become clear later, the main reason for this is that most of the numerical techniques used in previous investigations do not allow this viscous phase to be properly computed. In the present paper we employ DNS to focus on the transition from the slumping phase to the subsequent viscous regime(s) and on the characteristics of these viscous regimes themselves in the case of an instantaneous fixed-volume release in the lock-exchange configuration.

More precisely, we consider a semi-infinite horizontal channel of height H closed by a vertical wall located at x = -L 0 and containing a volume V 0 = L 0 H (in the two-dimensional case) of dense fluid restrained by a lock placed at x = 0. This lock is instantaneously opened at t = 0 so that a gravity current develops, the front position of which is located at x f (t). We assume molecular diffusion and surface tension effects to be negligible, so that the two fluids only mix owing to hydrodynamic effects, and no capillary regime has to be considered. In the whole paper, we only consider small density differences compatible with the Boussinesq approximation. All the computational results to be discussed are issued from two-dimensional simulations in plane channels except those of § 2 and § 3.3 in which intrinsically three-dimensional configurations are considered. In what follows we normalize lengths and time by H and √ H/g , respectively (g = g(ρ 2ρ 1 )/ρ 2 denoting the reduced gravity) and mark the resulting dimensionless quantities with a tilde. To clarify the raison d'être of some of the issues addressed in this paper, a brief discussion of the various successive regimes that can be expected in the above configuration is in order.

If viscosity is neglected, the solution at early times is that of the dam-break flow problem, namely xf = F r t, (1.1) where F r = (dx f /dt)/ √ g H is the Froude number. In a well-known paper [START_REF] Benjamin | Gravity currents and related phenomena[END_REF]) it was shown that F r = 0.5 for an energy-conserving current. However, high-Reynolds-number experiments rather indicate 0.45 < Fr < 0.5 (see e.g. [START_REF] Marino | The front condition for gravity currents[END_REF], a difference due to small wall-friction effects [START_REF] Bonometti | Wall effects in non-Boussinesq density currents[END_REF]. Equation (1.1) characterizes the slumping phase. If the receding front reaches the endwall during this phase, then xf = L0 for t = TS such that TS = L0 F r -1 .

(1.2) For much larger times ( xf L0 ) and provided viscous effects are still negligible, a self-similar regime during which buoyancy is balanced by inertia is approached. This is the inertial phase during which the current length in a plane channel obeys [START_REF] Huppert | The slumping of gravity currents[END_REF] lf = xf + L0 = 1.47 Ṽ 1/3 0 t2/3 .

(1.3)

In contrast, if inertia is negligible right after the opening of the gate or if it is long before the receding front reaches the endwall, the current is self-similar with [START_REF] Gratton | Self-similar viscous gravity currents: phase-plane formalism[END_REF] xf = αRe 1/2 t1/2 , (1.4) where the flow Reynolds number is defined as Re = g 1/2 H 3/2 /ν (ν denoting the kinematic viscosity) and α = 0.0922 in the two-dimensional case [START_REF] Hinch | Inclined to exchange: shocking gravity currents[END_REF] or α = 0.0739 in a cylindrical horizontal tube [START_REF] Séon | Transient buoyancydriven front dynamics in nearly horizontal tubes[END_REF]), H being the tube diameter in the latter case. This solution holds until the receding front reaches the endwall, say for t = TV such that TV = Re -1 ( L0 /α) 2 .

(1.5)

Then the solution ceases to be self-similar.

On the other hand, if inertia is negligible and xf L0 , a different self-similar asymptotic corresponding to the point-source approximation is approached in which lf = βRe 1/5 Ṽ 3/5 0 t1/5

(1.6) in a two-dimensional channel [START_REF] Huppert | The propagation of two-dimensional and axisymmetric viscous gravity currents over a rigid horizontal surface[END_REF] or

lf = β c Re 1/4 Ṽ 1/2 0 t1/4 (1.7)
in a horizontal cylindrical tube [START_REF] Takagi | The effect of confining boundaries on viscous gravity currents[END_REF], the prefactors β and β c being 1.133 and (12/35) 1/4 ≈ 0.765, respectively. For convenience, the viscous regimes with xf ∼ t1/2 and lf ∼ t1/5 or lf ∼ t1/4 will be referred to as the long reservoir (LR) and short reservoir (SR) regimes, respectively. Several of these predictions have been convincingly validated experimentally, while some others still deserve confirmation. Generally, experiments performed with axisymmetric currents [START_REF] Didden | The viscous spreading of plane and axisymmetric gravity currents[END_REF][START_REF] Huppert | The propagation of two-dimensional and axisymmetric viscous gravity currents over a rigid horizontal surface[END_REF] or with currents propagating in cylindrical tubes [START_REF] Séon | Transient buoyancydriven front dynamics in nearly horizontal tubes[END_REF][START_REF] Takagi | The effect of confining boundaries on viscous gravity currents[END_REF] have been more conclusive because influence of the boundaries is better controlled. Many other configurations involving viscous gravity currents have been considered as well but will not be discussed here. A comprehensive review of the field was given by [START_REF] Huppert | Geological fluid mechanics[END_REF].

In order to predict the succession of regimes experienced by a given current, the characteristic transition times to each regime have to be compared. In particular, for given flow conditions, it is first necessary to determine whether the inertial solution (1.3) or the SR viscous solution (1.6) occurs after the receding current hits the endwall. This can be determined by noting that the inertial phase is observed only if the dimensionless height of the resulting current goes below 0.075, which requires lf > lf I = 13.33 Ṽ0 according to the box-model prediction [START_REF] Huppert | The slumping of gravity currents[END_REF]. Moreover the transition length lf SR and time TSR at which the viscous force associated with the dam-break solution (1.1) starts to exceed that corresponding to the SR solution (1.6) after the current has been influenced by the endwall obey [START_REF] Huppert | The slumping of gravity currents[END_REF] lf SR = Re 2/7 Ṽ 5/7 0 , (1.8) TSR = Re 3/7 Ṽ 4/7 0 .

(1.9)

Hence the inertial phase can exist only if lf I < lf SR , which implies Re & 8650 Ṽ0 . Given the moderate Reynolds numbers used in the computations to be described below, this condition will never be met in this paper. Therefore the inertial regime (1.3) will be ignored in what follows, and we will only be concerned with situations in which the slumping phase solution (1.1) is followed either by the viscous LR solution (1.4) (which sets in for t = TLR ) or by the SR solution ((1.6) which sets in for t = TSR or (1.7)) or by both. More precisely, the succession of phases to be observed will depend on how TS , TSR , TV and TLR compare to each other.

Based on the above discussion, the structure of the present paper is as follows: Section 2 briefly describes the numerical techniques we use and shows that they allow a faithful reproduction of the characteristics of viscous gravity currents involving nearly immiscible fluids. Section 3 describes computational results which, to our knowledge, provide the first clear validation of the theoretical predictions for the viscous LR and SR regimes in two-dimensional channels, i.e. (1.4) and (1.6), respectively. We then discuss the consequences of the effects of an O(1) Schmidt number on previous computational predictions and examine the influence of sidewalls in experiments through three-dimensional computations and force-balance arguments. In § 4 we first make use of computational results to derive and validate a theoretical criterion predicting which type of viscous regime immediately succeeds the slumping phase. In the case in which the slumping phase is followed by the LR regime, we show that a secondary viscous transition from the LR regime to the SR regime later occurs and determine theoretically the corresponding transition time. Computations also help understand the gradual evolution of the current from the time the receding front hits the endwall to that at which the resulting disturbance reaches the front of the advancing current. In particular, they reveal that, under certain conditions, the LR regime can even set in a short time after the receding front has hit the endwall. We finally summarize our main conclusions in § 5 and present a unified diagram of the various evolutions a given viscous current may follow.

Numerical approach and validation

Most DNS of gravity currents reported in the literature made use of spectral or pseudo-spectral techniques because they offer excellent accuracy even with a moderate number of grid points (see [START_REF] Härtel | Analysis and direct numerical simulation of the flow at a gravity-current head. Part 1. Flow topology and front speed for slip and non-slip boundaries[END_REF][START_REF] Necker | High-resolution simulations of particledriven gravity currents[END_REF][START_REF] Necker | Mixing and dissipation in particle-driven gravity currents[END_REF][START_REF] Cantero | On the front velocity of density currents[END_REF]. However the numerical stability of these methods is compromised when sharp gradients have to be dealt with, such as the density gradient across the front of a high-Schmidt-number gravity current. (The Schmidt number is defined as the ratio of the kinematic viscosity to the binary molecular diffusion coefficient.) The usual remedy has been to stabilize the corresponding computations by considering a diffusion term in the density equation. Therefore the predictions made in the aforementionned studies correspond to a Schmidt number of order unity. This has almost no influence on the predictions during the slumping phase in which inertia effects dominate and is indeed representative of realistic conditions for viscous gazeous gravity currents. However most laboratory experiments make use of liquids (essentially fresh and salt water) for which the Schmidt number is rather O(10 3 ). To our knowledge, the only attempt to compute the viscous phase of a gravity current using a spectral approach has been that of [START_REF] Cantero | On the front velocity of density currents[END_REF] in which the viscous-spreading law was found to be closer to the l f ∼ t 3/8 prediction of [START_REF] Hoult | Oil spreading on the sea[END_REF] for a gravity current flowing on a fluid surface than to the expected result of [START_REF] Huppert | The propagation of two-dimensional and axisymmetric viscous gravity currents over a rigid horizontal surface[END_REF]. One can of course suspect the difference to be due to Schmidt number effects, as the theoretical model assumes that the released volume is strictly conserved and that the interface between the two fluids is kept sharp, which, at late times, is only true if molecular diffusion has negligible effects. In particular, a recent study by Bonometti & Balachandar (2008) points out that the front thickens tremendously with time in low-Reynolds-number gravity currents with an O(1) Schmidt number.

Based on these remarks, it is clear that a different numerical strategy is required to obtain reliable computational predictions of high-Schmidt-number viscous gravity currents which can be compared with theoretical predictions and most laboratory experiments. Obviously, the key point is to select a technique capable of handling stiff density profiles. To this end, we use the in-house JADIM code. This code was originally developed to compute two-phase incompressible flows of immiscible fluids with highly contrasted densities, especially bubbly flows (e.g. [START_REF] Bonometti | Transition from spherical caps to toroidal bubbles[END_REF], 2007). Here surface tension is set to zero, since we do not wish to consider the capillary regime of gravity currents. In the above finite-volume code, the transport equation for the local volume fraction of the heavy fluid (which is identical to that of the density field) is solved using a flux-corrected scheme. More precisely, the transport equation is split into successive one-dimensional sub-steps, one along each grid direction. Within each sub-step, the corresponding advective flux is discretized using a Zalesak scheme [START_REF] Zalesak | Fully multidimensional flux-corrected transport algorithms for fluids[END_REF]) which combines the use of a low-order (here first-order) and a high-order (here fourth-order) spatial discretization to preserve positivity and monotonicity of the volume fraction and to ensure second-order spatial accuracy. The momentum equations, in which no reference to the Boussinesq approximation is made, are discretized using second-order centred schemes. Time advancement is achieved through a Runge-Kutta (RK3)/Crank-Nicholson algorithm. The incompressibility condition is finally satisfied at the end of each time step, using a variable-density projection technique. Overall, the code is second-order accurate in space and firstorder accurate in time, owing to the Euler scheme employed to advance the density equation.

It has to be pointed out that while the set of equations used here corresponds to immiscible fluids, since molecular diffusivity is set to zero, the treatment of the density field on a discrete grid prevents the interface between the two fluids from remaining infinitely sharp. Actually, the density is found to vary from ρ 1 to ρ 2 over one to three grid cells, resulting in a large but finite effective Schmidt number. Bonometti & Balachandar (2008) employed the same code and showed the effective numerical Schmidt number to be O( 103 ), which is typical of experiments performed with liquids. This is essentially a qualitative estimate, since the numerical Schmidt number at a given point of the flow depends on the local grid size and orientation of the front with respect to the grid. Therefore a proper choice of the grid size has to be made to ensure that the front thickness remains much less than the viscous diffusion scale everywhere. Moreover, since this numerical diffusion spreads the interface over one to three cells, it would be impossible to obtain an accurate description of immiscible gravity currents at very late stages, when the current height is less than a few cells. We took care of these various aspects to ensure that the simulations presented in the paper are sufficiently resolved at all times. We refer the reader to [START_REF] Hallez | Effects of channel geometry on buoyancy-driven mixing[END_REF] for more details on the numerical technique as well as for validations with gravity-driven flows confined in tilted tubes and to Bonometti et al. (2008) for results obtained with the same code in non-Boussinesq density currents.

The computations to be described below were run in the lock-exchange configuration defined in § 1, either in two-dimensional or three-dimensional channels or in cylindrical tubes. The length-to-height (or length-to-diameter) ratio ranges from 32 to 102. The corresponding grids involve 80 nodes in the vertical direction for two-dimensional channels and 48 2 nodes in the cross-section for three-dimensional channels. In the streamwise direction, the number of grid points ranges from 512 to 3168 for all geometries, depending on the length-to-height ratio. The density ratio is set to ρ 2 /ρ 1 = 1.008. We assume the dynamic viscosity to be the same in both fluids, which implies that the kinematic viscosity is slightly larger in the light fluid.

Prior to discussing new findings provided by these computations, we find it useful to describe a validation test corresponding to a well-documented case, namely that of a current propagating in a long cylindrical tube with the lock located midway between the two endwalls, i.e. starting from an initial symmetric configuration. This

10 1 10 1 10 2 10 3 t t x f 1 xf ~ t 1/2 ĩf ~ t 1/4 ˜ Figure 1.
Evolution of the front position in a current resulting from a release in a cylindrical tube: -• -• -, computational result for a short release ( L0 = 1); --: computational result for a symmetric release (here L0 = 17); + + +, experimental data from [START_REF] Séon | Transient buoyancydriven front dynamics in nearly horizontal tubes[END_REF]; the thin dashed lines indicating the t 1/2 and t 1/4 regimes were plotted using the theoretical prefactors α = 0.0739 and β c = (12/35) 1/4 determined by [START_REF] Séon | Transient buoyancydriven front dynamics in nearly horizontal tubes[END_REF] and [START_REF] Takagi | The effect of confining boundaries on viscous gravity currents[END_REF], respectively.

problem was addressed, both theoretically and experimentally, by [START_REF] Séon | Transient buoyancydriven front dynamics in nearly horizontal tubes[END_REF].

We computed the same physical situation using the experimental set of parameters.

The tube was discretized with 32 nodes in the radial direction and 64 nodes in the azimuthal direction. The evolution of the front position recorded in the experiment performed under conditions Re = 790 and ρ 2 /ρ 1 = 1.008 with almost immiscible liquids (fresh and salt water) is plotted in figure 1, together with the corresponding computational prediction.

The theoretical prediction (1.4) obtained by the same authors for the LR regime is also reported for comparison. Computational results display an excellent agreement with experimental data all along the flow evolution which includes a slumping phase, a transitional stage and an LR viscous phase. Another computational result corresponding to a short release ( L0 = 1) is added with the corresponding theoretical prediction (1.7). In the latter case, the agreement on the time exponent is again excellent, and the numerical prefactor in (1.7) is recovered within 1 %.

From the above comparison, we conclude that our numerical technique accurately reproduces the experimental and theoretical results of both [START_REF] Séon | Transient buoyancydriven front dynamics in nearly horizontal tubes[END_REF] and [START_REF] Takagi | The effect of confining boundaries on viscous gravity currents[END_REF] for the LR and SR viscous regimes in a horizontal cylindrical tube, respectively. Therefore we are confident that our code will allow us to carry out numerical experiments to explore various aspects of viscous gravity currents of immiscible fluids.

Predictions and experiments in plane channels

3.1. Predictions in the LR regime In two dimensions, the theoretical prediction for the LR regime (1.4) (with α = 0.0922)

has not yet received any experimental validation. Figure 2(a) reports the results of several computations we ran to check this prediction, using large heavy fluid releases with L0 > 17 to ensure that the current becomes viscously dominated before it reaches the nearest endwall. In this figure, the front position xf and time t have been rescaled by the transition position xf LR and time TLR to the LR regime at the considered Reynolds number so that we define xf = xf /x f LR and t = t/ TLR . The transition takes place at the time at which the viscous force associated with the slumping solution starts to exceed that associated with the LR viscous solution. Equivalently, the transitional characteristics can be computed by equating the front positions predicted by (1.1) and (1.4), which yields

xf LR = α 2 ReF r -1 , (3.1) TLR = α 2 ReF r -2 . (3.2)
Introducing xf and t, (1.1) and (1.4) reduce to xf = t and xf = t1/2 respectively, whereas the transition obviously takes place at ( xf LR , tLR ) = (1, 1). With this renormalization, all computational curves collapse on a single master curve valid all along the current life. The difference between the time exponent (respectively the numerical prefactor) predicted by (1.4) and that found in computations is less than 1 % (respectively 2 %) for sufficiently well-established self-similar states. Hinch et al.

(2007) also showed the current height to be governed by

∂ h ∂ t = 1 3 Re ∂ ∂ x h3 (1 -h3 ) ∂ h ∂ x . (3.3)
The normalized current profile recorded at 45 successive times during the LR regime of a current corresponding to Re = 200 is plotted in figure 2(b). This profile is found to follow closely the theoretical prediction obtained by integrating (3.3). The two profiles corresponding to early stages of the flow indicate that the self-similar viscous profile is settled progressively, starting from the initial lock position. At early times, the local inertia/gravity balance results in a sharp profile with a prominent blob, and this picture holds until viscous effects dominate the front dynamics.

Predictions and experiments in the SR regime

The prediction (1.6) derived by [START_REF] Huppert | The propagation of two-dimensional and axisymmetric viscous gravity currents over a rigid horizontal surface[END_REF] for the SR viscous regime has been checked in several experiments, especially those of [START_REF] Rottman | Gravity currents produced by insantaneous releases of a heavy fluid in a rectangular channel[END_REF]. Nevertheless, the slopes inferred from the corresponding data reveal a significant variability. More precisely, the slopes deduced from figure 12 of the latter paper lead to time exponents between 0.23 and 0.27 approximately, which suggests a 1/4 exponent rather than the 1/5 exponent predicted by (1.6). A similar variability may be observed in the results of [START_REF] Marino | The front condition for gravity currents[END_REF] and is nicely illustrated in figure 9(d ) of the paper by [START_REF] Cantero | On the front velocity of density currents[END_REF]. (Note that there is a small misprint in this figure, the time exponent displayed by the solid line being actually -4/5). We also carried out some additional experiments in an horizontal 4 m long channel, 0.25 m wide and 0.5 m high. The two fluids were clear water and dyed salt water. In order to reach the SR regime, we performed lock exchanges with small volumes of salt water of typical density 1050 kg m -3 . The initial lengths L 0 of the released volumes were 0.044 m, 0.095 m or 0.15 m, and their initial heights h 0 were about 0.06 m, whereas different total heights H of fresh water were employed. The evolution of xf in the experiments of [START_REF] Rottman | Gravity currents produced by insantaneous releases of a heavy fluid in a rectangular channel[END_REF] and in our own experiments is reported in figure 3. The various runs correspond to different values of l f SR /L 0 (see (1.8) for the definition of l f SR ). The time exponent deduced from our experiments ranges from 0.15 to 0.27, these extreme two values being obtained with the same set of parameters (namely ρ 2 /ρ 1 , Re, L0 and h0 ). Finally, the evolution of the two currents of figure 3 that exhibit a 1/5 time exponent corresponds to an 'experimental' prefactor in (1.6) higher than its theoretical counterpart β by 18 % and 22 %, respectively. There may be multiple reasons for the large variability observed among experimental results. In particular, small perturbations generated by the opening of the gate may affect the free surface (an effect which does not exist in cylindrical tubes), and there may be small local variations in the slope of the channel. Also, even though the channels used in experiments have large aspect ratios, the flow is actually three-dimensional, owing to the no-slip condition imposed by the lateral walls. We will come back to this point later.

We ran several computations in the above situation with different Reynolds numbers and initial volumes. The corresponding results are displayed in figures 3 and 4. In particular, figure 4 shows the computed evolution of the front position for Re = 200 and two contrasted released volumes: a very short volume L0 = 1.5 and a much longer one L0 = 17. The run with L0 = 17 corresponds to the LR solution and is only reported here for comparison. In contrast, the run with L0 = 1.5 confirms that, after the slumping phase, a current generated by a sufficiently small release evolves in agreement with the SR self-similar solution (1.6). The numerical time exponent and prefactor are found to be 0.198 ± 0.001 and 0.93β, respectively. More generally, all computations leading to an SR viscous regime predict the 1/5 exponent with less than 2 % difference, despite the finite length of the domain. The prefactor is always found to lie between 0.86β and 1.02β. The above comparison indicates that present computations recover the theoretical self-similar prediction characterizing the SR viscous phase. [START_REF] Cantero | On the front velocity of density currents[END_REF] also briefly considered this regime but found the exponent of the front position to be closer to that corresponding to a gravity current flowing over another fluid [START_REF] Hoult | Oil spreading on the sea[END_REF], i.e. lf ∼ t3/8 , and the computational prefactor was more than 3β. We strongly suspect the discrepancy to be due to the O(1) Schmidt number used in these pseudo-spectral computations to ensure numerical stability. Indeed, such a large molecular diffusivity significantly reduces density gradients, which results in two effects. In the horizontal direction, the driving force is reduced, which in turn reduces the velocity of the current, at least for short times (Bonometti & Balachandar 2008). In the vertical direction, the current thickens dramatically, which reduces the velocity gradient at the bottom wall and thus the resistive viscous force. In particular, the characteristic length in the viscous stress is not the current height h(x) anymore but rather the penetration thickness √ νt. Equating the gravitational driving force with the viscous force then leads to the scaling lf ∼ t3/8 found by [START_REF] Hoult | Oil spreading on the sea[END_REF] instead of the t1/5 scaling, in line with the observations of [START_REF] Cantero | On the front velocity of density currents[END_REF]. Note that even though theoretical solutions for viscous currents with Sc = O(1) are not available, it is obvious that the self-similar behaviour (1.6) cannot subsist under such circumstances. This is made clear by considering the equation governing the reduced density (ρ(z)ρ 1 )/(ρ 2ρ 1 ) whose integral across the current height yields the governing equation for h(x, t). When molecular diffusion is absent, this equation is of the Reynolds type, leading to the self-similar solution (1.6). In contrast, when finite diffusion effects are present, it involves an additional diffusion flux at the interface, which leads to a modified evolution of h. Based on this remark, it is clear that SR viscous currents computed with spectral methods making use of an O(1) Schmidt number cannot obey the scaling (1.6).

Beyond such O(1) Schmidt numbers dictated by numerical constraints, an O(1) effective Schmidt number may also be encountered in the front region and along the interface in turbulent currents with immiscible fluids. In this context, one can also question the validity of (1.6) once such currents become viscously dominated, since turbulent diffusion may well have spread the interface between the two fluids in a way pretty similar to molecular diffusion. This could lead to results close to those obtained by considering laminar currents with an O(1) Schmidt number. If this happens to be the case, and we assume that finite-Schmidt-number effects result in a gradual shift in the time exponent of the evolution law from 1/5 to 3/8, it could explain the variability observed among the experiments, since they were generally carried out under turbulent conditions over a broad range of Reynolds number. Three-dimensional high-Schmidtnumber/high-Reynolds-number DNS could help check this conjecture.

Three-dimensional effects

To assess possible alterations of two-dimensional predictions by flow threedimensionality, we ran several three-dimensional computations in a channel of aspect ratio unity with either periodic conditions or no-slip conditions in the spanwise direction. A small disturbance was added to the initial density profile to trigger the transition to three-dimensionality in the former case. When the released volume is large enough for the LR regime to occur, both types of three-dimensional currents exhibit a t 1/2 evolution in the LR regime (figure 5a). Streamwise vortices are well visible on the body of the currents (figures 6a and 7a). Neither their spanwise spacing nor their typical strength seems to be significantly altered by the nature of the boundary condition on the sidewalls. More significant differences can be detected on the shape of the current fronts: as may be seen in figure 6(b), the no-slip condition on the sidewalls generates a marked rim at the extremities of the front, while the central part is pretty flat. The front is found to be steeper in presence of periodic lateral conditions with a clear presence of lobe and cleft structures (figure 7b). In presence of periodic boundary conditions, the front velocity is almost the same as that of the two-dimensional current, even though a tiny decrease of the distance travelled by the front may be detected in figure 5(a). This decrease certainly results from the small additional dissipation induced by the streamwise structures.

In contrast, compared to the two-dimensional reference, the transition from the slumping phase to the LR regime occurs significantly earlier for the current with no-slip conditions on the sidewalls. This is a direct consequence of the increased drag resulting from the presence of the lateral boundary layers. A rough estimate of the impact of this additional drag on the transition time to the LR regime may be achieved by evaluating the time TLR * at which the total viscous force starts to exceed the inertial force assuming that the dimensionless thickness of the lateral boundary layers is δ( t), while the bottom boundary layer is fully developed at that time and has a constant dimensionless thickness K/2. Equating the transition time TLR given by (3.2) with that provided by this momentum balance in the two-dimensional case yields K = (2α/F r) 2 . Then, in the three-dimensional case, the viscous force is increased by 1 + K/(2A δ( t)), where A denotes the channel aspect ratio, i.e. the width-to-height ratio. Provided the lateral boundary layers are still developing when t = TLR * , which is the case as long as TLR * < Re A 2 /4, one has δ( TLR * ) ≈ Re -1/2 T 1/2 LR * . Therefore the transition time now obeys TLR * (1 + K(Re/ TLR * ) 1/2 /2A) ≈ KRe/4.

(3.4) Equation (3.4) predicts that in the case of figure 5(a) (where Re = 660), the transition is advanced by about 43 % with respect to the two-dimensional situation. We determined the transition time graphically from this figure by intersecting the t 1 and t 1/2 lines tangent to the initial and final stages of the computational curves, respectively. This method provides TLR * = 27.1, while the two-dimensional result is TLR = 46.3, indicating that the transition is advanced by about 41.5 %, in good agreement with the above prediction. We also examined three-dimensional effects in currents resulting from much shorter releases leading to a direct transition from the slumping phase to the SR viscous regime. As shown in figure 5(b), the corresponding currents exhibit a t 1/5 evolution, whatever the boundary condition in the spanwise direction.

Possible three-dimensional effects in this regime were already considered by [START_REF] Didden | The viscous spreading of plane and axisymmetric gravity currents[END_REF] who studied currents with a constant inflow in channels of large aspect ratios. Using their estimates for all three forces in balance right at the time of the transition, it can be shown that the two-dimensional estimate (1.9) for this time is reduced by (1 -2Re -1 t/A δ( t)) 3/7 . Therefore, three-dimensional effects affect the corresponding transition time according to TSR *

(1 -2( TSR * /Re) 1/2 /A) 3/7 ≈ Re 3/7 Ṽ 4/7 0 .

(3.5)

In the case of the computations presented in figure 5(b) (where Re = 200), (3.5) predicts that three-dimensional effects decrease the transition time by about 23 %. Again, we also determined the transition time graphically from this figure by intersecting the corresponding t 1 and t 1/5 lines. This technique yields TSR * = 14.2, while the two-dimensional result is TSR = 19.0. The difference corresponds to a 25 % decrease, in good agreement with the above prediction.

To conclude, the three-dimensionality of the flow does not alter the two-dimensional self-similar scaling found for either LR or SR currents. The only discernible effect in both regimes is the decrease of the transition time (or position) due to the enhanced friction exerted on the currents by the lateral walls. Therefore it turns out that the variability of the evolution law of the front observed in the various experiments is not due to the three-dimensionality of the flow.

Prediction of successive regimes in a given current

4.1. Prediction of the earliest viscous regime The results presented in § 3 highlight the possibility to obtain either of the two scalings (1.4) and (1.6) in a given channel, depending on the released volume and flow Reynolds number. Therefore it is desirable to derive a criterion predicting whether the LR or SR viscous regime first succeeds the slumping phase. This goal may be achieved by comparing the transition times to the LR and SR viscous regimes, TLR and TSR . Given the direct transition from the slumping phase to a viscous regime in present computations, the transition position and time for the LR regime are provided by (3.1) and (3.2), respectively. The change of slope at point (1, 1) in figure 2 confirms these predictions. On the other hand, processing the data of figure 3 indicates that predictions (1.8) and (1.9) are fairly well recovered, both in experiments and in computations, when the slumping phase is succeeded by the SR regime. Therefore (1.9) appears as a relevant estimate of TSR . Then, expressing the condition TLR < TSR provides the minimum released volume Ṽ0 c required for the LR viscous regime to succeed the slumping phase. This condition is written as

Ṽ0 c Re = α F r 7/2 , ( 4.1)
where α = 0.0922 is the numerical prefactor involved in (1.4). (Note that since only two-dimensional currents are considered in this section, Ṽ0 may be replaced by L0 everywhere it appears.) We validated this criterion using 17 runs performed at Reynolds numbers ranging from 200 to 790 with various initial lengths L0 smaller or larger than L0c . The result of this validation is displayed in figure 8. Only in a narrow band located along the critical curve is there some doubt on the nature of the viscous regime that sets in after the slumping phase. Interestingly, F r is known to increase with Reynolds number towards [START_REF] Benjamin | Gravity currents and related phenomena[END_REF] prediction, F r = 0.5, and is roughly constant for Re > O(10 3 ), a condition achieved in most real flows. As a consequence, the relation (4.1) becomes quasi-linear at high Reynolds number and simplifies to Ṽ0 c = (3.3 ± 0.6)10 -3 Re for 0.45 < F r < 0.5. This condition can also be reversed to infer that given a released volume Ṽ0 the Reynolds number must be lower than O(300 Ṽ0 ) for the LR regime to occur.

The possible secondary viscous transition

As soon as the SR regime sets in, the current extends from one endwall to its head. Then the reservoir cannot be considered infinite any more, and the condition required to reach the LR regime can never be fulfilled at later time; i.e. the current stays in the SR regime whatever t > TSR . In contrast, if the slumping phase is succeeded by the LR regime, one can wonder whether the current can later experience a second transition from the LR regime to the SR regime. Referring to the characteristic times defined in the introduction, this could in principle be the case at late times if TLR is shorter than the time TS at which the receding current would hit the wall if it were staying in the slumping phase. Under this condition, the LR regime is expected to exist within a finite time interval, say TLR < t < TV , approximately, TV being the time at which the current hits the endwall during the LR regime. Then, for t TV , the conditions for the SR regime to occur are in principle fulfilled. To validate this scenario, we ran long computations with large enough released volumes. Figure 9 confirms that in such cases, the second viscous transition expected from the above argument indeed happens. More precisely, as predicted by (4.1), the current experiences a direct transition from the slumping phase towards the SR regime in cases in which L0 < L0c ( in figure 8, ᭛ and in figure 9). When L0 is close to L0c , the current reaches the LR regime but immediately departs from it towards the SR regime ( in figure 8). Finally, for larger releases, the current first reaches the LR regime and remains in this regime for a while. However, after the receding front has hit the wall, the viscous current eventually switches to the SR regime ( in figure 8, and O in figure 9).

10 100 10 t xf xf ~ t 1/2 xf ~ t ˜l f ~ t 1/5 ˜ Figure 9
. Evolution of the distance travelled by a two-dimensional numerical current corresponding to Re = 200 ( L0c = 2.87) for various released volumes: , L0 1.52; , L0 2.56; , L0 4.04; O, L0 8.04, -• -• -, slumping phase (1.1); --, LR regime (1.4); ---, SR regime (1.6); ᭹, transition from the slumping phase to the LR phase as predicted by (3.1) and (3.2); grey symbols: transition from the slumping phase to the SR regime as predicted by (1.8) and (1.9) for the largest two releases; black symbols: secondary viscous transition as predicted by (4.4) and (4.5) for the largest two releases. To compute these transitions, β was set to 0.986 to be consistent with the computational prefactor obtained in the SR regime.

The physical mechanism at work in this second transition may be better understood from figure 10 which compares the evolution of the current shape during a 'pure' LR regime without any endwall effect (say current 1) with that of a current during the transitional phase from the LR to the SR regime (say current 2). As can be seen in figure 10(b)-(e), a disturbance propagates downward after the receding front of current 2 hit the wall, even though the current dynamics are already fully controlled by viscous effects. Nevertheless, the velocity of the front of current 2 does not change immediately (figure 10b,c). Rather, the fronts of both currents still follow the same evolution until the disturbance reaches the front of current 2 (figure 10d ). At this point, the front of current 2 starts moving more slowly (figure 10e) and eventually joins the t1/5 evolution. It is of interest to notice that the above mechanism is similar to that observed by [START_REF] Rottman | Gravity currents produced by insantaneous releases of a heavy fluid in a rectangular channel[END_REF], even though in their case the current dynamics were still dominated by inertia at the time the receding current hit the wall. In particular they noticed that the advancing current did not 'feel' the endwall immediately but rather started to be affected only when the bore generated by the interaction of the receding current with the endwall reached the front. In the present computations, when the current is still dominated by inertial effects at the time of the reflection, the disturbance also takes the form of a bore propagating at the surface of the heavy fluid. This is easily discernible on space-time diagrams of the current height (see figure 6 of [START_REF] Rottman | Gravity currents produced by insantaneous releases of a heavy fluid in a rectangular channel[END_REF] for instance). In contrast, when the current is already in the LR viscous phase at the time of the reflection, as in figure 10, the disturbance merely consists in a local modification of the shape of the current that can only be detected through a comparison with the profile of an undisturbed current. The position of the front of the disturbance is easily detected in figure 10(b, c), whereas it is hard to infer from figure 10(d,e) the exact time at which the flow regime of current 2 changes. A much better tool for this is provided by space-time diagrams of the local height difference between current 1 and current 2 such as that of figure 11 to be discussed later. Using such diagrams, we found that the regime of the current always changes right at the time at which the disturbance catches the front. The secondary viscous transition is expected to take place at the time t = TLS at which the viscous force corresponding to the LR solution starts to exceed that corresponding to the SR solution. During the LR regime, the dimensionless force associated with the bottom viscous stress and the buoyancy force balance each other. The typical magnitude of the former force is Re -1 ũf xf / hLR , while that of the latter is h2

LR /2, where hLR and ũf denote the dimensionless averaged thickness and front velocity of the current, respectively. Using (1.4) and noting that ũf = dx f /d t, this first balance implies hLR = α 2/3 . (4.2)

Then, during the SR regime, the length lf of the current is given by (1.6). The above momentum balance then yields hSR ( t) = 2 5

1/3

β 2/3 Ṽ 2 0 Re 1/5 t-1/5 . (4.3)
The viscous forces in the two regimes must match at the transition time t = TLS , so that the balance of (4.2) and ( 4 Interestingly, setting xf = -L0 in (1.4) indicates that the receding current hits the endwall at t = TV = ( L0 /α) 2 /Re. Hence the prediction (4.4) implies that TLS is about 7.9 times larger than TV . Figure 9 qualitatively confirms the above prediction for TLS , even though the transition time is difficult to locate precisely, owing to the progressive nature of the transition.

4.3. Late occurrence of the LR regime An important ingredient in the derivation of (1.4) lies in the symmetry of the problem with respect to the lock position. This suggests that the LR regime can only be reached as long as the endwalls are far enough from the front of the advancing and receding currents. However our observations indicate that under certain conditions, a transitional LR evolution can arise even after the receding current has hit the endwall. Figure 11 illustrates the existence of such a late transition in the form of a space-time diagram representing the local height difference h( t, l) between current 1 with L0 1 (for which endwalls have no influence) and current 2 with L0 = 4.04 (which corresponds to the square symbols in figure 9). In the latter case, the receding front hits the endwall at t = TS ≈ 12, while the transition from the slumping phase to the LR regime is only expected to occur at t = TLR ≈ 18 according to (3.2). Indeed, for 0 6 t 6 TLR , the front position of both currents is found to evolve linearly (see the dashed line in figure 11), confirming that current 2 is still in the slumping phase when the receding front hits the endwall. As expected, when t > TLR , current 1 follows the theoretical prediction (1.4) for the LR regime (indicated by the solid line in figure 11). As soon as t > TS , the height difference h becomes non-zero in the upper part of the currents, say in the range 0 6 l 6 lD ( t), where lD ( t) denotes the downstream point reached by the reflected disturbance. As shown by the inclined line issued from the point ( t = TS , l = 0), this disturbance propagates with a constant speed and reaches the front of current 2 at t = TD ≈ 47. Within the time interval TS < t < TD , the evolution of the front position of both currents is found to be identical, as already observed in the example of figure 10. However the novel feature in figure 11 is that, even though current 2 was not yet in the viscous regime when the endwall started to influence its dynamics, its evolution closely follows the theoretical prediction for the LR regime for TLR < t < TD . This simply confirms and extends the conclusion of figure 10 that the evolution of the downward front is not influenced by the endwall as long as it has not been reached by the reflected disturbance. In the present case, the hierarchy of times is such that TS < TLR < TD , so that the front of the advancing current (more precisely, the downward part of the current corresponding to lD < l < lf ) switches to the LR regime for t ≈ TLR as if there were no endwall. Then for t > TD , current 2 slows down, and the position of its front departs gradually from the LR prediction until it eventually matches the SR prediction (1.6) for t ≈ 170. The situation of current 2 at t = TD is thus similar to that of any current experiencing a secondary viscous transition from the LR to the SR regime. Based on this remark, one expects TD to be identical to TLS . In the present example, (4.4) predicts TLS ≈ 45.9, which is indeed very close to the value of TD determined graphically from figure 11. Conversely, assuming TD = TLS provides a direct means to evaluate the speed of the reflected disturbance: at time TLS , (1.4) predicts lf ( TLS ) = L0 + αRe 1/2 T 1/2 LS . Hence the disturbance travels with a speed cD = lf ( TLS )/( TLS -TS ), where TS = L0 F r -1 , so that cD ≈ 0.40 in the above example.

Summary and concluding remarks

To obtain some new insight on the dynamics and evolution of viscous gravity currents, we carried out extensive Navier-Stokes computations in full-depth lockexchange configurations, for various released volumes and Reynolds numbers in the range 200 6 Re 6 800. The specificity of the numerical technique employed here lies in the treatment of the equation for the volume fraction of one of the fluids (or the density field) which allows us to deal with strong gradients. Hence, in contrast with pseudo-spectral approaches, the present numerical approach makes it possible to mimic the behaviour of almost immiscible fluids, which is crucial for the comparison with experimental results and theoretical predictions in the viscous regime. We first validated our code by predicting the evolution of currents propagating in horizontal cylindrical tubes, using the theoretical and experimental results of [START_REF] Séon | Transient buoyancydriven front dynamics in nearly horizontal tubes[END_REF] and [START_REF] Takagi | The effect of confining boundaries on viscous gravity currents[END_REF] in the LR and SR regimes, respectively. We then confirmed the prediction of [START_REF] Hinch | Inclined to exchange: shocking gravity currents[END_REF] for the LR regime of planar currents and that of [START_REF] Huppert | The propagation of two-dimensional and axisymmetric viscous gravity currents over a rigid horizontal surface[END_REF] for the corresponding SR regime, whereas available experiments including ours exhibit a significant variability of the time exponent in the latter regime. We also explored the influence of three-dimensional effects on the current dynamics. In particular, we showed that the boundary layers that develop along the sidewalls result in an earlier transition to the viscous regime. We obtained a theoretical estimate of the transition time to both viscous regimes through a force balance and showed that these predictions agree well with the computed evolutions.

We finally considered the successive transitions that may happen during the 'life' of a given current. First, we derived and validated a quantitative criterion predicting which minimum (respectively maximum) released volume Ṽ0 c is required for the LR (respectively SR) viscous regime to succeed directly the slumping phase. In cases in which the LR regime occurs, we then showed that the current later experiences a secondary viscous transition from the LR to the SR viscous regime and determined the corresponding transition time and position through a force-balance argument. Computations made it clear that in all cases, the current starts to obey the SR evolution (1.6) only at the time at which the reflected disturbance reaches the front of the advancing current, rather than at the time the receding current hits the endwall. In some cases, we observed that the LR regime sets in during the time interval separating these two events.

Keeping in mind that we focused on moderate Reynolds numbers for which the inertial regime governed by (1.3) does not occur, we may summarize our conclusions about the various 'routes' a given current may follow within a single diagram of states. This diagram is shown in figure 12. In brief, the SR regime directly succeeds the slumping phase only when the transition time TLR is larger than TSR . In contrast, an intermediate LR regime exists every time TLR < TSR . In such cases, the LR regime may set in before or after the receding front hits the endwall, depending on TLR being smaller or larger than TS . Then the current experiences a second transition to the SR regime at t ≈ TLS .

Figure 2 .

 2 Figure 2. (a) Evolution of the front position in a two-dimensional current reaching the LR viscous regime: computational results for Re = 200, 400, 600 and 790 (black lines) all collapse on the master curve defined by xf = t for t < 1 and xf = t1/2 for t > 1; (b) theoretical and computational predictions for the current profile corresponding to Re = 200; ---, theoretical prediction (3.3); , computational results taken at 45 successive times during the viscous phase; • • • •, two profiles obtained during the development of the viscous phase.

3 .

 3 Distance travelled by the front for various releases leading to the SR regime ( h0 = 1 except specified otherwise). Experiments of Rottman & Simpson (1983): ᭺, l f SR /L 0 = 6; ᭜, l f SR /L 0 = 10; M, l f SR /L 0 = 12; H, l f SR /L 0 = 18; , l f SR /L 0 = 25. Present experiments: I, l f SR /L 0 = 9.5; * , l f SR /L 0 = 10.8, h0 = 1/6; ×, l f SR /L 0 = 10.7; F, l f SR /L 0 = 10.8; +, l f SR /L 0 = 10.8, h0 = 1/6. Computational results: --, x f SR /L 0 = 3.04, Re = 200; -• -• -, x f SR /L 0 = 2.5, Re = 790; ---, lf ∼ t1/5 .

4 .

 4 Distance travelled by the front of a two-dimensional current as a function of time for Re = 200: --, LR regime reached with L0 = 17); -• -• -, SR regime reached with L0 = 1.5. The thin dashed (respectively dot-dashed) line corresponds to the prediction (1.4) (respectively (1.6) with a prefactor set to 0.93β).

5 .

 5 Comparison of the evolution of two-and three-dimensional currents at Re = 200 for a large release (LR) with L0 = 17 (a) and a small release (SR) with L0 = 1.72 (b) in a channel of aspect ratio A = 1: --, two-dimensional current; ---, three-dimensional current with periodic conditions in the spanwise direction; -• -• -, three-dimensional current with no-slip conditions on the sidewalls. The thin dashed line in figure 5(a) (respectively b) corresponds to the prediction (1.4) (respectively (1.6) with a prefactor set to 0.914β).

Figure 6 .

 6 Figure 6. Structure of a three-dimensional current in the LR regime with no-slip boundary conditions in the spanwise direction (Re = 660). (a) Vortical structures coloured by the strength W of the streamwise vorticity normalized by √ g /H ; (b) current surface plotted using the (ρρ 1 )/(ρ 2ρ 1 ) = 0.65 isosurface. The spatial scales are identical on both axes.

Figure 7 .

 7 Figure 7. Same as in figure 6 with periodic boundary conditions in the spanwise direction.

Figure 8 .

 8 Figure 8. Critical length L0c of the released volume required for the LR regime to succeed the slumping phase. Solid line: relation (4.1) evaluated using the Froude number obtained numerically. Each symbol corresponds to a two-dimensional computational run: , the LR regime succeeds the slumping phase; , the LR regime cannot occur and the slumping phase is directly succeeded by the SR regime; , intermediate cases for which the nature of the viscous regime that succeeds the slumping phase is difficult to determine unambiguously.

Figure 10 .

 10 Figure 10. Evolution of the profile of two currents corresponding to Re = 200: --, current 1 with L0 1 (the left wall is rejected far outside the figure); ---, current 2 with L0 = 8.04 (the left wall corresponds to the left side of the figure). Current 1 in the LR on snapshots. (a) Current 2 is in the LR regime. (b) and (c) The receding current has reached the endwall; the profile of current 2 is altered behind the front of the disturbance while its head is still in the LR regime. (d ) The disturbance just reaches the head of current 2. (e) Current 2 slows down and switches to the SR regime.

Figure 11 .

 11 Figure 11. Space-time diagram showing the height difference h(x, t) (in greyscale) between current 1 with L 0 1 and current 2 with L 0 = 4.04 (both currents correspond to Re = 200): ---, position of the front of current 2; --, LR solution (1.4) with α = 0.094; -• -• -, SR solution (1.6) with β = 0.974; , position of the transition from the slumping phase to the LR regime as predicted by (3.1) and (3.2); F, position of the transition from the slumping phase to the SR regime as predicted by (1.8) and (1.9).

  Figure 12. Diagram of the various possible evolutions of a given current according to the initial released volume L0 (or Ṽ0 ), Froude number F r and Reynolds number Re; Sl, LR and SR stand for slumping, LR and SR regimes, respectively. The numbers refer to the equation providing the expression of the corresponding transition time in terms of the above three parameters or to that providing the corresponding evolution law.

We thank Annaig Pedrono for her continuous and central support in the development of the parallel version of the code used in this research. The contributions of CALMIP, CINES and IDRIS, which provided most of the CPU time, are greatly appreciated.