S Ouchtout

B Rousseau

Y Favennec
email: yann.favennec@univ-nantes.fr

Finite element framework for modeling conducto-radiative transfers within heterogeneous media at both discrete and continuous scales

Keywords: conductive-radiative transfer, continuous-scale approach, discrete-scale approach, open-cell foams, Newton-Raphson, finite elements, parallel computing

The paper compares three distinct conducto-radiative models that can be used for modeling heat transfer at high temperature in view of being integrated within a topology optimization process. Finite elements are extensively used. Non-linearities are handled using a fixed-point method coupled with a Newton-Raphson linearization method. A three-dimensional square lattice is sandwiched between a hot and a cold plate. The discrete mesoscopic scale is the reference model. At all conduction-to-radiation regimes, and even when radiation predominates, it appears that the continuous-scale approach, which has been set-up thanks to proper homogenization laws, gives satisfactory results, which is not the case for the Rosseland approximation model. The homogenized continuous-scale model appears to be a good candidate for being used in the future in topology optimization problems.

Introduction

In the field of free-carbon heat generation and waste heat recovery at high temperatures (T > 1000 o C), there is a growing interest in the design of durable high temperature energy systems with high thermal efficiency such as gas-to-gas heat exchangers [START_REF] Pelanconi | Application of ceramic lattice structures to design compact, high temperature heat exchangers: material and architecture selection[END_REF], volumetric solar receivers [START_REF] Fend | Porous materials as open volumetric solar receivers: Experimental determination of thermophysical and heat transfer properties[END_REF][START_REF] Wu | Experimental investigations of porous materials in high temperature thermal energy storage systems[END_REF][START_REF] Capuano | Optimized volumetric solar receiver: Thermal performance prediction and experimental validation[END_REF], thermal protection systems [START_REF] Helber | Microstructure and gas-surface interaction studies of a lowdensity carbon-bonded carbon fiber composite in atmospheric entry plasmas[END_REF][START_REF] Borner | High temperature permeability of fibrous materials using direct simulation Monte Carlo[END_REF][START_REF] Banerji | Tomography-based radiative characterisation of decomposing carbonaceous heat shield materials[END_REF][START_REF] Paglia | Carbon-phenolic ablative materials for re-entry space vehicles: plasma wind tunnel test and finite element modeling[END_REF], and radiant tube inserts [START_REF] Pelanconi | Cellular ceramic architectures produced by hybrid additive manufacturing: A review on the evolution of their design[END_REF], among other interesting processes. These systems can be described as radiativeconvective heat exchangers [START_REF] Luque | Exploiting volumetric effects in novel additively manufactured open solar receivers[END_REF] containing porous structures which play a key role in promoting the volumetric propagation of both radiation and temperature advected by a fluid flow. In reason of their textural and chemical features, refractory macroporous ceramics such as open-cell foams [START_REF] Scheffler | Cellular ceramics: structure, manufacturing, properties and applications[END_REF] and media with porosity-gradient structures [START_REF] Pelanconi | Cellular ceramic architectures produced by hybrid additive manufacturing: A review on the evolution of their design[END_REF]) know an increasing concern. Indeed, for engineering applications such as heat exchangers working at high-temperature, they present interesting performances in terms of energy conversion. When convective transfer can be neglected, which is the case for stagnant fluids, a global comprehension of both conductive and radiative transfers taking place within the porous medium structure can allow one to design them smartly, in order to best satisfy a pre-established objective.

The main perspective induced by the presented work is the design of architectured materials, and, in particular, lattice materials. These materials are becoming increasingly popular in the context of additive manufacturing. To achieve this goal, we exploit a numerical tool we have developed these last years: a fast and accurate finite element solver dedicated to the solution of the radiative transfer equation. This solver, which is based on stabilized vectorial finite elements [START_REF] Hardy | Solution of the 2-d steady-state radiative transfer equation in participating media with specular reflections using supg and dg finite elements[END_REF][START_REF] Badri | Vectorial finite elements for solving the radiative transfer equation[END_REF], has been parallelized [START_REF] Badri | High performance computation of radiative transfer equation using the finite element method[END_REF] in such a way that the calculation can run on multiple processes. Further, its matrix-free version highly reduces the memory load on each process, so that complex geometries can be dealt with [START_REF] Jolivet | Deterministic radiative transfer equation solver on unstructured tetrahedral meshes: Efficient assembly and preconditioning[END_REF]. This solver, which uses efficient matrix solvers and preconditioners [START_REF] Badri | Preconditioned Krylov subspace methods for solving radiative transfer problems with scattering and reflection[END_REF], takes into account several kinds of radiative boundary conditions such as incoming collimated radiative intensity, emission, and specular reflexions [START_REF] Hardy | Specular reflection treatment for the 3d radiative transfer equation solved with the discrete ordinates method[END_REF]. The solver is also able to design ad hoc angular discretization [START_REF] Favennec | Ad hoc angular discretization of the radiative transfer equation[END_REF].

Macroporous ceramics have two phases, one is solid and the other is fluid. A fine modeling of heat transfers taking place in this two-phases medium should be made, ideally, at the scale of the discretized material, i.e. at the pore scale. If the solid phase is considered as an absorbing medium, the conductive problem is to be solved within the solid phase, the radiative problem is to be solved solved within the fluid phase, which can be considered as vacuum, and the energy transfer between the two problems takes place at the solid-void interface separating apart both phases. Such a discrete-scale approach, finely described in a previous work [START_REF] Badri | Conductive-radiative heat transfer within sic-based cellular ceramics at high-temperatures: A discrete-scale finite element analysis[END_REF], makes it possible to achieve a good level of details in terms of three-dimensional temperature fields, conductive flux fields, and radiative flux fields. This is important because, at high temperature, heat transfer modes are coupled in such a way that a given heat transfer mode is affected by the other modes, making their effective contribution non linear [START_REF] Perraudin | Numerical quantification of coupling effects for radiation-conduction heat transfer in participating macroporous media: Investigation of a model geometry[END_REF].

In order to apply the homogenization methodology, the study considered here deals with the computation of coupled conducto-radiative heat transfer problems involving high temperature in a three-dimensional square lattice, using a single cell as well as multiple cells, at both the discrete mesoscopic scale and at the related continuous macroscopic scale. Along with this simple geometry, a highly simple physical situation is dealt with, where the considered domain is sandwiched in between a hot plate and a cold plate. The main goal of the study is to evaluate, in such a simple situation, the amount of errors one gets, in terms of temperature, using the coupled model at the macroscopic continuous scale respectively to the reference discrete-scale model at the mesoscopic level. Added to the continuous-scale model, the Rosseland approximation, which assumes that the porous medium is considered optically thick [START_REF] Modest | Radiative heat transfer[END_REF], is also compared to the discrete-scale model to get a clear idea on approximated models which could be used in future optimization processes, knowing the bias and errors associated to them.

In this study, three models are presented for coupled conducto-radiative problems at mesoscopic and macroscopic scales. One is the discrete-scale coupled model, presented above. Another is the related continuous-scale coupled model. The last one is the Rosseland approximation model. The passage from all of them is detailed in the paper. Conducto-radiative problems being non linear, an iterative scheme is to be set-up in order to make the solution to converge. We propose an iterative scheme based on a weak coupling between both the radiative transfer equation and the heat conduction equation and, in order to improve drastically the convergence, the Newton-Raphson iterative method is derived on the non-linear heat conduction equation.

The paper is organized as follows. Section 2 details the three considered conducto-radiative models. Section 3 details the solution process thanks to finite elements for conduction problems, and vectorial finite elements for the radiative transfer equation. This section also details iterative schemes such as the fixed-point method and especially the Newton-Raphson linearization process. Section 4 is dedicated to the global set-up, with useful non-dimensional numbers, geometry definition, and other necessary tools used for validations, convergence analysis, and comparisons between models. Section 5 details and discusses the numerical results of the three presented models. It firstly presents the numerical results related to the continuous-scale coupled model, with validation and convergence analysis. This section then presents the numerical results related to the Rosseland model, with a comparison respectively to the continuous-scale coupled model, and a convergence analysis. It finally presents the numerical results related to the discrete-scale coupled model, with a comparison respectively to the continuous-scale coupled models, and, again, a convergence analysis. The conclusion which can be drawn from numerical results is twofold. At first, a numerical process has been designed to solve efficiently coupled conducto-radiative problems at high temperature, at both the discrete and continuous scales. The models have been validated thanks to data from the literature. Second, if effective properties are carefully calculated, the continuous-scale coupled model approximates the discrete-scale coupled model with moderate errors of maximum 5 %. The corollary of this is that the continuous-scale coupled model could be the one used in the future in optimization processes when coupled with the homogenization methodology.

Models of conducto-radiative transfers

Heat transfer involved in macroporous media at high temperature may include the three modes of heat transfer, namely conduction, convection, and radiation. The balance between these modes is influenced by operating conditions as well by textural parameters such as porosity or ligament size, and also physical properties, such as thermal conductivity, scattering coefficient, etc. At high temperature, heat transfer modes may be coupled so that, for example, the global heat flux is not as simple as the summation of individual fluxes [START_REF] Perraudin | Numerical quantification of coupling effects for radiation-conduction heat transfer in participating macroporous media: Investigation of a model geometry[END_REF][START_REF] Badri | Conductive-radiative heat transfer within sic-based cellular ceramics at high-temperatures: A discrete-scale finite element analysis[END_REF]. The knowledge of individual contributions as well as their coupling effect is a crucial route for optimizing processes involving high temperature, modifying adequately space-dependent physical properties at macro-scale, or unitary cells at meso-scale. Among the three heat transfer modes, only the conduction and radiation are considered in the study; fluid flows contribution shall be taken care of in a forthcoming study. Further, in this paper, three distinct conducto-radiative models are considered: i) the discrete-scale coupled model (DS-CM) in which two phases are considered, namely the void phase in which radiation occurs, and the solid phase in which conduction occurs. The transfer between these two physics takes place on the internal boundary between both domains; ii) the continuous-scale coupled model (CS-CM) that couples both physics of radiation and conduction within a single domain. The physical properties are effective ones calculated both from the ones of the discretescale model and from the geometry itself (porosity, pore size, etc.); iii) the related Rosseland approximation model (RM) that combines both physics of radiation and conduction in a single non-linear conduction equation.

