Inverse heat transfer analysis in a polymer melt flow within an extrusion die

Abstract

In forming processes like the extrusion of polymers, the temperature profile in the polymer melt flow through the die can be quite sharp, due to high viscous dissipation and very low heat conductivity. To predict accurately the temperature rise, the inlet temperature profile has to be taken into account. This profile is generally unknown because for creeping flow the temperature field is affected far downstream from the entrance of the die. This numerical study aims to restore the temperature field within the polymer from temperature measurements taken inside the die. The polymer flow is assumed to be an incompressible steady laminar flow of a Newtonian fluid. The numerical solution of the inverse problem is computed by using a classical conjugate gradient method. The analysis is important to decide the location of the thermocouples in an experimental die. The feasibility of restoring the temperature profile is illustrated for different thermal and flow conditions.

Keywords:

• Extrusion die
• Temperature profile
• Inverse method
• Conjugate gradient method

1. Introduction

The extrusion process is particularly useful when thermal and/or mechanical means are required to obtain a uniformly processed product in continuous operation. For plastic materials, it consists usually of the following sections: feeding, conveying and metering. In the metering section, the material is heated, pressured and homogenized [1]. The product is obtained by pushing the molten material through the extrusion die. In order to obtain the desired product quality and characteristics, knowledge of the inlet thermal state of the die is very important. So the die plays an important part in the extrusion. It is found that, because of this process and operating parameters, the thermodynamic properties of the product are strongly affected. In order to perform a numerical simulation of polymer flow and heat transfer in the extrusion dies, one has to specify inlet boundary conditions, which are generally unknown. Fortunately, for such a creeping flow the developing area for a velocity field is very short. The inlet velocity profile could be chosen as fully developed. On the contrary, the inlet temperature profile develops quite slowly and affects the temperature field far downstream. So, the problem is how to specify the inlet temperature profile. Hence the problem is referred to as a boundary inverse heat transfer problem: it belongs to the inverse heat conduction problems (IHCP), which were considered by Alifanov [2] and others. The IHCP have numerous important applications in various branches of engineering and science, including among others, estimation of unknown boundary heat flux [3,4], thermo-physical properties of materials [5,6], time-wise variation of the strength of energy sources located inside a medium [7,8] and estimation of initial conditions [9,10]. These problems are known to be ill-posed, in contrast to the direct heat conduction problems (DHCP), which are well-posed (i.e., the solution exists, it is unique and stable to small changes in the input data). A variety of numerical and analytical techniques have been proposed to approximate the solution of such IHCP [11,12]. The problem stated above is the problem of the initial temperature profile restoration. While the velocity field is known, one could expect the solution of the problem to be unique. There are several reports on this subject, for example: Liu, and Özisik [13] estimated the inlet temperature profile in the laminar duct flow. Subsequent investigations by Bokar and Özisik [14], Huang and Özisik [15], and Machado and Orlande [16] examined various aspects of this problem.

Recently, Hsu et al. [17] presented a two-dimensional inverse least square method to estimate both inlet temperature and wall heat flux of a steady laminar flow in a circular duct. Huang and Chen [18] have solved a non-stationary Navier–Stockes equation, to provide coefficients for the energy equation, but the velocity field itself does not depend on temperature. Among forced convection flows, the flow of a complex fluid in a narrow channel is obviously the case when coupling between the temperature field and the velocity field are essential, because its viscosity strongly depends on temperature. The inverse natural convection problem was also solved by Prud’homme and Nguyen [19]. Gejadze and Jarny [20] have presented a detailed analysis of IHCP coupled problems, when the velocity filed depends on the temperature field through viscosity and Navier–Stokes equations for a non-Newtonian fluid. Recently K.T. Nguyen and M. Prystay estimated the initial temperature profile and its evolution in polymer processing using surface temperature measurements and the conjugate gradient method (CGM) was employed to search for the minimum of the function [21].

