Manuel Girault
email: manuel.girault@ensma.fr

Julien Launay

Nadine Allanic

Pierre Mousseau

Rémi Deterre

Development of a thermal Reduced Order

Keywords: pseudoplastic fluid, elastomer, Galerkin projection, Modal Identification Method, optimization, intrusive sensor

Introduction

The diversity of natural resources and the increase of biosourced and recycled materials involve quality fluctuations of raw materials used in polymer forming process. The consequence is a variation of the polymer rheological behavior and thus difficulties during forming. Polymers, and especially elastomers, have a high pseudo-plastic behavior which concentrates the high shear rate near to the walls. Thus, the viscous dissipation, which highly depends on the material viscosity, is preponderant close to the wall. The resulting thermal gradient can be very high in the runners of industrial tools [START_REF] Launay | Intrusive measurement of polymer flow temperature[END_REF][START_REF] Wielhorski | Thermal balance between viscous heating and inlet thermal condition in non stationary polymer flow through a cylindrical die[END_REF][START_REF] Beaumont | Runner and gating design handbook[END_REF]. The knowledge of the flow rheology variations during forming process should be helpful to define production stability criteria.

Usually, the rheological apparatus used in production, and largely studied in literature, are based on the detection of pressure gap in derivation runners of rheometer arrangement [START_REF] Beaumont | Runner and gating design handbook[END_REF][START_REF] Agassant | Polymer processing: principles and modeling[END_REF].

Such apparatus are only devoted to mechanical parameters (flowrate, pressure) and allow measuring an average viscosity if the prerequisite condition of a stationary regime is met. Such average viscosity does not take into account the thermal heterogeneity of the material in the measurement area [START_REF] Laun | Capillary rheometry for polymer melts revisited[END_REF][START_REF] Laun | Pressure dependent viscosity and dissipative heating in capillary rheometry of polymer melts[END_REF]. Moreover, in most cases, the process conditions lead to unsteady flow or thermal regimes, which are hardly suitable for the required stationarity hypothesis previously mentioned.

Another possible way consists in exploiting the intrinsic thermal information of the flow (especially viscous heating) in order to estimate viscosity or slight change of viscosity. To study this possibility, we want to use the thermal measurements obtained with a cylindrical runner specifically equipped. This device has been developed in the context of former works [START_REF] Launay | Intrusive measurement of polymer flow temperature[END_REF] and is presented in Figure 1. It is thermally regulated and equipped with two intrusive measurement cells called TMC (Thermal Measurement Cell). These innovative sensors with an elliptic shape are introduced within the flow and contain 25 thermocouples of 50 µm diameter. Despite the perturbations of the flow induced by their invasive nature, these sensors allow to measure the temperature profile of elastomer or thermoplastic flows at the inlet and the outlet of the runner. The measurements obtained and numerically checked are sensitive to the process parameters and flow conditions, namely the rheological behavior of the material.

Through the association of these temperature measurements with a numerical procedure for parameter estimation, it is possible to consider an on-line indirect measurement of viscosity, which is a characteristic quantity of the flow rheological behavior.

A thermo-rheological model of the flow is required to solve the parameter estimation problem. In order to build a device able to perform on-line viscosity estimation, this model should be a Reduced Order Model (ROM) allowing fast resolution of the inverse problem.

Thus, the aim of this numerical preliminary work is to develop a thermo-rheological ROM for a steady incompressible viscous flow of a non-Newtonian pseudoplastic fluid.

A ROM is a model involving a small number of degrees of freedom (dof), which mimics the behavior of an actual system or a reference model of this system (model with a large number of dof sometimes called detailed model), whatever the time-varying boundary conditions and/or source terms, or for a range of values of some parameters. A ROM may be obtained by transformation of a detailed model or identified from data coming either from simulations of a detailed model or from measurements on the real system.

Among methods for building ROMs we may first cite those based on Proper Orthogonal Decomposition (POD) followed by a Galerkin projection. As an example, such an approach was used in [START_REF] Park | Karhunen-Loève Galerkin method with decimated sampling technique for the simulation of complex fluids defined in the phase space[END_REF] for building a ROM for the Fokker-Planck equation modelling dynamics of solute molecules in complex fluids, with an application to the flow of dilute polymeric liquids over a cylinder described by the FENE model. POD-Galerkin ROMs were also developed in [START_REF] Brigham | Inverse viscoelastic material characterization using pod reduced-order modeling in acoustic-structure interaction[END_REF] for a vibroacoustic system used for characterizing viscoelastic properties of solids by means of an inverse problem.

The Branch Eigenmodes Reduction Method (BERM) relies on a special spectral problem where the eigenvalue of each mode appears in a Steklov boundary condition. The resulting branch basis allows handling nonlinear problems. Reduction is then obtained by selecting or amalgamating the branch modes. The BERM was used for instance in [START_REF] Joly | Modal reduction of an advection-diffusion model using a branch basis[END_REF] for an advectiondiffusion problem featuring a rotating disk with variable velocity and multiple time-varying loads.

The Reduced Basis (RB) approximation method [START_REF] Quarteroni | Certified reduced basis approximation for parametrized partial differential equations and applications[END_REF] allows building ROMs for parametrized Partial Differential Equation (PDE). It uses a Greedy algorithm to compute the reduced basis thanks to solutions of a reference Finite Element model at optimally selected points in the parameter space. The ROM is then formed by Galerkin projection of PDE onto the reduced basis. The RB method provides errors bounds but in turn requires a reference FE model on which it relies. Therefore the construction of ROMs from experimental data does not appear to be possible with the RB approach. One should note that Greedy sampling methods are usually used in parameter domain whereas POD-based approaches are often used in temporal domain. Both techniques can be combined for parameter-time problems.

Both the goal oriented model-constrained approach developed in [START_REF] Bui-Thanh | Goal-oriented, modelconstrained optimization for reduction of large-scale systems[END_REF] and the Proper Generalized Decomposition (PGD) [START_REF] Chinesta | An overview of the proper generalized decomposition with applications in computational rheology[END_REF] allow the building of ROMs able to mimic dynamics of a system for both a time-varying load and a physical parameter range.

The former uses an optimization problem for computing the basis of space functions: an output error between solutions of the reference model and those of the ROM is minimized, assuming that each basis vector can be represented as a linear combination of snapshots. The approach has been applied to a compressible flow around a rotor blade, showing that for a small number of basis functions (≤10), stable ROMs were obtained whereas POD-Galerkin ROMs were unstable [START_REF] Bui-Thanh | Goal-oriented, modelconstrained optimization for reduction of large-scale systems[END_REF].

The latter (PGD) appears to provide a discrete solution in physical space, time and possibly parameters space rather than a continuous ROM. In [START_REF] Chinesta | An overview of the proper generalized decomposition with applications in computational rheology[END_REF] are gathered applications of PGD to some rheology-related problems where the interested reader will find additional references.

Applying PGD to our steady problem would consist to seek the temperature field as a sum of products of functions of physical space coordinates and functions of the rheological parameters defining the fluid viscosity. We would then obtain a discrete solution in the physical space as well as in the parameters space. Because the ROM we wish to build is intended to be used in an inverse problem for the estimation of rheological parameters, we have preferred to build a ROM depending explicitly on these parameters, thus allowing for direct differentiation with respect to these parameters in order to compute sensitivities.

The Modal Identification Method (MIM) allows the building of ROMs linking the system inputs (i.e. the quantities influencing its behaviorexternal loads and/or internal parameters)

to the outputs (either the whole fields of the problem variables, like temperature, velocity, etc., or their value in some specific locations). ROMs are identified from the knowledge of input-output data characterizing the system behavior, through the minimization of a quadratic functional based on the residues between reference output data on the one hand and solutions of the ROM on the other hand, corresponding to the same inputs. MIM hence uses optimization tools (here a Particle Swarm Optimization algorithm) for the identification of ROMs. So far, MIM has been mainly subject to works involving heat conduction [START_REF] Girault | Estimation of time-varying heat sources through inversion of a low order model built with the modal identification method from in-situ temperature measurements[END_REF][START_REF] Girault | Temperature regulation and tracking in a MIMO system with a mobile heat source by LQG control with a low order model[END_REF] or Newtonian fluid flows in forced convection [START_REF] Girault | Estimation of time-varying pollutant emission rates in a ventilated enclosure: inversion of a reduced model obtained by experimental application of the modal identification method[END_REF][START_REF] Rouizi | Model reduction by the Modal Identification Method in forced convection: Application to a heated flow over a backward-facing step[END_REF] or mixed convection [START_REF] Cordier | Reduced Order Modeling by Modal Identification Method and POD-Galerkin approach of the heated circular cylinder wake in mixed convection[END_REF] regimes. In particular, ROMs built by MIM have been used to solve efficiently inverse problems of heat conduction [START_REF] Girault | Estimation of time-varying heat sources through inversion of a low order model built with the modal identification method from in-situ temperature measurements[END_REF] and transport-diffusion [START_REF] Girault | Estimation of time-varying pollutant emission rates in a ventilated enclosure: inversion of a reduced model obtained by experimental application of the modal identification method[END_REF] and also thermal state-feedback control problems [START_REF] Girault | Temperature regulation and tracking in a MIMO system with a mobile heat source by LQG control with a low order model[END_REF]. Moreover, until now, inputs considered were mainly time-varying boundary conditions [START_REF] Girault | Temperature regulation and tracking in a MIMO system with a mobile heat source by LQG control with a low order model[END_REF][START_REF] Rouizi | Model reduction by the Modal Identification Method in forced convection: Application to a heated flow over a backward-facing step[END_REF] or internal sources [START_REF] Girault | Estimation of time-varying heat sources through inversion of a low order model built with the modal identification method from in-situ temperature measurements[END_REF][START_REF] Girault | Estimation of time-varying pollutant emission rates in a ventilated enclosure: inversion of a reduced model obtained by experimental application of the modal identification method[END_REF] rather than parameters. In the frame of heat conduction, an extension of MIM to the construction of a ROM having as inputs both a time-varying thermal load and a physical parameter in a predefined range (thermal conductivity) was addressed in [START_REF] Girault | Parametric low-order models in transient heat diffusion by MIM. Estimation of thermal conductivity in a 2D slab[END_REF].