The considered physics deals with a medium located between two plates on which a given temperature is prescribed; on all other boundaries, a symmetry is considered.

The discrete-scale coupled model (DS-CM)

In the discrete-scale coupled model, the solid phase domain is fully embedded within a box. This box, without the solid phase, is the void phase. The solid phase is denoted as Ω s ⊂ R 3 , and the void phase is denoted as Ω v ⊂ R 3 . The union of these two domains yields the whole domain treated in the previous section, i.e. Ω s ∪ Ω v := Ω. Further, the solid (resp. void) phase boundary is denoted ∂Ω s (resp. ∂Ω v ). The solid-void interface is denoted Γ = ∂Ω s ∩ ∂Ω v . Also, let n s (resp. n v ) denote the outward unit vector normal to the solid (resp. void) phase. Figure 1 schematically represents the two distinct domains as well as the immersion process of the solid domain within the void domain. The discrete-scale coupled model is quite similar to the one presented in [START_REF] Badri | Conductive-radiative heat transfer within sic-based cellular ceramics at high-temperatures: A discrete-scale finite element analysis[END_REF]; it is however, for the sake of clarity, given hereafter.

The steady-state heat conduction equation is solved in the solid phase only. Since this phase is assumed to be opaque, the energy transfer between both phases occurs on the boundary only and, consequently, no source term is involved in this partial differential equation. The conduction problem then consists in searching for a scalar-valued function T s : Ω s → R such that:

E DS-CM HCE := -∇ • λ s ∇T s = 0 ∀x ∈ Ω s ⊂ R 3 (1)
where λ s denotes the thermal conductivity of the solid phase. Similarly to both other models presented hereafter, temperature is prescribed at the two ends of the medium, so:

T s = T g- on ∂Ω s,D,-⊂ Γ; ∂Ω s,D,-= {x ∈ ∂Ω s , x < δ + min Ωs x }, (2)
T s = T g+ on ∂Ω s,D,+ ⊂ Γ; ∂Ω s,D,+ = {x ∈ ∂Ω s , x < -δ + max Ωs x }. (3)
In these relationships, the subscript 'D' stands for 'Dirichlet', and δ is a sufficiently small positive userdefined parameter so that the Dirichlet condition is applied on a boundary of sufficiently large enough area, i.e. |∂Ω s,D,± | > . T g+ and T g-are the given hot and cold Dirichlet temperatures, respectively.

The conduction equation in the solid phase is additionally supplied with the flux exchange condition which takes place on the boundary and which takes into account of both the outgoing and incoming radiation:

-λ s ∇T s • n s = ε s σ B n 2 s (T s ) 4 -ε s s•nv>0 I v s • n v ds (4)
In this equation, T s represents the solid-void interface temperature, ε s denotes the emissivity of the solid medium, σ B denotes the Stefan-Boltzmann's constant, n s is the refraction index of the solid phase. The first term in the right-hand side of the equality sign denotes the emission loss, while the second term denotes the incoming radiation flux.

We further use the steady-state radiative transfer equations for the void which is non absorbing and non scattering medium, and where no emission occurs:

E DS-CM RTE := s • ∇I v = 0 (5)
The non-zero temperature yielding emission on the common void-solid boundary Γ comes into play as the boundary-value condition (the hat symbol stands for the 'inlet' condition), and the specular reflection boundary condition is also prescribed (it is the same type as for the continuous-scale coupled model):

I v (x, s) = 1 π σ B ε s n 2 s (T s ) 4 + (1 -s )I v (x, s ) on Γ, s • n v < 0 (6)
with ρ s = 1 -ε s denotes the surface reflectivity. Eventually, on the two ends of the domain, hot and cold temperature are prescribed, yielding the emission condition:

I v (x, s) = 1 π σ B n 2 T 4 g∓ ∀x = {x min , x max } (7)
To sum-up on this configuration, a solid medium is immersed within a void region represented by the bounding box such that the solid medium does not touch any boundary of the bounding box. The radiation that goes out of these six boundaries of the bounding box reenters the bounding box due to specular reflection. On the junction between the solid and the void, there is an energy exchange: there is emission due to temperature (from solid to void), radiation absorption (from void to solid), and partial reflection of the radiation (from void to void). Additionally, on two extreme parts of the solid medium, the temperature is prescribed such that a temperature gradient is numerically created.

The continuous-scale coupled model (CS-CM)

In order to describe the mathematical model that governs the radiation physics coupled with the conduction one, in a continuous medium, let us denote the open bounded set Ω ⊂ R 3 that represents the spatial domain of interest. (In the present study, this domain is a parallelepiped, i.e. Ω = {(x, y, z

)|x ∈ [0, L x ], y ∈ [0, L y ], z ∈ [0, L z ], }.)
The boundary of the domain, denoted as ∂Ω, forms the union of Neumann and Dirichlet conditions, denoted as Γ N and Γ D , respectively, and Γ D itself is the union of Γ + D and Γ - D , such that different conditions can be applied on the boundary of the three-dimensional domain, such as prescribed temperature, adiabatic condition, reflection of radiative intensity, emission, etc. For the mathematical derivations, let n denote the outward unit vector normal to domain Ω.

In this work, we use the standard steady-state form of the radiative transfer equation which allows to describe the propagation of thermal radiation through a medium affected by absorption, emission and scattering processes. The mathematical statements, based on the principle of energy conservation, leading to this integrodifferential equation are given for example in [START_REF] Modest | Radiative heat transfer[END_REF][START_REF] Howell | Thermal radiation heat transfer[END_REF]. Its use requires to carefully check that the textural features of the porous medium allow to fulfill the assumptions of randomness, homogeneity, and continuity [START_REF] Baillis | Thermal radiation properties of dispersed media: theoretical prediction and experimental characterization[END_REF]. The condition of randomness is fulfilled because, as seen in next sections, the geometric clearance is higher than the radiative wavelengths, and the condition of homogeneity is fulfilled because the typical size of individual scatterers is smaller than the characteristic size of the medium. Eventually, homogeneity is essential for considering the porous medium as a continuous one. Also, though the generalized radiative transfer equation [START_REF] Taine | Generalized radiative transfer equation for porous medium upscaling: Application to the radiative fourier law[END_REF] is to be applied in low porosity media, the ordinary version of the radiative transfer equation is considered as valid in the high porosity media considered here:

E CS-CM RTE := s • ∇I(x, s) + βI(x, s) -σ 4π I(x, s )Φ(s, s ) ds -κ 1 π σ B n 2 T 4 = 0 (8)
In this equation, x = (x, y, z) ∈ Ω, s ∈ S 2 is the direction of propagation, β = κ + σ is the extinction coefficient, σ is the scattering coefficient, κ is the absorbtion coefficient, σ B is the Stephan-Boltzmann constant, n is the refraction index of the medium, and T is the temperature. Also, Φ(s, s ) is the phase function which can be determined by numerical means. For example, based on a statistical representation of the porous medium, the radiative distribution function identification (RDFI) method completely characterizes this scattering phase function [START_REF] Tancrez | Direct identification of absorption and scattering coefficients and phase function of a porous medium by a monte carlo technique[END_REF][START_REF] Bellet | Rdfi determination of anisotropic and scattering dependent radiative conductivity tensors in porous media: Application to rod bundles[END_REF]. In the present study, among other possible approximations, the Henyey-Greenstein function is used. For three-dimensional domains, the normalised Henyey-Greenstein function reads:

Φ(s, s ) = 1 4π 1 -g 2 (1 + g 2 -2gs • s ) 3 2 (9)
where g ∈] -1, 1[ is the anisotropy factor. The integro-differential eq. ( 8) for the transport of radiation comes with its boundary conditions. In order to mimic emission where temperature is prescribed, i.e. on x = x min and x = x max , and periodicity on other boundaries (with specular reflection), one lets, for s • n < 0:

I(x, s) = 1 π σ B n 2 T 4 ∀x = {x min , x max } (10)
I(x, s) = (1 -)I(x, s ) ∀x ∈ ]x min , x max [ ( 11
)
where the emissivity of the boundary, and s is deduced from s from the Householder rotation matrix R(n) = 2nn -I, following [START_REF] Comninos | Mathematical and Computer Programming Techniques for Computer Graphics[END_REF].

The steady-state heat conduction equation involving the radiation contribution consists in equaling both fluxes of conduction and radiation. The former is given by the Fourier law, and the latter is the radiative source term. One gets [START_REF] Modest | Radiative heat transfer[END_REF]:

E CS-CM HCE := ∇ • (-λ∇T ) + k4πI b -kG = 0 ( 12
)
with λ being the thermal conductivity. In this equilibrium equation, G denotes the radiative density and I b denotes black body source term: 13)

