Matthias Gille Levenson

Matthias Gille

Levenson Editing

Editing a 15th century political treatise using the computer: a back and forth between meaning and information

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

Editing a 15th century political treatise using the computer: a back and forth between meaning and information *

Matthias GILLE LEVENSON

Introduction: towards a comparative edition

I am working on the edition of part of the Regimiento de los prínçipes, a political treatise of the XIV th century, that is itself the translation of the famous De Regimine Principum of Gilles of Rome. I am focusing my study on the discourse of chivalry and nobility in this text. I am influenced by the so called "New Philology" 1 and by the works of Roger Chartier on the materiality of the texts: we cannot put aside this materiality when studying a text. As Roger Chartier points out, La sociologie des textes telle que l'a définie D.F McKenzie conduit à tenir chaque état d'une oeuvre comme l'un de ses états historiques, qu'il faut comprendre, respecter et, possiblement, éditer. Le concept d'un "ideal copy text", existant en deçà ou au-delà des différentes formes imprimées (ou manuscrites) d'une oeuvre, est une illusion que la critique textuelle doit abandonner au profit de l'analyse des effets produits sur le texte, ses lecteurs et, éventuellement, son auteur, par chacune de ses existences matérielles 2 .

My aim is to study the discourses that I have in the different manuscripts of my textual tradition, as well as its homogeneity: in other words, I want to evaluate how much of a text I have in the Regimiento. What can we say about the use of the different manuscripts that circulated in these 150 years ? What do they say, individually ? How were the texts read, received, or transformed by the actors that decided to re-produce them ? From an ecdotic point of view, I am not interested in the ideal text. I need to produce as many texts as I have witnesses. For this reason, I plan to do not a critical edition, but rather a comparative edition: I will not choose between "good" and "bad" variants in order to access and produce the original text. This presentation will focus on the methodology that is being used to produce this edition.

Data, information, meaning

Let's compare versions of an extract of the Regimiento. G, J, Z are the sigla of three witnesses, two manuscripts and an incunable. We are in the gloss of the chapter 21 of the last part of the text, that speaks of the funeral honours given by Alexander the Great to Darius' wife.

E llamó ý un muy grant sabio de los iudíos que dezían Apelles, que era de Jherusalém (J) E llamó ý un grand sabio de los iudíos que dezían Apelles, que era de Gerusalén (G) E llamó ý un muy grand sabio de los indios que dezían Apelles que era grand pintor (Z) One possible representation of the difference between these three sentences would be:

8 muy J Z | om. G 9 iudíos J G | indios Z 9 Jherusalém J G | grand pintor Z
The second apparatus entry means that in the witness Z, "iudío", was replaced by "indio". That is meaning, a translation made by the human mind of some information: a sequence of characters, a pipe "|" and another sequence of characters; this information is meaningful for those who know the conventions of the critical apparatus. This is one difference between us and the computer. The computer 'sees' the information, and may process it. From this, we get or we create a meaning by reading the information. In an article published in 2007, Chaims Zins analyses how searchers define data, information and meaning more than one hundred papers presented at an conference on Information Science. highlights that there are as many definitions of these concepts as people who try to define them3 . I decided to choose one definition from the 130 that the article presents.

Data are the basic individual items of numeric or other information, garnered through observation; but in themselves, without context, they are devoid of information. Information is that which is conveyed, and possibly amenable to analysis and interpretation, through data and the context in which the data are assembled 4 .

As for the term "meaning", the Cambridge Dictionary tells us:

The meaning 5 of something is what it expresses or represents 6 .

Meaning is information that makes sense. Meaning is the information when read by the human mind: this is why the same document can be either meaning or information, depending on what one does with it. This conceptual distinction is important to me, because it helps me to understand why and how I do the things I do to produce my texts. In this presentation I will talk of the dialectic between meaning and information, that is: what are the steps I have to take before the production of meaning, the publication of my edition ? How the conceptual distinction between 'information' and 'meaning' can help us (can help me) understand what I am doing ? How to describe my work in terms of back and forth between meaning and information ?

My use of the machine is based on the following assumption: computers can manipulate information better than we can, and they can help us in producing more homogeneous or consistent documents. But ultimately what is important for us as humans is to deal with meaning. How do we get from information to meaning? We will see that we have to complete successive tasks of translation between meaning and information to get the desired result, that is: a comparative edition for each of the transcribed witnesses, with a meaningful and readable apparatus 7 .

The workflow

The traditional way to edit a text is as follows: transcription of a witness, collation with the other witnesses, creation of the apparatus 8 . My own approach adds two steps (which are absolutely not original): the TEI encoding and a phase of addition of some grammatical information: lemmas 9 and Part of Speech 10 (PoS). If we try to think in terms of manipulation of the information, the steps I take are:

1. Acquisition of the information: the transcription of the texts 2. Structuring of the information: the TEI encoding 3. Enrichment of the information: lemmatisation and PoStagging 4. Production of more meaningful information: collation and creation of the apparatus 5. Translation of the information: production of readable documents: pdfs (and/or web-based interface)

The presentation will follow this workflow, step by step; I will try to show the links between all of them.

Acquiring the information: the transcription

If realized with the help of the computer, the transcription is an acquisition of the information, or a first translation between data and information: after that, we only have more or less structured information.

We use programs that are called OCR programms: optical character recognition, based on neural networks. I use Ocropy for recognizing the incunable. How does it work ? We train a model, and use this model to recognize new text. The model is trained by feeding it with veridic data, in our case text and image. This data is called ground truth. With this ground truth the programm will start creating a recognition model. With this first model, we will predict new text, and correct it to transform it into new ground truth, re-train the model, etc. Once a given error rate is reached, we can stop an use the model for transcribing the whole text. As surprising as it can seem, the ideal rate has not to be the lowest possible, because we don't want the model to be overfitting on the training data: we want it to be accurate on new data.

I will insist on the importance of the correction: a prediction using neural networks is based on probabilities, and thus is prone to errors 11 . 8 Depending on the methodology one chooses, we may add to this list the production of a stemma, and the selection of the "good" variants. I will not talk about these two steps in my presentation: I am not interested in re-creating a virtual original text; I am not advanced enough in my work to talk about the stemma codicum.

9 https://en.wikipedia.org/wiki/Lemma_(morphology). 10 https://en.wikipedia.org/wiki/Part_of_speech. 11 I have used the tool called Ocropus to predict the text of the incunable. The error rate would be too high for using it to predict the text of manuscripts, but more recent tools are being developped and are reaching an acceptable error rate with manuscripts: see Kraken (Benjamin KIESSLING, "Kraken -an Universal Text Recognizer for the Humanities", in: , DH2019 : Complexity, Utrecht, 2019) (open-source) or Transkribus (Philip KAHLE et al., "Transkribus-a service platform for transcription, recognition and retrieval of historical documents", in: 2017 14th IAPR International Conference on Document Analysis and Recognition (ICDAR), vol. 4, IEEE, 2017, pp. 19-24) (not fully open-source).

Structuring the information. Digital encoding and TEI

If not the most important, the XML-TEI encoding is one of the most important steps of my workflow, mainly because of the standard characteristic of the TEI12 . I will show in this part that the difference between information and meaning is not as obvious as it can seem first: a TEI conformant file, as the ideal description of a given text, is processable by the machine, but also readable and can be understood by the human being, thanks to its semantic nature. Why ? Because the TEI is a "conceptual model of textuality", as says Fabio Ciotti13 . It proposes an ontology on the structure of the text, not informatically speaking (the TEI is not what informaticians call an ontology), but rather in its philosophical meaning. The TEI data model doesn't explain or states on what is a text, but it does say what is inside a text, what can be the different components of a text: the different elements, but also their interelation. TEI gives us a semantism. It is human-readable, but also machine-readable: in other words, it is meaning, but it is also information that can be processed by the computer, to create multiple outputs, or indices, to extract some particular information, etc. Last but not least, it is a community driven standard: it evolves with the needs and the reflections of the community of users.

The paradox here is that the TEI encoding, for my project at least, is the step that can be the less automatized, precisely because of this semantic nature: we are here in the first phase of construction of the meaning, in the most important transition from information to meaning, a phase that has to be realized by the human being, by the editor.

Enrichment of the information: Lemmatisation and PoStagging

For my project, I need to be able to make efficient searchs within my text. For medieval languages, that can be tricky because the spelling of the words is not fixed yet: I need to add some information that will help me go over these issues. Grammatical annotations can be helpful here: I am talking about lemmatisation and PoStagging, that I mentionned before. The lemmatisation is the act of indicating for each word the canonical word it refers to; the PoStagging is about indicating the morphology of the word, with a code that indicate its nature, its gender, the mode if applied, the person, etc. This code, named tagset, is a convention: in my case I am using the EAGLES tagset14 . I will give an example: the legal aspects of my text have to be studied, as Jesús Velasco helped me understand. In Spanish, one particular mode and time is very often used in legal contexts: the future subjunctive. With a non PoStagged text, it would be really difficult to extract all of the occurrences. With a PoStagged corpus, if we want to find the verbs that are in future subjunctive, we will try to find all of the words whose PoS will match the code for future subjunctive. This is the step I call enrichment of the structured information: the interesting thing is that I can add this grammatical information to the TEI encoding: it will be really helpful in order to process the text, in particular to collate it, as I will show.

Comparison of the information. The collatio.

How do we get from n individual transcription to a compared text ? With this step I will be able to get to some meaning, and to answer to some my questions: how much of a text do I have ? What is the importance of the variation ? The collation is a way to add relevance to the information: from several witnesses, or, in my case, several transcription, I will produce a single document -per witness -that gathers all of the information about the variation. The process of collating several witnesses is the following: we want to compare the texts, so we need to find the equivalent texts (to align the texts) and to compare the content of each passage. These two steps are virtual: a human mind cannot dissociate them, but I need this distinction to simulate the process of collation that I will explain in detail.

The important idea is that here the words are represented by their attributes, form (the sequence of characters of the word), lemma and PoS. This deconstruction, as we will see, is a translation too, this time from meaning (the concept of "word") to information (the attributes that can identify a word). Another point is that we use here the computer to create more meaningful content, to clean up the apparatus from unuseful information: the computer deals with information to make the meaning emerge more easily. However, the created documents will not be critical editions, something that only the human mind can achieve; however this method could be used to create critical editions. It has issues: the main one here is related to what is called "overlapping", a situation that is not allowed by the specifications of the XML format, and that might be generated here.

Once the collation phase is done, the program will transform the result into XML, and inject this XML information to the individual TEI structured transcriptions. The interests of this method are multiple: it makes me gain a lot of time, and it allows me to have each particular manuscript as a base, with the apparatus as the result of the comparison between them.

Transformation of the information: the creation of the output

The final creation of a meaningful object: the edition. I use the typesetting language/program L A T E X for the pdf output, and I will show samples of the result. An illustration of what can be the web-based output of the edition15 , on the passage I've been talking about, can also be found here.

Conclusions

