Implementation of a parametric procedure allowing efficient positioning of heat sources: application to high-temperature composites thermoforming process

. This work describes the implementation of a simple procedure that helps to easily position the heating elements in press plates used in high-temperature composites thermoforming process. The developed method permits to obtain desired temperature pro ﬁ les on the surface of the press plates through two main steps. The ﬁ rst step consists in ﬁ nding out an appropriate parametric curve that de ﬁ nes the spatial location of the heating sources into the thickness of the press heating plates. The second step uses an inverse method that combines a stochastic optimization algorithm in conjunction with ﬁ nite element simulations. This second step serves for the adjustment of the position curve parameters to obtain a simulated temperature pro ﬁ le as close as possible to the expected one at the press plates surface. This easy-to-implement approach is shown to be very effective to rapidly obtain a suitable location of the heat sources that minimizes energy consumption.


Introduction
Improving energy efficiency of industrial systems is of utmost importance. Indeed, several production processes, like injection molding and composite thermoforming, use heating elements in order to heat a polymer-based material so that it takes a desired shape during the cooling stage. Such industrial processes are known to be highly energyconsuming because of the operating temperature levels and their mass production nature [1,2]. On the other hand, these processes are often not optimally designed leading to huge energy losses during their in-service life. Reducing the energy consumption can be achieved through several routes. Indeed, characterizing correctly material properties variation during its processing is mandatory as it helps to predict correctly heat exchanges during the molding and thus design optimal parts. In their work, Likozar and Krajnc [3] studied heat transfer in several elastomers using a pilot scale mold. They showed that a correct description of temperature dependent elastomer thermophysical properties leads to good agreement between measured and predicted temperatures throughout a temperature range between 20°C and 200°C. Le et al. [4] showed that by using an embedded thermocouple in the molding cavity, it is possible to characterize crystallization kinetics piezo dependence of an injected polypropylene under processing conditions leading thus to a better prediction of the molded part cooling. In their contribution, Cheheb et al. [5] investigated the thermal conductivity of a rubber dependency with respect to its state of curing using a hot disk technique. They showed that the thermal conductivity may vary as a function of vulcanization rate and hence better predictions of reaction process are possible if this variation is taken into account. In addition to classical characterization methods, thermophysical properties identification using inverse method is shown to be an interesting alternative. In their study, Kim et al. [6] showed that it is possible to successfully achieve valid estimation of the temperature-dependent of thermal conductivity and volumetric heat capacity by combining the temperature measured on the boundary of a thermoplastic prepreg and its numerical counterpart obtained by finite element analysis (FEA).
Another critical aspect consists in the optimization of heating elements by finding the right position and the right energy amount to be delivered. The difficulty behind such optimization is that the considered solution has to be thought from an industrial perspective. Hence, it should not by no means compromise the production time, reproducibility and quality of the processed products. Indeed, some reported works show that it is possible for existing processes to minimize energy consumption and enhance the molded parts quality at the same time [7]. In their contribution, Castro and Lee [8] studied using a finite element (FE) model the thermal design of a representative mold used in sheet molding compound (SMC) compression molding. They showed that steady state cure time for a given SMC material depends on the following key parameters: part thickness, the heating channel temperature, the initial SMC temperature, the distance between heating channels and the distance from the heating tube center to the mold surface. These authors highlighted in their study the pertinence of a local one-dimensional heat transfer for thermal analysis design of molds for thin SMC parts as long as the heating channels are not too far apart. This work was extended by Abrams and Castro [9] who developed charts that can be used to estimate the steady state cure time for new parts and thus to estimate the manufacturing cost. The analysis of the temperature uniformity in cavity surface of injection molds for filled thermosets was studied by Barone and Caulk [10]. They combined a quasi-steady heat conduction model with an optimization algorithm to determine the heating line positions and operating temperatures that minimize the spatial variation of temperature in cavity surface. They showed that optimal designs use considerably fewer heating lines in comparison with conventional heating system.
The design of efficient heating elements that permit to rapidly reach setting temperature is another topic that has drawn attention of some research teams. Indeed, Abdalrahman et al. [11] investigated the effect of the length, the profiles and the arrangement of thermal channels on the efficiency of the heat system by comparing the heating time and the temperature homogeneity at the surface. They concluded that the sectional shape of the heating channel had little influence on the thermal efficiency of the heating plates and that the most important parameter is their spatial distribution in the plates: the channels parallel arrangement is the model that allows the best thermal efficiency. Jeng et al. [12] have showed that in order to achieve rapid mold temperature control, a steam heating system combined with water cooling was more efficient than a water heating/water cooling system and helps to improve the part quality. In their study, Xiao and Huang [13] studied how to shorten the required heating time for the mold cavity surface to reach the setting temperature and to minimize the heterogeneity of the temperature field on the cavity surface. To do so, they developed an optimization tool by combining particle swarm optimization (PSO) algorithm with FEA to design the heating system of a rapid thermal cycling mold. Wang et al. [14] studied the rapid heat cycle molding process with steams and electric heating. The heating and the cooling channels were combined into one single system. The electric heaters were placed into annular channels filled with water. During the heating, the water placed into the annular gaps would avoid the overheating and improve the heat transfer. They have proven, with numerical simulation and experiments, the efficiency of this heating method by studying the influence of the gap size, the power density of the electric heaters and their layout on the heat transfer to the cavity surface. The work reported by Lin and Chou [15] was concerned by the improvement of the cooling system of an injection mold by dividing the problem into two major steps: (i) the optimization of the cooling time based on the energy balance in the mold and (ii) the optimization of the cooling channel disposition in the mold and the geometric shape factor.
The majority of the reported works in literature of which the authors are aware are basically concerned with heating and cooling systems embedded into mold and not with separated heating plates. Indeed, the challenge in the case of the heating plates is that these latter should be designed to be versatile with and an adaptive heat delivery capabilities according to the mold they meant to heat. In this work, we focus on the design of press plates heating system used in composite parts curing. Unlike conventional processes, the target temperature needed for the thermoforming of these materials is very high as it ranges from 180°C to 400°C. The present study describes a procedure developed to efficiently locate the heating element in the thickness of the plates.
After the context of this study is presented, Section 2 describes a simplified geometry of the studied system with its components. The proposed heat sources positioning technique is detailed in Section 3. Sections 4 and 5 present, respectively, the FE model used and the discussion of the obtained results.
2 Description of the studied heating plates geometry The heating press system serves in general to heat a mold containing polymer-like material sandwiched between the press plates as depicted schematically in Figure 1. The heat is transferred first by conduction from the heating elements to the plates. These latter, in turn, transfer a part of the received heat by conduction to the mold through their common contact surface. In this study, we focus on a 2D model of the heating press plate equipped with identical equidistance sources (along the x-axis) (Fig. 1). For sake of simplicity, each single heating element is considered as a black box that delivers the needed amount of energy regardless the used technology. These heating sources may be either electrical resistances or channels wherein a high temperature fluid flows.
In order to design such a heating system that is able to deliver the needed heat by the mold, one should deal first with several engineering issues: (1) What is the ideal energy amount to be delivered to the mold? (2) What are the characteristics of the heating elements (position, number)? and (3) What is the dynamics that should be set to monitor this complex system both in the heating and cooling phases?
Trial and errors techniques based on either experiment, analytical or numerical calculations are the three routes that are used to bring answers to the above-mentioned engineering issues concerning the heating system design. In most cases, the optimization from an engineering point of view relies on the use of simple analytical solutions of the heat equation considered mainly in steady state. While this approach allows obtaining a solution very quickly, it is unfortunately not suitable for complex systems particularly under unsteady state conditions. On the other hand, powerful tools such as FEA are extensively used as an alternative to simple analytical methods. Beside the increasing Central Processing Unit (CPU) capacity of computers, using FEA is still time consuming especially in large parameters set optimization. Thus, one has to develop efficient design procedures that achieve balance between the predictive capabilities of FEA of complex system behavior and computation time.