In the present work is presented an extension of the MIM aiming at the construction of a parametric thermal ROM for a steady incompressible flow of a non-Newtonian fluid in a cylindrical duct. The ROM is explicitly parametrized by two thermo-rheological parameters defining the dynamic viscosity of a pseudoplastic fluid: consistency index and power law index [START_REF] Agassant | Polymer processing: principles and modeling[END_REF]. The ROM allows the computation of the flow temperature profile at the channel outlet whatever the two input parameters in a predefined range of values. Due to the fact that an assumption of no thermodependence of viscosity is considered in this first approach, the developed model is more adapted to polymers having a low value of the activation energy of the viscosity [START_REF] Agassant | Polymer processing: principles and modeling[END_REF].

The paper is structured as follows. In section 2 are presented the geometry, equations and boundary conditions of the problem as well as the numerical reference model developed with a finite element software. In section 3 is presented an overview of the approach proposed for building a ROM of heat transfer parametrized by viscosity. The development of the ROM formulation is described in details in section 4. Only the energy balance equation is concerned as the velocity field is obtained by an analytical solution. A Galerkin projection is used for defining the ROM equations. Data generation and ROM identification are briefly presented in section 5. In section 6 the developed approach is applied to an elastomer cylindrical flow corresponding to an experiment performed with the device described in Figure 1. Several "viscosity-temperature" datasets are generated with the numerical reference model of the flow. Some of them are used to build a series of ROMs. The others are used for the ROMs validation: the outlet temperature profile computed with the identified ROMs is compared to the one computed with the reference model. Concluding remarks and prospects are given in section 7.

Pseudoplastic fluid flow in a cylindrical channel

Let us consider the steady incompressible viscous flow of a polymer melt in a channel. The polymer melt is considered as a non-Newtonian fluid. The complete domain Ω is composed of two sub-domains: Ω 𝑓 is the fluid domain inside the channel and Ω 𝑤 the solid domain in the channel wall.

General local governing equations

Because of the strong viscous nature of polymers, inertia terms are not taken into account and the flow is assumed to be steady. Buoyancy effects being also neglected [START_REF] Agassant | Polymer processing: principles and modeling[END_REF][START_REF] Byron | Transport Phenomena[END_REF], equations for conservation of mass, momentum and energy in the polymer flow hence write:

∇ ⃗ ⃗ . 𝑣 = 0 (1)
-∇ ⃗ ⃗ 𝑝 + ∇ ⃗ ⃗ 𝜏 = 0 ⃗ (2)
𝜌𝐶 𝑝 𝑣 . ∇ ⃗ ⃗ 𝑇 = ϕ 𝑣 + ∇ ⃗ ⃗ . (𝑘 𝑓 ∇ ⃗ ⃗ 𝑇) (3)
where 𝜌, 𝐶 𝑝 and 𝑘 𝑓 are respectively the density, specific heat capacity and thermal conductivity of the fluid. ϕ 𝑣 represents the power dissipated by internal forces. By considering an incompressible viscous fluid, the power dissipated by pure shear, can be expressed by:

ϕ 𝑣 = 𝜏 ∶ 𝜀̇ (4)
with 𝜀̇ the strain rate tensor.

The viscous stress tensor 𝜏 depends on the fluid viscosity. In the case of generalized Newtonian fluids, the dynamic viscosity 𝜂 depends on the generalized shear rate 𝛾̅ . This one is linked to the second invariant of strain tensor 𝐽 2 by [START_REF] Agassant | Polymer processing: principles and modeling[END_REF]:

𝛾̅ ̇= 2√𝐽 2 (𝜀) with 𝐽 2 (𝜀) = 1 2 ∑ 𝜀̇𝑖 𝑗 2 𝑖,𝑗 (5)
For polymer flows, the shear rate terms are much more important than the elongation ones in most of processes. This leads to the relation 𝛾̅ ̇= 𝛾̇ where 𝛾̇ is the shear rate (SI unit s -1 ). Thus, the viscous power can be expressed as:

ϕ 𝑣 = 𝜂𝛾̇2 (6)
Moreover in this work, we describe the viscosity of pseudoplastic fluids through a power law:

𝜂 = 𝐾𝛾̇𝑛 -1 (7)
where 𝐾 is the consistency index (SI units Pa.s n ) and 𝑛 < 1 is the power law index (dimensionless).

For polymers, the quality of viscosity curves fitting by power law depends on polymer molecular weight and molecular weight distribution [START_REF] Van Krevelen | Properties of polymers: their correlation with chemical structure; their numerical estimation and prediction from additive group contributions[END_REF] and also on the effect of fillers [START_REF] Leblanc | Rubber-filler interactions and rheological properties in filled compounds[END_REF].

Concerning the most important industrial processes (e.g. extrusion, injection) the shear rate range is mostly located in the pseudoplastic area of viscosity curves [START_REF] Agassant | Polymer processing: principles and modeling[END_REF]. Due to the high dependence of viscosity with shear rate (several orders of magnitude), the smaller variation of viscosity with temperature, currently expressed by an exponential law [START_REF] Launay | Effect of viscous dissipation in the prediction of thermal behavior of an elastomer cylindrical flow[END_REF] was neglected.

In the channel wall depth, diffusion is the lone heat transfer mode. The energy equation then comes down to the heat equation:

∇ ⃗ ⃗ . (𝑘 𝑤 ∇ ⃗ ⃗ 𝑇) = 0 (8)
where 𝑘 𝑤 is the channel wall thermal conductivity.

Established flow in a circular channel

Let us consider the special case of a circular channel as depicted in Figure 2.

Local governing equations

The flow is assumed to be established in the 𝑧-direction (channel axis) with a symmetry of the velocity field with respect to angular coordinate 𝜑. Hence the velocity field has a single nonzero component 𝑢 𝑧 depending on the radial coordinate 𝑟. The dynamic viscosity 𝜂 defined by Eq.( 7) therefore writes:

𝜂 = 𝐾 | 𝑑𝑢 𝑧 𝑑𝑟 | 𝑛-1 (9)
In the present work viscosity is assumed to be independent of temperature. Properties of the fluid and of the channel are also assumed to be temperature-independent. In addition buoyancy effects are neglected. As a consequence of these hypothesis, the velocity field does not depend on temperature.

Under the previously mentioned assumptions, mass conservation equation ( 1) and momentum conservation equation ( 2) have an analytical solution. The velocity field writes [START_REF] Agassant | Polymer processing: principles and modeling[END_REF]:

𝑢 𝑧 (𝑟) = ( 𝑛 𝑛 + 1 ) ( ∆𝑝 2𝐾𝐿 ) 1 𝑛 (𝑅 𝑛+1 𝑛 -𝑟 𝑛+1 𝑛 ) (10)
Where ∆𝑝 is the pressure difference between inlet and outlet sections of the channel of length 𝐿.

Taking into account the above mentioned assumptions, the energy equation (3) takes the following form describing steady heat transfer in the flow [START_REF] Agassant | Polymer processing: principles and modeling[END_REF][START_REF] Byron | Transport Phenomena[END_REF]:

𝜌𝐶 𝑝 𝑢 𝑧 𝜕𝑇 𝜕𝑧 = ∇ ⃗ ⃗ . (𝑘 𝑓 ∇ ⃗ ⃗ 𝑇) + 𝜂 ( 𝑑𝑢 𝑧 𝑑𝑟 ) 2 (11)
Using equations ( 9) and [START_REF] Quarteroni | Certified reduced basis approximation for parametrized partial differential equations and applications[END_REF], equation [START_REF] Bui-Thanh | Goal-oriented, modelconstrained optimization for reduction of large-scale systems[END_REF] writes:

𝜌𝐶 𝑝 ( 𝑛 𝑛 + 1 ) ( ∆𝑝 2𝐾𝐿 ) 1 𝑛 (𝑅 𝑛+1 𝑛 -𝑟 𝑛+1 𝑛 ) 𝜕𝑇 𝜕𝑧 = ∇ ⃗ ⃗ . (𝑘 𝑓 ∇ ⃗ ⃗ 𝑇) + 𝐾 ( ∆𝑝 2𝐾𝐿 𝑟) 𝑛+1 𝑛 (12)

Thermal boundary conditions

The boundary of the whole domain Ω is noted Γ. The fluid domain boundary is

Γ 𝑓 = Γ 𝑖𝑛 ∪ Γ 𝑓𝑤 ∪ Γ 𝑜𝑢𝑡 whereas the solid domain (wall) boundary is Γ 𝑤 = Γ 𝑤𝑑 ∪ Γ 𝑓𝑤 ∪ Γ 𝑤𝑛 .
Thermal boundary conditions are described in the following and shown in Figure 2.

At the channel inlet Γ 𝑖𝑛 , a radial temperature profile is assumed:

Γ 𝑖𝑛 ∶ 𝑇(𝑟, 𝑧 = 0) = 𝑇 𝑖𝑛 (𝑟) ; 0 ≤ 𝑟 ≤ 𝑅 (13)
At the channel outlet Γ 𝑜𝑢𝑡 , the flow is considered as thermally established, a homogeneous

Neumann condition is hence applied:

Γ 𝑜𝑢𝑡 ∶ 𝜕𝑇 𝜕𝑧 (𝑟, 𝑧 = 𝐿) = ∇ ⃗ ⃗ 𝑇(𝑟, 𝑧 = 𝐿). 𝑒 𝑧 ⃗⃗⃗ = 0; 0 ≤ 𝑟 ≤ 𝑅 (14)
On the channel wall external surface Γ 𝑤𝑑 , the temperature profile is assumed to be function of 𝑧:

Γ 𝑤𝑑 ∶ 𝑇(𝑟 = 𝑅 + 𝑒, 𝑧) = 𝑇 𝑤 (𝑧); 0 ≤ 𝑧 ≤ 𝐿 (15)
On the channel wall boundary Γ 𝑤𝑛 (𝑧 = 0 and 𝑧 = 𝐿), homogeneous Neumann conditions are assumed:

