S Nakouzi

J Pancrace

F M Schmidt

Y Le Maoult

F Berthet

Infrared curing simulations of liquid composites molding

Keywords: Composite, Infrared curing, Curing optimization, Radiation, Glass/Epoxy, Carbon/Epoxy. Pacs:

Infrared radiation is an effective energy source to cure thermosetting polymers. Its usage is expected to reduce curing time in comparison with thermal heating and mold thermally regulated. In addition, because of the polymerization mechanism and instant on-off control of this power, an improvement in the final properties of the material is also expected. In this paper, we studied the infrared interaction with carbon (or glass) fibers reinforced epoxy matrix, where Liquid resin infusion (LRI) is used to manufacture the composite. Temperature of the composite is a key parameter that affects its mechanical properties and is controlled by the infrared emitters and the exothermic heat released from the polymerization. Radiative heat flux is computed using the in-lab developed software RAYHEAT. Then, the heat flux (or absorbed energy for glass fibers) is exported to the finite element based program COMSOLMULTIPHYSICS where heat balance equation is solved. This equation is coupled with the exothermic heat released during the curing process in order to predict the composite temperature versus time and degree of cure. Numerical simulations will be performed on planar parts (sheet shape) as well as curvilinear shapes. Experimental validations of the infrared curing carbon (glass)epoxy composite system are presented in this paper Sheet surface temperature distribution are measured thanks to infrared camera. Kinetic parameters were estimated from differential scanning calorimeter (DSC) experimental data.

INTRODUCTION

A novel process of curing composites is the use of infrared heaters. Carbon reinforcement epoxy matrix is cured in an infrared oven (figure 1) composed of halogen lamps. This manufacturing process is flexible: Emitter's power which affects temperature distribution in the composite can be easily controlled. Different researchers 1 used IR radiation for curing thermosets, others 5 used it for the preheating process. In this study, Liquid Resin Infusion process is used to infuse the liquid resin (RTM6) through the reinforcement (figure 2).

IR OVEN