G(x) = 4π I(x, S) ds (
I b = I b (T ) := 1 π σ B n 2 T 4 (14)
The partial differential eq. ( 12) for the heat conduction comes with its boundary conditions. Temperature is prescribed on x = x min and x = x max , and adiabatic condition is prescribed everywhere else:

T = T - ∀x = {x min } (15)
T = T + ∀x = {x max } (16)
-∇T • n = T N = 0 x ∈ ]x min , x max [ (17)

The continuous-scale Rosseland model (CS-RM)

In the Rosseland approximation, the radiative source term is an explicit non-linear function of the temperature, according to, for instance, [START_REF] Modest | Radiative heat transfer[END_REF]:

E CS-RM := -∇ • (λ∇T ) -∇ • (αT 3 ∇T ) = 0 (18) with α = 16n 2 σ B 3β ( 19
)
where n is the refraction index, σ B the Stephan-Boltzmann constant, and β the extinction coefficient. Note that the Rosseland approximation is usually to be used only on optically thick media; the criterion βL > 5 is, for instance given in [START_REF] Sans | Caractérisation des propriétés thermophysiques de mousses céramiques à haute température[END_REF], with L = |x max -x min | the length towards the main direction of propagation. Added to this non-linear partial differential equation, the boundary conditions ( 15)-( 17) are to be applied.

Resolution in a finite element framework

Three distinct physical models have been presented in section 2 in order to solve the conducto-radiative coupled physics in a bounded domain, in both continuous and discrete scales.

The radiative transfer equation involving absorption, scattering and emission is involved in the continuousscale coupled model (CS-CM), see eq. ( 8), while it is mathematically simplified to a transport advection-type equation in the discrete-scale coupled model (DS-CM), see eq. ( 5). In the following, the solution process is given only for the full RTE, i.e. for E CS-CM RTE ; the process being easily transposed to the discrete-scale, i.e. for E DS-CM RTE . The solution of the RTE with vectorial finite elements is the subject of section 3.1; additional details can be found in [START_REF] Badri | Vectorial finite elements for solving the radiative transfer equation[END_REF] and accompanying papers.

Besides the radiative equation, a non-linear heat conduction equation is also to be solved. Due to its nonlinear character, the equation has to be linearized, and an iterative scheme is set-up to make the state converge to the right solution.

Section 3.1 deals with the finite element set-up and the variational formulation for solving the radiative transfer equation; section 3.2 deals with the variational formulation and iterative schemes for solving the nonlinear heat conduction equation; section 3.3 deals with the variational formulation and iterative schemes for solving the non-linear Rosseland approximation model; eventually, section 3.4 gives the global iterative process.

Vectorial finite elements and discrete ordinates for solving the radiative transfer equation

The discrete ordinate method is used to angularly semi-discretize the RTE. Choosing, in an appropriate way, N d couples of direction/weight, using, for example, the S N method, or any other method based on the discretization of the unit sphere (in the current study we use the uniformly refined octahedron, as presented in [START_REF] Favennec | Ad hoc angular discretization of the radiative transfer equation[END_REF]), the problem consists in solving, ∀m = 1, . . . , N d :

(s m • ∇ + β)I m (x) = σ N d n=1 ω m I n Φ m,n + κI b (20)
The vectorial finite element method consists in searching a vector of radiative intensities I (vectorial trial function) using a vectorial test function W which lies in the corresponding vectorial functional space V := Π N d i=1 V i . In order to derive the vectorial FEM variational formulations for equation [START_REF] Perraudin | Numerical quantification of coupling effects for radiation-conduction heat transfer in participating macroporous media: Investigation of a model geometry[END_REF], the first step is to convert coupled set of equations ( 20) into its equivalent vectorial form. To do that, let us first define the radiative intensities vector I, the collective direction vector S collecting all individual directions vectors S m , and a matrix Θ of size N d × N d that combines scattering coefficient σ s , weights ω m , and scattering phase function Φ m,n . These quantities are expressed as:

I =      I 1 I 2
. . .

I N d      , S =      s 1 s 2 . . . s N d      , Θ =      σ s ω 1 φ 1,1 σ s ω 2 φ 1,2 . . . σ s ω N d φ 1,N d σ s ω 1 φ 2,1 σ s ω 2 φ 2,2 . . . σ s ω N d φ 2,N d . . . . . . . . . . . . σ s ω 1 φ N d ,1 σ s ω 2 φ N d ,2 . . . σ s ω N d φ N d ,N d      (21)
Using these notations, the RTE in its discrete ordinates form, eq. ( 20), can now be reformulated in its equivalent vectorial form. The vectorial equation reads:

S • ∇I + βI -ΘI = κI b I d ( 22
)
with I d being the identity vector of same order as I. Here, S • I would give a vector, the ith component of which will be given by (S • I) i = s i • ∇I i . Let us also introduce the following notations:

A T B = i A i B i and (A : B) = A i B i .
The vectorial FEM weak formulation can now be built by multiplying the vectorial eq. ( 22) with the SUPG vectorial test function W + γS • ∇W and integrating over the domain of interest Ω. Based on that, the problem now reads: search I in V that satisfies:

Ω (S • ∇I + βI) T (W + γS • ∇W) dx - Ω (ΘI) T (W + γS • ∇W) dx = Ω (kI b,T I d ) T (W + γS • ∇W) dx ∀W ∈ V (23)
The Green theorem is then applied for introducing the inflow boundary conditions, so that the problem consists in searching I in V h that satisfies:

Ω (S • ∇W) T I dx + ∂Ω (S • n) : H [S•n>0] ) T W dx (24) + ∂Ω (S • n : H [S•n<0] : I in ) T W dx Ω (S • ∇I) T (γS • ∇W) dx + Ω (βI) T (W + γS • ∇W) dx - Ω (ΘI) T (W + γS • ∇W) dx = Ω (kI b I d ) T (W + γS • ∇W) dx W ∈ V
where H [S•n>0] is the vectorial indicator (Heaviside) function. This function results in zeros and ones depending on Boolean operations, e.g., H [S•n>0]i equals one if and only if si S • n > 0, and zero elsewhere.

In typical treated situations, hundreds of directions are needed in the discrete ordinate method, and a fine enough spatial discretization is also needed to ensure that the numerical solution is close enough to the real solution. This yields huge matrix systems after approximating [START_REF] Taine | Generalized radiative transfer equation for porous medium upscaling: Application to the radiative fourier law[END_REF]. To cope with such difficulties, domain decomposition approach is performed on this vectorial equation, see [START_REF] Badri | High performance computation of radiative transfer equation using the finite element method[END_REF], the use of efficient solvers based on Krylov-subspace methods are used, see [START_REF] Badri | Preconditioned Krylov subspace methods for solving radiative transfer problems with scattering and reflection[END_REF], and semi-matrix free approach is set-up on the top of [START_REF] Taine | Generalized radiative transfer equation for porous medium upscaling: Application to the radiative fourier law[END_REF], see [START_REF] Jolivet | Deterministic radiative transfer equation solver on unstructured tetrahedral meshes: Efficient assembly and preconditioning[END_REF]. Note that this solver handles both absorbing-scattering and non absorbing-non scattering media in the same way; just physical properties are set to zero for the latter kind of medium.

Newton-Raphson and finite elements for solving the heat conduction equation

In order to give the finite element discretization of E CS-CM HCE , we introduce the following Sobolev spaces [START_REF] Brezis | Functional analysis, Sobolev spaces and partial differential equations[END_REF]:

H 0 := {φ ∈ L 2 (Ω) / ∇φ ∈ (L 2 (Ω)) d }, (25)
H Γ -:= {φ ∈ H 1 (Ω), φ = 0 on Γ - D }, (26)
H Γ + := {φ ∈ H 1 (Ω), φ = 0 on Γ + D }. ( 27
)
To begin with, the partial differential equation eq. ( 12) is multiplied by a test function belonging to the above functional space and integrated over the whole domain of integration Ω. Then, in order to introduce the Neumman boundary conditions, we apply Green's formula on the first integrand so that, eventually, the weak formulation yields in finding

T ∈ (T -+ H Γ -) ∩ (T + + H Γ + ), such that Ω λ∇T • ∇v dx = Ω κ Gv dx - Ω (κ4σ B n 2 T 4 )v dx ∀v ∈ H 0 . (28)
Equation ( 28) is nonlinear. Therefore, starting from an initial guess condition, say T 0 , an iterative scheme is to be set-up coupled with a linearization of the equation. Three linearization methods are presented hereafter (these methods shall be compared in the numerical results dedicated section).

1. The ordinary fixed-point method consists in expressing the non-linear term with the state given at the previous iteration. Denoting k + 1 the current iteration, and k the previous iteration, this yields in searching

T k+1 ∈ (T -+ H Γ -) ∩ (T + + H Γ + ), k ∈ N + , that satisfies: Ω λ∇T k+1 • ∇v dx = Ω κ G dx - Ω κ4σ B n 2 T 4 k )v dx ∀v ∈ H 0 (29)
2. The improved fixed-point method consists in linearizing the non-linear term such that, at the current iteration, one writes T 4 = T k+1 T 3 k . This strategy yields in searching

T k+1 ∈ (T -+ H Γ -) ∩ (T + + H Γ + ) such that: Ω λ∇T k+1 • ∇v dx + Ω κ4σ B n 2 T 3 k T k+1 v dx = Ω κ Gv dx ∀v ∈ H 0 . (30)
3. the Newton-Raphson linearization method [START_REF] Whiteside | The mathematical papers of Isaac Newton[END_REF] directly applied on the weak formulation eq. ( 28) is here proposed as a third method. After some mathematical manipulations, the finite element discretization of the non-linear heat conduction problem, using the Newton-Raphson linearization method, consists in searching

T k+1 ∈ (T -+ H Γ -) ∩ (T + + H Γ + ) such that: Ω λ∇T k+1 • ∇v dx + Ω 4κ4π 1 π σ B n 2 T 3 k T k+1 v dx = Ω 3(κ4σ B n 2 T 4 k )v dx + Ω (κ G)v dx ∀v ∈ H 0 . (31)

Newton-Raphson and finite elements for solving the Rosseland approximation

The Rosseland approximation eq. ( 18) which couples both physics of conduction and radiation yields a nonlinear equation. To solve this equation, finite elements are used along with an iterative scheme based on a linearization of the equation. The same functional space as the one used in section 3.2 is introduced, so that the variational formulation reads:

Ω λ∇T • ∇ϕ dx + Ω αT 3 ∇T • ∇ϕ dx = 0 ∀ϕ ∈ H 0 . (32)
Equation ( 28) is nonlinear. Therefore, starting from an initial guess condition, say T 0 , an iterative scheme is to be set-up coupled with a linearization of the equation. Two linearization methods are presented hereafter (these methods shall be compared in the numerical results dedicated section).

1. The ordinary fixed-point method consists in expressing the non-linear term with the state given at the previous iteration. Denoting k + 1 the current iteration, and k the previous iteration, this yields in searching