In the present study, a two-dimensional conjugate heat transfer model has been developed for Newtonian materials being processed in an extrusion die. The aim of this work is to estimate the inlet temperature profile from the knowledge of temperature measurements inside the die. The governing equations for the fluid and wall regions are solved at the same time, by using a finite volume numerical method (Aquilon) [22]. An iterative numerical procedure is used to solve the direct, adjoint and sensitivity problems.

2. Direct problem equations

The steady laminar flow of an incompressible fluid is considered through an extrusion die ((---355----1) ), see figure 1. The velocity and the temperature fields are governed by the coupled equations of mass, momentum and energy, i.e., (---355--1) (---355--2) (---355--3) (---355--4) where, (---355--5) The following boundary conditions are considered for the velocity field:

a.	At the inlet of the channel: (---355--6-)
b.	At the outlet: (---355--6-)

At the internal surface, the usual no-slip boundary condition: (---355--6-) For the temperature field, the boundary conditions are taken as follows: at the entrance; (---355--7-) adiabatic conditions; (---355--7-) The outer surface temperature of the wall is fixed: (---355--7-)

Figure 1. Extrusion die.

The solution (T, u, w, p) of equations (1)(---355--1) –(7)(---355--7-) is denoted as the Direct Problem Solution. When the pressure drop is fixed, and the temperature profiles T₀(z) and T_p(x) are known, the Direct Problem Solution determines the temperature field within the extrusion die and the velocity field in the channel.

3. Inverse problem

3.1. Problem statement

For the inverse problem the inlet temperature profile T₀(z) is regarded as being unknown and is to be estimated. Additional data are required for the solution. They consist of temperature measurements (---355----2) given by Ns sensors located at appropriate positions (x_m, z_p), m = 1, … , Ns, p = 1, 2 inside the extrusion die. Such measurements may contain random errors, but all the other quantities appearing in the direct problem equations (1)(---355--1) –(7)(---355--7-) are considered to be known with a sufficient degree of accuracy.

The velocity field is established from the entrance of the channel, hence only the energy equation has to be inverted. The inverse problem is formulated by introducing the following least squares function: (---355--8) where T(x_m, z_p; T₀) are the temperatures computed from the Direct Problem Solution, at the measurement locations, by using an estimated inlet temperature profile T₀(z).

Minimization of the function (8)(---355--8) is achieved by constructing a sequence of approximations: (---355----3) according to the CGM such that; (---355--9) (---355--10) Where γⁿ is the descent length and pⁿ the descent direction, at iteration n.

Convergence of the above iterative procedure in minimizing quadratic functions like J(T₀) is guaranteed in [23]. To perform the iterations according to the equations (9)(---355--9) –(10)(---355--10) , the gradient of the function (---355----4) is needed.

3.2. Sensitivity problem

In order to develop the sensitivity problem equations, a variation εδT₀ of the inlet temperature profile is considered, the resulting temperature is denoted T⁺ = T(x, z; T₀ + εδT₀). The temperature sensitivity is then defined by: (---355--11)

Developing the direct problem equations (1)(---355--1) –(7)(---355--7-) for T(x, z; T₀) and T(x, z; T₀ + εδT₀) and then subtracting the resulting expressions, leads to (---355--12) (---355--13) (---355--14) (---355--15) (---355--16) (---355--17) where (---355----5) is the solution of the direct problem equations (1)(---355--1) –(7)(---355--7-) .

3.3. Adjoint problem

In a first step, the variation δJ(T₀) due to the variation δT₀ of the inlet temperature profile is developed according to the definition: (---355--18) From equation (8)(---355--8) , we have: (---355--19) and from equation (18)(---355--18) , it comes: (---355--20) where, (---355--21) By definition, the gradient (---355----6) satisfies: (---355--22)

In order to compute the gradient (---355----7) , the equation (20)(---355--20) has to be transformed under this new form. It is convenient to introduce an adjoint variable Ψ and the Lagrangian L(Ψ, T, T₀). (---355--23) and to take the variable Ψ solution of: (---355--24)