In order to use the computer, there is a series of translations to be made between meaning and information. One of them can be considered as a deconstruction of the aforementioned objects: the word, for example, will be represented by its attributes (string and grammatical information, lemma and PoS). The same can be said of philological methods: the deconstruction of the cognitive process of the collatio, for example, is what we have to look after if we want to simulate this ecdotic method. As we can see, the computer is here to accomplish the majority of the mechanisable tasks, and these tasks alone. The result is always to be controled by the editor: as I've said before, the correction phase is essential.

The very last translation operation would be, of course, the close reading of the text, in order to create the final meaning: I am talking about the comments on the edited sources, the construction of the discourse upon all of the texts. The apparatus shows a difference between "iudío" (jew) and "indio" (indian). How can I explain it ? The history of the incunable is important here: it was printed in 1494 by people that were somehow close to the Reyes Católicos, as I have demonstrated in my article; the edition of 1494 could be part of the cultural and political project of the Catholic Kings 16 . The change made by the incunable is not likely to be a correction; it seems to be a conscious and wanted modification, not an error or a misreading 17 : if there is no mention anywhere of the Greek painter Apelles being Jewish, neither is there any other mention of him being related to the Indians! This is why I think this modification can be interpreted as an act of censorship, the silencing of a minority that had been expulsed from the Peninsula two years before.

After the production: what comes next ?

I would like to highlight three concepts that are important to me, and to ask a question, to start the discussion.

Three ideas

Concerning the production of a scientific document as a digital edition can be considered, the most important points are, in my view, the: Citability of the information. A published and accessible web-based interface, as Elena Pierazzo calls them 18 is very often a work in progress: if we take this into account, that means there is a risk of inconsistency of the citations: how can one be sure that the text one cites will be the same in two months, or in a year ? This leads to the need to have consistent and fixed documents when they are published, or at least to record the different versions of the document. One solution is to use a version control software like git: each version is identified with the git hash, the "unique" identifier for the version of a project.

Perennity of the information. The lifespan of the web-based interfaces depends on funding a project can acquire. After a certain time, the funding ceases and we cannot be assured that the produced text will be accessible for the scientific community, either because the server doesn't exist anymore, or because the code is obsolete and cannot be read by modern browsers 19 . This is a real problem for the consistency of the scienfic works, a problem that takes its source in the separation between the information and its support: what is the value of an article based on a text that is not available anymore ? This is why the source (XML, in this case), as the bearer of the semantic value, should be available and stored to be kept for the long term, and why I think, at the very least, a pdf version of the edited text should be made available.

Accessibility of the information. The standardization of the information is an important step towards its processability; but if you don't know how to access this information, it can be considered as a lost/unuseful information. The International Image Interoperability Framework (IIIF 20 ) 16 Matthias GILLE LEVENSON, "L'évolution du Regimiento de los prínçipes (1345-1494), conditionnée par le pouvoir politique ?", in: Écritures du Pouvoir, ed. by Véronique LAMAZOU-DUPLAN, Ausonius Edition, vol. 2, Scripta Medievalia, 2019, pp. 137-148. 17 Regardless of the importance of the variant, we have here a good example of what the computer can and cannot do: when transcribing the text, I couldn't tell if the word was "iudío" or "indio". The context, and the indication that Apelles was from Jerusalem (information that doesn't appear in the incunable), helped me decide what was written in each witness. 18 Elena PIERAZZO, Digital Scholarly Editing: Theories, Models and Methods, Ashgate Publishing, 2015. 19 I don't think it is a good idea to rely on external initiatives like archive.org to make sure a scientific document is accessible. 20 https://iiif.io/ proposes standardized ways to retrieve an image and its metadata. The requests to access the images from a server are documented and are the same for each implementation. Distributed Text Services (DTS21 ) does the same, but with the text instead of the image: it provides guidelines to create a standard API (application programming interface) with a standard grammar to access a text. Hopefully, the use of this new standard will increase the interoperability of the texts, and will ease their reusability.

A question

Question: Where is the edition ? Is it in the TEI encoding, or in the output ?

Chaim ZINS, "Conceptual approaches for defining data, information, and knowledge", Journal of the American Society for InformationScience and Technology 58.4 (2007), pp.

479-493. 4 Ibid., p. 481. 5 I've decided to choose the term "meaning" because the conceptual opposition information vs meaning is far more productive that information vs knowledge: as far as we are concerned, the computer doesn't get the meaning of a text.6 https://dictionary.cambridge.org/dictionary/english/meaning7 The graphical variant will have to be excluded from the apparatus.

For Text Encoding Initiative, see here.

Fabio CIOTTI, "A Formal Ontology for the Text Encoding Initiative", Umanistica Digitale 2.3(Nov. 18, 2018).

See here. I have been using a dictionnary of forms, lemmas and PoS created for the project Freeling (Lluís PADRÓ, "Analizadores Multilingües en FreeLing", Linguamática 3.2(2011), pp. 13-20), that decided to follow this tagset.

The apparatus are still to be added.

https://w3id.org/dts


Emeric Tapachès

Richard Perez

Philippe Lauret

Marc Perez

Mathieu David

Mitigation of the Variability of a PV Fleet via Geographical Dispersion and Energy Storage Systems on the Reunion Island Non-Interconnected Grid

In this paper, the geographical dispersion of PV installations and Energy Storage System (ESS) are studied to smooth PV production on Reunion Island. Despite the small area of the Island, the results show that with dispersion at the total grid scale, the variability driven by weather conditions (unpredictable variability) is almost removed. Then a centralized EES management at grid scale and a few large EES units could mitigate the variability even further with a relative low cost.

I. INTRODUCTION

In 2017, the share of renewables in the electricity generation of Reunion Island was 32.4%. The rest comes from imported fossil energies (coal, fuel and diesel) which negatively impact both the local environment and the economy. In order to be no longer dependent on imported fossils fuels, the local authorities have the ambition to become self-sufficient for the production of electricity by 2030. To achieve this goal, a mix of available local renewable energies (i.e. biomass, hydro, solar, wind, geothermal, wave, OTEC) will probably be realized. Whatever the future energy mix, it will still have to provide inertia and frequency control reserves to the grid, and having a part of dispatchable generation in order to balance electricity supply and demand at any time.

The authors of this paper will focus on the photovoltaic (PV) and the penetration of this production into the island grid. PV has many advantages : more and more cost-effective, free energy source, easy to install, very low maintenance and operations, good life cycle, can be distributed on the grid and installed on roofs, parking, along roadsides, on greenhouse or over water bodies. This last point is attractive for Reunion Island because land pressure is high. On the other hand, a few drawbacks limit the high penetration of PV because the current production is: intermittent (driven by solar path), highly variable (affected by the weather), partially controllable (e.g. only curtailment is possible), difficult to forecast. PV changes the grid architecture and its regulation in comparison with the conventional generation means (e.g. nuclear, coal or hydro). As explained by [START_REF] Passey | The potential impacts of grid-connected distributed generation and how to address them: A review of technical and non-technical factors[END_REF] and [START_REF] Palmintier | On the Path to SunShot: Emerging Issues and Challenges in Integrating Solar with the Distribution System[END_REF], PV could negatively impact the grid balance/management. These impacts concern: planing the supply-demand balance (i.e. regulation of the grid frequency), voltage regulation, reverse power flow, inertia, frequency control reserve, protection of people and property, unintentional islanding and other. In Reunion Island, by experience and at the moment, the electric system operator cares mainly about grid frequency stability and therefore the aggregated PV production is a major concern. Concerning the local voltage fluctuations, they are contained by limiting the power factor of each installation having more than 36 kVA.

To secure the injection of PV on the grid, some operators limit the maximum instantaneous penetration rate of intermittent productions (e.g. about 30% for the non-interconnected grids in France) and/or limit the rate ramp of these productions (e.g. 10% of rated capacity per minute in Puerto Rico [START_REF] Vahan | Review of PREPA technical requirements for interconnecting wind and solar generation[END_REF]). Some approaches to mitigate PV production variability (i.e. to limit its ramp rate) exist: geographical dispersion, curtailment (only for positive ramp rate) or storage. Several papers have shown the smoothing effect of PV dispersed generation for different locations like Ontario [START_REF] Marcos | Smoothing of PV power fluctuations by geographical dispersion[END_REF], Andalusia [START_REF] Santiago | Influence of photovoltaic installation angles and geographical dispersion in the smoothing of photovoltaic fleet power fluctuations[END_REF], Hawaii [START_REF] Perez | Mitigating Short-Term PV Output Intermittency[END_REF] and Reunion [START_REF] David | Spatial and Temporal Variability of PV Output in an Insular Grid: Case of Reunion Island[END_REF]. Energy Storage System (ESS) to smooth PV production is also largely investigated in the literature. [START_REF] Sukumar | Ramp-rate control smoothing methods to control output power fluctuations from solar photovoltaic (PV) sources-A review[END_REF] represents a deep review of these approaches and classify three types of methods: (i) moving average and exponential smoothing based methods, (ii) filter based methods, and (iii) ramp-rate based algorithms. They found that moving average is a simple method largely used by researchers, however this method overbuild the capacity of the ESS and decrease its lifetime. Filter based methods and ramp-rate based algorithms have better performances. [START_REF] Gavriluta | Storage system requirements for grid supporting PV-power plants[END_REF] arrived at similar findings and also showed that ESS sized for power variability mitigation can also provide inertia emulation and primary control reserve, which increase the value of ESS for the grid. This paper presents a first investigation of possible solutions to mitigate aggregated PV production variability on Reunion Island in order to facilitate its large-scale integration. To this end, geographical dispersion and smoothing of aggregated PV output with ESS are applied at different spatial scale. The aim is to quantify the spatial footprint and the ESS size needed to limit the variability.

The paper is organized as follows: Section II presents the electrical grid of Reunion island while section III describes the methods and result of the mitigation of PV output variability. Section IV closes the paper with some conclusions, recommendations and perspectives.

II. CASE STUDY: REUNION ISLAND GRID

Reunion Island is a French overseas territory located in southern of the Indian Ocean. It has an area of 2,512 km 2 , its highest point is 3,070 meters above sea level and about 860,000 inhabitants live mainly on the coast. There is a unique Transmission/Distribution System Operator (TSO/DSO), EDF, and many energy suppliers (EDF, Albioma, Akuo, Quadran, individual producer, etc.). Fig. 1 shows the yearly global horizontal radiation and the 63 kV transmission grid which is connected to twenty-three 15 kV distribution grids. This island experiences a high spatial variability of solar resources and a very good potential on the coast. In 2017, the total electricity production was 2,985 GWh with the following sharing: (i) diesel/fuel/coal (67.6%), (ii) hydroelectricity (14.1%), (iii) bagasse (8.7%), (iv) photovoltaic (8.6%) and biogas (0.4%). In terms of power, the maximum instantaneous PV penetration was 37.8%. The cost of electricity in the Island is higher than metropolitan France. Hourly production was 201 €/MWh on average with a maximum of 454 €/MWh in 2017. In spite of these production costs, the price of electricity is around 140 €/MWh thanks to a national mechanism of tariff equalization.