Optimization method
Each individual heating source can be characterized by three main dimensional parameters, namely, its position coordinates in the plates thickness captured by the y i and x i coordinates and the diameter D i (Fig. 2). This latter parameter is generally correlated to the energy rate that can be generated by the heating source, be it a flowing hot fluid in a channel or an electrical resistor.
For simplicity, we consider in the following that the heating elements are identical delivering hence the same amount of energy. Also, the elements are considered equidistant in the x-axis direction, which means that the distance Dx = x i+1 À x i is constant for the studied configurations. The remaining unknown to be found consists in the vertical positions y i . With this choice, one still has n position parameters, corresponding to the heating sources number, to determine if the elements are considered individually. Obviously, this procedure is likely to be time consuming when used in conjunction of FE optimization loops. The proposed approach detailed below consists simply in defining a curve by m parameters along which the heat sources will be distributed. Hence, by taking m value smaller compared to n value one may drastically enhance the efficiency of heating element positioning step.

Selection of the parametric curve shape
A good starting point for defining the parametric curve is to have in mind that the desired temperature level is somewhat proportionally related to the in-depth position (distance from the plate surface) as well as to the energy rate generated by the heating element of course. This reasoning assumes obviously that heat is transferred mainly in one direction. As far as the studied problem is concerned, the heat transfer is close to a 1D problem thanks to the adjacency distribution of the sources. This adjacency condition can be obtained by maintaining moderate values for both distances Dx and Dy. It is worthy to mention that the distance Dx is inversely proportional to the heating elements number n. One should keep in mind that this assumption does not hold for the right and left edge neighboring sources. With all this said, a natural curve candidate is a parametric one defined by a homothetic transformation of the shape of the desired control temperature to be reached at the press plate surface. From this scaling, one defines the sources position y i along the ydirection knowing that x i was already fixed by the relation x i = Dx Â (i À 1) where i represents the heat source number. The distance d i introduced in Figure 2 can be easily calculated by the simple relation d i = H À y i where H denotes the plates thickness.