{ 𝜕𝑇 𝜕𝑧 (𝑟, 𝑧 = 0) = ∇ ⃗ ⃗ 𝑇(𝑟, 𝑧 = 0). 𝑒 𝑧 ⃗⃗⃗ = 0; 𝑅 ≤ 𝑟 ≤ 𝑅 + 𝑒 𝜕𝑇 𝜕𝑧 (𝑟, 𝑧 = 𝐿) = ∇ ⃗ ⃗ 𝑇(𝑟, 𝑧 = 𝐿). 𝑒 𝑧 ⃗⃗⃗ = 0; 𝑅 ≤ 𝑟 ≤ 𝑅 + 𝑒 (16)
On the channel wall internal surface Γ 𝑓𝑤 , i.e. the fluid/solid interface, continuity of both temperature and heat flux density (perfect contact) are assumed.

Numerical reference model

The Ansys Polyflow software [START_REF] Ghoreishy | Finite element analysis of a thermoplastic elastomer melt flow in the metering region of a single screw extruder[END_REF] has been used to model the system. This proprietary Finite

Elements software is particularly well suited for viscous fluids flows. This reference model will be called Detailed Model (DM) in the following.

The problem considered is a one-way coupled problem as temperature depends on velocity whereas velocity does not depend on temperature because viscosity is assumed to be independent of temperature. Hence for the fluid domain, only energy equation [START_REF] Chinesta | An overview of the proper generalized decomposition with applications in computational rheology[END_REF] has to be solved. However doing this with Ansys Polyflow is not convenient. Therefore the system of coupled equations ( 1) and ( 2) has first been solved for computing the velocity field, which is of course equal to the analytical profile [START_REF] Quarteroni | Certified reduced basis approximation for parametrized partial differential equations and applications[END_REF], and then equations ( 3) and ( 8) have been solved

to compute the temperature field in the fluid and solid domains respectively. The thermal boundary conditions ( 13)-( 16) are completed by following boundary conditions for velocity:

For both inlet and outlet sections, a homogeneous Neumann condition is applied:

Γ 𝑖𝑛 ∶ 𝜕𝑣 𝜕𝑧 (𝑟, 𝑧 = 0) = 0 ⃗ (17)
Γ 𝑜𝑢𝑡 ∶ 𝜕𝑣 𝜕𝑧 (𝑟, 𝑧 = 𝐿) = 0 ⃗ (18)
On the channel wall internal surface Γ 𝑓𝑤 , a no-slip condition is applied:

Γ 𝑓𝑤 ∶ 𝑣 (𝑟 = 𝑅, 𝑧) = 0 ⃗ (19)
Gaussian elimination has been used along with an iterative Picard scheme. A geometric mesh with 1.1 factor and minimal cell size 0.5 mm has been used to obtain convergence. Thus, the number of finite elements is 𝑁 𝑚𝑒𝑠ℎ =33322.

Heat transfer Reduced Order Model parametrized by viscosity: overview

The velocity field being analytically determined by Eq.( 10), model reduction is here concerned with energy equation only for both the fluid (Eq.( 12)) and the solid wall (Eq.( 8))

domains. Our goal is to build a Reduced Order Model with explicit dependence on the parameters governing the dynamic viscosity 𝜂 defined by Eq.( 9): consistency index 𝐾 and flow behavior index 𝑛.

The approach used in the present paper is based on the Modal Identification Method (MIM) [START_REF] Girault | Estimation of time-varying heat sources through inversion of a low order model built with the modal identification method from in-situ temperature measurements[END_REF][START_REF] Girault | Temperature regulation and tracking in a MIMO system with a mobile heat source by LQG control with a low order model[END_REF][START_REF] Girault | Estimation of time-varying pollutant emission rates in a ventilated enclosure: inversion of a reduced model obtained by experimental application of the modal identification method[END_REF][START_REF] Rouizi | Model reduction by the Modal Identification Method in forced convection: Application to a heated flow over a backward-facing step[END_REF][START_REF] Cordier | Reduced Order Modeling by Modal Identification Method and POD-Galerkin approach of the heated circular cylinder wake in mixed convection[END_REF][START_REF] Girault | Parametric low-order models in transient heat diffusion by MIM. Estimation of thermal conductivity in a 2D slab[END_REF]. The MIM consists in three main steps:

1) Defining the structure of the ROM equations able to adequately describe the involved physics (see section 4). This is the core of the paper;

2) Generating some input-output data representative of the system. Those data may come from in-situ measurements or, like in the present paper, from numerical simulations (see section 5.1 for a brief presentation and section 6.1 for the practical application);

3) Identifying the fixed constitutive parameters of the ROM equations through the minimization of a functional based on the quadratic residuals between the previously generated output data of the system on the one hand and the outputs of the ROM on the other hand, for the same input data (see section 5.2 for a brief presentation of technical aspects and section 6.2 for the practical application).

The MIM therefore aims to adjust the ROM constitutive parameters using optimization techniques, in order for the ROM to mimic the data characterizing the input-output responses of the system. Main conceptual differences with POD-G and PGD approaches are highlighted in section 5.3.

Reduced Order Model formulation

The general form of the ROM equations is given in section 4.1. Then, the Galerkin projection of the energy equation (Eq.( 12) in Ω 𝑓 and Eq.( 8) in Ω 𝑤 ), briefly presented in section 4.2, leads to a first structure for the ROM equations. An additional hypothesis allowing to separate the radial direction 𝑟, intrinsically linked to the flow behavior index 𝑛, from the 𝑧 direction, is made in section 4.3. A simplified treatment of diffusion terms is then done in section 4.4. The final form of the ROM is given in section 4.5 before some comments in section 4.6.

General form of ROM equations

For a steady-state problem, the ROM to be built takes the following general form:

𝐹 𝑘 (𝑋, 𝜃, 𝐾, 𝑛) = 0 ∀𝑘 ∈ [1; 𝑚] (20)
And:

𝑌 𝑗 𝑜𝑏𝑠 = 𝐺 𝑗 (𝑋, 𝜃) ∀𝑗 ∈ [1; 𝑁 𝑜𝑏𝑠 ] (21)
Equation ( 20), called ROM state equation, is a set of 𝑚 equations (𝑚 ≪ 𝑁 𝑚𝑒𝑠ℎ ), explicitly parametrized by 𝐾 and 𝑛, which allows computing a so-called state vector 𝑋 ∈ ℝ 𝑚 . This set of equations also depends on constitutive fixed parameters gathered in vector 𝜃(𝑚) ∈ ℝ 𝑁 𝑝𝑎𝑟𝑎𝑚 (𝑚) which have to be determined in order to build the ROM.

Equation [START_REF] Leblanc | Rubber-filler interactions and rheological properties in filled compounds[END_REF], called ROM observation equation, is a set of 𝑁 𝑜𝑏𝑠 equations allowing to compute the outputs of interest 𝑌 𝑗 𝑜𝑏𝑠 (some temperatures in our case) once 𝑋 computed.

𝐹 𝑘 , 𝑘 ∈ [1; 𝑚], 𝐺 𝑗 , 𝑗 ∈ [1; 𝑁 𝑜𝑏𝑠 ] as well as vector 𝜃(𝑚) depend on the considered problem and on additional hypotheses, as shown in the following sections.

Galerkin projection of energy equation

Energy equation (Eq.( 12) in Ω 𝑓 and Eq.( 8) in Ω 𝑤 ), parametrized by 𝐾 and 𝑛, is considered.

It is here assumed that the temperature field 𝑇(𝑟, 𝑧, 𝐾, 𝑛) can be written as a linear combination of some space functions 𝜙 𝑖 (𝑟, 𝑧), 𝑖 = 1, … , 𝑚:

𝑇(𝑟, 𝑧, 𝐾, 𝑛) ≈ ∑ 𝜙 𝑖 (𝑟, 𝑧)𝑋 𝑖 (𝐾, 𝑛) 𝑚 𝑖=1 (22)
where space functions 𝜙 𝑖 (𝑟, 𝑧), 𝑖 = 1, … , 𝑚 are a truncation of a basis of the Hilbert space formed by the ℒ 2 (Ω) space of square integrable functions equipped with the usual inner product 〈. , . 〉 Ω :

〈𝑢, 𝑣〉 Ω = ∫ 𝑢𝑣𝑑Ω Ω ( 23
)
where functions 𝑢 and 𝑣 are defined on Ω.

The 𝑋 𝑖 (𝐾, 𝑛), 𝑖 = 1, … , 𝑚 are the coefficients of the approximation.

Of course, our goal will be to find a ROM, hence corresponding to a small number 𝑚 of functions used in the temperature field decomposition.

Let us call ℛ the residual of energy equation (Eq.( 12) in Ω 𝑓 and Eq.( 8) in Ω 𝑤 ).

The Galerkin projection consists in forcing the residual ℛ, written with the approximation [START_REF] Launay | Effect of viscous dissipation in the prediction of thermal behavior of an elastomer cylindrical flow[END_REF], to be orthogonal to each 𝜙 𝑘 (𝑟, 𝑧), 𝑘 = 1, … , 𝑚, so that the projection of the residual onto the subspace of ℒ 2 (Ω) generated by the 𝜙 𝑘 would be null.

According to the inner product [START_REF] Ghoreishy | Finite element analysis of a thermoplastic elastomer melt flow in the metering region of a single screw extruder[END_REF], the Galerkin projection of the energy equation writes:

〈ℛ, 𝜙 𝑘 〉 Ω = ∫ ℛ𝜙 𝑘 (𝑟, 𝑧)𝑑Ω Ω = 0 ∀𝑘 ∈ [1; 𝑚] (24)
Taking into account the respective forms of energy equation on the fluid domain Ω 𝑓 and on the solid domain Ω 𝑤 , Eq.( 24) writes:

〈ℛ, 𝜙 𝑘 〉 Ω = 𝜌𝐶 𝑝 ( 𝑛 𝑛 + 1 ) ( ∆𝑝 2𝐾𝐿 ) 1 𝑛 ∫ (𝑅 𝑛+1 𝑛 -𝑟 𝑛+1 𝑛 ) 𝜕𝑇 𝜕𝑧 𝜙 𝑘 𝑑Ω Ω 𝑓 ⏟ 𝒯 -𝐾 ( ∆𝑝 2𝐾𝐿 ) 𝑛+1 𝑛 ∫ 𝑟 𝑛+1 𝑛 𝜙 𝑘 𝑑Ω Ω 𝑓 ⏟ 𝒫 -∫ ∇ ⃗ ⃗ . (𝑘 𝑓 ∇ ⃗ ⃗ 𝑇)𝜙 𝑘 𝑑Ω Ω 𝑓 -∫ ∇ ⃗ ⃗ . (𝑘 𝑤 ∇ ⃗ ⃗ 𝑇)𝜙 𝑘 𝑑Ω Ω 𝑤 ⏟ 𝒟 = 0 ∀𝑘 ∈ [1; 𝑚] (25)
𝒯, 𝒫 and 𝒟 are respectively the transport, production and diffusion terms.