T k+1 ∈ (T -+ H Γ -) ∩ (T + + H Γ + ), k ∈ N + , that satisfies: Ω λ∇T k+1 • ∇ϕ dx + Ω αT 3 k ∇T k+1 • ∇ϕ dx = 0 ∀ϕ ∈ H 0 . ( 33
)
where H 0 , H Γ + and H Γ -are respectively the Sobolev spaces given in previous section, ϕ is the test function and T 0 is initial guess condition. 2. the Newton-Raphson linearization method [START_REF] Whiteside | The mathematical papers of Isaac Newton[END_REF] directly applied on the weak formulation eq. ( 32) is here proposed as a second method. After some mathematical manipulations, the finite element discretization of the non-linear heat conduction problem, using the Newton-Raphson linearization method, consists in searching

T k+1 ∈ (T -+ H Γ -) ∩ (T + + H Γ + ) such that: Ω λ∇T k+1 •∇ϕ dx+ Ω 3αT 2 k T k+1 ∇T k •∇ϕ dx+ Ω αT 3 k ∇T k+1 •∇ϕ dx = Ω 3αT 3 k ∇T k •∇ϕ dx ∀ϕ ∈ H 0 (34)

Numerical algorithm for solving the coupled problem

Added to non-linear conduction equations, the coupling itself between both radiation and conduction physics makes the problem to be non-linear. It is suggested the two physics are solved iteratively, one after the other. The following algorithm 1 has been set-up for the continuous-scale coupled model.

A linear steady-state heat conduction equation is firstly solved in order to get the initial guess (see line 5), which is needed to start the iterations. The outer loop (see line 7) is used to solve successively the radiative and then the conduction equation until stabilization is reached. The stabilization criterion is based on the infinite-based norm (in space) on temperature difference between the current iteration and the previous one, within the outer loop, but out of the next inner loop. The inner loop (see 12) is used for the non-linear heat conduction equation, until a stabilization criterion is satisfied; this one is based on the spatial infinite-based norm of the difference between the current iteration and the previous one, within the inner loop.

The set-up for the DS-CM can be easily derived from this one.

Algorithm 1 Iterative algorithm for solving the continuous-scale coupled conducto-radiative problem. T (k,l) ← T in 10:

calculate

I b = I b (T (k,l) ) 11:
solve RTE return I h 12:

for l ← 1 to l max do 13:

solve Conduction using I h return T (k,l)

14:

e (l) ← T (k,l) -T (k,l-1) ∞ 15:
if e (l) < ε l then

T (k,l-1) ← T (k,l)
19:

end for 20:

l(k) ← l 21:
T k ← T (k, l)

22:

E (k) ← T (k) -T (k-1) ∞ 23: if E (k) < ε k then 24:
break 25:

end if

26:

T in ← T (k) 27: end for

Set-up numerical validations and comparisons

This preliminary section details the set-up of all necessary tools used for validations, convergence analysis, and comparisons between models, all these being given in the next section.

Non-dimensional numbers

Hereafter are given the necessary non-dimensional numbers used in the next section, for characterizing tests, plot figures, and do proper comparisons:

A dimensionless temperature Θ(x) as well as the reference dimensionless temperature Θ 0 are defined using temperature of both the hot and the cold plates:

Θ(x) = T (x) T max ; Θ 0 = T min T max (35)
The optical thickness parameter τ (x) and the optical thickness of the medium, τ 0 , are defined by, respectively:

τ (x) = x 0 β dx = βx; τ 0 = βL (36)
In the next section, most plots will consider dimensionless temperature Θ as a function of the optical thickness τ (x).

The albedo of single scattering, also known as the Schuster number [START_REF] Viskanta | Heat transfer by simultaneous conduction and radiation in an absorbing medium[END_REF], is defined as the ratio between the scattering coefficient and the absorption coefficient. Both these relationships will be used in th next section:

ω = σ β ; κ β = 1 -ω (37)
Note that, at continuous scale, the albedo of single scattering is related to the emissivity, through:

ω = 1 -ε (38)
The conduction-to-radiation parameter [START_REF] Modest | Radiative heat transfer[END_REF][START_REF] Viskanta | Heat transfer by simultaneous conduction and radiation in an absorbing medium[END_REF], also known as the Plank number, or the Stark number, determines the relative role of both the conduction and radiative physics: conduction predominates for large N, while radiation predominates for low N.

N = λβ 4σ B T 3 max (39)

Geometry definition

The three-dimensional geometry used for continuous-scale models is a rectangular parallelepiped (note that this geometry has been prefered to the ordinary cube in order to save computational resources):

Ω = - L 2 , L 2 × - L 6 , L 6 × - L 6 , L 6 (40)
The chosen related discrete-scale geometry is the one defined by Perrausin and Haussener [START_REF] Perraudin | Numerical quantification of coupling effects for radiation-conduction heat transfer in participating macroporous media: Investigation of a model geometry[END_REF], and shown schematically in fig. 2-left. It is composed of squared-section bars. Symbols L, , and d nom denote the full length (it is the one used in the continuous-scale geometry), the squared-section size, and the nominal diameter of the pore, respectively. This geometry in fact is a single cell of a much bigger open-cell foam which is composed of these unit cells, see fig. 2-right. From geometrical considerations, both the porosity and the pore nominal diameter can be expressed explicitly as a function of the pore nominal diameter:

p( ) = L 3 -12 2 L + 16 3 L 3 ; d( ) = L -2 2 (41)
Besides, following [START_REF] Guévelou | A simple expression for the normal spectral emittance of open-cell foams composed of optically thick and smooth struts[END_REF][START_REF] Kribus | The promise and challenge of solar volumetric absorbers[END_REF], the effective extinction coefficient relative to a discrete-scale media is an explicit function of both the porosity and the pore nominal diameter, see eq. ( 45). So, for the given geometry composed of squared bars as presented in fig. 2, for a given full length L, and for a given effective extinction coefficient β, the bar width is to be solution of:

3 - 288L 384 2 - 5βL 3 384 + 5βL 4 768 = 0 (42)
This third-degree equation has three solutions. With the given values of L and β, the unique positive real solution is given by [START_REF] Murray | mathematical handbook of formulas and tabel[END_REF]:

= 2 -Q cos 1 3 ϕ - 1 3 a 1 ( 43
)
with

a 1 = -288L 384 ; a 2 = -5βL 3 384 ; a 3 = 5βL 4 768 ; Q = 3a 2 -a 2 1 9 R = 9a1a2-27a3-2a 3 1 54 ; ϕ = arccos R/ -Q 3
As reported in table 2, for exemple for L = 3.8 × 10 -4 m, and β = 2624.67, the porosity p has been found to be approximately equal to 0.9127, the nominal diameter of the pore d nom ≈ 1.5963 × 10 -4 m, and the squared section width ≈ 3.0869 × 10 -5 m. 4.3. Set-up physical properties from discrete to continuous scale and vice-versa Effective physical properties of a porous medium are functions of properties of the different phases at the discrete scale as well as on geometric characteristics. The determination of such effective properties by upscaling-type approaches is an active research subject as been denoted very recently [START_REF] Le | Homogenized models of conduction-advection-radiation heat transfer in porous media[END_REF]. In the radiative transfer community, usual routes are based on numerical tools such as Monte-Carlo or ray tracing simulations on representative elementary volumes [START_REF] Hendricks | Absorption/scattering coefficients and scattering phase functions in reticulated porous ceramics[END_REF][START_REF] Guévelou | Representative elementary volumes required to characterize the normal spectral emittance of silicon carbide foams used as volumetric solar absorbers[END_REF]. Following for example [START_REF] Guévelou | A simple expression for the normal spectral emittance of open-cell foams composed of optically thick and smooth struts[END_REF][START_REF] Kribus | The promise and challenge of solar volumetric absorbers[END_REF][START_REF] Kamiuto | Modeling of composite heat transfer in open-cellular porous materials at high temperatures[END_REF], the effective conductivity λ, the effective extinction coefficient β, and the effective emissivity ε are given by:

L d nom

λ = 1 3 (1 -p)λ s (44) β = 4.8 1 -p d nom (45) ε = 1 -(1 -ε s ) × (1 -0.8p) ( 46
)
where λ s , and ε s are the thermal conductivity, the emissivity of the solid phase involved in the discrete-scale model, and p and d nom denote the porosity and the nominal diameter of the pore, respectively. Note that the single scattering albedo and the emissivity in the continuous-scale case are denoted by ω and ε, respectively, while, however, in the discrete-scale case, they are denoted by ω s and ε s , respectively. In order to perform the passage from the discrete scale to the continuous scale, and vice-versa, let us first point out that the continuous-scale coupled model is fully defined by these four parameters: i) the effective thermal conductivity λ, ii) the effective absorption coefficient κ, iii) the effective scattering coefficient σ, and iv) the effective emissivity of the medium . Besides, the discrete-scale coupled model is fully defined by these four parameters: i) the solid-phase thermal conductivity λ s , ii) the solid-phase emissivity ε, iii) the porosity p, and iv) the nominal diameter of the pore d nom .

The knowledge of effective properties at the continuous scale, for a given geometry (so for a given porosity and nominal pore diameter) gives directly, using eqs. ( 44) and ( 46), both the thermal conductivity and the emissivity of the solid phase:

λ s = 3 λ 1 -p (47)
ε s = 1 - 1 -ε 1 -4 5 p (48)
At the discrete scale, the knowledge of the solid-phase properties λ s and ε s along with the porosity p and the pore nominal diameter d nom , determines the continuous-scale model, as follows. Firstly, eqs. ( 44) and (46) give directly the effective thermal conductivity λ and the effective emissivity ε. Then, combining eqs. [START_REF] Guévelou | Representative elementary volumes required to characterize the normal spectral emittance of silicon carbide foams used as volumetric solar absorbers[END_REF] to [START_REF] Viskanta | Heat transfer by conduction and radiation in absorbing and scattering materials[END_REF] and ( 44) to (46) yields the Stark number N as an explicit function of λ, p, and d nom , see eq. ( 49). The effective absorption coefficient κ is then expressed as an explicit function of the Stark number N, the solid-phase emissivity ε s , the porosity p, and the effective thermal conductivity λ, see eq. (50). Eventually, the effective scattering coefficient σ is expressed as an explicit function of solid-phase emissivity ε s , the porosity p, and the effective absorption coefficient κ, see eq. (51).