According to this stationary condition, the following adjoint problem equations are obtained: (---355--25) (---355--26) (---355--27) (---355--28) (---355--29) (---355--30)

Together with the gradient equation: (---355--31)

3.4. Conjugate gradient algorithm

The minimization algorithm by conjugate gradient (CG) involves the resolution of the three problems (direct, adjoint, and sensitivity) at each step of the iterative process. The overall CG algorithm may be summarized as follows:

1.	Choose an initial guess T0(z);
2.	Solve the direct problem given by equations (1)(---355--1) –(7)(---355--7-) , to obtain T(x, z; T0);
3.	Knowing the computed T(xm, zp; T0) and the measured temperatures (---355----8) , at the sensor location compute the least square function J(T0);
4.	Solve the adjoint problem;
5.	Compute the gradient (---355----9) ;
6.	Compute the direction of descent pn;  If n = 0, (---355----10) .  Otherwise (---355----11) with: (---355--32)
7.	Solve the sensitivity problem with: δT0 = pn;
8.	Compute the step size: (---355--33)
9.	Correct the inlet temperature profile: (---355--34)
10.	Terminate the iteration when the convergence criterion is satisfied. Otherwise,
11.	Go to step 2.

4. Results and discussion

In order to examine the feasibility and then the accuracy of the inverse analysis for estimating the unknown inlet distribution temperature T₀(z), by using the CGM. Two kinds of distribution are considered: a constant function and a smooth function. The effects of flow conditions and measurement errors are also analyzed.

4.1. Numerical solution of the direct problem

The computational technique used to solve the direct problem is based on a finite volume method (FVM) with discrete components of the governing mass, momentum, and energy equations. The staggered grid for Marker and Cell method (MAC) was employed to separate the governing equations. An augmented Lagrangian method [22] is used to compute the solution of the coupling velocity–pressure equations (1(---355--1) –4)(---355--4) . The grids are non-uniformly spaced in the z direction. The mesh spacing is taken smaller near the wall. The grids in the principal direction of flow are uniformly spaced (200 × 160). This distribution was dictated by the desire to capture the details of the viscous dissipation, because the temperature rise occurs in the region very near the wall and remains very low in the center of the channel die. The accuracy of the numerical solutions of the direct problem was studied in a previous article for a Newtonian and a pseudo-plastic fluid [24]. These solutions were compared to analytical solutions; the comparison confirms that the selected grids are adequate. The same solver and grids are used to compute the adjoint and the sensitivity variables, because in the cases which are considered the highest variations of the sensitivity and the adjoint variables are also located very near the wall–polymer interface.

The problem is solved for a geometry illustrated in figure 1. The flow of a Newtonian polymer is studied in a parallelipedic channel. The height, the length, and the width of the channel are respectively h = 2 × 10⁻³ m, L = 0.2 m, and l = 20 × 10⁻³ m, and the solid wall thickness is H = 15 × 10⁻³ m. A coordinate system (see figure 1) is considered such that Ox is the principal direction of the flow, z is the thickness coordinate, and y the width coordinate. It is assumed that the channel ratio l/h is sufficiently large such that a two-dimensional solution T(x, z) is valid at the medium plane (y = 0).

The properties of the polymer melt are taken as ρ_f  = 1000 kg m⁻³, λ_f  = 0.2 W m⁻¹K⁻¹ η = 1000 Pa s and C_pf  = 2000 J kg⁻¹ K⁻¹. The pressure drop in the channel die is taken as Δp = 30 MPa and the specified temperature at the outer wall of the die is T_p(x) = 473 K. The thermal conductivity, density and heat capacity of the metallic wall are taken for this numerical study, as λ_w = 10 W m⁻¹K⁻¹, ρ_w = 7800 kg m⁻³ and C_pw = 400 J kg⁻¹ K⁻¹.