Four steps are successively performed, for which the detailed developments are not given here.

Noting that in cylindrical coordinates

𝜕𝑇 𝜕𝑧 = ∇ ⃗ ⃗ 𝑇. 𝑒 𝑧 ⃗⃗⃗ in the transport term, both the transport term 𝒯 and the diffusion term 𝒟 are integrated by parts using the Green formula:

∫ 𝑢 ⃗ . ∇ ⃗ ⃗ 𝑓𝑑Ω Ω = ∫ 𝑓𝑢 ⃗ . 𝑛 ⃗ 𝑑Γ Γ -∫ 𝑓∇ ⃗ ⃗ . 𝑢 ⃗ 𝑑Ω Ω (26)
This makes the boundary terms appear (Dirichlet-type for transport term and Neumanntype for diffusion term). Further treatment is required in order to make Dirichlet-type boundary conditions on Γ 𝑖𝑛 and Γ 𝑤𝑑 appear in the diffusion boundary terms. The relation ∇ ⃗ ⃗ (𝑓𝑔) = 𝑓∇ ⃗ ⃗ 𝑔 + 𝑔∇ ⃗ ⃗ 𝑓 is used for that matter.

2. The known thermal boundary conditions ( 13) to ( 16) are introduced in terms where they explicitly appear;

3. The temperature field approximation [START_REF] Launay | Effect of viscous dissipation in the prediction of thermal behavior of an elastomer cylindrical flow[END_REF] is inserted in integral terms on domain Ω and in remaining boundary terms;

4. Finally the resulting equation is written in compact form.

The Galerkin projection [START_REF] Clerc | Particle swarm optimization[END_REF] can finally be written as follows:

∑ 𝒜 𝑘𝑖 (𝐾, 𝑛)𝑋 𝑖 (𝐾, 𝑛) 𝑚 𝑖=1 = ℬ 𝑘 (𝐾, 𝑛) ∀𝑘 ∈ [1; 𝑚] (27)
Or under matrix-vector form:

𝒜(𝐾, 𝑛)𝑋(𝐾, 𝑛) = ℬ(𝐾, 𝑛) (28)
With:

𝒜 𝑘𝑖 (𝐾, 𝑛) = 𝜌𝐶 𝑝 ( 𝑛 𝑛 + 1 ) ( ∆𝑝 2𝐾𝐿 ) 1 𝑛 ( ∫ (𝑅 𝑛+1 𝑛 -𝑟 𝑛+1 𝑛 ) 𝜙 𝑖 𝜙 𝑘 𝑑Γ Γ 𝑜𝑢𝑡 ⏟ [𝒜 1 ] 𝑘𝑖 -∫ 𝜙 𝑖 ∇ ⃗ ⃗ . [(𝑅 𝑛+1 𝑛 -𝑟 𝑛+1 𝑛 ) 𝜙 𝑘 𝑒 𝑧 ⃗⃗⃗ ] 𝑑Ω Ω 𝑓 ⏟ [𝒜 2 ] 𝑘𝑖 ) + 𝑘 𝑓 ∫ ∇ ⃗ ⃗ 𝜙 𝑖 . ∇ ⃗ ⃗ 𝜙 𝑘 𝑑Ω Ω 𝑓 + 𝑘 𝑤 ∫ ∇ ⃗ ⃗ 𝜙 𝑖 . ∇ ⃗ ⃗ 𝜙 𝑘 𝑑Ω Ω 𝑤 ⏟ [𝒜 3 ] 𝑘𝑖 -𝑘 𝑓 ∫ ∇ ⃗ ⃗ (𝜙 𝑖 𝜙 𝑘 ). 𝑛 ⃗ 𝑑Γ Γ 𝑖𝑛 -𝑘 𝑤 ∫ ∇ ⃗ ⃗ (𝜙 𝑖 𝜙 𝑘 ). 𝑛 ⃗ 𝑑Γ Γ 𝑤𝑑 ⏟ [𝒜 4 ] 𝑘𝑖 ∀𝑘 ∈ [1; 𝑚], ∀𝑖 ∈ [1; 𝑚] (29)
And:

ℬ 𝑘 (𝐾, 𝑛) = 𝜌𝐶 𝑝 ( 𝑛 𝑛 + 1 ) ( ∆𝑝 2𝐾𝐿 ) 1 𝑛 ∫ 𝑇 𝑖𝑛 (𝑟) (𝑅 𝑛+1 𝑛 -𝑟 𝑛+1 𝑛 ) 𝜙 𝑘 𝑑Γ Γ 𝑖𝑛 ⏟ [ℬ 1 ] 𝑘 -𝑘 𝑓 ∫ 𝑇 𝑖𝑛 (𝑟)∇ ⃗ ⃗ 𝜙 𝑘 . 𝑛 ⃗ 𝑑Γ Γ 𝑖𝑛 ⏟ [ℬ 2 ] 𝑘 -𝑘 𝑤 ∫ 𝑇 𝑤 (𝑧)∇ ⃗ ⃗ 𝜙 𝑘 . 𝑛 ⃗ 𝑑Γ Γ 𝑤𝑑 ⏟ [ℬ 3 ] 𝑘 + 𝐾 ( ∆𝑝 2𝐾𝐿 ) 𝑛+1 𝑛 ∫ 𝑟 𝑛+1 𝑛 𝜙 𝑘 𝑑Ω Ω 𝑓 ⏟ [ℬ 4 ] 𝑘 ∀𝑘 ∈ [1; 𝑚] (30)
Equation ( 27) is a system of 𝑚 algebraic equations for 𝑚 unknowns that are the coefficients 𝑋 𝑖 , 𝑖 = 1, … , 𝑚, of the temperature field approximation [START_REF] Launay | Effect of viscous dissipation in the prediction of thermal behavior of an elastomer cylindrical flow[END_REF].

As shown by equations ( 29) and (30), 𝒜(𝐾, 𝑛) ∈ ℝ 𝑚×𝑚 and ℬ(𝐾, 𝑛) ∈ ℝ 𝑚 in equation ( 27)

both depend on space functions 𝜙 𝑖 (𝑟, 𝑧), 𝑖 = 1, … , 𝑚. These functions could be defined by different means: they might be arbitrarily chosen, obtained by POD or, in the frame of an identification method as in the present work, obtained via an optimization algorithm.

However, in order to simplify the ROM identification procedure, we are going to make an additional hypothesis about the form of space functions 𝜙 𝑖 (𝑟, 𝑧), 𝑖 = 1, … , 𝑚.

Additional hypothesis: separation of r and z dependencies

In [START_REF] Girault | Parametric low-order models in transient heat diffusion by MIM. Estimation of thermal conductivity in a 2D slab[END_REF] a heat diffusion ROM explicitly parametrized by thermal conductivity was developed.

As thermal conductivity was assumed to be uniform throughout the whole domain, it was possible to extract it from the integral diffusion term arising in the Galerkin projection. As a consequence, all matrices and vectors defined by integrals depending on space functions in the Galerkin projection were independent of thermal conductivity. In the frame of the Modal Identification Method, the space functions themselves were not identified. The special integral form of matrices and vectors arising from the Galerkin projection was not taken into account and the components of matrices and vectors in the ROM were identified through an optimization problem.

For the problem studied in the present paper, the flow behavior index 𝑛 of the pseudoplastic fluid is intrinsically associated with radial coordinate 𝑟, as shown by heat transport and production terms in local energy equation [START_REF] Chinesta | An overview of the proper generalized decomposition with applications in computational rheology[END_REF] and consequently in equations ( 29) and (30)

defining 𝒜(𝐾, 𝑛) and ℬ(𝐾, 𝑛) appearing in equation [START_REF] Rabinowitsch | Über die viskosität und elastizität von solen[END_REF]. The flow behavior index 𝑛 cannot therefore be extracted from integral terms [𝒜 1 ] 𝑘𝑖 and [𝒜 2 ] 𝑘𝑖 in eq.( 29) and [ℬ 1 ] 𝑘 and [ℬ 4 ] 𝑘 in eq.( 30). In order to make further simplifications in the ROM formulation, an additional hypothesis is thus made for space functions 𝜙 𝑘 (𝑟, 𝑧), 𝑘 = 1, … , 𝑚, which are now assumed to be written as:

𝜙 𝑘 (𝑟, 𝑧) = 𝜓 𝑘 (𝑟)𝜉 𝑘 (𝑧), 𝑘 = 1, … , 𝑚 (31)
This hypothesis allows to separate the radial direction 𝑟, intrinsically linked to the flow behavior index 𝑛, from the 𝑧 direction.

Using (31) the temperature field approximation defined by eq.( 22) is written as:

𝑇(𝑟, 𝑧, 𝐾, 𝑛) ≈ ∑ 𝜓 𝑖 (𝑟)𝜉 𝑖 (𝑧)𝑋 𝑖 (𝐾, 𝑛) 𝑚 𝑖=1 (32)
Injecting additional hypothesis (31) into equations ( 29) and ( 30) and taking into account the expressions of the divergence of a vector and of the gradient of a scalar in cylindrical coordinates, yields the following expressions for 𝒜(𝐾, 𝑛) and ℬ(𝐾, 𝑛):

𝒜 𝑘𝑖 (𝐾, 𝑛) = 𝜌𝐶 𝑝 ( 𝑛 𝑛 + 1

) ( ∆𝑝 2𝐾𝐿 )

Where 𝑉 𝑇,𝑖𝑛 , 𝑉 𝐷,𝑖𝑛 , 𝑊 𝑤 and 𝑊 are vectors of size 𝑚 which depend on space functions 𝜉 𝑖 (𝑧), 𝑖 = 1, … , 𝑚 and/or on their partial derivatives with respect to 𝑧. Vector 𝑊 𝑤 also depends on the wall temperature 𝑇 𝑤 (𝑧) of Dirichlet boundary condition [START_REF] Girault | Estimation of time-varying pollutant emission rates in a ventilated enclosure: inversion of a reduced model obtained by experimental application of the modal identification method[END_REF]. Their expressions are given in Appendix B.

Matrices 𝑀 𝑇 , 𝑀 𝐷1 , 𝑀 𝐷2 , 𝑀 𝐷3 and vectors 𝑉 𝑇,𝑖𝑛 , 𝑉 𝐷,𝑖𝑛 , 𝑊 𝑤 , 𝑊 are independent of parameters 𝐾 and 𝑛. In the frame of an approach such as the MIM based on the identification of a ROM, their components are identified using optimization techniques. In fact functions 𝜉 𝑖 (𝑧), 𝑖 = 1, … , 𝑚 are not computed, they are embedded in these matrices and vectors to be identified. A space discretization in direction 𝑧 is hence avoided and a reduction in this direction is obtained.