Composite is exposed to radiation in the infrared oven. This one is composed of nine halogen lamps provided by TOSHIBA LIGHTING company. Short wave halogen lamps (in the range 0.6 -6µm and is the wavelength at the maximum lamp's relative intensity) made by TOSHIBA LIGHTING company were employed for perform the curing process. The lamps have variable power levels and a nominal power of 1000W. Tungsten filament is the infrared emitter. It is surrounded by quartz glass and a kind of halogen gas-filled lamp.

Air aspiration

The infrared radiation is absorbed by the carbon ply present on the top surface of the composite, in addition this surface is a convective boundary surrounded by the air present in the oven T air (t). Thermal insulation (Kerlane) is in contact with the frame of the composite, the edge effect becomes negligibly small, that is we can ignore the heat exchanges all arround its contour. The bottom of the composite is in contact with a 12 mm thick PTFE block, this material was taken into account in the numerical simulations. we chose this in order to limit heat exanges with the aluminium mold and study the IR interactions. Bottom surface of the PTFE is in contact with the aluminium mold where its experimental temperature was introduced as a boundary conditions in the finite element software COMSOL. During curing, Incident radiation is absorbed by the top surface of the carbon/epoxy composite. Then the composite is thermally heated in the thickness via conduction from its top surface. Figure 4 presents the IR density distribution as well as the boundary conditions.

RAY TRACING METHOD

Thermal radiant heat flux upon the top surface of the composite is computed in MATLAB using the in-lab developed algorithms 7. Theses algorithms are based on ray tracing method. In ray tracing, we model radiation as rays of light. The propagation of the rays is only based on the laws of geometrical optics. The general idea consists to simulate the interactions between the ray of light emitted by a source (in our case halogen lamps) and every object present in the scene (Infrared oven): multiple lamps, reflectors, composite. In the ray tracing software, coiled tungsten is the infrared emitter. It is assumed to be a Lambertian grey body and it is modeled as a cylinder with equivalent diameter 8. That is, we define a number N of rays emitted from the lamp; this number is chosen in a way to have a good agreement between the calculations costs and the needed precision of the results. Cosson and al. [7] calculated the convergence of the radiative heat flux to the analytical solution [START_REF] Leuenberger | ASME Annual Meeting[END_REF], as a function of the number N of the rays. Then we start from a random point P on the infrared emitter and trace the ray emitted in a random direction. This direction is defined by two parameters and . These parameters are defined by stochastic variables:

where and are independent uniform stochastic variables 12 in the range . Each ray emitted by the source is defined by its origin, direction, and spectral intensity. Then the ray is tested against all objects in the scene to determine what it hits and if it intersects the composite (fig. 5). The different interactions between the ray and the objects present in the scene change the ray properties. In our ray tracing software, assumptions found in the literature are made for the different optical properties of lamps, reflectors, recovery surface and the carbon/epoxy system. Ray tracing method was validated in an anterior work in the laboratory with an analytical solution given by view factor 13.

RADIATIVE PROPERTIES

A well known of the IR heat flux absorption on the top surface necessitate determination of the radiative properties for the sequence of plies: recovery surface film and de-molding fabric. Using an IR spectrometer, we measured radiative properties and determined the fraction of IR rays absorbed by the composite. Recovery surface transmit 85% of the infrared radiation to the demolding fabric (figure 6), than this one has high absorption spectral bands; 0.04% of the IR rays is transmitted to the composite (figure 7). In a second hand, spectral measurements in transmission for 1mm thick carbon fiber reinforcement sheet show that carbon reinforcement is optically opaque to the IR radiation (figure 8). These results were also found in the literature. To resume our spectral measurement, surface radiation transport occur on the top surface of the composite, our objective is to determine the ray fraction that is reflected when rays impact this surface. Spectral measurements in reflection for the sequence of the recovery surface film, the demolding fabric and the composite, compacted together, are plotted in figure 9. In the spectral range of the short wave halogen lamps, we assume a constant reflectance and we take the mean value in the spectral band 0.6-6 μm (R = 12%). Glass reinforcement as it is presented in the literature is optically transparent to the infrared radiation; this is not the case of our measurement in the spectral transmission band 1.3-25µm. As a first approach for the Glass reinforcement is opaque to the infrared radiation, in addition, the demolding fabric that is in contact with the fibers, transmit 0.04% of the infrared radiation to the reinforcement, that is, first approach that we can treat this case as a surface radiation transport.

MODEL DEFINITION

The composite measures 160 mm-by-140 mm and is 6 mm thick. The physics inside the composite14 is described by a heat balance (1) coupled with the exothermic heat (2) released during polymerization of epoxy matrix.

(1)

In the previous equation, denotes the thermal conductivity of the composite, equals the temperature, is the ultimate heat of reaction of the resin, is the degree of curing. The material property parameters required in this analysis are listed in table 1.

(2)

(3) Where is a diffusion factor, is diffusion rate constant with an Arrhenius type of dependence, is a constant having the order of one, the glass transition temperature, is a free volume fraction of the polymer, is a free volume fraction for a temperature equal to , is the thermal expansion coefficient. For an amorphous polymer which is the case of our composite, and .

Composite

Transverse thermal conductivity.

[18] Longitudinal thermal conductivity.

[18]

EXPERIMENTAL

Experimental measurement was used to validate simulations of the whole Infrared curing composites process.

FIGURE 11

Schematic presentation of thermocouples at section A-A of the composite.

FIGURE 10 Specific heat capacity of the carbon fiber reinforcement

In the experimental set up, carbon fibers placed in the mold were preheated by the infrared lamps present in the oven. Then the infrared heaters are turned off and the RTM6 resin at 80°C is infused through the reinforcement using the Liquid Resin Infusion process. In order to refer to the same operating conditions between experiments and simulations (we have not addressed in our simulations the filling step during which resin flows through the fiber reinforcement), infrared heaters are turned on after the infusion takes place and stabilization of the resin in the mold. IR curing the composite is flexible and controllable; we controlled the heat flux in order to optimize a curing cycle. In another word, we imposed a high heat flux to initiate and accelerate reaction in the epoxy resin, and then the heat flux is turned off, every time when we reached 195°C (figure 15) in temperature (because the degradation temperature of the matrix is 220°C). Figure 12 shows the heat flux value imposed on the top surface of the composite (at the center).

The heat transfer coefficient, h, is obtained using a semi-empirical expression (from COMSOL library) as a function of temperature and distance from the leading edge of the composite (figure 14). The bottom of the composite is in contact with a 12 mm thick Teflon block, which is in contact with the aluminium mold. Aluminium temperature was imposed as a boundary condition in the finite element software COMSOL MULTIPHYSICS.

Curing process is achieved with 47% of the nominal power of the infrared heaters. Because the degradation temperature of the matrix is 220°C, we turned off the lamps when we reached 195°C in temperature. Then the lamps were turned on when we decreased to 160°C in temperature. We introduced the experimental operating conditions as input data in our simulations: oven and composite geometry, air temperature, initial temperature of the composite, mold temperature.

SENSITIVITY ANALYSIS

Because of carbon fibers embedded in the epoxy matrix are very strong IR absorbers at all wavelengths oh the IR spectrum , heat flux on the top surface of the composite is absorbed or reflected (figure 5), IR interaction with the composite is a surface radiation transport. The ray fraction absorbed by the composite is: . Where is the mean reflectance, is the heat flux that impacts the top surface of the composite and is the heat flux absorbed by the composite. The heat flux absorbed by the composite is proportional to the reflectance of the top surface of the composite. Sensitivity analysis were performed in the time range 0<t[s]<358 (the heat flux is turned off at t=358s, figure 12) in order to determine temperature sensitivity to the inward heat flux received on the top surface of the composite. Sensitivity of the composite temperature to the heat flux received on its top surface is calculated using equation 7: (7) From figure 19, we can deduce that composite temperature is very sensitive to the inward heat flux (which is proportional to the reflectance of the composite top surface) and is monitored by the infrared heaters.

CONCLUSIONS

Until now, we simulated the whole curing process of the carbon/Epoxy system, and we validated experimentally theses simulations, these results were presented previously. Numerical simulations for fast curing the thermoset polymers and their experiment validation results are presented in this paper. We presented also radiative properties measurements for the carbon fiber and glass fiber reinforcement, for the recovery surface and for the de-molding fabric.

In our future work, we will simulate the curing process of composites with an industrial geometry; a specific mold will be used in the Liquid Resin Infusion process.

FIGURE 18

Heat flux sensitivity on composite temperature

FIGURE 1 FIGURE 3

 13 FIGURE 1 Schematic presentation of the Liquid Resin Infusion process

FIGURE 2

 2 FIGURE 2 Infrared oven

FIGURE 4 FIGURE 5

 45 FIGURE 4 Boundary conditions

FIGURE 6 FIGURE 7 FIGURE 9 FIGURE 8

 6798 FIGURE 6 Spectral transmission of the recovery surface

TABLE 1

 1 Thermo physical properties used in numerical simulations

	Property	Material	Expression

ACKNOWLEDGEMENT

Financial support for this work was provided by TOSHIBA LIGHTING (France) Company. The authors are grateful to Mr. D. Goyot for his support and his attention to this work.


Vignes Jérémie
email: jeremie.vignes@mines-albi.fr

Dusserre Gilles
email: cgilles.dusserre@mines-albi.fr

D E Almeida

Jean Dalmasso

Frédéric

Numerical simulations of clay tiles compression

Keywords: Clay, free compression test, rheological behaviours, friction law, inverse analyse, forming simulation

During the pressing step, the clay tiles undergo stresses which result in the appearance of defects. A rheological study, based on free compression tests, allowed to confirm the Elasto-viscoplastic behaviour of the clay. The different constitutive parameters were estimated by fitting the force-displacement experimental curves using the optimisation algorithm (ES Metamodel) implanted in the commercial software Forge 2009®. The influence of the tribological parameters was studied using squeezing numerical simulations of a full tile. The numerical model was validated with experimental squeezing test of technological specimen with a tile lug. Then, we have compared experimental force with the numerical one and deduced that the clay/tool interface is not perfectly sliding. A friction Tresca's law was used to model the clay/tool interface. Numerical results showed that the actual geometry of tile lug didn't allow to form correctly the tile. Several areas undergo tensile stress, air traps ,... A new geometry of tile lug was proposed in order to limit this phenomenon. Using a simplified defect criteria (Latham and Cockroft), the numerical model allowed to locate the areas where there is a risk of crack.

Introduction

During the shaping of a clay tile by pressing, the clay undergo stresses which can have for consequence the appearance of cracks, pulling out, and distortion. The manufacturing process of clay tiles can be divided into five steps: preparation of the raw material, extrusion, pressing,and finally drying and cooking. It is necessary to limit the appearance of defects in order to improve the yield of the process. In this paper, we will present a study which aims to determinate the triborheological behaviour of the clay used to produce tiles Numerical simulation is performed to understand the phenomenon leading to the appearance of defects. The study focuses on the understanding of phenomenon involved in the third step. There are different types of tests to identify the behaviour of the clay. There are the methods originate of the ground mechanical [START_REF] Zenzri | Encadrement du moment extrême dans l'essai au scissomètre[END_REF] and the fluid mechanical [START_REF] Mehmet | Simulation of lubricated squeezing flow of a Herschel-Bulkley fluid under constant force[END_REF]. Like explained in the thesis of E. Foudrinier [START_REF] Foudrinier | Etude numérique et expérimentale du procédé d'extrusion de pâtes argileuses[END_REF], it is essential to choose a test involving a state of stress representative of the actual geometry to determinate the behaviour of the clay. For example, a capillary rheometer was used to study the extrusion step. We decided to use a squeeze flow test [START_REF] Leider | Squeezing flow between parallel disks: I. Theoretical analysis[END_REF][START_REF] Covey | Use of the parallel-plate plastometer for the characterisation of viscous fluids with a yield stress[END_REF] on a clay cylindrical sample to find the rheological behaviour. This test can be easily modelled by a 2D axisymetric model. The analytical solutions given by Engmann [START_REF] Engmann | Squeeze flow theory and applications to rheometry: A review[END_REF] allowed to calculate the pressing strength according to the different parameters of the behaviour laws. The SE metamodel algorithm [START_REF] Emmerich | Meta-model Assisted Evolution Srategies[END_REF] available in the optimization module of the commercial software Forge 2009 ® allowed to identify the parameters of the constitutive law. The element deletion method [START_REF] Lee | Application of element deletion method for numericalanalyse of cracking[END_REF] implemented in the software Forge 2009 ® is based on a critical value of the Latham & Cockcroft [START_REF] Cockroft | Ductility and the Workability of Metals[END_REF] criteria. The simulation performed by using the identify behaviour and damage models.

Identification of the constitutive law

A squeeze flow test was used to determine the parameters of the constitutive law. This test involved a state of stress similar to those involving in the pressing step. The identification of the rheological parameters is complex because it requires to identify the friction law at the same time. To get around this problem, the pressing plates were modified in order to ensure a totally sticky contact. The tests were performed on a electro-mechanical testing machine, with a land cell 30 kN. To analyse a wide range of speed strain rates, we used several crosshead displacement rates. The clay samples were cylindrical. The dimension can't be given in this article for confidential reason. As it can be seen on the figure 1, the tested clay shows an Elasto-visco-plastic behaviour.The instantaneous answer of the material is linear and is supposed to follow the generalized Hooke's law (eq 1) until the flow stress is reached.

(1)

The clay follows then a visco-plastic behaviour. The evolution law of the flow stress is given by (eq 2).

(2)

An inverse analyse was carried out to obtain the different parameters of the constitutive law using the optimisation module of the commercial software Forge 2009 ® . The optimisation algorithm is based on a metamodel assisted evolution strategy. The rheological parameters identified were the Young's modulus, the consistency, the sensitivity with the strain rates and the sensibility with the strain hardening. The Poisson's ratio was fixed at 0,45, like it is commonly given in the literature. Asset of parameters was identified for each crosshead displacement rate. The mean value of each evaluation was finally calculated for each parameters in order to take into account the whole experimented strain rate range. For confidential reason the final parameters can't be given. Strain rate and hardening sensitivity parameters are low. However, if these parameters are neglected, the numerical strength / displacement curve diverges from the experimental curve. of the test is quite good. The quality of the model given by the optimisation algorithm can be controlled with the calculation of the quadratic gap between the two curves. For the model presented, the mean quadratic gap was about three percent. The Figure 2 shows the evolution of a cylinder clay sample during a numerical and experimental free compression test. In this case the friction between the sample and the pressing plates was greased. The numerical model allows to find the evolution of the convex of the sample.

Modelling of the damage during the free compression test

To model the damage observed during the pressing test and introduce into the numerical simulation an element deletion method, available in Forge 2009 ® , was used. This technique compares the value of a damage criterion calculated in each element to a limit value corresponding to the occurrence of the failure. If the criterion value is higher than this limit value, the element is automatically deleted. We chose the Latham and Cockroft criterion C defined by (eq 3) where is the plastic strain, are the principal stresses, and are respectively equivalent stress and strain.

.

For each test, the sample height Hd, where the material begin to damage has been measured. The critical value defined by the value of the criterion is reached in the simulation when the height of the sample is Hd. Whereas damage begin to occur on material heterogeneities, the model describes an heterogeneous criterion field due to finite element discretization. The result is then very sensitive to the mesh size, and a very fine mesh is necessary to provide a crack pattern similar to the experimental one. Moreover, the element deletion method is not well suited to model the crack propagation because of a bad evaluation of the stress concentration at the head of the crack. It results in a fracture pattern not so sharp as observed experimentally.

Determination of the friction parameter

It is quite difficult to identify the friction parameter from a free compression test. The analyse of the evolution of the sample profile could be used to determinate the friction parameter but this method is difficult to set up. We decided to use an inverse analyse on a new compression test. For this test, we pressed a small clay sheet in a prototype reproducing the strains involved in the forming of a real tile. The geometry of the tools allows to stamp a tile lug in the sample. The tools and the clay sample are presented in the figure 3.

The inverse analyse allowed to identify a Tresca parameter corresponding to an imperfectly sliding contact. The figure 4 show the strength / displacement curves obtained for the new test. We can observe that the numerical strength curve was positioned in the dispersion range of the experimental strength curve. The numerical model allow to reproduce different steps of the pressing test

Validation of the numerical model with a specific specimen

We have identified the different parameters of the rheological behaviours and damage laws (with the free compression test) and the friction law (with the prototype pressing test). Now, we can use the numerical model to reproduce the prototype pressing test by taking into account the damage. The model reproduces correctly the geometry of sample. This helps us to validate the different parameters determined. The Mises equivalent stress field (figure 3) helps to define the critical areas on the tile. The figure 5 illustrate that the sharp edge, in the lower face of the lug, will be difficult to form without damage. The damage begins in the same time on the numerical and the experimental sample. The crack appears on the top of the lug. The thinning areas are favourable to tears. Nevertheless, the numerical model does not reproduce perfectly the actual behaviour of the material. The deletion element method creates a lost of material during the pressing. The strength obtained with the numerical model, when the damage is taken into account, is lower than the actual strength. Furthermore, it is necessary to mesh the part with very fine elements. This is why we have used an automatic remeshing during the pressing. When the criterion of Latham and Cockroft, on an element, was near the critic value, this element was divided into ten small elements. This allowed to limit the lost of volume. The major problem of this method is the dramatic increase of the computation time.

Conclusions and perspectives

This study allowed to develop a methodology to identify the different parameters of the constitutive law of the clay. Estimating a single friction parameter from a complex test allows to obtain a mean value taking several contact conditions involved in the actual process. Moreover, it allows to discriminate visually and classify the different clay available according to their capacity to result to a plastic forming. The influence of the tribo-rheological parameters was put in evidence with the numerical simulation of a clay tile forming. The friction law plays a lot in the forming quality. When the clay slide in the tool, the tile is correctly formed but when there are some adhesion areas, the clay is cut through. This leads to the creation of damage areas. The numerical simulation will be a new solution to improve the actual geometry of the tiles. The commercial software Forge 2009 ® can supply a map representive of the state of internal stress. This will allow to develop new tools and new geometries to increase the yield of the process. Indeed, the analyse of the numerical stress fields helped to divided by two the internal stress at the end of the pressing step.

Figure 1 :

 1 Figure 1: Experimental and numerical strength / displacement curves of a free compression test on a clay sample

Figure 2 :

 2 Figure 2: Experimental and numerical compression test on a clay sample

Figure 3 :

 3 Figure 3: Pressing tools and sample used to identify the friction law

Figure 4 :

 4 Figure 4: Experimental and numerical strength / displacement curves of a free compression test on a clay sample

Figure 5 :

 5 Figure 5: Mises stress field on the clay sample at the end of the test

Acknowledgements

The present study supported by Terreal company. The authors wish to acknowledge this industrial partner for this support and his help.


Myriam Dauphin
email: myriam.dauphin@valeo.com

Stéphane Albin

Christine Roucoules

Yannick Le Maoult

Fabrice Schmidt²

Mouna El Hafi

Modélisation des échanges radiatifs dans un système d'éclairage automobile par la méthode de Monte-Carlo

Keywords:

La géométrie complexe des optiques des systèmes d'éclairages automobiles peut générer des points de focalisation du rayonnement sur les parties plastiques du système. Cet article porte sur la modélisation des transferts radiatifs entre une lampe à incandescence et un réflecteur dont le revêtement est en aluminium. L'objectif a été d'évaluer le modèle de lampe adapté pour le calcul de la température d'une plaque opaque, irradiée par un système « lampe-réflecteur ». Pour cela, nous avons comparé trois modèles de lampes : une source ponctuelle, une source cylindrique et une source issue de mesures goniométriques. Le calcul de rayonnement a été mis oeuvre au moyen d'une méthode de Monte Carlo associée à des techniques de suivi de rayons. On a démontré que la manière de modéliser la lampe a un impact sur la température de la plaque.

Introduction

Les simulations thermiques effectuées au sein de Valeo Lighting Systems permettent de prédire les températures des matériaux des systèmes d'éclairages automobiles en vue d'éviter leur dégradation. Des simulations utilisent la méthode des ordonnées discrètes et permettent d'estimer les points chauds avec une précision de +-5°C sur les températures [START_REF] Chenevier | Thermal simulation in lighting systems -5 days/5degrees[END_REF]. Cependant, cette méthode ne permet pas de traiter précisément les réflexions spéculaires et n'est pas adaptée aux phénomènes de focalisation du rayonnement rencontrés dans certains cas. C'est pourquoi, des méthodes de « raytracing » sont appliquées aux cas automobiles ( [START_REF] Moore | Thermal Analysis of automotive lamps using the ADINA-F coupled specular radiation and natural convection model[END_REF], [START_REF] Filipuzzi | Thermal Behavior of Bulb and LED based Automotive Rear Lamps[END_REF]). Dans l'article de Filippuzzi [START_REF] Filipuzzi | Thermal Behavior of Bulb and LED based Automotive Rear Lamps[END_REF], un cas simple d'une lampe irradiant une plaque valide la distribution de température sur le bulbe de la lampe 1 par le biais d'un couplage rayonnement-CFD du système lampe-plaque.

Dans ce travail, nous faisons le choix de découpler le calcul du rayonnement du calcul de CFD. Notre objectif est de valider le modèle d'émission d'une lampe à incandescence (figure 1) dans le cadre d'une méthode numérique de Monte-Carlo permettant de prédire les températures atteintes par une plaque opaque, irradiée par un système « lampe-réflecteur » (figure 2). La distribution de température, calculée en face arrière de la plaque, est comparée aux mesures thermographiques. Ce cas d'étude permet de valider indirectement le modèle d'émission de la lampe et la modélisation des surfaces réfléchissantes. Dans un premier temps, on s'intéressera aux propriétés radiatives de la lampe. Puis, nous détaillerons les simulations effectuées pour trois modèles de lampes et nous comparerons nos résultats avec des mesures par thermographie infrarouge.

Propriétés radiatives de la lampe

Le rayonnement thermique du filament, donné par la relation [START_REF] Chenevier | Thermal simulation in lighting systems -5 days/5degrees[END_REF], fait intervenir l'émissivité du tungstène et la température, Tf. L'émissivité, ελ(Tf), dépend de la température et varie fortement avec la longueur d'onde. Une précédente étude [START_REF] Dauphin | Towards thermal model of automotive lamps[END_REF] a permis d'interpoler les données de la littérature sur la bande spectrale [0,25μm; 28μm] [START_REF] Devignes | Propriétés des rayonnements optiques[END_REF][START_REF] Ln | Monochromatic emissivity of tungsten in the temperature range 1200-2600K and in the wavelength range 0.4-4μm[END_REF][START_REF] Fu | Physical Modeling of filament light sources[END_REF].

      f f f T L T T L       (1)
Des mesures en transmission du bulbe ont également été réalisées [START_REF] Dauphin | Towards thermal model of automotive lamps[END_REF] au moyen d'un spectromètre IRTF. La synthèse de ces données est présentée figure 4. Le bulbe transmet 93% du rayonnement émis par le filament, 93 . 0  b  , dans la bande [0,3µm ; 5,5 µm] (Figure 3). Au-delà de 5,5 µm, le bulbe est opaque affectant l'énergie partant du filament : 4).

        f f b f T L T T L         (Figure

Modélisation du système « lampe-réflecteur-plaque »

Nous utilisons la méthode numérique de Monte-Carlo comme une méthode statistique de calcul intégral pour estimer les distributions de densités de flux. Elle permet de tirer les rayons selon des lois de densités de probabilité et d'introduire aisément des effets radiatifs dans le processus de simulation [START_REF] Howell | The Monte Carlo method in radiative transfer », Hemisphere[END_REF]. Cette méthode, implémentée dans l'environnement, EDStar, développé par la communauté scientifique StarWest [9], utilise la bibliothèque de synthèse d'images PBRT (Physically Based Rendering Tool) [START_REF] Pharr | Physically Based Rendering: from theory to implementation[END_REF].

Modèle d'émission de la lampe

Pour cette étude, nous proposons de comparer trois modèles de sources radiatives de puissance équivalentes, 25W, dans le but d'identifier un modèle de lampe adapté et simple (Tableau 1). Les sources se distinguent par la distribution spatiale des rayons. La première source est ponctuelle. La deuxième source est modélisée par un filament cylindrique d'où sont tirés les rayons selon une pdf « Lambertienne » 2 . L'énergie de la lampe est uniformément répartie sur la totalité des rayons. En chaque point source, les vecteurs directions des rayons sont uniformément distribués dans l'hémisphère sortant selon la loi de Lambert. Enfin, la troisième source, « source goniométrique », est issue des mesures de l'intensité radiative de la lampe pour différents angles sur la bande spectrale [0,4; 1,1 µm]. En pratique, il s'agit d'un fichier de rayons, chaque rayon étant affecté de la même valeur unitaire. La polaire d'émission de la lampe est donc recréée par les variations de densités de rayons sur la surface du bulbe (Figure 5). Nous utilisons seulement les données portant sur les vecteurs directions. Ainsi, la source goniométrique est vue comme une source de recalage spatial pour notre modèle de lampe recréant fidèlement le trajet optique des rayons.

Modèle de réflexion du réflecteur

Notre algorithme, se basant sur les propriétés statistiques de la méthode de Monte-Carlo, consiste à effectuer un test sur une valeur seuil pour lequel un nombre aléatoire sur l'intervalle unité est tiré de manière uniforme. Ce nombre, selon qu'il soit inférieur ou supérieur à cette valeur seuil, détermine le futur du rayon (absorbé ou réfléchi). Un grand nombre de rayons permet de prendre en compte la fraction réelle d'énergie absorbée ou réfléchie par le réflecteur. On suppose une réflexion spéculaire idéale avec une valeur seuil constante ρ R =0,8 de sorte que 80% du rayonnement est réfléchi. L'algorithme permet ainsi de ne pas modifier l'énergie du rayon au cours des évènements qu'il subit avant d'impacter la plaque.

Validation du modèle par thermographie infrarouge

Banc expérimental et simulation thermique

Deux expériences permettent de valider le modèle radiatif d'un système lampe-réflecteur. Dans l'expérience 1, l'irradiation d'une plaque opaque de type ABS-N est faite avec lampe seule (Figure 6). En simulation, seules les sources ponctuelles et goniométrique sont testées pour ce cas. Dans l'expérience 2, l'irradiation de la plaque est faite avec le système lamperéflecteur (Figure 7). La distribution de densité de flux sur la plaque irradiée, calculée dans EDStar, est appliquée comme condition aux limites sur la face avant, afin de réaliser une simulation thermique 3D dans le logiciel commercial Comsol Multiphysics®. Le résultat à l'issu de chaque configuration est ensuite comparé aux mesures par caméra infrarouge-FLIR SC325 réalisées dans la même configuration d'irradiation. La somme des écarts quadratiques des profils thermiques issus de l'expérience 1 (Figure 8, Figure 9) s'élève à 3 à 4 % pour la source goniométrique, et 5 à 6 % pour le cas ponctuel. Par ailleurs, la courbe relative à la source goniométrique surestime le profil thermographique de 2 °C au maximum. A l'inverse, la courbe relative à la source ponctuelle sous-estime ce profil thermographique de 2 °C. Pour l'expérience 2 (Figure 10, Figure 11), les écarts quadratiques s'accentuent. Ils varient entre 10 et 14,5 % pour la source ponctuelle, 8 et 14,3 % pour la source goniométrique et la source cylindrique. Cette fois ci, la source ponctuelle surestime la thermographie tandis que la source cylindrique et la source goniométrique sous estiment la thermographie. Par rapport à l'expérience 1, il y a donc eu une inversion de la tendance liée à la présence du réflecteur. Le modèle de réflexion étant le même pour les trois sources, on peut d'ores et déjà dire que cette inversion est liée aux différences de distribution spatiale du rayonnement des sources. Dans l'expérience 2, on constate que la source cylindrique approche mieux le profil goniométrique. Cette dernière étant notre source de recalage spatial, nous pouvons conclure qu'il n'est pas nécessaire d'affiner le chemin optique des rayons et de prendre en compte les effets de réfraction du bulbe.

Résultats

Analyse des distributions de températures sur la plaque

Analyse des distributions de densités de flux sur la plaque

Le tableau 4 présente la répartition de la puissance radiative dans le système (expérience 2). Les puissances obtenues sur la plaque et le réflecteur sont comparables. La puissance absorbée par la plaque dans le cas ponctuel est même inférieure de 1% à celle obtenue avec la source goniométrique. Les écarts constatés sur les profils thermiques ne peuvent s'expliquer de manière globale, il faut s'intéresser à la distribution des densités de flux pour détecter des éventuels points de focalisation.

 

sur la plaque (2) afin de comparer les différences de répartition des de flux entre la source goniométrique et la source ponctuelle. L'expérience 1 (Figure12) montre que la source goniométrique est plus directive vers le centre de la plaque avec une irradiation plus importante d'environ 170 W/m² par rapport à la source ponctuelle. Figure 13, on observe l'effet inverse : la source ponctuelle irradie plus au centre de la plaque avec un écart de densité de flux significatif jusqu'à 600W/m² (zones gris médium) et très localement 1000 W/m² (zones noires). L'inversion de la tendance dans les profils thermiques est bien liée à un problème de dispersion des rayons.

Figure 1 :Figure 2 :

 12 Figure 1 : Lampe P21W-Philips -filament (L= 6mm, d= 0,75mm)

Figure 3 : 4 :

 34 Figure 3 : Transmission du bulbe Figure 4 : Influence du bulbe sur le rayonnement transmis

Figure 5 : 1 :

 51 Figure 5 : Source goniométrique : représentation 3D d'un échantillon de 1000 rayons-Données Valeo Lighting Systems

Figure 6 :Figure 7 :Tableau 3 :

 673 Figure 6 : Expérience 1 : lampe seule Distance filament-plaque = 7cm

Figure 8 :Figure 9 :Figure 10 :

 8910 Figure 8 : Expérience 1-Profil XX'

5 :

 5 Figure 12Figure13

Figure 12 :

 12 Figure 12 : Expérience 1 -

Bulbe de la lampe : enveloppe en verre protégeant le filament de l'oxydation (figure1).

Conclusions et perspectives

L'introduction du réflecteur dans le système radiatif montre que la source ponctuelle conduit à un résultat numérique surestimant la température de 15°C par rapport au résultat de la simulation réalisée avec la source goniométrique. Nous avons ainsi démontré que le modèle de source ponctuelle n'est pas adapté à notre système.

Le calcul fait avec le modèle de source cylindrique, quant à lui, se recale thermiquement sur celui de la source « goniométrique ». Cela confirme qu'il n'est pas nécessaire de prendre en compte les effets de réfraction du bulbe.

En revanche, les modèles retenus, goniométriques et cylindriques, induisent des résultats qui sous-estiment la mesure thermographique de 6°C au maximum de température. En perspective, nos efforts se focaliseront donc sur l'amélioration du modèle de réflexion pour approcher les mesures expérimentales.

Remerciements

Les auteurs souhaitent remercier Valeo Lighting Systems pour leur soutien.


Sawsane Nakouzi
email: sawsane.nakouzi@altran.com

Florentin Berthet
email: fberthet@mines-albi.fr

Yannick Le Maoult
email: yannick.lemaoult@mines-albi.fr

Fabrice Schmidt
email: fabrice.schmidt@mines-albi.fr

Simulations of an Infrared Composite Curing Process

Keywords: Curing composite, infrared oven, Radiation, Optimization, Epoxy resin, Carbon fibers

A 3D composite geometry (carbon fiber reinforced epoxy matrix) has been optimally cured using a well-designed infrared oven. The main purpose of this study is to validate a FE thermal model based upon commercial software COMSOL Multiphysics™ developed for the prediction of the infrared incident heat flux on the top surface of the composite during the curing process. Through the FE thermal model, an optimization study on the percentage power of each infrared heater was performed in order to optimize the incident IR heat flux uniformity on the composite. In order to achieve validations, an experimental set-up using an IR camera and some thin thermocouples has been developed in order to achieve temperatures measurement. The recorded temperatures at different positions on the 3D composite part will be compared to computed temperatures.

Introduction

Epoxy resins have several applications in the aerospace and automobile industry. Because of their good adhesive properties, superior mechanical, chemical and thermal properties, and resistance to fatigue and micro cracking, they produce high performance composites. In a second hand rising energy costs have also promoted an interest in using infrared energy to manufacture these composites [START_REF] Karkanas | Cure modelling and monitoring of epoxy/amine resin Systems II. Network formation and chemoviscosity modelling[END_REF][START_REF] Chern | On-line Processing of Unidirectional Fiber Composites Using Radiative Heating: I. Model[END_REF]. In the technology presented here, the composite is cured using an IR oven which includes halogen lamps (fig. 1.a). The liquid resin infusion (LRI) process is used to manufacture the composite, whereby liquid resin is infused through a fiber reinforcement previously laid up in a onesided mold (fig. 1.a). The 3D geometry used to validate numerical simulations is a simplified underbody car and its dimensions are sketched in Figure 1 The objective of this paper is to compare the temperature measurements of the 3D composite geometry realized using an infrared camera and thin thermocouples (fig. 1) versus numerical simulations.

Experimental set-up

In order to measure surface temperatures, we have used a microbolometer camera equipped with a lens's focal length mm 18 F  . Its resolution is 320 X 240 pixels. Its spectral range corresponds to the average infrared  

m 13 3 . 7  
. The image frequency is 60 Hz. In Fig. 2, the position of the camera is sketched as well as the corresponding measured surfaces   3 2 S , S .

Fig.2 -Positions of IR camera and measured surfaces

The infrared camera is located outside the oven. The axis of the camera ' 3 n  has been adjusted such as to be parallel to the normal 3 n  of the surface 3 S . For the temperature measured on surface 2 S , the angle  between the axis of camera and normal of surface has been adjusted to  60 . This value guarantees the accuracy of the measures and especially that the emissivity of such dielectric materials is not varying too much up to an angle of 60°. In order to ensure this assumption, spectrometric measurements of emissivity at different angles were performed but not addressed in the present paper. The apparent emissivity of the composite in the spectral range of the camera 95 . 0 13 7

  

has been calculated thanks to the apparent reflectivity 13 7 R  [3].

13 7 13 7 R 1      (1)
The apparent reflectivity has been deduced from spectral reflectivity measurements (fig. 3), performed thanks to an infrared spectrometer Bruker Vertex 70. It is to be noted that this spectrometer is equipped of two detectors allowing to cover all the IR spectral range   m 25 to 8 . 0

 . An important detail is that the resulting measurement includes all the different films used for resin infusion (recovery surface, demoulding fabric) as well as carbon fibers.

In addition, four K-type thin thermocouples of 80 m diameter have been placed at different positions inside the composite geometry (numbers 1 to 4, fig. 1.b), and especially where the IR camera has not the capability to achive measurements.

Numerical modeling

In previous papers [START_REF] Nakouzi | Optimization of the incident IR heat flux upon a 3D geometry composite part (Carbon/Epoxy)[END_REF][START_REF] Nakouzi | Simulations of an Infrared Composite Curing Process[END_REF], a FE thermal model based upon software COMSOL Multiphysics™ was fully detailed, for the prediction of the infrared incident heat flux on the top surface of the composite during heating. Although we have applied ray tracing method previously [START_REF] Nakouzi | Curing simulation of composites coupled with infrared heating[END_REF], it is not convenient to use it for an industrial application because of the computational expense when solving the radiation problems and time-consuming resulting from coupling ray tracing method with the optimization algorithms. The method employed by the FE heat transfer solver to compute the radiation heat transfer between opaque diffuse surfaces is based upon radiosity which is defined as the total radiative flux leaving a surface. The radiosity method subdivides a scene into surface patches that possess Lambertian optical characteristics and subsequently determine a condition in which the light transfer between the patches constitutes an energetic equilibrium [START_REF] Siegel | Thermal Radiation Heat Transfer[END_REF]. The fraction of radiative energy absorbed by the composite has been calculated via view factor computations thanks to the hemicube method [START_REF] Cohen | The hemicube: a radiosity approach for complex environments[END_REF]. In order to achieve numerical simulations, the radiative power of each halogen lamps is computed using the following relationship:

  t 4 f f f SB L S T T P    (2)
Where SB  is the Stefan-Boltzmann coefficient,

  f f T 
the integrated emissivity of tungsten filament, f T the temperature of tungsten filament and t S the equivalent surface of filament. It is noticeable that the emissivity is strongly dependent of the temperature of filament. This last parameter has been previously measured by comparing different technics such as color temperature, the resistivity method and a semi-analytical model based upon flux balance [START_REF] Monteix | Recent Issues In Preform Radiative Heating Modeling[END_REF].

The composite was meshed into 8248 tetrahedral elements. A very dense mesh is used for the volume material discretization along and near the edges highlighted in fig. 4.a, as higher temperature gradients will occur here, while a coarser mesh is used in the remaining volume. The resulting temperature distribution is represented in fig4.b. Due to optimization procedure [START_REF] Nakouzi | Optimization of the incident IR heat flux upon a 3D geometry composite part (Carbon/Epoxy)[END_REF], the maximum temperature difference throughout the composite geometry is less than 5°C. Again, the positions and associated numbers of thermocouples and surfaces are sketched in fig.

Discussion

In fig. 5, measured temperatures issued from both thermocouples and IR camera, and computed ones are plotted versus time. The three-dimensional numerical model is developed in the finite element software COMSOL Multiphysics™, where the heat balance equation is coupled with the cure kinetic model of the resin [START_REF] Ryan | Kinetics of Epoxy Cure: A Rapid Technique for Kinetic Parameter Estimation[END_REF]. The computed radiative heat flux is exported to COMSOL Multiphysics™, and imposed as a boundary condition on the top surface of the composite. This numerical model allows the calculation of the temperature distribution in the composite during curing as well as the degree of curing. From a global point of view, the agreement is unexpectedly good considering all the input data involved in the numerical model. At position 1, the measured temperature using thermocouple (TC1) and computed one are in good agreement, especially during temperature rise. However, a temperature difference of 5°C is calculated during the steady state. At position 2, the surface temperature measured by IR camera is slightly higher than those measured by thermocouple (TC2). Nevertheless, the measured and computed temperatures exhibit again a good agreement.

Fig.5 -Comparison between measured and computed temperatures

At position 3, again the best agreement between experimental and numerical temperatures occurs during temperature rise (transient state). A temperature difference of 8°C is observed during the steady state, both with IR camera and thermocouple measurement. Indeed, the proximity of this thermocouple of a vacuum port leads to underestimate the measured temperature. From these results, we can conclude that the measured temperatures on horizontal surfaces (S 2 and S 4 ) have the best agreement with computed temperatures during the steady state. In addition, the measured temperatures on inclined surfaces (S 1 and S 3 ) have the best agreement with computed temperatures during the transient state. These differences may be attributed to the fluctuation of air temperature inside the oven. Only one air temperature measurement has been recorded using a thermocouple located at the center of the oven (above the S 2 surface).The recorded air temperature is plotted versus time in fig. 6. 5), it is noticeable that slight fluctuations are observed around t=5000s. This fluctuation can be correlated with those observed on air temperature. This fluctuation is due to the variation of room temperature where curing process has been achieved (air flow, opening and closing of doors …). In addition, the root mean square error Err including the four thermocouples measurements is calculated thanks to the following equation and plotted in fig. 7:

    100 T T T Err 4 1 i 2 i exp 4 1 i 2 i num i exp        (3)