The polymer melt enters the extrusion die at x = 0 with a temperature profile T₀(z). Two specific examples are considered: a non-uniform profile (figure 2) in the first example and a uniform one in the second example. For both cases, the direct problem is solved to generate the velocity profile in the channel (figure 3) and the measured temperature (---355----12) (figure 4). These data (---355----13) are then used as ‘experimental’ data. In a first step, the data are considered without noise, and then a zero mean Gaussian noise is added to the data. For both experiments, the temperature measurements are given by 40 sensors located at z_p, p = 1, 2 in the die wall; see figure 1. In practice, the sensor number will be less than 40, and then the inlet temperature profile T₀(z) has to be parameterized in the following form: (---355--35) where the set of np basis functions ω_k(z) is a priori given, with the parameter number less than the sensor number, (---355----14) . This parameterization is not considered in this article.

Figure 2. Inlet temperature profile (example 1): Δp = 30 MPa and T_p(x) = 473 K.

Figure 3. The x velocity profile (examples 1 and 2): Δp = 30 MPa and T_p(x) = 473 K.

Figure 4. Measured temperatures: Δp = 30 MPa and T_p(x) = 473 K.

4.2. Numerical results of the inverse problem

4.2.1. Example 1 – non-uniform inlet temperature profile-data without noise

The inlet temperature profile (figure 2) is chosen as the developed solution of a Newtonian fluid between two parallel isothermal walls. It is characterized by two peaks close to the die walls, which are due to high viscous dissipation [16]. In order to examine the accuracy of the inverse analysis by using the CGM, two different temperature guesses (---355----16) are studied: a uniform profile and a parabolic one.

1.

Uniform guess The computation starts from the uniform guess (---355----17) . Figure 5 illustrates the effects of the number of iterations on the accuracy of the estimation. It is observed that when the number of iterations increases, the final value of the least square criterion decreases and a reasonably inlet temperature solution is obtained after 40 iterations. Although the optimal number of iterations is difficult to be determined a priori, the plot of the iterative process, figure 6, confirms the efficiency of the CGM. The temperature residuals, equation (36)(---355--36) , are plotted on the figure 7 for different iteration numbers: n = 20, 40, and 100. (---355--36) Their magnitude remains small compared to that of (---355----18) . The estimated temperature profile at 100 iterations is within ±0.2 K of the exact inlet profile, except near the wall where it is within ±0.5 K of the correct profile.

2.

Parabolic guess Figure 8 shows the comparison between the estimated inlet temperature profile and the exact one. For this case, the temperature inlet profile is represented quite well, as well as for the peak and near the wall, after 150 iterations.

Figure 5. Exact and estimated inlet temperature profiles for different iterations. Example 1 – data without noise.

Figure 6. Least square criterion J(T₀(z)). Example 1 – data without noise.

Figure 7. Temperature residuals E ⁿ(x_m, z_p). Example 1 – data without noise.

Figure 8. Exact and the estimate inlet temperature profile for 150 iterations. Example 1 – data without noise.

The norm J_e of the deviation between the estimated inlet temperature (---355----19) and the exact temperature profile (---355----20) is computed at the inlet boundary of the channel die. The variability of J_e and (---355----21) versus the iteration number can be seen on figure 9. After 100 iterations, no significant corrections to the estimated profile are observed (J_e ≈ constant), while the criterion (---355----22) is still decreasing.

Figure 9. Evolution of least square criterions J(T₀(z)) and J_e. Example 1 – data without noise.

4.2.2. Example 2 – uniform inlet temperature profile-data without noise

In polymer extrusion, the material is plasticated in an extruder before passing through the die for shaping purposes. The material is normally assumed to have a uniform temperature at the entrance of the die and the die is usually controlled at a given temperature. As such, it is hoped that the extrudate will have a uniform temperature. For this reason, we try here to estimate a uniform temperature profile. In order to check the algorithm, three typical initial guesses are considered:

1.	a uniform initial guess,
2.	a fourth-order polynomial with a maximum at the centerline, and
3.	a non-uniform profile with two peaks close to the die walls.

Figures 10, 11, and 12 show the reconstructed inlet temperature profiles. For the three cases, the estimated solutions oscillate around the exact one, and there are little differences between the exact and the estimated temperature at the centerline of the channel die, except near the walls where the viscous dissipation effects are important.

Figure 10. Exact and estimated inlet temperature profiles for different iterations: case (1) – Example 2.

Figure 11. Exact and estimated inlet temperature profile: case (2) – Example 2.

Figure 12. Exact and estimated inlet temperature profile for different iterations: case (3) – Example 2.

Figures 13 and 14 show the effects of the initial guess on the convergence of the conjugate gradient algorithm. It is observed that the reconstruction of a uniform inlet temperature profile is more difficult than the previous non-uniform one. Due to the thermal history of the material in the extruder, a uniform temperature profile at the die entrance is not realistic.

Figure 13. Evolution of the criterion J(T₀(z)). Influence of the initial guess – Example 2.

Figure 14. Evolution of the criterion J_e. Influence of the initial guess – Example 2.

It must be noted that for both the aforementioned examples, the residual principle is not applicable to stop the CG algorithm, because the data measurements are considered without noise.

4.2.3. Example 1 – non uniform inlet temperature profile-with noisy data

Since all experimental data are corrupted by noise, the challenge is to develop a stable algorithm, less susceptible to noise. To assess the effect of noisy data on the reconstructed inlet temperature profile, a zero mean Gaussian noise (σ = 0.1 K) was added to the exact data at sensor locations. The exact and the noisy data are shown in figures 15 and 16.

Figure 15. Temperature measurements z₁ = 0.01495 m – Example 1.

Figure 16. Temperature measurements ; z₂ = 0.0170 m – Example 1.

Figure 17 shows the difference between the reconstructed inlet temperature profiles obtained using exact and noisy data. The temperature peak was represented quite well by the inverse solution after 150 iterations. It can be seen that the estimated profile using noisy data is not symmetric. This confirms that the level of the noise is not exactly the same in the upper and lower wall (see figures 15 and 16). We remark that, due to the existence of noise in data, it is meaningless to proceed with the iterations to a very small criterion value (figures 18 and 19), so for this case, the residual principle is applicable. After n = 150 iterations, the r.m.s of the residuals reaches the final value (---355----25) , this value is slightly less than the standard deviation of the noise. The optimal number of iterations given by the residuals principal (---355----26) , is n = 90, as shown on figures 18 and 19.

Figure 17. Exact and estimated temperature profiles with exact and noisy data – Example 1.

Figure 18. The criterion J_e with exact and noisy data – Example 1.

Figure 19. The output criterion J(T₀(z)) with exact and noisy data – Example 1.

4.2.4. Sensitivity study

The effect of the flow conditions on the reconstructed inlet temperature profiles is observed by solving the inverse problem for different values of the pressure drop Δp: (case (1) Δp = 3 MPa, case (2) Δp = 10 MPa, case (3) Δp = 30 MPa, and case (4) Δp = 300 MPa). The inlet temperature profile to reconstruct is symmetric (like for the case studied in section 4.2.1). The initial guesses (---355----27) for each case is taken uniform, equal to 473 K.

Figure 20 shows the exact and estimated inlet temperature profile T₀(z) and the variations of the output criterion J(T₀(z)) are plotted on figure 21. The best results are obtained for the cases (2) and (3); the estimated profile is in excellent agreement with the exact solution except near the wall. In case (1), there are small deviations between the exact and the estimation at the centerline of the die channel. There are significant differences for the case (4), especially at the centerline of the channel. These results are confirmed by a sensitivity analysis.

Figure 20. Exact and reconstructed temperature profiles for each pressure drops – Sensitivity study.

Figure 21. Output criterion J(T₀(z)) for each pressure drops – Sensitivity study.