On the contrary, the 𝑟-dependence needs to be retained in the model because of terms 𝑟 𝑛+1 𝑛 in transport and production terms of equations ( 33) and (34). Considering the physics involved and our goal to build a ROM parametrized by the flow behavior index 𝑛, the radial direction is in fact "irreducible" due to its intrinsic link to 𝑛.

Further simplification of the diffusion terms is going to be made in the following.

Simplified form of diffusion terms

First derivatives of functions 𝜓 𝑘 (𝑟), 𝑘 = 1, … , 𝑚 which appear in some diffusion terms in equation (33) are sources of numerical errors. Moreover they should be calculated at each iteration step of the particle swarm optimization algorithm used for the ROM parameters identification and for each particle of the swarm. Thus their calculation may turn out to be time consuming, especially if the derivation is performed with a compact scheme on a fine radial mesh.

In addition, diffusion terms in equation (33) require the identification of three 𝑚 by 𝑚 matrices 𝑀 𝐷1 , 𝑀 𝐷2 , 𝑀 𝐷3 to build the ROM.

In order to reduce the number of parameters to be identified and to avoid numerical difficulties involved by the computation of first derivatives of functions 𝜓 𝑘 (𝑟), we have

chosen
Therefore, instead of using equations (33) and (34) for which matrices 𝑀 𝐷1 , 𝑀 𝐷2 , 𝑀 𝐷3 ∈ ℝ 𝑚×𝑚 and vector 𝑊 𝑤 ∈ ℝ 𝑚 need to be identified and first derivatives of space functions 𝜓 𝑘 (𝑟), 𝑘 = 1, … , 𝑚, have to be computed, we have chosen to use equations ( 35) and (36) for which matrix 𝑀 𝐷 ∈ ℝ 𝑚×𝑚 and vector 𝑉 𝐷,𝑤 ∈ ℝ 𝑚 are going to be identified.

Another consequence of these simplifications of diffusion terms is that we do not need to evaluate functions 𝜓 𝑘 (𝑟), 𝑘 = 1, … , 𝑚 in the wall depth (𝑟 ∈ ]𝑅; 𝑅 + 𝑒]) anymore, but only inside the channel (𝑟 ∈ [0; 𝑅]).

Final form of the Reduced Order Model

Usually the flow rate 𝑄 in the channel is given rather than the pressure difference ∆𝑝. The flow rate is defined as:

𝑄 = ∫ 𝑣 . 𝑛 ⃗ 𝑑𝑆 𝑆 = ∫ 𝑢 𝑧 (𝑟)2𝜋𝑟𝑑𝑟 𝑟=𝑅 𝑟=0 (37)
According to the velocity profile [START_REF] Quarteroni | Certified reduced basis approximation for parametrized partial differential equations and applications[END_REF] in the channel, the flow rate (37) writes:

𝑄 = 𝜋 ( 𝑛 3𝑛 + 1 ) ( ∆𝑝 2𝐾𝐿 ) 1 𝑛 𝑅 3𝑛+1 𝑛 (38)
The Reduced Order Model is composed of two sets of equations.

The first one is a set of 𝑚 algebraic equations:

∑
The second set of equations allows the computation of some temperatures of interest once vector 𝑋 has been computed. Here we focus on a particular observable set consisting in 𝑁 𝑟 temperatures referred to as 𝑌 𝑗 , 𝑗 = 1, … , 𝑁 𝑟 , along the radial direction at a given position 𝑧 * .

Using eq.( 32) and setting 𝛼 𝑖 = 𝜉 𝑖 (𝑧 * ), 𝑖 = 1, … , 𝑚, one obtains:

𝑌 𝑗 = 𝑇(𝑟 𝑗 , 𝑧 * ) = ∑ 𝜓 𝑖 (𝑟 𝑗 )𝜉 𝑖 (𝑧 * )𝑋 𝑖 𝑚 𝑖=1 = ∑ 𝜓 𝑖 (𝑟 𝑗 )𝛼 𝑖 𝑋 𝑖 𝑚 𝑖=1 𝑗 = 1, … , 𝑁 𝑟 (42)
The ROM finally consists in equations ( 39) and (42) (respectively corresponding to general forms [START_REF] Van Krevelen | Properties of polymers: their correlation with chemical structure; their numerical estimation and prediction from additive group contributions[END_REF] and ( 21)), where 𝒜 and ℬ are defined by ( 40 The discrete values 𝜓 𝑖 (𝑟 𝑗 ) of space functions 𝜓 𝑖 (𝑟), 𝑖 = 1, … , 𝑚 are computed using (43) and used to compute integral terms in ( 40) and ( 41) and are also used in (42). Building a ROM valid for any value of a third parameter (say 𝑄 for instance) in addition to 𝐾 and 𝑛 is possible but requires to use temperature datasets for different values of these three parameters.

Comments on the Reduced Order Model

Inlet temperature profile

As mentioned in section 4.3, because our goal is to build a ROM parametrized by the flow behavior index 𝑛, the radial dependence has to be kept in the model because of terms 𝑟 𝑛+1 𝑛 .

Hence the radial distribution 𝑇 𝑖𝑛 (𝑟) of the inlet temperature does not involve additional difficulties. In this work the developed ROM is valid only for a specific given profile 𝑇 𝑖𝑛 (𝑟) which is injected in the two first terms of the equation (41) during the ROM construction.

Building a ROM valid for any distribution 𝑇 𝑖𝑛 (𝑟) could be possible, for instance by using an approximation 𝑇 𝑖𝑛 (𝑟) = ∑ 𝑐 𝑘 𝑔 𝑘 (𝑟) 𝑞 𝑘=1

on a basis of known functions 𝑔 𝑘 (𝑟) with the coefficients 𝑐 𝑘 as additional explicit parameters in the ROM. In such case the ROM identification would require to use temperature datasets corresponding to several couples (𝐾, 𝑛) and several sets (𝑐 𝑘 , 𝑘 = 1, 𝑞).

Wall temperature profile

As shown in Appendix B, the wall temperature distribution 𝑇 𝑤 (𝑧) appears in vector 𝑊 𝑤 of equation ( 34):

[𝑊 𝑤 ] 𝑘 = ∫ 𝑇 𝑤 (𝑧)𝜉 𝑘 (𝑧)𝑑𝑧 𝑧=𝐿 𝑧=0 , 𝑘 = 1, … , 𝑚
The particular form of vector 𝑊 𝑤 is not taken into account in the ROM. In fact the term 34) is even replaced by 𝑘 𝑤 [𝑉 𝐷,𝑤 ] 𝑘 to get equation ( 36)