Step 1: Scale normalization of the control temperature curve
The idea here is to position the heat sources in such a way to replicate isotherms with a shape similar to the one of the desired control temperature. This shape can be obtained by  first applying scale normalization to the control temperature curve (Eq. (1)).
T con ðxÞ ¼ T con ðxÞ À T con min T con max À T con min ; ð1Þ where T con min and T con max represent respectively the minimal and maximal control temperature values.
The obtained curveT con ðxÞ has the same shape as T con (x) but evolves between 0 and 1. A graphical representation of this scaling is shown in Figure 3a and b.

3.3
Step 2: Definition of the position curve y(x, a, b) of the heat sources From the normalized curveT con ðxÞ, a new parametric curve that fits inside the press plates dimensions (L Â H) can be introduced by equation (2). yðx; a; bÞ ¼ a ÂT con ðxÞ þ b: Each parameters set a and b defines a unique positioning curve. The result obtained thanks to this second step is illustrated in Figure 3c. At this stage, one needs only to define the desired number of heat sources to be placed along the position curve y(x, a, b) making the model ready to be tested using FE (Fig. 3d and e).
By varying these parameters, one can sweep a wide range of positioning curves that are likely to generate isotherms of similar shape to the control temperature demonstrating the strength of this technique (Fig. 4).
From a practical point of view, a minimum safety thickness should be taken into account in the definition of the curve y(x, a, b) in a way that the sources do not touch the edges of the plate. This thickness can be determined based on material strength considerations such as maximum stress/strain constraints.

Optimization procedure of the position curve parameters
It is worthy to mention that the determination of T con (x) is based on processed part quality criteria and obeys itself to an optimization procedure that is beyond the scope of this work. Hence, we assume in the following that T con (x) is already known and the main goal is to find out the optimal parameters of the curve y(x, a, b) that minimize the mean square difference DT(a, b) between FE simulated temperature T sim (x) and the target control temperature at the surface of the press plate T con (x). This difference, known also as the objective function, is given by equation (3).
where n is the number of evaluation points that correspond to the chosen number of heat sources in the plate. The extreme points located near plate right and left edges are not taken into account in the calculation of DT(a, b) since they do not meet the 1D heat transfer assumption as already discussed. The overall optimization procedure for the determination of the best heat sources location is illustrated by the chart in Figure 5.
The loop starts by an initial guess parameters set that defines the first position curve. T sim (x) obtained by FE simulations is then compared to T con (x) and the difference is captured by DT(a, b). The stochastic algorithm "pattern search" implemented in Matlab software [6,7] drives the evolution of the a and b parameters set through an optimization loop until the minimum value of DT(a, b) is reached.