These observations led us to investigate the influence of the flow and thermal conditions on the estimated inlet profiles. This study shows that the Cameron number, (---355----28) , can be used to characterize these conditions. The sensitivity of the measured temperature (---355----29) to the inlet profile T₀(z) is analyzed by considering the less favorable conditions where the variation δT₀(z) is taken equal to zero everywhere except in the interval Δz located at the centerline of the channel (z = 0.016 m). The resulting sensitivity field θ(x, z), solution of equations (12(---355--12) –17)(---355--17) is calculated for each pressure drops Δp. Figures 22, 23 and 24 show respectively the computed reduced values θ(x, z)/max_x,_z|θ(x, z)| within the solid wall. It is observed that the maximum of the sensitivity θ(x, z) is always reached at the interface polymer/wall, and that its location on x-axis depends on the Cameron number.

Figure 22. Sensitivity θ(x, z)/max_x,_z|θ(x, z)| in the upper wall : Δp = 3 MPa and Ca = 0.01.

Figure 23. Sensitivity θ(x, z)/max_x,_z|θ(x, z)| in the upper wall : Δp = 30 MPa and Ca = 0.1.

Figure 24. Sensitivity θ(x, z)/max_x,_z|θ(x, z)| in the upper wall : Δp = 300 MPa and Ca = 1.019.

The results obtained in figure 22 correspond to ‘adiabatic’ flow conditions and Ca = 0.01. Only the sensors located at the exit of the channel (x ≥ 0.16 m) are sensitive to the variation δT₀(z). Then the length of the channel is too short to estimate accurately the inlet temperature in the centerline of the channel.

Figure 23 corresponds to ‘transient’ flow conditions and Ca = 0.1. The effect of the variation δT₀(z) appears clearly along the whole channel. The sensitivity is great enough to estimate the inlet temperature profile at the centerline of the channel.

Figure 24 corresponds to the ‘established’ flow conditions and Ca = 1.019. The average temperature of the flow remains constant, the sensitivity along the flow axis is close to zero as soon as x ≥ 0.06 m.

5. Conclusions

The classical conjugate gradient algorithm was developed for computing the numerical solution of the inverse heat transfer problem, which consists in the determination of the unknown inlet temperature profile of the melted polymer at the entrance of an extrusion die. The additional temperature data are measured in the solid wall of the die. The method was validated using exact and noisy data. The sensitivity analysis has shown that the most favorable conditions are characterized by the Cameron number value Ca = 0.1. These numerical results are encouraging and the applications of the method to the actual experimental data for a pseudo-plastic polymer are in progress.

Nomenclature

1.	Δp = Pressure drop, Δp = P0 − PL
2.	cp = Heat capacity
3.	η = Dynamic viscosity
4.	F = Source term of energy
5.	T0 = Inlet temperature
6.	pn = Direction of descent at iteration n
7.	J = Functional to be minimized
8.	(---355----30)  = Measured temperatures
9.	E = Temperature residuals

Subscripts

1.	m = Sensor location
2.	p = Sensor location
3.	e = Inlet die
4.	w = Wall
5.	f = Fluid

Superscripts

1.	n = Iteration number
2.	p = Pressure
3.	P0 = Inlet pressure
4.	PL = Outlet pressure
5.	(---355----31) = Fluid velocity
6.	(---355----32)  = Mean velocity
7.	T = Temperature
8.	n = Vector normal
9.	Ns = Number of sensors

Greek symbols

1.	γ = Shear rate
2.	γn = Step size at iteration n
3.	βn = Conjugate coefficient at iteration n
4.	Ψ = Adjoint variable
5.	θ = Sensitivity variable
6.	ρ = Mass density
7.	δ = Dirac delta function
8.	Ω = Spatial domain

Acknowledgments

The numerical computation were supported by the laboratory TREFLE, UMR 8508, Laboratoire interétablissements CNRS, ENSAM, ENSCPB, Université Bordeaux 1.