Except at the beginning and the end of the heating where the error is greater, it is noticeable that the average value is less than 4%. Again, this is a good result considering the great number of parameters involves in this kind of numerical simulations.

Future works

Different temperature measurements thanks to IR camera and thermocouples allow validating computed temperatures. In order to simulate the curing process, 3D radiation numerical simulations include also a kinetic model based upon modified Kamal-Sourour equation [START_REF] Nakouzi | Simulations of an Infrared Composite Curing Process[END_REF][START_REF] Ryan | Kinetics of Epoxy Cure: A Rapid Technique for Kinetic Parameter Estimation[END_REF]. For estimating the influence of IR radiation upon mechanical properties of the final composite piece, it would be important to achieve mechanical tests (DMA, 3-points bending, shear test …). For example, comparisons may be performed between samples cured by conduction mode (or autoclave), and in other hand by infrared oven. In addition, careful examinations of MEB images would lead to better understand the microstructure of composite piece obtained using IR curing process.

Fig. 1 -

 1 Infrared Composite Curing Process

Fig. 3 -

 3 Fig.3 -Spectral reflection  R on the top surface of composite (Recovery surface + demoulding fabric+ carbon fibers)

Fig. 4 -

 4 3D mesh of the composite geometry and computed temperature distribution

Fig. 6 -

 6 Fig.6 -Recorded air temperature versus time

Fig. 7 -

 7 Fig.7 -Root mean square error versus time