N = 4.8λ(1 -p) 4σ B T 3 max d nom (49) κ = 4σ B T 3 max N 1 -(1 -ε s )(1 -4 5 p) λ (50) σ = 1 1 -(1 -ε s )(1 -4 5 p) -1 k (51)
Note that the anisotropy coefficient involved in the Henyey-Greenstein approximation function eq. ( 9) cannot be determined directly by this upscaling-type methodology. It could however be determined by the radiative distribution function identification (RDFI) method, for example, which completely characterizes, among other quantities, the scattering phase function [START_REF] Tancrez | Direct identification of absorption and scattering coefficients and phase function of a porous medium by a monte carlo technique[END_REF][START_REF] Bellet | Rdfi determination of anisotropic and scattering dependent radiative conductivity tensors in porous media: Application to rod bundles[END_REF]. In the current study, the anisotropy coefficient has been set to zero yielding isotropic scattering. This can be justified numerically, see section 5.1.3, because, in the configuration dealt with in this study, output temperature is little sensitive to this parameter.

Some more quantities

Some useful quantities are now defined; these are used in order to plot one-dimensional non-dimensional averaged temperature out of results obtained in three-dimensional domains, as well as to quantify errors. The non-dimensional x-averaged temperature T (x) is the non-dimensional averaged temperature on a given yzplane located at a specific x = x (this plane is denoted as A x). For the discrete scale (for the continuous-scale problems, the indicator function is taken equal to one everywhere), it is:

T (x) := 1 T max A x T (x)1 [x∈Ωs] dy dz A x 1 [x∈Ωs] dy dz . (52)
Further, in order to compare the three distinct models to each other, we denote by e α-β the difference between the computational results obtained from models α and β:

e α-β (x) = T α (x) -T β (x). ( 53
)
and these two ordinary norms are used (the first is an averaged quadratic error norm while the second gives the maximum error along the main axis of propagation of the heat flux):

e 2 (u) = 1 L L 0 u 2 (x) dx; e ∞ (u) = max x∈[0,L] |u(x)| (54)

Numerical results

This section details and discusses the numerical results of the three presented models. Section 5.1 firstly presents the numerical results related to the continuous-scale coupled model, with validation and convergence analysis. Section 5.2 then presents the numerical results related to the Rosseland model, with a comparison respectively to the CS-CM, and a convergence analysis. Section 5.2 finally presents the numerical results related to the DS-CM model, with a comparison respectively to the CS-CM, and a convergence analysis.

Continuous-scale coupled model (CS-CM) 5.1.1. Validation

The validation of the three-dimensional continuous-scale coupled model is performed against the onedimensional pioneering work of Viskanta et al. [START_REF] Viskanta | Heat transfer by simultaneous conduction and radiation in an absorbing medium[END_REF][START_REF] Viskanta | Heat transfer by conduction and radiation in absorbing and scattering materials[END_REF]. Based on these references, the characteristic length L is chosen to be 3.8 × 10 -4 m. The mesh needed for the approximation of the functional space is chosen to be sufficiently fine enough, and quasi-uniform, with mean size h = L/80. Such a choice yields N e = 281, 247 tetrahedral elements and N v = 47, 068 vertices. Further, in order to have sufficiently fine enough angular discretization, the octahedron has been refined twice, so that N d = 128 directions are considered.

Following [START_REF] Viskanta | Heat transfer by simultaneous conduction and radiation in an absorbing medium[END_REF], the first validation is performed on an absorbing non scattering medium. Let us point out that our model is three-dimensional. Further, it involves emission on both hot and cold boundaries, it involves body emission, and it involves specular reflection on other boundaries. Specular reflections, in the present case, follow the 1-DP method, see [START_REF] Hardy | Specular reflection treatment for the 3d radiative transfer equation solved with the discrete ordinates method[END_REF]. As such, this first validation step implicitly validates all the cited mathematical operators. Following rigorously [START_REF] Viskanta | Heat transfer by simultaneous conduction and radiation in an absorbing medium[END_REF], the dimensionless temperature Θ 0 is set to 0.5, with T min = 1, 111 K and T max = 2, 222 K, the optical thickness of the medium τ 0 is set to 1 by setting the absorption coefficient κ to 2624.67 cm -1 .

The results are plots of adimensionalized (averaged) temperature along the main direction of propagation, x, for Stark numbers equal to 10, 1, 0.1, and 0.01 (the Stark number is modified with the value of the thermal conductivity parameter, following eq. ( 39)). Figure 3 presents such plots for both our results, and those from Viskanta et al. [START_REF] Viskanta | Heat transfer by simultaneous conduction and radiation in an absorbing medium[END_REF] (fig. 2, page 68). This figure shows that the line plots (which represent our three-dimensional results) coincide with the empty marks which represent the results obtained by Viskanta et al. in 1D. Notice however that the non-dimensional temperature for the three-diensional case is an averaged one, performed over the transverse section of the main direction of propagation of the flux, i.e. e x , see eq. ( 52).

We notice that in the case of N=10 the temperature behavior is almost linear; this is explained by the fact that the conduction physics is dominant over the radiation physics. While the Stark number decreases, the radiation physics is no more negligible and, consequently, the obtained temperature profile moves away from the linear curve, see for example orange and red curves in fig. 3, for N = 0.1 and N = 0.01, respectively.

The properties needed to run this first test case are summarized in table 2 (Test A). Next, fig. 4 presents the three-dimensional temperature fields for both extreme cases, for N = 10, i.e. when the conduction physics predominates, and for N = 0.01, i.e. when radiation predominates, respectively. It can be seen from this figure how temperature is somehow homogenized due to radiation, as compared to the case where there is conduction only.

Following [START_REF] Viskanta | Heat transfer by conduction and radiation in absorbing and scattering materials[END_REF], the second validation test case considers also the scattering operator, added to all others already considered in the previous validation test case. Our numerical results are plotted and compared against those cited, setting the Stark parameter N to the value 0.1, the extinction coefficient β = 2624.67 m -1 , and the optical thickness to the value τ 0 = 1 m. The validation is performed for several albedos, i.e. ω equal to 0, 0.5, and 1, respectively, and for two temperature ratio θ 0 equal to 0.5 and 0.1, respectively. All the results are summarized in table 2. Figure 5 presents all the corresponding plots in terms of variation of the dimensionless temperature Θ along the main direction of propagation, for different values of ω and Θ 0 . In can be seen from presented curves that our threedimensional results very well coincide with those one-dimensional of Viskanta; this validates completely the numerical modeling of coupled conducto-radiative problem in three-dimensions, with operators of absorption, scattering, emission, and reflection. Please note that other strategies have been developped for validation of the radiative transfer equation only; see our accompanying papers in the references section.

Convergence analysis

This section presents and discusses convergence results for the three schemes presented in section 3.2. As a first test case, a scattering medium is considered, with temperature ratio θ 0 = 0.1, albedo ω = 0.5, and Stark number N = 0.1. Figure 6 (top-left) shows, for this particular case, the error evolution with respect to the iterations of the inner loop, for the very first iteration of the outer loop (see algorithm 1, line 14). For this particular case, the three iterative schemes are convergent. However, the improved fixed-point method (in orange) converges faster than the ordinary fixed point method (in red), and the Newton-Raphson method (in green) outperforms both fixed-point-type methods. Please note (it is not reported here in plots for the sake of conciseness) that, in all considered cases, when the scheme was convergent, the stabilization was reached after three iterations of the outer loop, at most. We also indicate that the numerical resolution was carried out with the tolerance condition, for both the outer and inner loops, ε k and ε l , equal to 10 -4 (cf. algorithm 1).

As a second test case, a non scattering medium is considered, i.e. σ = 0, with the Stark number N = 0.1. Figure 6 (top-right) shows, for this particular case, the error evolution with respect to the iterations of the inner loop, for the very first iteration of the outer loop. It can be seen that both fixed-point-type methods diverge. The Newton-Raphson method is the only one that converges, and so, very efficiently.

As a third test case, the Rosseland model RM is considered along with ω = 0, Θ 0 = 0.5 and N = 0.1. Figure 6 (Bottom-left) shows that the Newton-Raphson method converges much faster than the ordinary fixed point method. Note also that the fixed-point method did not converge for the particular case N=0.01.

As a fourth test case, the discrete-scale coupled model is considered with ω s = 0, Θ 0 = 0.5 and N = 0.01. Figure 6 (Bottom-right) shows, for this particular case, the error evolution with respect to the iterations of the inner loop, for the very first iteration of the outer loop. It can be seen that both fixed-point-type methods diverge; the Newton-Raphson method is the only one that converges, and so, very efficiently.

With the numerical test that we performed, we can also conclude that the computational schemes converge quickly when conduction predominates (the physics is close to being linear), but, when radiation predominates, the order of convergence is lower (see Figure 7), and the Newton-Raphson method is the only one that gives satisfactory results (see Table 1).

Effect of the anisotropy coefficient

The knowledge at the discrete scale of both the solid-phase conductivity and emissivity, for a given porosity and pore nominal diameter enables us to evaluate both absorption and scattering coefficients of the related effective medium. However, the phase function cannot be determined in such a straightforward way, except using numerical tools such as, for example, the radiation distribution function identification method.

Figure 8 presents on the left the averaged temperature in the continuous-scale effective medium for two very different anisotropy coefficients involved in the Henyey-Greenstein approximation function: g = -0.9 for modeling strong backward scattering, and g = +0.9 for modeling strong forward scattering. Such curve has been built for two radiation-to-conduction regimes, with the Stark number equal to 0.01, and 0.0001. Figure 8 presents on the right the temperature difference between these both extreme cases and the isotropic scattering case. From this figure, it is seen that the averaged temperature is little sensitive to the anisotropy coefficient (errors are less than 2 %), at least for this considered test case where the effective medium is sandwiched between two given temperatures.