-𝑘 𝑤 [𝑊 𝑤 ] 𝑘 𝑑𝜓 𝑘 𝑑𝑟 (𝑅 + 𝑒) in equation (
and finally equation ( 41). Thus, in this work, the developed ROM is valid for a specific given wall temperature distribution 𝑇 𝑤 (𝑧) which does not appear explicitly in the ROM equations.

𝑇 𝑤 (𝑧) is in fact implicitly taken into account in the temperature data used for the ROM identification. Building a model valid for any distribution 𝑇 𝑤 (𝑧) would be conceivable by considering approximations of 𝑇 𝑤 (𝑧) and functions 𝜉 𝑘 (𝑧) on some bases of known functions, at the expense of adding many explicit parameters to the ROM.

Taking into account the dependence on angular coordinate

The ROM formulation has been developed for the axisymmetric case with no dependence on angular coordinate 𝜑. Should the wall temperature distribution be 𝑇 𝑤 (𝜑, 𝑧), vector 𝑊 𝑤 in equation ( 34) would be such as

[𝑊 𝑤 ] 𝑘 = ∫ ∫ 𝑇 𝑤 (𝜑, 𝑧)𝜉 𝑘 (𝜑, 𝑧)𝑑𝜑 𝜑=2𝜋 𝜑=0 𝑑𝑧 𝑧=𝐿 𝑧=0
, 𝑘 = 1, … , 𝑚.

However, 𝜑 would be handled just as 𝑧. Therefore, except for equation (42) that could take into account several angular positions, the final ROM formulation (equations (39), (40) and (41)) would be the same as in the axisymmetric case, with functions 𝜉 𝑖 (𝜑, 𝑧), 𝑖 = 1, … , 𝑚 embedded in the matrices and vectors to be identified. The ROM constitutive parameters would be identified in order to fit the particular non-axisymmetric temperature datasets.

Data generation and Reduced Order Model identification

Data generation

The ROM parameters identification requires some input-output data, that are a set of 𝑁 𝑐 𝑖𝑑 couples (𝐾, 𝑛) 𝑖 𝑑𝑎𝑡𝑎 , 𝑖 = 1, … , 𝑁 𝑐 𝑖𝑑 in a chosen range of values and the corresponding temperatures 𝑌 𝑗 𝑑𝑎𝑡𝑎 ((𝐾, 𝑛) 𝑖 𝑑𝑎𝑡𝑎 ), 𝑖 = 1, … , 𝑁 𝑐 𝑖𝑑 , 𝑗 = 1, … , 𝑁 𝑟 along the radial direction at the chosen position 𝑧 * . In the present work, the temperature datasets are computed with the Finite Elements detailed model briefly presented in section 2.3.

Reduced Order Model identification

In order to effectively build the ROM consisting of equations ( 39) and (42) with 𝒜 and ℬ defined by ( 40) and ( 41) respectively and space functions 𝜓 𝑖 (𝑟), 𝑖 = 1, … , 𝑚 defined by equation ( 43), components of matrices 𝑀 𝑇 , 𝑀 𝐷 ∈ ℝ 𝑚×𝑚 , vectors 𝑉 𝑇,𝑖𝑛 , 𝑉 𝐷,𝑖𝑛 , 𝑉 𝐷,𝑤 , 𝑊 ∈ ℝ 𝑚 , matrix 𝛽 ∈ ℝ 𝑚×𝑁 𝑏 and vector 𝛼 ∈ ℝ 𝑚 , need to be identified. Matrix 𝑀 𝐷 which corresponds to diffusion terms in eq.( 40) is symmetric (see term 𝒜 3 + 𝒜 4 in eq.( 29) or see eq.( 33) along with Appendix A), so only 𝑚(𝑚 + 1) 2 ⁄ components are identified to build 𝑀 𝐷 .

The ROM construction is recast into a parameter estimation problem. For a given model order 𝑚, all parameters to be identified, whose total number is 𝑁 𝑝𝑎𝑟𝑎𝑚 (𝑚) = 𝑚(𝑚 + (𝑚 + 1) 2 ⁄ + 5 + 𝑁 𝑏 ) with 𝑁 𝑏 = 15, are gathered in vector 𝜃(𝑚) ∈ ℝ 𝑁 𝑝𝑎𝑟𝑎𝑚 (𝑚) and identified through the minimization of a functional 𝒥 𝑖𝑑 (𝜃(𝑚) ) based on the quadratic deviation between temperature data computed with the detailed model (see previous section 5.1) on the one hand and corresponding temperatures computed with the ROM on the other hand, for all couples (𝐾, 𝑛) 𝑖 𝑑𝑎𝑡𝑎 , 𝑖 = 1, … , 𝑁 𝑐 𝑖𝑑 :

𝒥 𝑖𝑑 (𝜃(𝑚)) = ∑ ∑ (𝑌 𝑗 (𝜃(𝑚), (𝐾, 𝑛) 𝑖 𝑑𝑎𝑡𝑎 ) -𝑌 𝑗 𝑑𝑎𝑡𝑎 ((𝐾, 𝑛) 𝑖 𝑑𝑎𝑡𝑎 )) 2 𝑁 𝑟 𝑗=1 𝑁 𝑐 𝑖𝑑 𝑖=1 (44)
In order to quantify the global quality of the ROMs identification process, we use the mean quadratic deviation 𝜎 𝑖𝑑 (𝑚) (or root mean square of the residues) between reference data computed by the DM and corresponding ROM outputs, defined as:

𝜎 𝑖𝑑 (𝑚) = √𝒥 𝑖𝑑 (𝜃(𝑚)) (𝑁 𝑐 𝑖𝑑 × 𝑁 𝑟 ) ⁄ (45)
The optimization procedure for a given model order 𝑚 is shown on Figure 3. The minimization of 𝒥 𝑖𝑑 (𝜃(𝑚)) is performed with a Particle Swarm Optimization (PSO) algorithm [START_REF] Clerc | Particle swarm optimization[END_REF]. The stopping criterion is based on the lowest absolute value 𝜎 𝑖𝑑 (𝑚) (𝜃 𝑏𝑒𝑠𝑡,𝑎 (𝑚)) obtained. The process is terminated if either this quantity is not improved for a pre-defined number of successive iterations or if the accuracy wished by the user is reached: 𝜎 𝑖𝑑 (𝑚) (𝜃 𝑏𝑒𝑠𝑡,𝑎 (𝑚)) < 𝜀 𝑢𝑠𝑒𝑟 . A model of order 𝑚 = 1 is first of all identified, meaning that a vector 𝜃(1) is obtained by minimizing 𝒥 𝑖𝑑 (𝜃(1)). Then, a model of order 𝑚 = 2 is built by minimizing 𝒥 𝑖𝑑 (𝜃(2)). Models of higher orders are obtained by repeating the building process until the accuracy wished by the user is reached or until 𝜎 𝑖𝑑 (𝑚+1) ≈ 𝜎 𝑖𝑑 (𝑚) .

Vector 𝜃(𝑚) of structural parameters of the order 𝑚 model are used as initial guesses for the corresponding parameters in vector 𝜃(𝑚 + 1) to be estimated in the order 𝑚 + 1 model.

Main conceptual differences with POD-G and PGD approaches

The POD-G approach [START_REF] Park | Karhunen-Loève Galerkin method with decimated sampling technique for the simulation of complex fluids defined in the phase space[END_REF][START_REF] Brigham | Inverse viscoelastic material characterization using pod reduced-order modeling in acoustic-structure interaction[END_REF] is a two-step procedure. First space functions 𝜙 𝑘 (𝑟, 𝑧), 𝑘 = 1, … , 𝑚, would be computed by performing a Proper Orthogonal Decomposition (POD) on some temperature datasets corresponding to a set of (𝐾, 𝑛) couples. Functions 𝜙 𝑘 (𝑟, 𝑧) would then be plugged in the Galerkin projection corresponding to equations ( 27), ( 29) and (30) to form the ROM. In a different way, for the MIM, no POD is performed and the ROM parameters (including space functions) are identified simultaneously through an optimization procedure. Unlike MIM and POD-G, Proper Generalized Decomposition (PGD) [START_REF] Chinesta | An overview of the proper generalized decomposition with applications in computational rheology[END_REF] does not rely on datasets coming either from measurements or simulations to build the ROM: in PGD, space functions and potentially functions of parameters of interest, are found a priori.

Applying PGD to our steady problem would consist to seek the temperature field as a sum of products of functions of physical space coordinates and functions of the rheological parameters defining the fluid viscosity. A complete decomposition would write 𝑇(𝑟, 𝑧, 𝐾, 𝑛) ≈ ∑ 𝜓 𝑖 (𝑟)𝜉 𝑖 (𝑧)𝐾 𝑖 (𝐾)𝑁 𝑖 (𝑛) 𝑚 𝑖=1

but a decomposition of the form (32) is also possible, among others. Whichever decomposition is used, PGD would lead to a discrete solution in the physical space as well as in the parameters space. In the present work, MIM aims at building equation (39) for the computation of 𝑋 ∈ ℝ 𝑚 as a continuous function of 𝐾 and 𝑛, temperatures being then computed using (42).

Results and discussion

In this section, the feasibility of the method is presented by considering an elastomer flow experimentally tested with the experimental set-up previously described (Figure 1). The flow length is 𝐿 = 205 mm and the radius 𝑅 = 5 mm. The die is made of stainless steel, (𝑘 𝑤 = 49 W.m -1 .K -1 ), its external radius is 15.5 mm thus corresponding to 𝑒 = 10.5 mm. An experiment performed with a rotation speed of the extrusion screw of N = 70 RPM (corresponding to a mass flowrate of 𝑄 𝑚 = 41.2 kg.h -1 ) is retained as an example. The temperature profile measured at the outlet section is given in Figure 4. It is compared to the numerical outlet temperature profile, obtained with the physical model described in section 2.1. The inlet temperature profile is also shown in Figure 4. Pa.s n and 𝑛 = 0.529. The temperature profiles imposed on the die wall and the channel inlet correspond to the measured temperature profiles [START_REF] Launay | Analysis of the thermal and kinetic homogenization of an elastomer flow[END_REF]. As shown in Figure 4, the physical model fits well the experimental results, with a root mean square of 𝜎 = 1.4 K. Such numerical validation was also obtained with lower flowrates [START_REF] Launay | Analysis of the thermal and kinetic homogenization of an elastomer flow[END_REF].

As mentioned in section 5.1, the identification of the ROM parameters requires some inputoutput data. Thus, in the following work, the physical model is used by keeping this experimental configuration (thermal properties and boundary conditions). Only values of 𝐾 and 𝑛 are modified to provide the datasets necessary to the ROM building.

Input-output data for ROM identification and validation

Instead of a Cartesian regular mesh of this (𝐾, 𝑛) space, an Improved Hypercube Sampling (IHS) approach [START_REF] Beachkofski | Improved distributed hypercube sampling[END_REF] has been used in order to cover the space with a limited number of (𝐾, 𝑛)

couples. Using IHS the coordinates in the parameters space are regularly spaced out but the set is formed so that any two distinct (𝐾, 𝑛) couples do not share a common 𝐾 or 𝑛 value.

Identification of a series of ROMs

The radial temperature profiles computed by the Detailed Model (DM) for 𝑁 𝑐 𝑖𝑑 = 60 (𝐾, 𝑛) couples (blue diamond dots cloud in Figure 5) are shown in Figure 6. For a better understanding of how the temperature profile depends on parameters 𝐾 and 𝑛, temperature profiles for 6 (𝐾, 𝑛) couples among the 60 used for identification are shown in Figure 7.

Using the 60 input-output data a series of ROMs of order 𝑚 = 1 to 5 has been built by minimizing the functional 𝒥 𝑖𝑑 (𝜃(𝑚)) defined by equation (44).

The second line of Table 1 shows the value of 𝜎 𝑖𝑑 (𝑚) defined by equation ( 45) as a function of the ROM order 𝑚. As expected, 𝜎 𝑖𝑑 (𝑚) decreases with respect to order 𝑚, from 1.5878 °C for the order 1 ROM to 9.5 10 -3 °C for the order 5 ROM.

The mean quadratic error 𝜎 𝑖𝑑,𝑔𝑙𝑜𝑏(𝐾,𝑛) (𝑚) (𝑟 𝑗 ), 𝑗 = 1, … , 𝑁 𝑟 , for which the mean is made only on the 𝑁 𝑐 𝑖𝑑 = 60 (𝐾, 𝑛) couples, is also defined:

𝜎 𝑖𝑑,𝑔𝑙𝑜𝑏(𝐾,𝑛) (𝑚) (𝑟 𝑗 ) = √ ∑ (𝑌 𝑗 (𝜃(𝑚), (𝐾, 𝑛) 𝑖 𝑑𝑎𝑡𝑎 ) -𝑌 𝑗 𝑑𝑎𝑡𝑎 ((𝐾, 𝑛) 𝑖 𝑑𝑎𝑡𝑎 )) 2 𝑁 𝑐 𝑖𝑑 𝑖=1 𝑁 𝑐 𝑖𝑑 ⁄ (46)
Figure 8 shows the value of 𝜎 𝑖𝑑,𝑔𝑙𝑜𝑏(𝐾,𝑛)

(𝑚)
as a function of 𝑟, for model order 𝑚=1 to 𝑚=5. It can be observed that:

• For a given model order 𝑚, the error is of the same order of magnitude along the radius.

• The error tends however to be higher close to the duct wall, which is logical since temperature gradient is much higher close to the wall than elsewhere in the duct.

• For almost all positions, the error decreases with model order 𝑚, however slightly from 𝑚=3 to 𝑚=5. This is consistent with the global error 𝜎 𝑖𝑑 (𝑚) shown in Table 1.

Whereas 𝜎 𝑖𝑑 (𝑚) is a global quantity, Figure 9 shows the difference between radial temperature profiles computed by the DM on the one hand and by the order 5 ROM on the other hand, for each of the 60 (𝐾, 𝑛) couples: for all couples, residues are comprised between ± 0.025°C for 0 ≤ 𝑟 ≤ 4 mm and almost never exceed ± 0.05°C except in the very near wall zone where shear is maximal and temperature gradients are the largest.

Numerical validation of the identified ROMs

Once identified, the ROMs are tested with the 𝑁 𝑐 𝑣𝑎𝑙 = 30 (𝐾, 𝑛) couples (red square dots cloud in Figure 5) as inputs. The radial temperature profiles computed by the DM for these 30 (𝐾, 𝑛) couples are shown in Figure 10. As examples, temperature profiles for 6 (𝐾, 𝑛) couples among them are also shown in Figure 11.

Similarly to 𝜎 𝑖𝑑 (𝑚) (𝑟 𝑗 ), 𝑗 = 1, … , 𝑁 𝑟 between temperatures computed by the DM and corresponding ROM outputs are respectively defined for these 30 couples. The last line of Table 1 shows the value of 𝜎 𝑣𝑎𝑙 (𝑚) as a function of ROM order 𝑚. Just as 𝜎 𝑖𝑑 (𝑚) , 𝜎 𝑣𝑎𝑙 (𝑚) decreases with respect to order 𝑚, from 1.7155 °C for the order 1 ROM to 9.7 10 -3 °C for the order 5 ROM. As expected for this set of (𝐾, 𝑛) couples which had not been used for the construction of the ROMs, 𝜎 𝑣𝑎𝑙 (𝑚) is slightly higher than 𝜎 𝑖𝑑 (𝑚) for a given order 𝑚, except for 𝑚 = 2.

Figure 12 shows the value of 𝜎 𝑣𝑎𝑙,𝑔𝑙𝑜𝑏(𝐾,𝑛)

(𝑚)
as a function of 𝑟, for model order 𝑚=1 to 𝑚=5.

The observations corresponding to Figure 8 can also be made for .

In order to add insight to the global quantity 𝜎 𝑣𝑎𝑙 (𝑚) , Figure 13 shows the difference between radial temperature profiles computed by the DM on the one hand and by the order 5 ROM on the other hand, for each of the 30 (𝐾, 𝑛) couples: as for the models identification, residues for all couples are comprised between ± 0.025°C for 0 ≤ 𝑟 ≤ 4 mm and almost never exceed ± 0.05°C except in the very near wall zone.

The results of these validation tests allow us to conclude that the identified ROMs may be

Asset of using ROM in parameter estimation problem

Assumed values for the real flow 𝐾 = 11618 Pa.s n and 𝑛 = 0.529, which were deduced from measurement of pressure drop at several flowrates with the experimental device, have been used to compute the outlet temperature profile with both DM and order 5 ROM. The root mean square of the residues between both models is 𝜎 = 9.1 10 -3 °C for these parameters,

showing again that the ROM provides results very close to those of the reference model.

However Figure 4 shows that the temperature profile given by DM (and consequently ROM)

is not in full agreement with the profile measured on the experiment, especially at the location where temperature is maximal. Considering the simple geometry of the system and all precautions taken to properly handle experimental boundary conditions, we may assume that the DM built with the Ansys Polyflow software should be a good model of the experiment.

The discrepancy between simulated and measured outlet temperature profiles could thus be due to the fact that actual values of 𝐾 and 𝑛 for the real flow are different of the estimated values 𝐾 = 11618 Pa.s n and 𝑛 = 0.529. Better estimates for 𝐾 and 𝑛 could possibly be obtained from temperature measurements through an inverse problem, which would be solved efficiently thanks to the ROM, for which the computing time is less than 1 s compared to 390 s for the DM.

Conclusion

An approach based on the Modal Identification Method has been developed to build a thermorheological Reduced Order Model (ROM) for a steady incompressible flow of a pseudo-plastic fluid in a circular runner. The dynamic viscosity has thus been described by a shear rate power law defined by two parameters: consistency index 𝐾 and flow behavior index 𝑛.

Dependence of viscosity on temperature has not been taken into account and the analytical expression of the velocity field has been used. The aim was to build a ROM able to compute the radial temperature profile at the channel outlet for any (𝐾, 𝑛) couple in a predefined range.

Thus the ROM had to be explicitly parametrized by 𝐾 and 𝑛.

An approximation of the temperature field on a reduced set of 𝑚 space functions has been considered. In a first step a Galerkin projection of the energy equation has led to a first formulation of the ROM equations. In a second step an additional hypothesis about the form of space functions, consisting in separating the radial direction 𝑟, intrinsically linked to the flow behavior index 𝑛, from the 𝑧 direction, has allowed simplifying the ROM formulation. In a third step, a simplification of some diffusive terms has been made, allowing to reduce the number of fixed parameters to be identified during the ROM construction.

A numerical example has then been proposed. The fixed constitutive parameters of the ROM have been identified through an optimization procedure using temperature data for a set of 60

Figure 1 :

 1 Figure 1: Experimental device constituted of a cylindrical runner instrumented with two TMC (Thermal Measurement Cells).

Figure 2 :

 2 Figure 2 : Flow in circular channel. Schematic view including boundary conditions

  ) and (41) respectively. It is explicitly parameterized by 𝐾 and 𝑛 which define the dynamic viscosity. Parameters 𝜌, 𝐶 𝑝 , 𝑄, 𝑅, 𝑘 𝑓 , 𝑘 𝑤 in (40) and (41) are considered as fixed and known in the model. The inlet temperature profile 𝑇 𝑖𝑛 (𝑟) is also assumed to be known. Matrices 𝑀 𝑇 , 𝑀 𝐷 ∈ ℝ 𝑚×𝑚 and vectors 𝑉 𝑇,𝑖𝑛 , 𝑉 𝐷,𝑖𝑛 , 𝑉 𝐷,𝑤 , 𝑊 ∈ ℝ 𝑚 need to be identified to build the ROM, as well as coefficients 𝛼 𝑖 and space functions 𝜓 𝑖 (𝑟), 𝑖 = 1, … , 𝑚 .In order to ease the estimation of functions 𝜓 𝑖 (𝑟), 𝑖 = 1, … , 𝑚 , these are written as the linear combination of chosen functions which are here Bernstein basis polynomials[START_REF] Joy | Bernstein Polynomials[END_REF] 𝐵 𝑘,𝑁 𝑏 (𝑟 * ), 𝑘 = 1, … , 𝑁 𝑏 with 𝑁 𝑏 = 15 and 𝑟 * = 𝑟 𝑅 ⁄ ∈ [0; 1]: 𝜓 𝑖 (𝑟) = ∑ 𝛽 𝑖𝑘 𝐵 𝑘,𝑁 𝑏 (𝑟 * ) 𝛽 𝑖𝑘 , 𝑖 = 1, … , 𝑚, 𝑗 = 1, … , 𝑁 𝑏 with 𝑁 𝑏 = 15 are thus estimated instead of the 𝜓 𝑖 (𝑟 𝑗 ), 𝑖 = 1, … , 𝑚, 𝑗 = 1, … , 𝑁 𝑟 with 𝑁 𝑟 = 100.