Modeling of heat exchange in the press plates
The heat transfer within the press plates is governed by the heat equation [16] with a source term: where r, C p and l are respectively the density, the heat capacity and the thermal conductivity of the press plates material. Q represents the energy rate generated by the  heat sources considered here as punctual sources located at distance y(x, a, b) from the bottom edge of the press plate.
The boundary conditions are shown schematically in Figure 6. On the bottom edge, an adiabatic condition is applied in order to take into account the isolating material that is generally used to reduce energy loses. On the side and top edges, the boundary conditions are considered to be convective. The convection coefficient h 1 is equal to 5 W·m À2 ·K À1 on the side faces. On the top edge where the mold is meant to be placed, a high value convection coefficient h 2 = 100 W·m À2 ·K À1 is used. The room temperature is set to 20°C.
In 2D, the press plate has a rectangular shape (0.15 m Â 0.5 m) made of a mold steel with the following typical thermophysical characteristics: r = 7874 kg·m À3 , C p = 475 J·kg À1 ·K À1 and l = 45 W·m À1 ·K À1 .
The FE software COMSOL © is used to solve the equations governing the heat transfer in the press plates for a heating period fixed to 300 s and a value of Q set to 28.5 kW·m À1 in order to reach desired temperature levels (200-450°C).
The selected mesh consists in several thousand of 6nodes triangular elements. The mesh is refined near the heat sources and the edges of the plate as shown in Figure 6.   Profile 1 : T con ðxÞ ¼ T con min þ ðT con max À T con min Þ⋅ Profile 2 : T con ðxÞ ¼ T con min þ ðT con max À T con min Þ Ã e À70⋅ xÀ L 4 ð Þ 2 þ ðT con max À T con min Þ Ã e À70⋅ xÀ 3L where L represents the length of the plate along the x-axis (Fig. 4). While T con min was kept constant, several values of T con max were tested (Tab. 1). This gives a total of five target control temperatures for each profile.
With only two parameters to optimize and a light weight mesh FE heat transfer model, the calculation is quite quick as it takes a few minutes only under an HP-Zbook PC equipped with an Intel i7-4800 QM, processor and a 24 GB of RAM. One can see from Figure 7 that 50 iterations are necessary to converge to the solution with the default setting of the "pattern search" Matlab routine. For practical reasons, we introduce hereã andb that simply represent normalized values of a and b parameters introduced in equation (2).
The starting guess fixed arbitrarily to (ã ¼ 0:1 and b ¼ 0:1) leads to huge DT ðã;bÞ difference between simulated and target control temperatures. Despite this unsuitable initial guess, the optimization routine is capable of finding out, more or less quickly,ã andb values that minimize the objective function DT ðã;bÞ.
The obtained results are given in Table 2 and the simulated and target control temperatures curves are depicted in Figure 8.
For readability concerns, only runs #1, #3 and #5 are plotted in Figure 8. One can see that a good agreement between the simulation and the target control temperature is obtained in the case of the linear profile with a maximum deviation ranging between 2.61°C and 3.40°C. However, in the case of the profile #2 (bi-bump shaped function) the mismatch starts to be noticeable from run #3 corresponding to T con max ≥ 350°C. Figures 8 and 9 show that this discrepancy is mainly due to the difference in the maximal and minimal temperature region between simulated and target control temperatures. In addition to this, a greater difference is also observed at the extreme points where the temperature difference can be easily greater than 5°C because of the boundary effect evidenced by Figure 10. This figure also shows that the used number of heating elements seems to be sufficient to respect the 1D heat transfer assumption.
The obtained location of the heat sources for the tested profiles (a) profile #1 and (b) profile #2 for the runs #1, #3 and #5 is depicted in Figure 11 where the thickness H of the plate is amplified by a factor of 2.8 for readability purpose.   These positions are coherent with the adopted strategy. Indeed, as the T con min is kept constant and T con max is increased, the shape parameterã variation is sufficient to minimize DT ðã;bÞ. This shape variation works very well especially for moderate T con max À T con min values as they meet the 1D heat transfer assumption.

Application to the reduction of energy consumption
After validation of the feasibility of the proposed heat sources positioning approach, it is interesting to analyze the influence of energy consumption reduction on the evolution of the temperature difference and the obtained sources position. To do so, several calculations are performed using different heat energy coefficient applied to run #4 of the profile #1 (Tab. 1) which is taken as a reference. Table 3 summarizes the obtained results. One can see that better matching between simulated and target control temperature is obtained until a coefficient of 0.7. The difference starts to increase bellow a coefficient of 0.6 to attain as much as 20.23°C. Figure 12 shows the obtained temperature profiles compared to the target one.
The position of the heat sources depicted in Figure 12 shows that the reduction of energy delivered by the heat sources leads to a vertical translation of the latter until they are too close to press plate surface.
These last results show clearly that an energy criterion can be directly incorporated into the optimization procedure permitting to find out the best location that meets the desired temperature profile while minimizing energy consumption at the same time.

Summary
In this study, a simple yet efficient method for heat sources positioning within press plates used for high temperature thermoforming process is presented. The position of the heating sources is defined by a parametric function originated by a two scaling steps of the desired target control temperature. The basic idea is to reproduce isotherms in the press plates thickness having similar shape as the one of the control temperature. Then, an optimization loop is used to drive the determination of optimal parameters set of the position parametric curve that defines the best location of the heat sources.
Two target quite different control temperature profiles were tested with five sets of maximum and minimum temperatures values each. The results show that the Fig. 9. Comparison between absolute temperature differences obtained for the tested profiles (a) profile #1 and (b) profile #2. The bars from left to right correspond to run #1, #3 and #5 respectively.  proposed technique is efficient provided the 1D heat transfer assumption is valid. In addition to this, it was shown that this technique can be used in order to reduce energy consumption. This work will be extended using a FE model with a more realistic heat source geometry that takes into account the whole system configuration (press heating plates with the mold containing a polymer like material). A validation with regard to experimental data is of great interest and will hence be also considered.