In the whole rest of the study, the anisotropy coefficient has been set to zero yielding isotropic scattering. Left: temperature profiles with N ∈ {0.01; 0.001} and g ∈ {0.9; 0; -0.9}. Right: temperature difference for the same values of N .

Rosseland

The methodology for solving the Rosseland approximation based on finite elements coupled with the Newton-Raphson method is first validated. The method of manufactured solution is used, see for example [START_REF] Gao | A fast-forward solver of radiative transfer equation[END_REF][START_REF] Hardy | Solution of the 2-d steady-state radiative transfer equation in participating media with specular reflections using supg and dg finite elements[END_REF] for other examples in the field of radiation. The manufactured solution is built such that the same Neumann and Dirichlet conditions as that involved in the CS-CM are used. Further, a cubic interpolation has been performed on the CS-CM solution for the case N = 1, ω = 0, Θ 0 = 0.5 and τ 0 = 10, yielding a reference solution. All operators involved in the Rosseland model are applied on this reference solution, yielding a source term. Solving the Rosseland model with this source term gives a numerical solution which is to be compared with its reference solution. The comparison is presented in fig. 9-bottom in brown color. As can be seen, the error between the manufactured solution and its reference solution is negligible (in this case the quadratic error is around 0.1 %). This validates numerically the developed methodology based on finite elements and linearization for solving the Rosseland approximation model.

The Rosseland approximation is then compared with respect to the continuous-scale coupled method, in terms of averaged temperature. Three test cases, all recorded in table 2, are used for the comparison. Figure 9 plots the averaged temperature along the e x axis, as well as the difference with temperature obtained by the CS-CM. Furthermore, table 2 gives errors ê2 , and ê∞ between both the Rosseland model and the CS-CM.

A first test case deals with a non-scattering medium, for maximum temperature ratio Θ 0 equal to 0.5, an optical thickness of 1, and for Stark numbers equal to 10, 1, 0.1, and 0.01 (see fig. 9 and table 2, test A). As can be seen, the maximum error is about 5 % at maximum for the case exhibiting high radiation. Other tests, where radiation is less significant, give errors around 1 %.

The following test case (B) deals with errors evolution for different albedos and maximum temperature ratio. As can be seen from table 2 and fig. 9, maximum error is about 10 %. The higher the albedo, the higher the error. Comparing test B-1 and B-2, we conclude that the lower the maximum temperature ratio, the higher the error.

We deduce from these comparisons that the temperature profiles obtained with the Rosseland model tend towards the ones obtained with the continuous-scale couple model when conduction dominates; otherwise, errors become large.

Convergence analysis is also performed on the Rosseland model. The error evolution in terms of iterations is plotted in fig. 6 (bottom-right) for both the ordinary fixed-point method and the Newton-Raphson method, following eq. ( 33) and eq. ( 34), respectively. From this figure, both methods converge to the solution, for this particular considered case, and the Newton-Raphson method is more efficient in terms of convergence speed.

discrete-scale coupled model (DS-CM)

This section presents some coupled conducto-radiative problems similar to those of previous sections, but this time at a discrete scale. The physical properties as well as the discrete geometry have been chosen accordingly to fit with physical properties of the continuous-scale coupled model, following the derivation given in section 4.

In order to perform the immersion of the solid phase within a bounded convex domain representing the void, we consider a cube of size just slightly bigger than L Ωs .

Firstly, a single cell of the open-cell foam is considered. Then, three and a half cells are considered. Temperature in solid parts as well as radiative density are plotted for both cases, see fig. 10.

For the single-cell case, the geometry as well as physical properties have been set up in such a way that the averaged temperature along the main direction of propagation could be compared with those related to the effective medium. The way this has been done has been explained before in section 4. A quasi-uniform mesh is built for both computational geometries. The mesh size requirement has been set to h i = L Ωs /80 for the spatial discretization. This lead to numbers for degrees of freedom equal to (N e = 263, 931 tetrahedral elements ; N v = 48, 731 vertices) and (N e = 1, 426, 938 tetrahedral elements ; N v = 235, 729 vertices) for the solid and the void, respectively. The number of directions has been set to N d = 128 in order to reach the convergence relatively to the angular discretization. To run such test, 62 processes have been used. Plots of averaged temperature are recorded in fig. 11 for test A (i.e. modifying the Stark number N while keeping null albedo and maximum temperature ratio equal to one half.) and for test B (i.e. modifying the albedo while keeping other quantities of interest constant). From these plots, it is seen that, both CS-CM and DS-CM follow, roughly speaking, the same trend. The main differences between the two models appear close to boundaries where Dirichlet conditions are prescribed. In theses regions, the CS-CM overestimates the temperature gradient. However, the maximum error of the continuous model is about 5 %, and the quadratic error is lower than 1 %.

For the three-cells case, similarly as above, the geometry as well as physical properties have been set up in such a way that the averaged temperature along the main direction of propagation could be compared with those related to the effective medium. A quasi-uniform mesh is built for both computational geometries. The mesh size requirement has been set to h i = L Ωs /80 for the spatial discretization. This lead to numbers of degrees of freedom equal to (N e = 331, 861 tetrahedral elements ; N v = 66, 489 vertices) and (N e = 1, 753, 681 tetrahedral elements ; N v = 289, 968 vertices) for the solid and the void, respectively. The number of directions has been set to N d = 128 in order to reach the convergence relatively to the angular discretization. To run such test, 62 processes have been used. Plots of averaged temperature are recorded in fig. 12 for test A (i.e. modifying the Stark number N while keeping null albedo and maximum temperature ratio equal to one half.), both for the CS-CM and the DS-CM. From these plots, it is seen that, both CS-CM and DS-CM follow the same trend. The main difference between the two models appear close to boundaries where the Dirichlet condition is prescribed (through eqs. (3) and ( 16)). In the hottest region, the CS-CM slightly overestimates the temperature decrease; the maximum error of the continuous model is about 1 % only.

The same comparison has been performed between the constinuous-scale coupled model and the Rosseland approximation model, see fig. 13. From this plot, it is seen that, both CS-RM and DS-CM roughly speaking follow the same trend, but with large differences. Overall, the temperature is over-estimatated, except in the colder region. Errors are about 5 %, where there were only about 1 % for the continuous-scale coupled model.

Computation times

Table 1 gives computation costs related to the three models (the discrete-scale couple model, the continuous scale coupled model, and the approximation Rosseland model), for the three numerical schemes (fixed-point method, improved fixed-point method, and Newton-Raphson), for four radiation-to-conduction regimes (for N=10, N=1, N=0.1 nd N=0.01). The computation time takes into account of spatial mesh reading, mesh partitioning, RTE solving, non-linear heat conduction solving, solution writing, etc. From this table, the following conclusions are drawn.

For any Stark number and for any linearization numerical scheme: the Rosseland model is much cheaper than the continuous-scale coupled model, and the continuous-scale coupled model is much cheaper than the discrete-scale coupled model, but this conclusion is to be balanced with the accuracy of the obtained results discussed above.

For any Stark number, and for any physical model, when convergence is reached, the Newton-Raphson iterative scheme is more efficient than other fixed-point type methods in terms of computation time. However, the difference time is not much because, for all treated cases, only two outer iterations were necessary (the RTE was solved twice for all cases); while the time difference comes from the number of times the linearized heat conduction equation is solved (see fig. 6), having in mind that the solution of a single iteration of the non-linear heat equation is very small. This table eventually shows that, for small Stark numbers (N<0.1), fixed-point type methods may diverge (computation times are set to infinity in Table 1) while the Newton-Raphson linearization scheme is always convergent.

Conclusion

The paper compared three distinct conducto-radiative models that can be used for modeling heat transfer at high temperature in view of being integrated within a topology optimization process. Finite elements are extensively used for all models, either using standard Galerkin for the Rosseland approximation as well as for other conduction equations, or in the vectorial form based on the stabilized version of finite elements for the radiative transfer equation. The non-linearity is handled using an iterative scheme. Experiments given here show that the Newton-Raphson strategy converge much faster than other fixed-point type methods.

The considered test case was a three-dimensional square lattice sandwiched between a hot and a cold plate. The model at the discrete level, at the pore scale, was the reference model. For small Stark numbers, when the 1: Comparison of the computation costs related to the three models (the discrete-scale couple model -DS-CM, the continuous scale coupled model -CS-CM, and the approximation Rosseland model -RM), for the three numerical schemes (fixed-point method -F-P, improved fixed-point method -I F-P, and Newton-Raphson -N-R), for four radiation-to-conduction regimes (for N=10, N=1, N=0.1 and N=0.01). Values of the Test A of Table 2 have been used, with 62 processes.

conduction-to-radiation parameter is small, i.e. when radiation predominates, it appears that the continuousscale approach, which has been set-up thanks to proper homogenization laws, gives satisfactory results so that, at most, the error found was about 5 % for a test case involving a single cell. For the case involving 3 × 3 × 3 cells, the maximum error was less, around 1 %. The Rosseland model, because it is very simple and easy to be integrated within an optimizer, has also been tested. This one, however, gives much larger errors: the Rosseland model over-estimates temperature almost everywhere, except in the colder region. The Rosseland model cannot be used confidently within a topology optimization process. The homogenized continuous-scale coupled problem, however, gives satisfactory results and can be used confidently in topology optimization problems.

In the future, some other more complex 3D architectures will be considered. Doing so, relationships of homogenized physical properties, as functions of the textural parameters taken at meso-scale will be consolidated. This will be highly fruitful information for the topology optimization of processes involving conducto-radiative transfers. The main application induced by the presented work is the optimization of architectured materials, and, in particular, lattice materials.

Characteristics

Figure 1 :

 1 Figure 1: Computational domains used in the discrete-scale approach. Left: the solid phase Ωs only; middle: the void phase Ωv only; right: immersion of Ωs into Ωv.