In non-interconnected French territories, to avoid the risk of failure in the electricity supply, the grid operator (i.e. EDF) can legally disconnect PV systems higher than 3 kVA when the instantaneous PV penetration on the grid exceeds a regulatory threshold. Recently EDF has increases this threshold from 30 to 35% by increasing the resilience of the grid ( i.e. increasing the primary frequency control reserve and modifying the offloading plan). Mitigation of PV output ramps will be another way to increase this threshold.

III. MITIGATION OF PV OUTPUT VARIABILITY

A. Input data

For this study, we used gridded hourly global irradiance spanning the year 2011:

• Hourly Global Horizontal Irradiance (GHI) ground measurement from 34 pyranometers distributed throughout the Island • Hourly GHI satellite estimations using the SUNY model [START_REF] Perez | Validation of the SUNY Model in a Meteosat Environment[END_REF] from Meteosat-7 images (pixel: 5.4 km by 5.4 km) PV output was derived from these data with the System Advisor Model (SAM) software [START_REF] Gilman | SAM Photovoltaic Model Technical Reference Update[END_REF] of the US National Renewable Energy Laboratory (NREL).

B. Geographical dispersion

It is well established that geographical dispersion of PV generation connected to a common grid smoothes the aggregated PV output signal. However, is this "smoothing" effect significant on Reunion Island where the largest dimension does not exceed 70 km? Assuming a uniform spatial distribution of the PV fleet, the effect of geographical dispersion is equivalent to the spatial averaging of the GHI. Thus, GHI data have been averaged considering different footprints:

• A single point (ground measurement stations),

• 1, 4, 9, 16 and 25 pixels covering the island,

• All pixels covering the electric grid (55 pixels),

• All pixels covering the island (84 pixels), Amongst theses footprints, four representative scales are highlighted in this paper: (i) a single point, (ii) the largest distribution grid area (i.e. 9 pixels and 270 km 2 ), (iiii) the total grid area (about 1,650 km 2 ) and (iv) the area of the island (2,520 km 2 ).

Solar variability was evaluated thanks to the standard deviation of the hourly ramps of the GHI ( ). The ramp is 2 defined as follow:

GHI where i ε Ω (1) ΔGHI i 1hr = GHI i - i-1
where is the set of diurnal hours of the year and N the Ω number of elements of . It is also interesting to evaluate the Ω "smoothing" effect for the extreme ramps because the power size of ESS depends on these rare events. Thus, we also computed the maximum (3) and the 99th percentile ( 4 It seems that despite the small area of the island, the total PV production can benefit from geographical dispersion of PV installations. This is related to the high diversity of microclimates due to the specific relief of the island. The variability (i.e. standard deviation) and the extreme ramps decrease rapidly with the footprint. For a footprint of 270 km 2 (i.e. the largest distribution grid area of the island), the measured ramps are close to those observed under clear sky conditions which are totally predictable. Unpredictable ramps are almost removed except a low percentage of extreme events which occur when large clouds passe over the island.

C. Energy Storage System (ESS)

Some ESS with short time response (~milliseconds) like flywheels or Li-ion batteries can be used to smooth PV output ramps. Trivially, more important the "smoothing" is, more important the size and the cost of batteries are. Thus, as previously, it is interesting to study the size/cost of ESS as a function of the footprint of dispersed PV fleet (still under uniform spatial distribution hypothesis). To this end, the approach used by Perez et al. [START_REF] Perez | The Cost of Mitigating Short-term PV Output Variability[END_REF] was implemented in this work. Briefly, they calculate the PV clearness index defined as the ratio of actual PV output to PV output under clear sky conditions. Then, simple moving average (SMA) method is applied on the PV clearness index for a given time window. The result is multiplied by the PV output under clear sky conditions to obtain the "smoothed" PV output. This loop is repeated by increasing the time windows until the largest ramp of the "smoothed" PV output is lower than a given ramp limitation. Therefore, the mitigation is applied only to the unpredictable ramps. Namely, ramps driven by weather conditions. At the end, the power that passes through the ESS (P ESS ) is calculated as the difference between the actual and the "smoothed" PV outputs divided by the square root of the ESS efficiency. The ESS power and energy capacity are chosen respectively as is the maximum absolute value of P ESS and the maximum daily amplitude of the cumulative energy in the ESS. A control of the ESS considering ideal perfect PV output forecast was also considered in this work. To this end, centered moving average method is applied instead of the SMA method in the algorithm described above. Indeed, anticipation of future ramps allows reducing the ESS size. Fig. 3 shows the result of both approaches (i.e. with and without forecasts) and for different footprints. The PV outputs were normalized to the installed PV capacities, the ESS efficiency was set to 90% and the hourly ramp limitation was fixed at 25% of the installed PV capacity. As ESS investment cost is almost proportional to the energy capacity, the cost is also evaluated assuming Li-ion batteries with a unit cost of 250 €/kWh. Fig. 3: ESS capacity and cost to limit the hourly ramps of PV generation to 25% of installed power considering 4 different footprints.

Economically, under our hypothesis and for our purpose the best choice is to limit the ramps of the aggregated PV production at the largest footprint. In another words, mitigating the ramps of a single PV installation with individual batteries is not cost effective. Moreover, mitigation of the aggregated PV production at grid scale offers the possibility to use one or a few large ESS units connected to the transmission grid (i.e. centralized ESS). This configuration has the advantage of reducing the investment, operation and maintenance costs. If batteries are used, this configuration also has the advantage of making it easier to apply the safety standards of the batteries and allowing optimized operation. However, in order to better value the storage, other services (reserve, voltage regulation, peak shaving, etc.) should be taken into account to formulate relevant conclusion on EES location on the grid [START_REF] Delille | Participation of Energy Storage in the Advanced Management of Power Systems : organizational, Technical and Economic Approaches in Distribution Grids[END_REF]. Fig. 3 shows also that the PV forecast should provide an additional margin to lower ESS costs, as well as facilitating PV integration into the generation mix.

These results can be compared with those of [START_REF] Perez | Mitigating Short-Term PV Output Intermittency[END_REF]. Indeed, for mitigation of hourly ramps at 25% of installed PV capacity and for a footprint of 2500 km 2 , they found the results presented in Table I. In this case, the ESS energy capacity is similar to the values found for California and Hawaii. The value of the ramp limitation depends on some grid characteristics. In order to have a wider view, Fig. 4 shows the storage size as a function of the ramp limitation per installed PV capacity at total grid scale. The maximum calculated PV output ramp over one hour is about 30% of the total PV installed capacity. This value is also observed by EDF on real production. Thus, beyond this value ESS is not necessary. At the opposite, under 22%, the limitation is too restrictive and ESS size/cost becomes prohibitive. Between 23% and 30%, the ESS cost is lower than 300 EUR per kW of installed PV capacity, namely 10-25% of current PV installation cost [START_REF] Fu | U.S. Solar Photovoltaic System Cost Benchmark: Q1 2018[END_REF].

IV. CONCLUSION

The results show that to increase the penetration of PV on the insular grid of La Reunion island, the solution could be to disperse first, quite uniformly, PV installations on the grid area. Secondly, a few ESS units with a centralized management by the grid operator would allow a further mitigation of the variability. That would lead to almost remove the unpredictable ramps (driven by weather). In addition, the forecast of the total PV production could both reduce ESS cost and anticipate the grid management by the operator. This last point shows that the challenge of forecasting of PV connected to the grid is not at a single plant scale, but at grid scale. Moreover, forecast at grid scale should be more reliable.

Based on this work, further improvements will be carried out. Firstly, time and space resolution of GHI data will be improved. Indeed, in March 2017 the satellite Meteosat-7 was replaced by Meteosat-8. This satellite offers a time and space resolution respectively of 15 min and 1 km. Thus, ramp mitigation over shorter time step will be studied. Secondly, regarding the geographical dispersion of PV generation, non-uniform spatial distribution will be investigated. The challenge is to optimize the distribution of PV systems in order to both minimize the variability and maximize the total production. Thirdly, regarding the smoothing of PV output with ESS, as explained in the introduction, other methods, like filter based methods, will be implemented. These improvements should reduce further unpredictable variability and the ESS cost.

Fig. 1 :

 1 Fig. 1: Global horizontal yearly radiation map and the transmission grid (63 kV) The penetration of PV into the grid depends on the existing grid configuration, in particular inertial response and available frequency control reserves. The following points present a few characteristics of the studied electric grid: • Installed power capacity: (i) diesel generators (211 MW), (ii) coal/bagasse-fired plants (210 MW), (iii) PV (187 MW), (iv) hydroelectricity (137 MW), (v) fuel turbines (80 MW), (vi) wind turbines (15 MW), (vii) biogas-fired plants (4 MW) and demand-side response capacities (15 MW). • Start-up time of the controllable facilities: (i) hydroelectricity (3 min), (ii) fuel turbine (15 min), (iii) diesel engine (1 hr) and (iv) coal/bagasse-fired plant (12 hr)

Fig. 2 .

 2 Fig. 2. presents the maximum values of these three metrics for each footprint, and also the theoretical ramps of GHI under clear sky conditions (). As shown by Fig.2, ΔGHI cs,1hr geographical dispersion has no effect on ramps observed under clear sky.

Fig. 2 :

 2 Fig. 2: Mitigation effect of spatial averaging on maximum hourly GHI ramps, associated P99 and standard deviation of hourly GHI ramps, 2011, Reunion Island.

Fig. 4 :

 4 Fig. 4: Storage size (power and energy) at total grid scale as a function of the ramp limitation

Table I

 I ESS energy capacity to limit the hourly ramps of PV generation to 25% of installed power for 4 different regions

	Location	EES energy capacity in
		kWh storage /kW PV
		with forecast	without forecast
	Mississippi	1.29	2.71
	California	0.42	1.17
	Hawaii	0.40	1.01
	Reunion Island	0.3	0.8

ACKNOWLEDGMENTS

This work was supported by the European Regional Development Fund, the Regional Council of Reunion Island, EDF-SEI and ADEME through the PEPS project.


Anne Gégout

Aurélie Gueudin

Clémence Karmann

Modelling Estimation process Theoritical results Simulation studies Network inference for truncated gaussian data

Let us look at the likelihood L jk (Σ * jk , y ) of the pairs (Y j , Y k ), j < k where :

L jk (σ, y ) = 1 a,b=0
φ ab,jk (σ, y j , y k )n ab (y j , y k ), with:

Anne Gégout, Aurélie Gueudin and Clémence Karmann

Network inference for truncated gaussian data Modelling Estimation process Theoritical results Simulation studies

Theorotical tools -Pairs likelihood

Let us look at the likelihood L jk (Σ * jk , y ) of the pairs (Y j , Y k ), j < k where :

L jk (σ, y ) = 1 a,b=0
φ ab,jk (σ, y j , y k )n ab (y j , y k ), with:

• n ab (y j , y k ) = 1 ζ j =a,ζ k =b où ζ l = 1 si y l ∈ [a l , b l ] \ {0}, 0 si y l = 0.
Let us look at the likelihood L jk (Σ * jk , y ) of the pairs (Y j , Y k ), j < k where :

L jk (σ, y ) = 1 a,b=0
φ ab,jk (σ, y j , y k )n ab (y j , y k ), with:

• n ab (y j , y k ) = 1 ζ j =a,ζ k =b où ζ l = 1 si y l ∈ [a l , b l ] \ {0}, 0 si y l = 0.

• φ 11,jk (σ, y j , y k ) = f (y j , y k , σ)

Anne Gégout, Aurélie Gueudin and Clémence Karmann Network inference for truncated gaussian data Modelling Estimation process Theoritical results Simulation studies

Theorotical tools -Pairs likelihood

Let us look at the likelihood L jk (Σ * jk , y ) of the pairs (Y j , Y k ), j < k where :

L jk (σ, y ) = 1 a,b=0
φ ab,jk (σ, y j , y k )n ab (y j , y k ), with:

• n ab (y j , y k ) = 1 ζ j =a,ζ k =b où ζ l = 1 si y l ∈ [a l , b l ] \ {0}, 0 si y l = 0.
• φ 11,jk (σ, y j , y k ) = f (y j , y k , σ)

• φ 01,jk (σ, y j , y k ) = φ 01,jk (σ, y k ) = [a j ,b j ] c f (x, y k , σ)dx • φ 10,jk (σ, y j , y k ) = φ 10,jk (σ, y j ) = [a k ,b k ] c f (y j , y , σ)dy

Theorotical tools -Pairs likelihood

Let us look at the likelihood L jk (Σ * jk , y ) of the pairs (Y j , Y k ), j < k where :

L jk (σ, y ) = 1 a,b=0
φ ab,jk (σ, y j , y k )n ab (y j , y k ), with:

• n ab (y j , y k ) = 1 ζ j =a,ζ k =b où ζ l = 1 si y l ∈ [a l , b l ] \ {0}, 0 si y l = 0.
• φ 11,jk (σ, y j , y k ) = f (y j , y k , σ)

• φ 01,jk (σ, y j , y k ) = φ 01,jk (σ, y k ) = [a j ,b j ] c f (x, y k , σ)dx • φ 10,jk (σ, y j , y k ) = φ 10,jk (σ, y j ) = [a k ,b k ] c f (y j , y , σ)dy • φ 00,jk (σ, y j , y k ) = φ 00,jk (σ) = [a j ,b j ] c ×[a k ,b k ] c f (x,

Context

Let Y := (Y (1) , . . . , Y (n) ) be an i.i.d. n-sample du vecteur Y . For every j < k, let us denote:

• L (n) jk (Σ * jk , y) the likelihood of the pairs sample (Y (i) j , Y (i) k ) 1≤i≤n
given by:

L (n) jk (σ, y) = n i=1 L jk (σ, y (i) ). • L (n) jk (Σ * jk , y) the log-likelihood L (n) jk (Σ * jk , y) = log(L (
(i) j , Y (i) k ) i=1,...,n .
Estimator of the covariance matrix Estimator Σ (n) of Σ * is defined by:

Σ (n) jk = argmax |σ|≤1 L (n) jk (σ, y).
(1)

Step 2: precision matrix Θ * Idea : estimate Θ * thanks to the graphical Lasso (Friedman et al. ( 2007)) for the GGM is to maximize the penalized log-likelihood of the Gaussian model.

Θ(n) = argmax Θ 0 log det(Θ) -trace(Θ Σ(n) ) -λ n ||Θ|| 1,off . (2)
with Σ(n) being the empirical covariance matrix of X Step 2: precision matrix Θ * Idea : estimate Θ * thanks to the graphical Lasso (Friedman et al. ( 2007)) for the GGM is to maximize the penalized log-likelihood of the Gaussian model.

Θ(n) = argmax Θ 0 log det(Θ) -trace(Θ Σ(n) ) -λ n ||Θ|| 1,off . ( 2
)
with Σ(n) being the empirical covariance matrix of X Problem:

Σ (n) unachievable! replaced by Σ (n)
Modelling Estimation process Theoritical results Simulation studies

Covariance matrix estimation Precision matrix estimation

Step 2: precision matrix Θ * Idea : estimate Θ * thanks to the graphical Lasso (Friedman et al. ( 2007)) for the GGM is to maximize the penalized log-likelihood of the Gaussian model.

Θ(n) = argmax Θ 0 log det(Θ) -trace(Θ Σ(n) ) -λ n ||Θ|| 1,off . ( 2
)
with Σ(n) being the empirical covariance matrix of X Problem:

Σ (n) unachievable! replaced by Σ (n)
Estimator of the precision matrix Estimator Θ (n) of Θ * is defined by:

Θ (n) = argmax Θ 0 log det(Θ) -trace(Θ Σ (n) ) -λ n ||Θ|| 1,off .
(3)

Modelling

Estimation process Covariance matrix estimation Precision matrix estimation

Theoritical results

Simulation studies

Anne Gégout, Aurélie Gueudin and Clémence Karmann Network inference for truncated gaussian data

Modelling Estimation process Theoritical results

Simulation studies

Tools

Step 1: [START_REF] Mei | The landscape of empirical risk for non-convex losses[END_REF], under suitable assumptions, good properties of the population risk can be carried to the empirical risk, even in an non-convex case control of

Σ (n) j,k -Σ * j,k
Anne Gégout, Aurélie Gueudin and Clémence Karmann Network inference for truncated gaussian data

Modelling Estimation process Theoritical results

Simulation studies

Tools

Step 1: [START_REF] Mei | The landscape of empirical risk for non-convex losses[END_REF], control of Σ

(n) j,k -Σ * j,k
Step 2: Ravikumar & al (2011), estimating the concentration matrix under sparsity conditions without specific distributional assumptions, but rather analyze the estimator in terms of the tail behavior of max j,k Σ (H2) Let j < k. Let g be a function:

(n) j,k -Σ * j,k Anne
σ ∈ [-1 + δ, 1 -δ] → E L (n) jk (σ, y) . • (-1 + δ) et (1 -δ) are not critical points of g ,
• g admits a finite number of critical points,

• all the critical points of g , different from Σ * jk , are non-degenerated, that is:

σ = Σ * jk , g (σ) = 0 ⇒ g (σ) = 0.
Anne Gégout, Aurélie Gueudin and Clémence Karmann Network inference for truncated gaussian data (H2) Let j < k. Let g be a function:

σ ∈ [-1 + δ, 1 -δ] → E L (n) jk (σ, y) . • (-1 + δ) et (1 -δ) are not critical points of g ,
• g admits a finite number of critical points,

• all the critical points of g , different from Σ * jk , are non-degenerated, that is:

σ = Σ * jk , g (σ) = 0 ⇒ g (σ) = 0.
(H3) Technical (not written). The underlying intuition is to limit the influence of "non-edged" terms on edge terms.

Modelling Estimation process Theoritical results

Simulation studies

Intermediate result on the estimated covariance matrix Σ (n) Proposition Under (H1) and, (H2), let 0 < ρ < 1. There exist known constants B,

C and D such that if n satisfies n log n ≥ C log B ρ
, then the estimated covariance matrix Σ (n) satisfies:

P Σ (n) -Σ * ∞ ≥ D log n n log B ρ ≤ p(p -1) 2 ρ,
where ||A|| ∞ = max j,k∈{1,...,p}

|A jk | is the infinite norm of matrix A ∈ R p 2 .

Mei et al. (2017)

Anne Gégout, Aurélie Gueudin and Clémence Karmann Network inference for truncated gaussian data

Modelling Estimation process Theoritical results

Simulation studies

Final result about Θ (n) Theorem Under (H1), (H2) et (H3). Let c > 2. There exist known constants B, C and D such that for all n satisfying n log n > f (B, C , D, c, Σ * ) and

λ n ∝ log n n log Bp c , we have with probability 1 - 1 p c-2 : (a) || Θ (n) -Θ * || ∞ ≤ D log n n log Bp c . (b) E ( Θ (n) ) ⊂ E (Θ * )

Simulations design

• p = 100 covariates, n = 500 observations

• Gaussian X with a chain structure (99 edges) :

X 1 ←→ X 2 ←→ . . . ←→ X 100

• p = 100 covariates, n = 500 observations

• Gaussian X with a chain structure (99 edges) :

X 1 ←→ X 2 ←→ . . . ←→ X 100
• Different truncation thresholds: Identical thresholds : a = -0.5, b = 2 Non identical thresholds a = seq(-1, 0), b = seq(0.5, 3)

• p = 100 covariates, n = 500 observations

• Gaussian X with a chain structure (99 edges) :

X 1 ←→ X 2 ←→ . . . ←→ X 100
• Different truncation thresholds: Identical thresholds : a = -0.5, b = 2 Non identical thresholds a = seq(-1, 0), b = seq(0.5, 3)

• Detection rate over 50 repetitions and comparaison with naive graphical Lasso directly on truncated data Y.

  covariance matrix Σ * • empirical covariance matrix from Y: poor performance • Likelihood Maximisation of Y : 2 p terms... covariance matrix Σ * • empirical covariance matrix from Y: poor performance • Likelihood Maximisation of Y : 2 p terms... estimator of Σ * term to term from the pairs (Y j , Y k ). We estimate each coefficients Σ (n) jk by maximising the log-likelihood of the sample (Y covariance matrix Σ * • empirical covariance matrix from Y: poor performance • Likelihood Maximisation of Y : 2 p terms... estimator of Σ * term to term from the pairs (Y j , Y k ). We estimate each coefficients Σ

  jk by maximising the log-likelihood of the sample (Y

  For all j < k, |Σ * jk | = 1. (H1) For all j < k, |Σ * jk | = 1.

  For all j < k, |Σ * jk | = 1.

  Structure "hub".

FigureFigure : Figure :

 :: Figure: Représentation graphique des trois graphes utilisés dans cette Anne Gégout, Aurélie Gueudin and Clémence Karmann Network inference for truncated gaussian data

Figure :

 : Figure:Detection rates obtained by our method (estimated truncation points) and by the graphical Lasso applied directly to the truncated data. The configurations of the truncation points are "identical only". Detection rates are obtained over 50 independent repetitions for n = 500 observations of p = 100 variables. "hub" structure; red = true edges.

  y , σ)dxdy .

	Modelling Modelling
	Estimation process Estimation process	Covariance matrix estimation Covariance matrix estimation
	Theoritical results Theoritical results	Precision matrix estimation Precision matrix estimation
	Simulation studies Simulation studies
	Modelling
	Estimation process
	Covariance matrix estimation
	Precision matrix estimation
	Theoritical results
	Simulation studies
	Anne Gégout, Aurélie Gueudin and Clémence Karmann Anne Gégout, Aurélie Gueudin and Clémence Karmann	Network inference for truncated gaussian data Network inference for truncated gaussian data

Anne Gégout, Aurélie Gueudin and Clémence KarmannNetwork inference for truncated gaussian data


domain: phys math

text:
Lorenzo Campana
email: lorenzo.campana@inria.fr

Mireille Bossy
email: mireille.bossy@inria.fr

Jean Pierre Minier
email: jean-pierre.minier@edf.fr

A Lagrangian stochastic model for rod orientation in turbulent flows

Keywords: Turbulent flow, anisotropic particles, stochastic modelling

Suspension of anisotropic particles can be found in various applications, e.g. industrial manufacturing processes or natural phenomena (micro-organism locomotion, ice crystal formation in clouds). Microscopic ellipsoidal bodies suspended in a turbulent fluid flow rotate in response to the velocity gradient of the flow. Understanding their orientation is important since it can affect the optical or rheological properties of the suspension (e.g. polymeric fluids). In this work, the orientation dynamics of rod-like tracer particles, i.e. long ellipsoidal particles (in the limit to infinity of the aspect-ratio) is studied. The size of the rod is assumed smaller than the Kolmogorov length scale but sufficiently large that its Brownian motion need not be considered. As a result, the local flow around a particle can be considered as inertia-free and Stokes flow solutions can be used to relate particle rotational dynamics to the local velocity gradient tensor A ij = ∂u i ∂x j . The orientation of a rod is described as the normalized solution of the linear ordinary differential equation for the separation vector R 12 between two fluid tracers. Separation evolves under the action of the velocity gradient tensor. Simultaneously, a re-normalization procedure R 12 R 12 is introduced to obtain the unit-vector p aligned with the rod. In this frame, the rod orientation is described by a Lagrangian stochastic model, assuming that cumulative effects of the velocity gradient tensor on the observation time interval fluctuate with a Gaussian distribution. Indeed, cumulative velocity gradient fluctuations are here represented by a white-noise tensor such that it preserves the incompressibility condition. Large observation timescale (overall objective of the work) justifies the Gaussian distribution hypotheses, with a decorrelation timescale equal to the Kolmogorov one τ η . Finally, the Lagrangian stochastic model is tested in the case of homogeneous isotropic turbulence.

Introduction

Turbulent flows with suspended particles of non spherical shape are a common occurrence in many industrial and natural processes. Industrial processes include pulp making and papermaking [START_REF] Lundell | Fluid mechanics of papermaking[END_REF], as well as soot emission from combustion processes [START_REF] Ryan | In-situ measurements of the mixing state and optical properties of soot with implications for radiative forcing estimates[END_REF]. Natural processes include the dispersion of pollen species in the atmosphere [START_REF] Sabban | Measurements of pollen grain dispersal in still air and stationary, near homogeneous, isotropic turbulence[END_REF], the dynamics of icy clouds [START_REF] Andrew | Precipitation development in stratiform ice clouds: A microphysical and dynamical study[END_REF], and the cycle of plankton such as diatoms [START_REF] Musielak | Nutrient transport and acquisition by diatom chains in a moving fluid[END_REF]. In many of these applications, the flow is highly turbulent and has an effect on the rotational dynamics, alignment trends and correlations of anisotropic particles (such as fibers, discs or more general shapes). In a turbulent flow, a small ellipsoidal particle rotates in response to the velocity gradients tensor along its Lagrangian trajectory. Because these Lagrangian velocity gradients are controlled by the small scales, they are similar in many different turbulent flows and have been the focus of extensive study [START_REF] Meneveau | Lagrangian dynamics and models of the velocity gradient tensor in turbulent flows[END_REF].

The earliest investigation into the motion of non-spherical particles in a carrier fluid is that of [START_REF] Barker | The motion of ellipsoidal particles immersed in a viscous fluid[END_REF]. For the special case of an axis-symmetric ellipsoid, he derived the evolution equation for the orientation vector as function of the local velocity gradient tensor. In turbulent flows the velocity gradient tensor A ij = ∂u i ∂x j fluctuates and is dominated by small scale motions of the order of Kolmogorov length scale l η . Much work has focused in rod-like particles whose size is smaller than the l η . Studies of the orientation dynamics of such particles in turbulent flows have included those of [START_REF] Shin | Rotational and translational dispersion of fibres in isotropic turbulent flows[END_REF] and [START_REF] Pumir | Orientation statistics of small particles in turbulence[END_REF] using isotropic turbulence data from direct numerical simulation (DNS), those of [START_REF] Zhang | Ellipsoidal particles transport and deposition in turbulent channel flows[END_REF] and [START_REF] Ph Mortensen | Dynamics of prolate ellipsoidal particles in a turbulent channel flow[END_REF] for particles in channel flow turbulence using DNS. In many numerical studies, Lagrangian tracking is most often used to determine the particle trajectories and simultaneous time integration of the Jeffery equation along the trajectory leads to predictions of the particles' orientation dynamics. Generic properties of the orientation dynamics, such as the variance of the fluctuating orientation vector or its alignment trends may also be studied by making certain assumptions about the Lagrangian evolution of the carrier fluid's velocity gradient, in particular about its symmetric and skew-symmetric parts, the strain tensor S ij = (A ij + A ji ) 2 and the rotation tensor A number of theoretical studies have been based on the assumption that these flow variables obey isotropic Gaussian statistics, e.g. are the result of linear Ornstein-Uhlenbeck processes (e.g. see [START_REF] Brett K Brunk | Observations of coagulation in isotropic turbulence[END_REF]; [START_REF] Pumir | Orientation statistics of small particles in turbulence[END_REF]; [START_REF] Wilkinson | A model for alignment between microscopic rods and vorticity[END_REF] and [START_REF] Vincenzi | Orientation of non-spherical particles in an axisymmetric random flow[END_REF]).

Ω ij = (A ij -A ji ) 2.
For small tracers particles whose size is smaller than the Kolmogorov scale, the local flow around the particle can be considered to be inertia-free and Stokes flow solutions can be used to relate the rotational dynamics of the particles to the local velocity gradient. To understand the dynamics of ellipsoidal particles in turbulence, there are needs to extend the understanding of the Lagrangian statistics of the velocity gradient tensor, and to include the orientational dynamics, resulting from the integrating Jeffery's equation Jeffery (1922) along the particle trajectory. This is a challenging problem, both because of the complexity of statistically quantifying the particle orientation with respect to the velocity gradient tensor. In this frame, the rod orientation is described by a Lagrangian stochastic model, assuming that cumulative effect of the velocity gradient tensor on the observation time interval fluctuate with a Gaussian distribution. Indeed, cumulative velocity gradient fluctuations are here represented by a white-noise tensor such that it preserves the incompressibility condition.

The Lagrangian stochastic model involves three timescales: the Kolmogorov timescale τ η , the Lagrangian integral timescale of the fluid T L and the integration timescale ∆t. The first two timescales are physical characteristic timescales. The third one, ∆t, represents the 'observation' time-scale. Large observation timescale, such as pump clogging (objective of the modeling) justifies the Gaussian distribution hypotheses, with a decorrelation timescale equal to the Kolmogorov one τ η . Besides, the development of both model and suitable numerical scheme for a large integration time steps, for the stochastic differential equation associated, is a difficult task to address for the separation equation. In this context the focus on the orientation information only is a crucial point. This paper presents some advances in modeling of rod orientation in the context of large observation time-scale simulations and results are reproduced in the context of isotropic homogeneous turbulence.

Equations of motions

Generally speaking, the equation of motion for a rod-like particle (ellipsoid in the limit of infinity aspect-ratio see Fig. 1) in a turbulent velocity field u(r, t) is considered. The rod is assumed to be naturally buoyant and smaller than the smallest length scale characterizing fluid motion, i.e. the Kolmogorov length l η , but sufficiently large that their Brownian motion do not need be considered. The motion of the rod or "tumbling rate" is determined by the particle orientation and the velocity gradient tensor (as a particular case of the general Jeffery's equation (Jeffery 1922)):

ṗi = Ω ij p j + S ij p j -p i p k S kl p l (1)
where p is a unit director along the symmetric axis of the particle, and S ij and Ω ij are the symmetric and skew-symmetric parts of the velocity gradient tensor, respectively. However, this approach leads to a non-linear equation that in the framework of stochastic modeling is difficult to deal with. The same dynamics can be solved by considering a linear equation for the separation of two fluid particles, particularly, the equation of motion has been derived considering the dynamics of two fluid tracers described by the vector R 12 and eventually normalizing the solution.

In more detail, the separation vector R 12 = R 1 -R 2 between two fluid particles with trajectories R n (t) = R(t; r n ) passing at t = 0 through the points r n satisfies the equation

Ṙ12 (t) = v(R 1 , t) -v(R 2 , t).
(2)

In the rest of the paper the subscript R 12 is replaced by R.

Considering an incompressible fluid flow where the particles generally separate. In smooth velocity field, for separation R much smaller than the viscous scale of turbulence l ≪ l η , i.e. the so-called Batchelor regime (Batchelor 1959), the velocity difference between two fluid particles can be expressed as

v(R 1 , t) -v(R 2 , t) ≈ A ij (t)R with the Lagrangian strain matrix A ij = ∂u i ∂x j R 2 (t), t .
In this regime, Eq. ( 2) simplifies to the ordinary differential equation

Ṙ(t) = A ij (t)R(t) (3)
leading to the linear solution

R(t) = W ij (t)R(0) (4)
where W ij is evolution matrix or the deformation gradient tensor that characterizes the distortion withstand to the fluid element. The evolution matrix provides a Lagrangian description of the fluid stretching. In other words, for threedimensional flows, this process can be visualized by considering a sphere that is distorted into a tri-axial ellipsoid ad it is stretched by the flow (see Fig. 1). This means that passive vectors along with thin rods-like become preferentially aligned with the longest principal axis of the ellipsoid, and at long times approach (after 10τ η ) perfect alignment with the eigenvector corresponding to the maximum stretching has been studied by [START_REF] Ni | Alignment of vorticity and rods with lagrangian fluid stretching in turbulence[END_REF].

The general solution for the three dimensional case is determined by products of random matrices. The evolution matrix W ij may be written as

W ij (t) = T exp t 0 A ij (s)ds = ∞ n=0 t 0 A ij (s n )ds n ⋯ s3 0 A ij (s 3 )ds 2 s2 0 A ij (s 2 )ds 1 (5)
This time-order exponential form (T ), in general, is not very useful for the direct computation except the particular case of a short correlated strain (see sec. Stochastic model for the orientation). The basic idea of the present approach, relies on the result that, in almost realization of the strain gradient tensor, the matrix 1 t ln W ij W ji stabilizes as t → ∞ (see [START_REF] Falkovich | Particles and fields in fluid turbulence[END_REF] for details). To give some intuitive insight, as already mentioned, considering some fluid volume, like a sphere, which evolves into an elongated ellipsoid at later time. As time increases, the ellipsoid is more and more elongated; furthermore the long time evolution of Eq. ( 3) is characterized by the following result

λ L = lim t→∞ 1 t ln R(t) R(0) (6)
which for large t tends with probability one to the largest Lyapunov exponent λ L (finite value), governing the chaotic properties of the particles trajectories in the turbulent flow. In this way it is possible to solve the rod-like orientation through a linear ordinary differential equation Eq. ( 3); where the initial conditions are W ij (0) = I (identity matrix) and R(0) = p(0). The solution for rod orientation in Eq. ( 1) is obtained by normalizing the solution of Eq. (3):

p(t) = R(t) R(t) . ( 7
)

Stochastic model for the orientation

The evolution of p depends upon the velocity gradient tensor A ij (t), therefore a stochastic model for it is presented. Other Lagrangian stochastic model are present in literature, but so far this description has been always investigated in the framework of direct numerical simulations (DNS). Gaussian processes have been proposed for the velocity gradient statistics. For instance, [START_REF] Pumir | Orientation statistics of small particles in turbulence[END_REF] and [START_REF] Vincenzi | Orientation of non-spherical particles in an axisymmetric random flow[END_REF] have considered an Ornstein-Uhlenbeck process for A ij with a different correlation time scale for the symmetric and skew-symmetric parts. This is more realistic, since it is known that in turbulence the correlation time scale for the rotation rate is significantly longer than that the strain rate. Otherwise, more refined model for the velocity gradient has been obtained by [START_REF] Campana | Numerical scheme for a lagrangian stochastic model describing rods orientation[END_REF], introducing the RFDA model, which overcomes the Gaussian description, predicting a variety of local, statistical, geometric and anomalous scaling properties of 3-D turbulence. It is important to mention as a reference, the work of [START_REF] Chevillard | Orientation dynamics of small, triaxial-ellipsoidal particles in isotropic turbulence[END_REF]; which examines in very detail the orientation dynamics of the anisotropic particle in isotropic homogeneous turbulence, using both DNS and different stochastic models. In this framework it is important to underline that these models are tested in a DNS context, or rather, they have been developed for an observation time of the order of the Kolmogorov timescale (τ ≃ τ η ).

The solution of the differential equation such as Eq. ( 3) with a stochastic input A ij is given by Eq. ( 5) with the matrix W ij (t) involving stochastic integrands over the time. The case of a short correlated velocity gradient tensor (Kraichnan 1968) allows for a complete solution. More specifically, when the observation time of the dynamics is much larger than the correlation time of the velocity gradient tensor, i.e. t ≫ τ η , W ij may be viewed as a continuous product of independent random matrix. Under this assumption, the instantaneous velocity gradient tensor can be decomposed in two contributions: mean field and fluctuations

A ij X(t, ω), t = ⟨∂ i u j X(t, ω), t ⟩ + ξ ij (t, ω) (8)
where the ⟨⋅⟩ means the ensemble average and ξ ij (t, ω) is a white-noise that can be rewritten as

dF ij (t, ω) = ξ ij (t, ω)dt.
This stochastic forcing (fluctuation) is of the form dF ij = b ijkl dW kl , where dW kl represent a Wiener process, i.e ⟨W ij ⟩ = 0 and ⟨dW ij dW kl ⟩ = δ ik δ jl dt leading to,

⟨dF ij dF kl ⟩ = b ijmn b mnkl dt (9)
Incompressibility, isotropy and parity invariance impose the form of the fourth-order tensor b ijkl ,

b ijkl = b 1 δ ij δ kl + b 2 δ ik δ jl + b 3 δ il δ jk . (10)
with the coefficients b

1 = (- √ 3 3) τ 2 η , b 2 = 0.5( √ 3 + √ 5) τ 2 η and b 3 = 0.5( √ 3 - √ 5) τ 2
η which depend exclusively on the Kolmogorov time scale τ η (see [START_REF] Perry | A closure for lagrangian velocity gradient evolution in turbulence using recent-deformation mapping of initially gaussian fields[END_REF] for details). Now, replacing the decomposition for A ij into the rods' orientation Eq. (3), it is possible to re-formulate the stochastic differential equation (SDE)

dR i (t) = ⟨A ij ⟩ R j dt + dF ij ○ R j (11)
interpreted in the Stratonovich form (○). The natural form for a direct calculation of the SDE is the Itô form, so applying the transformation from Stratonovich to Itô the SDE in the index notation becomes

dR i = āii + 1 2 ((b 1 + b 2 + b 3 ) 2 + 2b 2 1 ) R i dt + 3 j=1 āij R j (1 -δ ij )dt + (b 2 + b 3 )dW ii + b 1 3 j=1 dW jj R i + (b 2 + b 3 ) 3 j=1 R j (1 -δ ij )dW ij (12)
where āij are the elements of ⟨A ij ⟩ and dW ij is a matrix composed by nine independent Brownian motion. Here, the repeated index is not an implicit summation and δ ij is the Kronecker's symbol. To attain the complete description of the rod's orientation p the SDE is coupled with the normalization procedure, as outlined in the Eq. ( 7).

It is worth to notice that, when the observation time-step ∆t becomes smaller with respect to the Kolmogorov inner time scale τ η the model for the orientation is no longer valid. This constraint reflects the argument that on a large enough observation time scale (meaning precisely ∆t ≫ τ η ), the stochastic model for the orientation is coherent with the physical description of the Lagrangian stochastic model. In particular, it should be stressed that there is a strong interplay between the physical aspects of the model (which leads to a formulation in terms of SDE) and the numerical aspects of the practical simulations as also has been remarked by [START_REF] Minier | Probabilistic approach to turbulent two-phase flows modelling and simulation: theoretical and numerical issues[END_REF].

Results and Discussion

This section presents some results obtained from the numerical simulations of the Lagrangian stochastic model described by Eq. ( 12) and Eq. ( 7). The model has been performed using a numerical integration scheme based on splitting algorithm, which solves the symmetric and skewsymmetric part of the velocity gradient tensor separately (Campana et al. N.D.). The validation has been restricted to the homogeneous isotropic turbulent (HIT) case, so that the mean velocity gradient tensor is assumed to be zero ⟨A ij ⟩ = 0 . As discussed above, the validation of the model has to be performed for an integration time-step bigger than the Kolmogorov one. Three simulations have been performed: varying the value of the ratio between the observation timescale and the Kolmogorov time scale (the setup for three different cases are detailed in Tab. 1). Independent sample of rods have been initialized at time zero, imposing to each of them a uniform distribution on the unit sphere (U(S 2 )) for the three cases. In the context of HIT, due to the hypotheses of homogeneity and isotropy of the fluid field, the rod orientation vector should remain uniformly distributed on the unit sphere. Firstly, results are presented in the Cartesian coordinate system. The whole sample of orientation vectors p remain, for long time simulation, distributed on the unit sphere as showed in Fig. 2. Regarding the first and second moments of the orientation vector (p 1 , p 2 , p 3 ), the Tab. 2 presents measures of the error of mean and variance of the stochastic process evaluating at final time. The empirical mean of any generic moments is denoted as

⟨f (p)⟩ M C = 1 N p sample f (p simualtions ) (13)
and the error is usually decomposed in two contributions: the bias part due to the integration scheme and the Monte Carlo error due to the size of sample as, Table 1: Summary of the simulations performed for the Lagrangian stochastic model. For each cases the parameters are referred to: the time-step (∆), the ratio between the observation and Kolmogorov timescales (τ τ k ), the numbers of particles used to performed the Monte Carlo (MC) simulations (N p ), the final time (T ) and the initial distribution imposed (I.d).

E E[f (p)] -⟨f (p)⟩ M C 2 = O ∆t 2α + O Var(f (p)) N p (14)
The values of the errors are listed in Tab. 2; these are small and of the same order for the three values of the integration time-step ∆t. The explanation of this last result relies into the fact that, analyzing the two sources of error, it appears that the total error is dominated by the Monte Carlo approximation, even for the large integration time-step. In Fig. 3 is presented the time evolution of the mean and variance of the rod orientation vector. Both of the mean (see Fig. 3a) and variance (see Fig. 3b) account for the ergodicity of the simulated process. This fact is of main importance regarding the construction of the model (Eq. 6), as discussed above. In fact, the mean oscillates around zero (which is the exact value) for the three component of the orientation vector p (see Fig. 3a) with an error (Tab. 2) that is produced only by the MC simulations for the three cases. Similar remarks can be done looking the variance (see Fig. 3b), in this case the value of the components of p fluctuating around the exact value of 1 3 for each of them (Tab. 2).

To complete the description of the numerical results and to have, at the same time, a more physical characterization of the model itself, the results on the empirical probability distribution function (p.d.f ) are presented. To analyze the orientation dynamics on a non-spherical particle, it is convenient to move from Cartesian coordinates (p 1 , p 2 , p 3 ) to spherical coordinates (r, ϕ, θ) according to the usual transformations

r = p 2 1 + p 2 2 + p 2 3 , ϕ = arctan p 2 1 + p 2 2 p 3 , θ = arctan(p 2 p 1 ) (15)
with r ⩾ 0, 0 ⩽ ϕ ⩽ π and 0 ⩽ θ ⩽ 2π. In Fig. 4 are presented the marginal distributions for the two spherical angles: the Fig. 4a and Fig. 4b) show respectively the empirical pdf (ϕ) and pdf (θ) for the three cases. The three corresponding curves are all in good agreement with the theoretical marginals p.d.f :

p.d.f (θ) th = 1 4π π 0 sin ϕdϕ = 1 2π p.d.f (ϕ) th = 1 4π 2π 0 sin ϕdθ = 1 2 sin ϕ. (16)
Finally, on account of the unitary length of the orientation vector p (r = 1), the empirical joint probability density function of the orientation must take the form of p.d.f (ϕ, θ) in spherical coordinates. It is showed in Fig. 5 and its twodimensional projection Fig. 6. This empirical p.d.f (ϕ, θ) has been obtained through the use a Gaussian smoothing kernel that wrings the surface near the boundaries of the plots. It appears that the empirical joint p.d.f is in agreement with the theoretical one, p.d.f (ϕ, θ) th = 1 4π sin ϕ.

(17)

Conclusions

In this study, the rod orientation is described by a Lagrangian stochastic model using dedicated splitting scheme for time integration along large observation time-scale. The cumulative effects of the velocity gradient tensor are assumed fluctuating with a Gaussian distribution with short-time correlations. A new stochastic differential equation (SDE), which describes the separation between two fluid particles, has been proposed. This SDE coupled with a renormalization procedure provides the description of the rod's orientation. The model have been tested numerically in the HIT case and three reference time-scale orders have been examined, showing a good agreement between theoretical and empirical, both for the marginal and joint probability density function expressed in spherical coordinates. This preliminary study serves as a base to go further both on the modeling and the numerical development.

Figure 1 :

 1 Figure 1: (Colour online) Typical ellipsoidal shape of a fluid element along the turbulent trajectory. The fluid element center of mass X(t) is supposed to follow the evolution of the fluid particle and the deformation is governed by the statistics of the fluid velocity gradients along the trajectory, ∂ i u j (X(t), t). The fluid element is always assumed of ellipsoidal shape with the three semi-axes ordered as d 1 ⩾ d 2 ⩾ d 3 . The rod-like particle is described by taking d 1 → ∞. This figure is adapted from Biferale et al. (2014).

Figure 2 :

 2 Figure2: (Colour online) Sample orientation vector evaluated at finale time (T = 10 3 ) distributed on the sphere with a number of iterations 10 5 for the case 3.

Figure 3 :

 3 Figure 3: (Colour online) Mean (a) and variance (b) for the three component of the orientation vector p as a function of time, obtained for the case 3.

10 thFigure 4 :

 104 Figure 4: (Colour online) P.d.f of the angle ϕ in (a) and θ in (b) in the spherical coordinates system. Results from case 1 (blue solid line), case 2 (green dashed line) and case 3 (red dot-dashed line).

  Figure 5: (Colour Online) Empirical joint P.d.f in spherical coordinates produces with the case 3.

Figure 6 :

 6 Figure 6: (Colour online) Projected empirical joint P.d.f in spherical coordinates produces with the case 3.

  10 th International Conference on Multiphase Flow, ICMF 2019, Rio de Janeiro, Brazil, May 19 -24, 2019

  10 th International Conference on Multiphase Flow, ICMF 2019, Rio de Janeiro, Brazil, May 19 -24, 2019

	(a)	0.02
	< p 1 >	0.00
		-0.02 0.02
	< p 2 >	0.00
		-0.02 0.02
	< p 3 >	0.00
		-0.02	0	200	400	600	800	1000
					t
	(b)
	1 ]	0.34
	Var[p	0.33
	2 ]	0.34
	Var[p	0.33
	3 ]	0.34
	Var[p	0.33
		0		200	400	600	800	1000
					t

Table 2 :

 2 Estimation of the weak error for the first, second moment and the function g = g(p 1 , p 2 , p 3 ) = p 1 p 2 p 3 . The p

	Error	Case 1	Case 2	Case 3
		0.00157 0.00360 0.00089
	E[p i ] -⟨p i ⟩ M C	0.00113 0.00005 0.00046
		0.00028 0.00019 0.00172
		0.00060 0.00019 0.00092
	E[p 2 i ] -⟨p 2 i ⟩ M C	0.00079 0.00021 0.00117
		0.00138 0.00005 0.00211
	E[g] -⟨g⟩ M C	0.00005 0.00023 0.00023

i indicates the exact value, the computed one. The table is organized as following: the error for the mean and variance is computed for each coordinates in order i = 1, 2, 3.

th International Conference on Multiphase Flow,ICMF 2019, Rio de Janeiro, Brazil, May 19 -24, 2019

Acknowledgments

This work has been supported by EDF R&D (projects PTHL of MFEE and VERONA of LNHE) and by the French government, through the Investments for the Future project UCA JEDI ANR-15-IDEX-01 managed by the Agence Nationale de la Recherche


Yunfei Ma

Jiahuan Cui

Nagabhushana Rao Vadlamani

Paul Tucker

O Tatin

Effect of Fan on Inlet distortion: A Mixed fidelity Approach

Keywords: inlet distortion, fan effects, mixed fidelity

Inlet distortion is typically encountered during the off-design conditions on civil aircra s and in S-ducts in military aircra s. It is known to severely deteriorate the performance of a gas-turbine engine. As the intakes get shorter, there is an increased interaction between the inlet distortion and the downstream fan. Previous studies in the literature use low-fidelity methods (RANS or URANS) to address this unsteady interaction, due the substantial computational cost associated with the high fidelity methods like LES/DNS. On the other hand, it is well known that the low order methods like RANS have limitations to accurately represent the distorted flows. In this paper, we propose a mixed-fidelity approach and employ it to study the intake-fan interaction at an affordable computational cost. e results demonstrate that there are two ways through which the fan affects the separated flow. Firstly, the suction effect of the fan alleviates the undesired distortion by 'directly' changing the streamline curvature, intensifying the turbulence transport and closing the recirculation bubble much earlier. Secondly, the enhanced turbulence in the vicinity of the fan feeds back into the initial growth of the shear layer by means of the recirculating flow. is 'indirect' feedback is found to increase the turbulence production during the initial stages of the shear layer. Both the direct and indirect effects of the fan significantly suppress the inlet distortion.

INTRODUCTION

Flow distortion is typically encountered on the engine intakes and in the duct flows under the off-design conditions.

e flow over an intake lip separates specifically during the take-off at high angles of incidence and under severe crosswinds during taxing or landing of an aircra . e distorted flow at the inlet convects downstream and deteriorates the performance of the fan. In some extreme cases, it can even lead to the catastrophic events like stall and surge. Interestingly, as the inlets become shorter, the proximity of the fan to the source of the distortion is shown to affect the distortion recovery ([ ]). In order to numerically investigate this interaction, it is crucial to accurately resolve the following flow regimes: (a) the separated flow domain and (b) the influence of a fan.

Most of the research in the literature addressed the fandistortion study using low fidelity methods such as URANS and RANS [ , , ]. It is well known that the predictive capability of these low fidelity approaches are satisfactory at the design point. However, they still suffer from severe limitations under off-design conditions, involving flow separation, distortion and unsteadiness [ , , ]. In particular, when predicting flows with large-scale, low-frequency turbulence, the results may vary significantly among the different RANS models [ ]. For these flows, eddy resolving simulations such as DNS/LES and hybrid LES/RANS are demonstrated to yield much more promising results [ , , ].

On the other hand, the fan influence can be represented either by resolving all fan blades or with a body force. While the former approach is computationally expensive [ ], the Body Force Method (BFM) is much more economical[ ] and practically feasible. BFM can be further classified into the smeared and standard methods, based on the extent to which the geometry has been modelled. e Immersed Boundary Method with Smeared Geometry (IBMSG) is typically used to model the force field generated by a set of rotating blades. is approach assumes an infinite number of blades in the circumferential direction; which avoids the computational overhead of capturing the detailed blade geometry and/or incorporating moving boundaries to represent the fan. In this approach, an inviscid (or wall-normal) force and a viscous (or parallel) force due to the rotating blades are added to the circumferentially averaged Navier-Stokes equations. e inviscid component of the force, formulated by Marble [ ], ensures that the flow follows the blade metal angle. e viscous component, introduced by Xu [ ], accounts for the losses encountered over the blades using a force-velocity relation. Xu [ ] has successfully employed this approach to investigate the distortion transfer in a high pressure turbine.

On the other hand, the standard Immersed Boundary Method (IBM) is generally used to resolve the features of a real geometry. Peskin [ ] proposed IBM based on Eulerian and Lagrangian variables, linked by the interaction equations involving smooth appximation of the Dirac-delta function.

is approach has been successfully applied by Fadlun [ ] to study the flow past a backward facing step and by Defoe [ ], who modelled a rotor blade and predicted the noise propagation under the influence of an inlet flow distortion.

In this paper, a mixed-fidelity approach is applied to accurately investigate the influence of the fan on the inlet distortion at an affordable computational expense. e standard IBM proposed by Peskin [ ] is used to model the distortion generator upstream of the fan. e unsteady separated regime downstream of the distortion generator is captured using a high-fidelity eddy resolving approach. Further downstream, IBMSG is used to model the force field due to a rotating fan. Subsequently, the recovery of distortion in the presence of a fan is examined using both the mean and turbulent characteristics.

. NUMERICAL FRAMEWORK

Figure illustrates the computational domain and the boundary conditions considered in the current study. is simplified setup is motivated by the experimental studies on the Darmstadt Rotor [ , , , , ], albeit under different operating conditions. In these studies, the distortion generators were designed to reproduce the flow conditions in a real engine within the laboratory. Measurements by Lieser [ ] and Bi er [ ] show that the compressor performance is largely sensitive to the distortions encountered at the tip. Hence, a periodic distortion generator is placed upstream of the tip of the fan in order to reproduce the distortion encountered over the intake lip at high angles of a ack. e test case employs the original duct and the rotating fan from the Bi er's test rig [ , ], with the periodic beam installed upstream of the fan. It features a 30 deg sector duct with a distortion generator ('beam') of height 'H' and length '1.5H' placed at an axial distance of '12.5H' from the inlet. e fan is positioned at a streamwise distance of . H from the beam. All the spatial quantities mentioned in the following sections are normalised by the beam height H. e velocity is normalised by the velocity, u ∞ , measured at the maximum height of the beam (Fig. ), which corresponds to the outer edge velocity of a separated shear layer. e primary objective of the current study is to capture the distortion generated on the lower wall. Hence, an invscid boundary condition is imposed on the upper wall, which ensures that the pressure distribution due to the spinner is well represented at a reduced computational cost. International Standard Metric Conditions (P 0 = Pa and T 0 = . K) are applied at the inflow and the massflow rate is fixed at . kg/s. is massflow rate corresponds to the peak efficiency point at % rotational speed ( . rad/s). A radial equilibrium boundary condition is imposed at the outflow. Periodicity is imposed in the circumferential direction. e extent of the sector (30 • corresponds to 5H) is sufficient enough to ensure that the structres are decorrelated in the circumferential direction. In a mixed fidelity framework, the beam is represented using the conventional IBM [ ]. e separated flow downstream of the beam is captured using the eddy resolving approach, while the force field of a rotating fan is replicated using IBMSG [ ]). is approach avoids the need to capture detailed blade geometry or incorporate mov-ing boundaries. Instead, it represents the force field generated by the rotating fan blades and captures the suction effect of the fan, thereby substantially reducing computational cost.

. METHOD e present simulations are carried out using a Rolls-Royce's in-house CFD code, HYDRA [ , ]. Both IBM and IBMSG are implemented into the numerical framework and are thoroughly validated [ ]. For the IBM part, the distortion generator is modelled using feedback forcing, proposed by Goldstein et al. [ ]. It has the following expression,

f (x, t) = α ∫ t 0 ∆udt + β∆u, ∆u = u(x 0 , t) -u 0 (x 0 , t) ( )
where the subscript 0 represents the solid boundary. e coefficients α and β are negative constants because essentially they function as a proportional-integral (PI) feedback controller evaluated by the velocity difference. Eventually, the controller enforces u = u 0 on the immersed boundary. From the viewpoint of control theory, the value of the force is independent of the two coefficients α and β, under the condition that the computation is convergent. To ensure the convergence, the two coefficients should be carefully selected. In fact, they are associated with two important parameters: the

< -β -β 2 -2αk α , ( )
where the order of k is . e fan modelled with IBMSG can be regarded as a set of an infinite number of infinitesimally thin blades. e forces are circumferential averaged into every cell in the blade region and turn the flow to the local desired direction. Cao et al. [ ] introduced this method when studying intake separation and find this model is capable to capture the key flow features. Hence, we use this normal force model

f n = - u l • n ∆t -RH S • n, ( )
where n is the normal vector of blade surface, and RH S is the right hand side of Navier-Stokes equation. e parallel force from Watson et al. [ ] is also used, which is

f p = -K(r)ρu 2 rel , K(r) = 4k 1 s 2 + k 1 , ( )
where s is the fraction of span and the coefficient k 1 need to be calibrated by the existing performance data. e blockage λ = 1 -1 2

(t 1 + t 2 ) S(r) cos β , ( )
where t 1 and t 2 are the blade thickness for the pressure part and suction part, and S(r) is the surface pitch, as is shown in Figure

. Mechanism

To further investigate how fan interacts the separation bubble through turbulence, this section delves into the mechanism by which the fan has reduced the recirculation region. In the previous section, Figure b and c demonstrates a substantial increase in both the production and convection. However, the location of such increase varies, and consequently may present different traits. To examine this in more details, the individual contributions from the dominant terms (in order): u ′ u ′ ∂U/∂ x, v ′ v ′ ∂V/∂ y and u ′ v ′ ∂U/∂ y to the TKE production are analysed.

e results demonstrate that there are two distinct ways through which the fan effects the separated flow: shear flow dominant (Zone ) and strain rate dominant (Zone ), shown in Figures a andb. In Zone , the primary effect of the fan is observable between x = 2H -4H where the recirculation region decreases due to the change of the streamline curvature. Both u ′ u ′ and v ′ v ′ increase significantly within this region by around % and % respectively. e increase in v ′ v ′ is much more pronounced than u ′ u ′ , in- we consider this increase in Reynolds stresses to be a 'direct effect' which is a ributed to the additional strain rate caused by the curvature change.

In zone , Figures a and b also demonstrate an increase in both the streamwise and wall normal fluctuations within the shear layer at the edge of the beam between x = 0 -2H. Clearly, the spreading rate of the shear layer is much more predominant in the presence of the fan. is is a ributed to the 'indirect effect' of the fan, where the enhanced turbulence generated in the vicinity of the fan feeds back into the origin of the shear layer by means of the recirculating flow. is is also evident from the Figures a andc where an increase in the TKE and TKE convection is observable within the recicurlating region. Consequently, this 'indirect' feedback is found to further intensify the turbulence production within the shear layer (Figure b) increasing the Reynolds stresses

( u ′ u ′ , v ′ v ′ ).
As expected, u ′ v ′ ∂U/∂ y is found to be the major contributor to the total production of TKE in Zone within a spreading shear layer. e increased TKE in Zone is subsequently convected downstream into the Zone through the mainstream flow and resulting in the cyclic feedback between both the zones. Hence, although the flow in the Zone is not affected directly by the fan, it is still a key source of turbulence due to the feedback from the reverse flow. It supplements the turbulence generated in the vicinity of the fan in the Zone and also contributes to the earlier rea achment of the separation bubble. In this paper, we propose a mixed-fidelity approach and employ it to study the intake-fan interaction at an affordable computational cost. e results demonstrate that there are two ways through which the fan effects the separated flow: Firstly, the suction effect of the fan (effective upto almost half the chord length upstream of the fan) alleviates the undesired distortion by 'directly' changing the streamline curvature, intensifying the strain rate and turbulence transport thereby closing the recirculation bubble much earlier. Secondly, the enhanced turbulence in the vicinity of the fan feeds back into the initial growth of the shear layer by means of the recirculating flow. is 'indirect' feedback is found to increase the turbulence production and spreading rate of the shear layer. Both these direct and indirect effects of the fan significantly suppress the inlet distortion.

From an engineering perspective, the current study shows that the loading of the fan blade can be conveninently altered to maximise the stagnation pressure recovery by suppressing the distortion. A tip loaded fan can increase the streamline curvature driving early rea achment of the flow. Further studies will be carried out to address these aspects and explore the possibility of some potential flow control methods to alleviate the inlet distortion.

[ ] V Jerez Fidalgo, CA Hall, and Y Colin. A study of fandistortion interaction within the nasa rotor transonic stage. Journal of Turbomachinery, ( ): , .

[ ] Sebastian Barthmes, Jakob P Haug, Andreas Lesser, and Reinhard Niehuis. Unsteady cfd simulation of transonic axial compressor stages with distorted inflow. In Advances in Simulation of Wing and Nacelle Stall, pages -. Springer, .

[ ] John Seddon and Eustace Laurence Goldsmith. Intake aerodynamics. Amer Inst of Aeronautics &, .

[ ] Paul Tucker and Yan Liu. Turbulence modeling for flows around convex features. In th AIAA Aerospace Sciences Meeting and Exhibit, page , .

[ ] Nicolas Gourdain, Laurent YM Gicquel, and Elena Collado. Comparison of rans and les for prediction of wall heat transfer in a highly loaded turbine guide vane. Journal of propulsion and power, ( ): -, .

[ ] Yangwei Liu, Hao Yan, Lipeng Lu, and Qiushi Li. Investigation of vortical structures and turbulence characteristics in corner separation in a linear compressor cascade using ddes. Journal of Fluids Engineering, ( ): , .

[ ] PR Spalart and DR Bogue.

e role of cfd in aerodynamics, off-design.

e Aeronautical Journal ( ), ( ): -, .

[ ] Feng Gao, Wei Ma, Jinjing Sun, Jérôme Boudet, Xavier O avy, Yangwei Liu, Lipeng Lu, and Liang Shao. Parameter study on numerical simulation of corner separation in lmfa-naca linear compressor cascade. Chinese Journal of Aeronautics, ( ): -, .

[ ] Jeffrey Slotnick, Abdollah Khodadoust, Juan Alonso, David Darmofal, William Gropp, Elizabeth Lurie, and Dimitri Mavriplis. Cfd vision study: a path to revolutionary computational aerosciences.

.

[ ] Roberto Verzicco, Jamaludin Mohd-Yusof, Paolo Orlandi, and Daniel Haworth. Large eddy simulation in complex geometric configurations using boundary body forces.

AIAA journal, ( ): -, .

[ ] FRANK E Marble. ree-dimensional flow in turbomachines. High Speed Aerodynamics and Jet Propulsion, : -, .

[ ] L Xu. Assessing viscous body forces for unsteady calculations. In ASME Turbo Expo : Power for Land, Sea, and Air, pages -. American Society of Mechanical Engineers, .

[ ] Charles S Peskin. e immersed boundary method. Acta numerica, : -, .

[ ] EA Fadlun, R Verzicco, Paolo Orlandi, and J Mohd-Yusof. Combined immersed-boundary finite-difference methods for three-dimensional complex flow simulations.

Journal of computational physics, ( ): -, .

[ ] Jeffrey J Defoe and Zoltán S Spakovszky.

α

  must be large enough to keep the force frequency 1 2π |α| much larger than any other frequencies in the flow [ ]. is means that the feedback force should change rapidly and fit the flow to the desired direction within a unsteady period. Also, to ensure the numerical stability, the time step should also meet the requirement ∆t

Figure .

 . Figure . Experiment se ings

Figure .

 . Figure . Sketch of blade blockage

.

  Time averaged flow fieldA er flushing out the initial transience, statistics are collected for around 94H/U ∞ .e maximal deviation of u ′ u ′ extracted from two different times (with the interval of 50H/U ∞ ) is only %, indicating the convergence is acceptable. Figure compares the mean velocity profiles at different streamwise locations on a carpet plot. A line joining the locus of inflectional points of the velocity profiles is also overlaid. e recirculation region can be also seen in Figure.As noted from the instantaneous flow, the extent of the recirculation zone is significantly reduced due to the fan. e flow reattaches at an axial location which is more than a beam height upstream of the fan leading edge.

  Figure a shows the contours of the time-averaged turbulent kinetic energy.e TKE in the shear layer and in the rea aching regime has increased by around -% in the axial direction in the presence of fan.is is partially caused by a higher production in the shear layer (Fig.b) and partially by a stronger convection in the middle-rear part of the recirculation region (Fig.c).

Figure .

 . Figure . Total pressure distribution

Figure .

 . Figure . Recirculation region

  Figure compares the velocity gradients ∂U/∂ y to the contours of the Reynolds stress u ′ v ′ (Fig. c) without/with fan. It is evident that the contribution to the TKE production is primarily affected by an increase in the Reynolds stress u ′ v ′ . On the other hand, the change in the velocity gradients is marginal.

  (a) Contours of TKE (b) Contours of TKE Production (c) Contours of TKE convection

Figure .

 . Figure . Turbulent Kinetic Energy Statistics

  Figure . Main Reynolds stresses

  Martin Bi er, Fabian Wartzek, Simon Übelacker, Heinz-Peter Kähler Schiffer, and Christian J. Characterization of a distorted transonic compressor flow using dualluminophore pressure-sensitive paint. In th Pacific Symposium on Flow Visualization and Image Processing (PSFVIP-), fedOA (Federico II Open Archive), Naples, Italy, . [ ] Fabian Wartzek, Felix Holzinger, Christoph Brandste er, and Heinz-Peter Schiffer. Realistic inlet distortion patterns interacting with a transonic compressor stage. In Advances in Simulation of Wing and Nacelle Stall, pages -. Springer, . [ ] Teng Cao, Paul Hield, and Paul G Tucker. Hierarchical immersed boundary method with smeared geometry. In th AIAA Aerospace Sciences Meeting, page , . [ ] PI Crumpton, P Moinier, and MB Giles. An unstructured algorithm for high reynolds number flows on highly stretched grids. Numerical methods in laminar and turbulent flow, pages -, . [ ] L Lapworth. Hydra-cfd: a framework for collaborative cfd development. In International Conference on Scientific and Engineering Computation (IC-SEC), Singapore, June, volume , . [ ] D Goldstein, R Handler, and L Sirovich. Modeling a noslip flow boundary with an external force field. Journal of Computational Physics, ( ): -, . [ ] Rob Watson, Jiahuan Cui, Yunfei Ma, and Paul Hield. Improved hierarchical modelling for aerodynamically coupled systems. In ASME Turbo Expo : Turbine Technical Conference and Exposition. American Society of Mechanical Engineers, . [ ] Yangwei Liu, Hao Yan, Yingjie Liu, Lipeng Lu, and Qiushi Li. Numerical study of corner separation in a linear compressor cascade using various turbulence models. Chinese Journal of Aeronautics, ( ): -, . [ ] Peter Bradshaw. Effects of streamline curvature on turbulent flow. Technical report, DTIC Document, .

				Effects
	of boundary-layer ingestion on the aero-acoustics of
	transonic fan rotors. Journal of Turbomachinery,
	( ):	,	.
	[ ] JA Lieser, C Biela, CT Pixberg, H-P Schiffer, S Schulze,
	A Lesser, CJ Kähler, and Niehuis R. Compressor rig test
	with distorted inflow using distortion generators.	.
	[ ] R Niehuis, A Lesser, Axel Probst, Rolf Radespiel, Sonja
	Schulze, Christian Kähler, Frank Spiering, and Norbert
	Kroll. Simulation of nacelle stall and engine response.
	.
	[ ] Simon Übelacker, Rainer Hain, and Christian J Kähler.
	Flow investigations in a stalling nacelle inlet under dis-
	turbed inflow. In Advances in Simulation of Wing and
	Nacelle Stall, pages	-. Springer,	.
	[ ]

ACKNOWLEDGMENTS

e authors gratefully acknowledge the UK Turbulence Consortium grant (EPSRC grant EP/L / ), the financial support from the Chinese Scholarship Council, and the technical support from Rolls-Royce.