4. 6 . 1

 61 Temperature outputs at other locations Equation (42) corresponds to a series of temperatures along the radial direction at a given position 𝑧 * . Adding temperatures 𝑌 𝑗 ′ , 𝑗 = 1, … , 𝑁 𝑟 at another position 𝑧 ′ is very simple. Using eq.(32) and setting 𝛾 𝑖 = 𝜉 𝑖 (𝑧 ′ ), 𝑖 = 1, … , 𝑚, one would have: 𝑌 𝑗 ′ = 𝑇(𝑟 𝑗 , 𝑧 ′ ) = ∑ 𝜓 𝑖 (𝑟 𝑗 )𝜉 𝑖 (𝑧 ′ )𝑋 𝑖 𝑚 𝑖=1 = ∑ 𝜓 𝑖 (𝑟 𝑗 )𝛾 𝑖 𝑋 𝑖 𝑚 𝑖=1 𝑗 = 1, … , 𝑁 𝑟 It would simply add 𝑚 parameters 𝛾 𝑖 , 𝑖 = 1, … , 𝑚 to be identified in the ROM. 4.6.2 Known fixed physical and geometrical parameters Density 𝜌, specific heat capacity 𝐶 𝑝 and thermal conductivity 𝑘 𝑓 of the fluid, wall thermal conductivity 𝑘 𝑤 , mass flow rate 𝑄, channel internal radius 𝑅, channel length 𝐿 and wall depth 𝑒 are all known fixed parameters in the ROM which is then valid only for the values used for all simulations of the DM providing the reference temperature data. 𝜌 and 𝐶 𝑝 could have been embedded in matrix 𝑀 𝑇 of (40) and in vector 𝑉 𝑇,𝑖𝑛 of (41). 𝑘 𝑓 and 𝑘 𝑤 could have been respectively embedded in vectors 𝑉 𝐷,𝑖𝑛 and 𝑉 𝐷,𝑤 of (41). However they have been retained explicitly in the ROM equations in order to provide a better understanding of the different transport and diffusion terms. Due to the simplifications made in sections 4.3 and 4.4, 𝐿 and 𝑒 do not appear explicitly in the ROM equations. Both are in fact implicitly taken into account in the temperature data used for the ROM identification. If one wants to build a ROM for other specific values of any of these fixed parameters, the identification procedure has to be performed again with temperature data corresponding to the new specific parameter values.

Figure 3 :

 3 Figure 3: Sketch of optimization procedure for ROM identification, for a given model order 𝑚.

Figure 4 :

 4 Figure 4: Comparison of experimental and numerical outlet temperature profiles of an elastomer flow. The inlet temperature profile is drawn in dashed line.

  Chosen ranges of parameter values are 𝐾 ∈ [5000; 20000] and 𝑛 ∈ [0.3; 0.6]. The input data for the identification of a series of ROMs are the set of 𝑁 𝑐 𝑖𝑑 = 60 couples (𝐾, 𝑛) depicted by the blue diamond dots in Figure 5. The output data are the corresponding 60 temperature profiles along the radial direction (𝑁 𝑟 = 100) at the channel outlet (𝑧 * = 205 mm), computed with the detailed model. Once the ROMs identified, some test cases are conducted for the ROMs validation: a second set of 𝑁 𝑐 𝑣𝑎𝑙 = 30 couples (𝐾, 𝑛) depicted by the red square dots in Figure 5 is used as input data for the ROMs and the resulting 30 outlet temperature profiles are compared with the profiles computed with the detailed model for the same couples (𝐾, 𝑛).