Figure 2 :

 2 Figure 2: Schematic representation of the discrete-scale geometry composed of crossed bars.

Figure 3 :Figure 4 :

 34 Figure 3: Comparison of our computational results with the results from Viskanta et al. [31] (Fig. 2 Page 68) for different conduction-to-radiation parameters N, (N=0.01, 0.1, 1, 10 ) with Θ 0 = T min /Tmax = 0.5 (Tmax = 2222.22 K, τ 0 = κL = 1.0 (κ = 2624.6719 m -1 ) and σ = 0.

τFigure 5 :

 5 Figure 5: Comparison of our computational results with the results from Viskanta et al. [39] (Fig. 2 Page 146). Effect of scattering albedo on temperature distribution. Conduction-radiation parameter N = 0.1 and τ 0 = 1.0, for Θ 0 = 0.5 and Θ 0 = 0.1.

)Figure 6 :

 6 Figure6: Convergence plots for the three schemes: "ordinary" fixed-point method, "improved" fixed-point method, and Newton-Raphson. Top-left: continuous-scale coupled model CS-CM on a scattering medium (ω = 0.5, N = 0.1 and Θ 0 = 0.1). Top-right: continuous-scale coupled model CS-CM on a non-scattering medium (ω = 0, N = 0.1 and Θ 0 = 0.5). Bottom-left: Rosseland model RM on a non-scattering medium (ω = 0, N = 0.1 and Θ 0 = 0.5). Bottom-right: discrete-scale coupled model DS-CM on a non-scattering medium (ωs = 0, N = 0.01 and Θ 0 = 0.5).

Figure 7 :N 2 • 10 - 2 τ

 72102 Figure7: Convergence history of the non-linear heat conduction computational scheme using the Newton-Raphson method for different values of Stark numbers N (from N=10, i.e. when radiation predominates, down to N=0.01, i.e. when conduction predominates) using the CS-CM model and with properties taken from Test A (see Table2).

Figure 8 :

 8 Figure 8: Effect of the anisotropy coefficient on temperature distributions for the continuous-scale coupled model CS-CM using different values of anisotropy coefficients and Stark numbers, with {σ = 1312, κ = 1312, θ 0 = 0.1, Tmax = 2222, and T min = 222}. Left: temperature profiles with N ∈ {0.01; 0.001} and g ∈ {0.9; 0; -0.9}. Right: temperature difference for the same values of N .

1 ; 1 0Figure 9 :

 119 Figure9: Rosseland model. Left: averaged adimensioned temperature T along the axis of propagation x. Right: difference ê(x) between these temperatures between the Rosseland model and the related continuous-scale coupled model. Top: test A (N = {10, 1, 0.1, 0.01}, ω = 0, Θ 0 = 0.5, τ 0 = 1). Middle: test B-1 (N = 0.1, ω = {0, 0.5, 1}, Θ 0 = {0.5, 0.1}, τ 0 = 1). Bottom: additional test (N = {10, 1, 0.1, 0.01}, ω = 0, Θ 0 = 0.5, τ 0 = 10); the brown curves in "line" and "dotted" are respectively the manufactured "exact" and approximate solutions with a non-nul source term g.

Figure 10 :

 10 Figure 10: Radiative density field visualized for the subdomain void part of Ωv and temperature field visualized for all the solid part Ωs, for one and three representative elementary volumes, with N = 0.01.

Figure 12 :Figure 13 :

 1213 Figure12: Results of the discrete-scale approach "DM" on the left and its comparison with the continuum-scale approach "CM" on the right (V al CM -V al DM ); in the porous medium with 3 VER representative elementary volume in each direction and with a porosity of 91 %.

Table 2 :

 2 Parameters and results of different errors obtained in the comparison of the models CM, RM and DM. Note that: "C", "R" and "D" mean respectively the continuous model, the Rosseland model and the discrete model.

1 :

 1 input parameters for the geometry; physical properties 2: load spatial and angular discretizations 3: define finite element spaces 4: k ← 0 5: solve -λ∆ T = 0 return T

	8:	l ← 0
	9:

6: T in ← T 7: for k ← 1 to k max do

Table

	CM	324 s 336 s 371 s 325 s 339 s 400 s	338	+∞	+∞	347	+∞	+∞
	CS-CM	69 s	70 s	70 s	70 s	77 s	78 s	74 s	+∞	+∞	77 s	+∞	+∞
	RM	9 s	-	10 s	10 s	-	13 s	11 s	-	17 s	11 s	-	+∞

Acknowledgements

The authors would like to thank the FreeFEM team for providing the very efficient finite element library [41]. which has been extensively used in this research. Values used in [31,39] have been converted into SI units using the online unit converter "https://converter.eu" and, in order to do the comparison, plots have been imported thanks to the open source library PlotDigitizer "http://plotdigitizer.sourceforge.net" [42]. Salih Ouchtout thanks the FEDER and the I-SITE NExT Innovez program for funding its post-doctoral support.


E Hazemann

Y Zhang

K Inal

Ch A Gandin

M Long

N T Niane

V Maguin

M Leduc

V Jaquet

M Bellet

Anisotropic Mechanical Behaviour of Nickel Base Single Crystal Superalloy CMSX-4 for Prediction of Grain Recrystallization during Homogenization Heat-treatment Following Solidification

In the present work, the methodology to identify parameters of an elastic viscoplastic behaviour law is presented. Orthotropic elasticity and Hill48 plastic anisotropy for as-cast CMSX-4 are considered. In order to obtain the corresponding parameters by inverse analysis, an experimental database is built by tensile/relaxation tests at constant temperatures between 800°C and 1200°C on specimens oriented along <001>, <110> and <111> directions at different strain rates. Tests are performed with a resistive heating machine, using infrared thermography and digital image correlation. The obtained mechanical behaviour is very dependent on the orientations of the specimens at all tested temperatures, and especially for the deformation mechanisms. Recrystallization behaviour is also very relevant to the specimens' orientation.

Introduction

In jet engines, high pressure turbine blades undergo the most severe operating conditions. To face both extreme temperature and stress, single crystal nickel base superalloys are extensively used as they show higher resistance to creep and thermal fatigue compared to conventional microstructures. Their composition is characterized by the absence of grain boundaries strengthening elements. As a result, the presence of unexpected grain boundaries, caused by defects such as recrystallized grains (RGs), is seriously detrimental to the mechanical performances of the final part. Recrystallization is related to the accumulated plastic strain introduced during solidification and cooling of the molten alloy when casting a part. Being a thermally activated phenomenon, recrystallization may occur during the homogenization heat treatment that follows solidification. In order to predict accurately RGs occurrence in a turbine blade, a representative descriptionthrough finite element calculationof stresses and strains formed during the casting process is required. Therefore, an anisotropic elastic-viscoplastic behaviour law of the single crystal nickel-based superalloy in the as-cast state is necessary. Different works have studied the conditions of recrystallization occurrence with different thermomechanical process such as indentation, compression or traction at different temperatures [Zam2007, Cox2003, Pan2013, Li2015]. The effects of deformation temperature, level of plastic strain, annealing temperature have been established in these studies. However, stresses were studied at a constant temperature, which is not representative of the mechanical load undergone by a turbine blade during solidification. These studies are focused on the conditions of recrystallization occurrence, not aiming at the understanding of the elementary occurrence mechanism. The present paper gives a methodology to identify parameters of a behaviour law for as-cast CMSX-4 with elastic and plastic anisotropy considerations. Thermomechanical tensile tests with constant cooling rate and loading ramp are presented for three different orientations of specimens to reproduce solidification and cooling paths of CMSX-4 in a turbine blade. These specimens are subjected to subsequent standard homogenization heat treatment to reveal recrystallization occurrence.

Materials and Methods

Tensile samples oriented along <001>, <110> and <111> directions

Second generation nickel-based superalloy CMSX-4 in the as-cast state was used for this investigation. Its composition is given in table 1. Several plates of 50 x10 x145 mm 3 and 80 x10 x 165 mm 3 were cast in Bridgman type furnace using respectively grain selector and single crystal seeds. Casting process followed industrial procedure at Safran Advanced Turbine Airfoils. Tensile specimens were machined by Electro Discharge Machining (EDM) within plates of 50 x10 x145 mm 3 for samples with tensile direction aligned with <100> crystallographic orientation. Samples with tensile direction aligned with <110> or <111> crystallographic orientation were obtained within plates of 80 x10 x 165 mm 3 by EDM. Samples are cut off by 1 mm from each surface to avoid surface defects. The exact orientation of tensile specimens was measured by X-Ray Diffraction (XRD). Only specimens with less than 10° of misorientation were used in this study. Standard electrochemical etching technique was performed in Safran Advanced Turbine Airfoils to ensure that neither stray grains, multi-crystals nor freckles were present in gauge part of tensile specimen. Additional etching using 5 % HF, 27.6 % HNO3, 67.4 % H2O for 60 minutes removed residual surface affected by EDM. As-cast microstructure cartographies of one specimen oriented along the <001> direction were obtained by Scanning Electron Microscopy (SEM) with Zeiss Supra 40 apparatus [Zei2022]. Local composition was measured on two specimens oriented along <001> by Energy Dispersive Spectroscopy (EDS) at different heights to account for macrosegregation with Brücker system and Maya3 Tescan apparatus [Bru2022, Pol2022]. For SEM imaging and EDS analyses, samples were cut on a XY plane (figure 1) and mechanically polished to a final step of colloidal silica. A high voltage of 20 kV was used for the imaging.

Table 1. Nominal composition of CMSX-4.

Element Al Co Cr Hf Mo Ni

Re Ta Ti W wt% 5.6 9.7 6.4 0.1 0.6 Balance 3.0 6.5 1.0 6.4

Tensile tests at high temperature for behaviour law identification