Figure 5 :

 5 Figure 5: Sets of (𝐾, 𝑛) couples generated by Improved Hypercube Sampling. The set of blue diamond dots is used as input data for ROMs identification. The set of red square dots is used as input data for ROMs validation.

Figure 6 :

 6 Figure 6: Outlet temperature profile computed by DM for 60 (𝐾, 𝑛) couples used for ROMs identification.

Figure 7 :

 7 Figure 7: Outlet temperature profile computed by DM for 6 (𝐾, 𝑛) couples among the 60 used for ROMs identification.

Figure 8 :

 8 Figure 8: ROMs identification. Mean quadratic error between ROM and DM as a function of 𝑟 (mean on 60 data (𝐾, 𝑛) couples), for model order 𝑚=1 to 𝑚=5.

Figure 9 :

 9 Figure 9: ROMs identification. Difference between outlet temperature profiles computed by DM and order 5 ROM for 60 (𝐾, 𝑛) couples.

Figure 12 .

 12 Values of error 𝜎 𝑣𝑎𝑙,𝑔𝑙𝑜𝑏(𝐾,𝑛) (𝑚) (𝑟) are also very close to those of 𝜎 𝑖𝑑,𝑔𝑙𝑜𝑏(𝐾,𝑛) (𝑚)

  used for any (𝐾, 𝑛) couple in the construction ranges 𝐾 ∈ [5000; 20000] and 𝑛 ∈ [0.3; 0.6]with a good knowledge about their accuracy with respect to the DM. In particular, the order 5 ROM is able to reproduce the temperatures computed by the DM within ± 0.05°C maximum except in the very near wall zone.

Figure 10 :

 10 Figure 10: Outlet temperature profile computed by DM for 30 (𝐾, 𝑛) couples used for ROMs validation.

Figure 11 :

 11 Figure 11: Outlet temperature profile computed by DM for 6 (𝐾, 𝑛) couples among the 30 used for ROMs validation.

Figure 12 :

 12 Figure 12: ROMs validation. Mean quadratic error between ROM and DM as a function of 𝑟 (mean on 30 test (𝐾, 𝑛) couples), for model order 𝑚=1 to 𝑚=5.

Figure 13 :

 13 Figure 13: ROMs validation. Difference between outlet temperature profiles computed by DM and order 5 ROM for 30 (𝐾, 𝑛) couples.

(

  𝐾, 𝑛) couples in the construction ranges 𝐾 ∈ [5000; 20000] and 𝑛 ∈ [0.3; 0.6]. These data were computed by a Finite Elements reference model experimentally validated on an instrumented apparatus. A series of reduced models of order 𝑚 = 1 to 5 has been identified and then tested for a second set of 30 (𝐾, 𝑛) couples. The order 5 ROM has been able to reproduce the radial outlet temperature profile computed by the reference model with a r.m.s. error of about 10 -2 °C. The good quality of the obtained results allows us to expect that such ROMs can be used for estimating 𝐾 and 𝑛 from temperature measurements. Moreover, thanks to their small number of degrees of freedom, ROMs allow fast resolution of the inverse problem for parameters estimation and hence can be used for on-line estimation of viscosity. An intermediate step could consist in improving developments by taking into account the thermodependence of the viscosity in the establishment of the ROM model. Appendices Appendix A: matrices in Eq.(33) Matrices in Eq.(33) have the following expressions: ∀𝑘 ∈ [1; 𝑚], ∀𝑖 ∈ [1; 𝑚]:

  Where 𝑀 𝑇 , 𝑀 𝐷1 , 𝑀 𝐷2 and 𝑀 𝐷3 are 𝑚 by 𝑚 matrices which depend on space functions 𝜉 𝑖 (𝑧), 𝑖 = 1, … , 𝑚 and/or on their partial derivatives with respect to 𝑧. Their expressions are given in

	Appendix A.
	ℬ 𝑘 (𝐾, 𝑛) = 𝜌𝐶 𝑝 ( 𝑛 + 1 𝑛	) (	∆𝑝 2𝐾𝐿	)	1 𝑛	[𝑉 𝑇,𝑖𝑛 ] 𝑘 ∫ 𝑇 𝑖𝑛 (𝑟) (𝑅 𝑟=𝑅 𝑟=0	𝑛+1 𝑛 -𝑟	𝑛+1 𝑛 ) 𝜓 𝑘 (𝑟)𝑑𝑟
	+ 𝑘 𝑓 [𝑉 𝐷,𝑖𝑛 ] 𝑘 ∫ 𝑇 𝑖𝑛 (𝑟)𝜓 𝑘 (𝑟)𝑑𝑟 𝑟=𝑅 𝑟=0	-𝑘 𝑤 [𝑊 𝑤 ] 𝑘	𝑑𝜓 𝑘 𝑑𝑟	(𝑅 + 𝑒)
	+ 𝐾 (	∆𝑝 2𝐾𝐿	)	𝑛+1 𝑛	𝑟=0 [𝑊] 𝑘 ∫ 𝑟 𝑟=𝑅	𝑛+1 𝑛 𝜓 𝑘 (𝑟)𝑑𝑟	∀𝑘 ∈ [1; 𝑚]
										1 𝑛	[𝑀 𝑇 ] 𝑘𝑖 ( ∫ (𝑅 𝑟=𝑅	𝑛+1 𝑛 -𝑟	𝑛+1 𝑛 ) 𝜓 𝑖 (𝑟)𝜓 𝑘 (𝑟)𝑑𝑟	)
										𝑟=0
	+ [𝑀 𝐷1 ] 𝑘𝑖 (𝑘 𝑓 ∫ 𝑟=𝑅 𝑟=0	𝑑𝜓 𝑖 𝑑𝑟	(𝑟)	𝑑𝜓 𝑘 𝑑𝑟	(𝑟)𝑑𝑟
	+ 𝑘 𝑤 ∫ 𝑟=𝑅+𝑒 𝑟=𝑅	𝑑𝜓 𝑖 𝑑𝑟	(𝑟)	𝑑𝜓 𝑘 𝑑𝑟	(𝑟)𝑑𝑟	-𝑘 𝑤	𝑑(𝜓 𝑖 (𝑟)𝜓 𝑘 (𝑟)) 𝑑𝑟	(𝑟 = 𝑅 + 𝑒))	(33)
										𝑟=𝑅
	+ ([𝑀 𝐷2 ] 𝑘𝑖 + [𝑀 𝐷3 ] 𝑘𝑖 ) 𝑘 𝑓 ∫ 𝜓 𝑖 (𝑟)𝜓 𝑘 (𝑟)𝑑𝑟
										𝑟=0
							𝑟=𝑅+𝑒
	+ [𝑀 𝐷2 ] 𝑘𝑖 𝑘 𝑤 ∫ 𝜓 𝑖 (𝑟)𝜓 𝑘 (𝑟)𝑑𝑟	∀𝑘 ∈ [1; 𝑚], ∀𝑖 ∈ [1; 𝑚]
							𝑟=𝑅

  to replace diffusion terms in equation (33) by a unique term [𝑀 𝐷 ] 𝑘𝑖 without taking into account the particular form of diffusion terms. This is possible since parameters 𝐾 and 𝑛 do not appear in the diffusion terms of equation (33) which is thus replaced by:

	𝒜 𝑘𝑖 (𝐾, 𝑛) = 𝜌𝐶 𝑝 ( 𝑛 + 1 𝑛	) (	∆𝑝 2𝐾𝐿	)	1 𝑛	[𝑀 𝑇 ] 𝑘𝑖 ( ∫ (𝑅 𝑟=𝑅 𝑟=0	𝑛+1 𝑛 -𝑟	𝑛+1 𝑛 ) 𝜓 𝑖 (𝑟)𝜓 𝑘 (𝑟)𝑑𝑟	)	(35)
		+ [𝑀 𝐷 ] 𝑘𝑖					∀𝑘 ∈ [1; 𝑚], ∀𝑖 ∈ [1; 𝑚]
	In a similar manner, the diffusion term inside the wall depth in equation (34), that is
	-𝑘 𝑤 [𝑊 𝑤 ] 𝑘	𝑑𝜓 𝑘
	ℬ 𝑘 (𝐾, 𝑛) = 𝜌𝐶 𝑝 ( 𝑛 + 1 𝑛	) (	∆𝑝 2𝐾𝐿	)	1 𝑛	[𝑉 𝑇,𝑖𝑛 ] 𝑘 ∫ 𝑇 𝑖𝑛 (𝑟) (𝑅 𝑟=𝑅 𝑟=0	𝑛+1 𝑛 -𝑟	𝑛+1 𝑛 ) 𝜓 𝑘 (𝑟)𝑑𝑟
								𝑟=𝑅
		+ 𝑘 𝑓 [𝑉 𝐷,𝑖𝑛 ] 𝑘 ∫ 𝑇 𝑖𝑛 (𝑟)𝜓 𝑘 (𝑟)𝑑𝑟	+ 𝑘 𝑤 [𝑉 𝐷,𝑤 ] 𝑘
								𝑟=0
		+ 𝐾 (	∆𝑝 2𝐾𝐿	)	𝑛+1 𝑛	𝑟=0 [𝑊] 𝑘 ∫ 𝑟 𝑟=𝑅	𝑛+1 𝑛 𝜓 𝑘 (𝑟)𝑑𝑟	∀𝑘 ∈ [1; 𝑚]

𝑑𝑟

(𝑅 + 𝑒), is replaced by 𝑘 𝑤 [𝑉 𝐷,𝑤 ] 𝑘 . Equation (34) is thus replaced by:

Table 1 :

 1

Identification and validation of ROMs. Values of 𝜎 𝑖𝑑

(𝑚)

and 𝜎 𝑣𝑎𝑙

(𝑚)

as functions of ROM order 𝑚.

  and 𝜎 𝑖𝑑,𝑔𝑙𝑜𝑏(𝐾,𝑛)

	(𝑚)	(𝑟 𝑗 ), 𝑗 = 1, … , 𝑁 𝑟 , mean quadratic deviations 𝜎 𝑣𝑎𝑙 (𝑚) and
	𝜎 𝑣𝑎𝑙,𝑔𝑙𝑜𝑏(𝐾,𝑛) (𝑚)