Tests are conducted in tensile machine DEDIMET developed at CEMEF [Pig2018, Gao2022]. A schematic view of the sample in the machine is presented in figure 1. The sample is heated by Joule effect. All the tests are carried out under primary vacuum conditions to limit oxide formation on sample surface. Tests are performed at a constant nominal temperature between 1000 and 1200 °C. Setpoint temperature is monitored by a type K thermocouple welded at mid hight of the gauge length. A bichromatic pyrometer and an Infra-Red (IR) camera measure maximum temperature and temperature field, respectively, within the gauge length. The upper jaw is fixed, setpoint velocity is imposed on the lower jaw. VIC-3D is used to analyse images to obtain displacement and deformation of the speckles in the gauge length with respect to a reference image [Sut2009]. For each sample, 3 cycles of tensile/relaxation phases are performed with respectively 0.01, 0.003 and 0.008 mm. s -1 set point velocities for the three tensile cycles. Duration of each tensile phase depends on set point temperature and is adapted to introduce up to 20% of deformation within the gauge length at the end of the test. Durations were calculated using preliminary tensile tests with one tensile/relaxation phase. Heating rate and cooling rate were 6 °C. s -1 and -7.5 °C. s -1 , respectively. A temperature stabilisation step of 100 s precedes the first tensile/relaxation cycle.

Tensile tests at high temperature, constant cooling rate and constant load rate for recrystallization criterion identification

Additional thermomechanical tests with constant cooling rate and load ramp are performed in DEDIMET to study recrystallization occurrence. This type of test is presented in figure 2. The tests parameters are gathered in table 2. Three different tests with various final loading 𝐹 𝑒𝑛𝑑 were carried out on specimens oriented along the <100>, <110> and <111> directions. Set point temperature was monitored with a bichromatic pyrometer.

Homogenization heat treatment and subsequent etching to highlight recrystallization occurrence

Specimens used for non-isothermal tests were subjected to the standard homogenization heat treatment in a primary vacuum heat treatment furnace Carbolite STF16-180 [Dat2022]. Specimens were attacked for 2 hours with Kalling etchant to highlight RGs when recrystallization occurred [Zam2007].

Mathematical model

Behaviour law

Only main equations are presented in this work. More details can be found elsewhere [Gao2022]. The mathematical model is implemented in a library developed at CEMEF (Cimlib).

For small strains, one can assume an additive decomposition for the strain rate: 𝜺̇= 𝜺̇𝒆 𝒍 + 𝜺̇𝒗 𝒑 + 𝜺̇𝒕 𝒉 (1) With 𝜺̇𝒆 𝒍 , 𝜺̇𝒗 𝒑 and 𝜺̇𝒕 𝒉 respectively the elastic, viscoplastic and thermal contributions to the strain rate tensor. The elastic strain rate is adapted for orthotropic properties. To consider anisotropic plasticity, Hill48 criterion adapted for orthotropic material is chosen [Hill1948]. The developed expression for the behaviour law is as follows:

σ ̅ = (σ Y0,T ref + Q T ref * (1 -exp(-bε ̅))) * exp [ Q 1 R ( 1 T - 1 T ref )] + k T ref ε ̅ ̇m * exp [ Q 2 R ( 1 T - 1 T ref )] (2)
With Q 1 and Q 2 activation energies and T ref the reference temperature.

Boundary conditions of the model

Only the gauge length is considered in the model. Temperature profile is directly imposed and velocity along tensile direction Z is considered as the boundary conditions during the simulation. They both are adapted from experimental measurements. More detailed methodology is presented in [Gao2022].

Identification strategy

An inverse analysis software developed at CEMEF is used to identify the different parameters of the behaviour law. Optimization strategy is similar to the one presented in [Gao2022]. Cost function expression is as follows:

𝐹 𝑐 = √ ∑ (𝐹

𝑇 ̇𝑖𝑛𝑖𝑡

Considering the high total number of parameters used in the model, elastic properties were chosen as inputs, which are tabulated with temperature. Values are adapted from [Her1996]. First step of identification strategy was to identify all the parameters for the orientation <001>. After identifying all the parameters for <001> orientation, grouped optimizations are performed for one temperature and three different orientations. All the parameters are identified in this step, including the plastic anisotropy factor.

Results and Discussion

Model development

The results of grouped optimization for specimens oriented along <001> direction at 1021°C, 1116°C and 1201°C are presented in figure 3. The 3 tensile/relaxation cycles are indicated by the vertical dashed lines for the specimen deformed at 1021°C. One can see very good agreement between simulated and experimental loads is obtained. Even though elastic properties are given as inputs of the model, simulated elastic part agrees well with experimental data. One can conclude that the model can predict accurately the behaviour of CMSX-4 between 1021°C and 1201°C for the orientation <001>. Grouped optimization for <001> orientation and more temperatures are ongoing. Model inputs for anisotropy consideration are Euler angles. Results obtained for elastic part between <001> and <111> directions at 1000°C have been compared with analytical calculations to validate the model. Grouped optimizations for one temperature and three orientations are ongoing to identify plastic anisotropy factor.

Recrystallization occurrence

Optical micrographs of specimens oriented along <111> direction after non-isothermal load ramps to 7000 N and 9600 N, homogenization heat treatment and attacked by Kalling etchant are presented in figures 4 and 5, respectively. No recrystallization is visible on gauge length for the specimen loaded up to 7000 N, while numerous RGs are visible on almost all gauge length for specimen loaded up to 9600 N. Additionally, one can observe strong necking for the latter, indicating large plastic strain compared to the sample loaded to 7000 N. Interestingly, most of the recrystallized grain boundaries appear to follow the direction of the primary dendritic arms. Recrystallized grain size is quite large, indicating that only few orientations were favoured to grow. The results for recrystallization occurrence of three different orientations and two different load ramps are gathered in table 3. The absence of recrystallization is indicated with a -mark. Recrystallization occurrence is evaluated qualitatively between + and +++ marks, with + the smallest recrystallized zone and +++ the largest. For ramp to 9600 N, one can conclude that specimen orientation has significant effect on recrystallization occurrence: recrystallized area is larger for specimen oriented along <111> direction than <100> direction, while no recrystallization is visible for specimen oriented along <110> direction. A recrystallization criterion can be identified in term of load ramp at constant cooling rate between 7000 and 9600 N for <111> oriented specimens. For <110> orientation, a load ramp up to a higher value than 9600 N is necessary. For <100> orientation, two more tests are required up to 5000 N and 7000 N to identify a recrystallization criterion. -+++ Following steps include the simulation of non-isothermal load ramps for the different orientations with the previously presented behaviour law to identify a recrystallization criterion among all simulated quantities. These simulations are ongoing to choose the quantity that fits the best recrystallization behaviour with orientation consideration.

Conclusions

A methodology to identify parameters of anisotropic elastic viscoplastic behaviour law for as-cast CMSX-4 has been presented. Results show that parameters identified for <100> orientation and temperature between 1021°C and 1201°C agrees well with experimental results. Developments to include more temperatures for parameters identification are ongoing. Plastic and elastic anisotropies are considered in the model, identification of plastic anisotropy factor for different temperatures is ongoing. Elastic anisotropy has been validated with analytic calculations. Orientation appears to have a strong effect on recrystallization behaviour for specimens tested on nonisothermal load ramps with subsequent homogenization heat treatment. Ongoing works include direct simulation on thermomechanical paths to identify a recrystallization criterion. Additionally, several experimental tests are still required to identify a restrained window of recrystallization occurrence for <001> and <110> oriented specimens.

  The machine is monitored by a Labview code [Lab2022]. During heating and cooling, a zeroset point load is fixed to prevent stresses caused by thermal dilatation/contraction. With Joule effect heating, strain is localized due to high vertical temperature gradient. Strain and displacement are tracked locally within gauge length with Digital Image Correlation (DIC). Temperature resistant speckles made of Aremco Pyro-Paint 634-ALP are applied on samples surface [Are2022]. VIC-Snap is used for images acquisition [Sut2009]. Two cameras Point Grey Grasshopper2 are used for the tests in the same vertical axis with an angle inferior to 30° between them [Fli2022]. They are combined with Schneider-kreuznach xenoplan 2.0/280901 15218330 optics [Sch2022]. Acquisition frequency varies between 5 and 10 Hz.

Figure 1 .

 1 Figure 1. Schematic view of tensile testing machine DEDIMET.

Figure 3 .

 3 Figure 3. Grouped optimization for <001> oriented specimens at 1021, 1116 and 1201 °C.

Figure 4 .

 4 Figure 4. <111> specimen after non-isothermal loading up to 7000 N, standard homogenization heat treatment and subsequent etching.

Figure 4 .

 4 Figure 4. <111> specimen after non-isothermal loading up to 9600 N, standard homogenization heat treatment and subsequent etching.

Table 2 .

 2 Testing parameters of non-isothermal tests.

	𝑇 𝑖𝑛𝑖𝑡	𝑇 ̇𝑖𝑛𝑖𝑡	𝑇 ̇𝑡𝑒𝑛𝑠𝑖𝑙𝑒	𝑇 ̇𝑒𝑛𝑑	𝐹 𝑖𝑛𝑖𝑡	𝐹 𝑒𝑛𝑑	𝑡 𝑡𝑒𝑛𝑠𝑖𝑙𝑒	𝑡 𝑡𝑖𝑛𝑖𝑡.	𝑡 𝑠𝑡𝑎𝑏.
	(°C)	(°C.s -1 )	(°C.min -1 )	(°C.s -1 )	(N)	(N)	(s)	(s)	(s)
	1200	6	-10.6	-7.5 1000 7000	2300 500 60
						9600

Z X Y Figure 2. Schematic representation of non-isothermal tests.

Table 3 .

 3 Recrystallization results for non-isothermal tests of three different orientations.

	𝐹 𝑒𝑛𝑑 (N) 7000 9600
	<001>		++
	<110>	-	-
	<111>

Acknowledgments

The authors are thankful to all the CEMEF staff that helped for preparation of specimens, microscopic observations and mechanical discussion: Cyrille Collin, Gilbert Fiorucci, Christophe Pradille, Salim Kraria, Guillaume Corvec, Arnaud Pignolet and Francis Fournier.