Int J Mater Form DOI 10.1007/s12289-015-1234-y
ORIGINAL RESEARCH
Numerical simulation and optimization of the injection blow molding of polypropylene bottles - a single stage process J. Biglione1 · Y. B´ereaux1 · J.-Y. Charmeau2 · J. Balcaen2 · S. Chhay2
Received: 11 December 2014 / Accepted: 16 March 2015 © Springer-Verlag France 2015
Abstract Single stage injection blow molding process, without preform storage and reheat, could be run on a standard injection molding machine, with the aim of producing short series of specific hollow parts. The polypropylene bottles are blow molded right after being injected. This implies that the preform has to remain sufficiently malleable to be blown while being viscous enough to avoid being pierced during the blow molding stage. These constraints lead to a small processing window, and so the process takes place between the melting temperature and the crystallization temperature, where the polypropylene is in its molten state but cool enough to enhance its viscosity without crystallizing. This single stage process introduces temperature gradients, high stretch rate and high cooling rate. Melt rheometry tests were performed to characterize the polymer behavior in the temperature range of the process, as well as Differential Scanning Calorimetry. A viscous Cross model is used with the thermal dependence assumed by an Arrhenius law. The process is simulated through a finite element code (POLYFLOW) in the ANSYS Workbench framework. The geometry allows an axisymmetric approach. The transient simulation is run under anisothermal conditions and viscous heating is taken into account. Sensitivity studies
J. Biglione
[email protected] Y. B´ereaux
[email protected] 1
INSA-Lyon LaMCoS UMR5259, Site de Plasturgie, B.P.807, 01108 Oyonnax, France
2
INSA-Lyon IMP-LMM UMR5223, Site de Plasturgie, Oyonnax, France
are carried out and reveal the influence of process parameters such as the material behavior, the blowing pressure and the initial temperature field. Thickness measurements using image analysis are performed and the simulation results are compared to the experimental ones. The simulation shows broad agreements with the experimental results. An optimization loop is run to determine the optimal initial thickness repartition. Design points are defined along the preform and the optimization modifies the thickness at these locations. Keywords Injection blow moulding · Finite element method · Polypropylene · Oscillatory rheometry · Dynamical mechanical analysis · Optimization · Thickness measurement · Random copolymer
Introduction Numerical simulation in polymer processing has to aim towards more accuracy and more predictability to efficiently assist the design of polymer parts. Indeed numerical simulation can cut down the number of trials required to set the process, thus saving time and money. Moreover, numerical simulation should not only be used to predict the process but also to optimize it. The notion of optimization encompasses process parameters and geometry of the part. To do so, the physics of the process have to be well understood and the polymer characterized in the temperature and strain range of the process. The injection blow molding single stage process is a variant of the conventional injection stretch blow molding process. Here the preform is blow molded right after being injected, without storage and reheat. This single stage process does not involve a rod-stretching stage.
Int J Mater Form
An overview of the numerical simulation of blowing processes is proposed by [2] In the case of the conventional process, the polymer used is most of the time the PolyEthylene Terephtalate (PET). The particularity of the PET is its strain hardening behavior, allowing it to resist during the process and not fail. Polypropylene (PP) is much more frequently used in injection molding than PET, but much less so in injection blow molding. The polymer used here is a PP random copolymerized ethylene. General properties of PP are presented by [7] where copolymer PP are said to be more transparent and more resilient to shocks than the homopolymers PP but their Young modulus, melting point and density are smaller. Numerical simulation of the blow molding stage of the process requires a choice of constitutive law of the material consistent with the type of polymer and the type of strain it is submitted to. Fatt and Ouyang [8, 15, 22, 23] and [9] simulated the injection stretch-blow molding of PET parts using a viscoplastic model. Huang et al. [11] studied this same process through the visualization of the preform growth. Bellet et al. [2] explain that in the case of a molten polymer, a generalized Newtonian or a viscoelastic law is legit, as it is the case in the extrusion blow molding process. Hopmann et al. [10] use an integrative method to simulate this process. A first step concerns the heating stage and a second step is dedicated to the stretch blow molding stage. The structural response of the bottle to a given load is then studied using the commercial finite-element software ABAQUS. The modelisation of the material behavior is clearly of importance, but, nonetheless, it has to be run under the rights process conditions. The knowledge of the temperatures of the process, as well as the actual pressure imposed to the preform is determinant. Fatt and Ouyang [8, 15, 23] and [9] represent the air evolution in the hollow part during the blow molding stage. A thermodynamical law links the volume evolution to the establishment of the pressure along the preform surface. More recently, [14] studied a free blow and simulated it using the commercial finite element software ABAQUS. They show that it is more accurate to simulate the evolution of the air flow than imposing a constant normal force only. Schmidt et al. [17] study the pressure evolution in a hollow part and highlight the need to take it into account to simulate a blow molding process. Bellet et al. [2] made the same observation in their review of the numerical simulation of blowing processes. Bordival et al. [4] used a thermodynamical model to simulate the air pressure, taking into account the impact of the volume variation in the evolution of the pressure in the simulation of the stretch-blow molding of PET bottles. Once the material behavior is identified and the simulation is set up, it is of interest to optimize the process to help to the design of real parts. Thibault et al. [21] propose an optimization of the preform geometry and of the process
parameters in the case of injection stretch-blow molding of PET bottles. A mathematical description of the preform is used to optimize the geometry of the preform. Bordival et al. [4] optimize the initial temperature repartition in injection stretch-blow molding of PET bottle. The optimization is done on the initial temperature repartition to obtain a final homogeneous thickness repartition using a Globalized Nelder-Mead algorithm proposed by [13], which consists in probabilistic restarts of the classic Nelder-Mead method, allowing further searches for the minimum value of the cost function. In [12], prediction of the optimal thickness distribution in extruded blow molded preform is carried out. A FEM code is used to simulate the process. A NewtonRaphson scheme is used, with under-relaxation factor to facilitate the convergence and to avoid oscillations. In the first part of this work, the material characterization of the polymer in the molten state is presented, as well as the constitutive equation chosen. The process conditions are then presented and the principle of the numerical simulation follows. Sensitivity studies are carried out and reveal the influence of process parameters such as the material behavior, the blowing pressure and the initial temperature field. The results are then discussed and compared to the effective thickness repartition measured by image analysis. It is sometimes claimed that the blow molding process is too constrained by itself to allow for meaningful differences to be observed in the thickness repartition. Our results prove this is not quite the case. Finally, optimization of the initial geometry is presented.
Characterization Random copolymers are often used in blow molding processes. Random copolymer PP consists of chains mostly composed of propylene blocks and of a few ethylene units. The presence of ethylene in the formulation grants a better flexibility to the polypropylene. This flexibility is suited to processes involving high stretches. It is a semi-crystalline thermoplastics. Considering semicrystalline polymers implies that the shaping process takes place at temperature higher than the crystallization temperature. The glass transition temperature of the PP is in the range of −13 ◦ C and it’s melting temperature is above 150 ◦ C. This particularity enables the process to take place at lower temperature than the conventional process. The PP viscosity has to be high enough so that even under strain the viscosity is sufficiently high to resist to air pressure. This is why an injection blow molding grade was chosen in this work: PPR 3221 provided by Total Petrochemical. We first investigate the melting temperature and solidification temperature of the polymer through differential scanning calorimetry (DSC).
Int J Mater Form
Differential scanning calorimetry
Dynamical rheometry
The Fig. 1 represents the results of the DSC test run on a TA instrument Q20 machine at 10 ◦ C.min−1 . The first heating cycle goes from room temperature to 250 ◦ C, the cooling cycle from 250 ◦ C to -50 ◦ C and the second heating cycle from -50 ◦ C to 250 ◦ C. The two peaks during the fusion in the second heating are due to the presence of two types of structure in the morphology of the PP formed under quiescent crystallization during the cooling phase. This was also seen in [5]. The two types of morphology correspond to the radial lamellas and to the tangential ones. The tangential lamellas correspond to the peak at the highest temperature, thus the radial lamellas correspond to the peak at the lowest temperature. A fusion enthalpy of 71.75J.g−1 is measured in the first heating. Given the reference enthalpy is equal to 165J.g−1 as chosen by [6], a crystallinity of 44 % is deduced. The analysis of those results leads to the properties presented in Table 1. The crystallization is a phenomena which strongly depends on the cooling rate. The question of the crystallization occurring in the present process is of order. To bring some hints to this problem, DSC trials have been run at different cooling rates: after heating the material up to 200 ◦ C and waiting for the stabilization of the temperature, a cooling ramp is imposed from 200 ◦ C to room temperature. The Fig. 2 shows the DSC curves of these experiments. One can see that the faster the cooling rate, the lower the temperature at which the crystallization occurs. The same observation has been made by [18, 20] and [19]. Those works highlight the effect of cooling rate on the evolution of the crystallization of isotactics PP.
The rheological behavior is investigated under shear in molten state, using oscillatory parallel plates, on a TA instrument ARES 2k FRT. Frequency sweeps are performed at different temperatures. The viscosity is then deduced and the effects of shear and of temperature are highlighted. The temperature of 180 ◦ C is chosen as the reference temperature. The shift factors aT , allowing for time temperature superpositions, are determined and the master curve can be obtained (Fig. 3). From this master curve, the shear thinning effect is studied and a shear thinning law is fitted. The temperature dependence is deduced from the curve log(aT ) = f ( T1 ) represented in Fig. 4. The Cross model is used to represent the shear thinning behavior:
Fig. 1 Diagram of the DSC test at 10 ◦ Cmin−1
η=
a T η0 1 + (aT λγ˙ )1−n
(1)
With η0 the zero-shear viscosity, λ the characteristic time and n the shear thinning exponent. The effect of the temperature is assumed by the use of an Arrhenius law: aT = e
Ea R
1 1 T +TK − Tref +TK
(2)
with R the perfect gas constant, Ea the material activation energy, TK equal to 273.15 ◦ C and Tref the reference temperature. The values for the parameter of the viscosity model are presented in Table 2.
2
1st heating cooling 2nd heating
1.5
Heat flow [W.g- 1 ]
1
0.5
0
-0.5
-1
-1.5 -100
-50
0
50
100
150 o
Temperature [ C]
200
250
300
Int J Mater Form Table 1 DSC Results at 10 ◦ C.min−1
Glass transition [◦ C]
Crystallization point [◦ C]
Melting point [◦ C]
Crystallinity [%]
-17
113
149
44
Numerical simulation
The mold opening phase lasts 6s and is also a transient heat conduction problem but with a heat flux boundary condition imposed at the outside of the preform instead of the previous mold temperature. The air temperature is set at 35 ◦ C and the convection coefficient equal to 10W.m2 .s−1 , as in a free convection problem. Finally, the blow molding stage is simulated. The preform is molded while the polymer is in its molten state and the high cooling rate involved in the process does not promote the crystallization. As seen previously, the faster the cooling, the lower the crystallization temperature. Plus, in the beginning of the project, as the final set of parameters were not already defined, trials with the broach temperature imposed at 90 ◦ C have been carried out. The blow molding of parts was a failure and the final parts exhibited defects similar to cold drawing. Those defects were obviously due to the presence of crystalline structure. On the contrary, with a broach temperature fixed at 150 ◦ C, the preform was not viscous enough to withstand the air pressure and the final parts were pierced. The temperature of 135 ◦ C for the broach and a temperature of 135 ◦ C for the mold have proven to be the only working set of parameters for the given geometry at the time. Both DSC experiments at different cooling rates and actual observation tend to legitimate the assumption that crystallization happens after the blow molding stage, and that during the process the material remains in its molten state.
Process conditions The process allows to produce hollow parts on a standard injection molding machine. The preform is blow molded without prior stretching by a rod. This is a sequential process and a special mold had to be designed with internal movements enabling the sequence of the following steps: 1. injection molding of the preform; 2. holding and cooling of the preform; 3. opening of the mold and transfer from the injection molding location to the blow molding location; 4. blow molding of the preform; 5. opening of the mold and ejection of the part; In this work, the injection molding is not taken into account but can be simulated using a FEM software (Moldflow). We suppose here that no residual stresses are left in the preform after the injection and that the preform is considered homogeneous with an initial temperature equal to the injection temperature (240 ◦ C). The holding and cooling phase, during 10s, are simulated through a transient heat conduction calculation with a temperature of 135 ◦ C imposed both on the inside of the preform (broach temperature) and on the outside of the preform (mold temperature).
Fig. 2 DSC curves for trials at several cooling rates
3
1oC.min-1 5oC.min-1 10oC.min-1 15oC.min-1 30oC.min-1 50oC.min-1
Heat flow [W.g-1]
2.5
2
1.5
1
0.5
0 20
40
60
80
100 o
Temperature [ C]
120
140
160
Int J Mater Form Fig. 3 Master curve at a reference temperature of 180 ◦ C and the fitted Cross law
105
Master curve points Cross law
η [Pa.s]
104
3
10
102 -2 10
10-1
100
-1
101
102
103
ω [rad.s ]
–
Geometry The Fig. 5 represents the geometry of the problem with the different domains and boundaries, along with the conditions imposed on them in the numerical model. Two domains are defined: the mold and the preform. In the mold, two boundaries are set up: the contact side, which correspond to the inside of the mold, and the other sides which are not used in the calculation. In the preform, four boundaries are defined: –
the exterior boundary, which will be taken into account in the contact detection with the mold;
Fig. 4 Shift factors and the fitted Arrhenius law
–
–
the interior boundary, where the pressure is imposed constant along the boundary, with its value following a transient evolution; the upper boundary, where a condition of zero normal velocity and zero tangential force is imposed (perfect slip); the remaining boundary is the axis of symmetry;
Meshing of the preform The mesh is composed of quadrilateral elements only in the preform. The side of the mold where the contact takes place is refined so that the detection of the preform is enhanced.
0.4
0.2
log(aT)
0
-0.2
-0.4
-0.6
-0.8
shift factors Arrhenius law
2
2.05
2.1
2.15
2.2
1/T [1E-3 K-1]
2.25
2.3
2.35
Int J Mater Form Table 2 Parameters of the viscosity model η0 [Pa.s]
λ [s]
n
Ea [J.mol−1 ]
Tref [◦ C]
2.6E4
1.37
0.35
34 600
180
This leads to a total of 6827 nodes for 6084 faces. The Fig. 6 shows the total meshed geometry, as well as a zoom on the lower part of the preform. Finite elements approach Polyflow is used to simulate the problem. This a FiniteElement Method (FEM) Computational Fluid Dynamics (CFD) program dedicated to the resolution of fluid problems involving viscous and viscoelastic fluid flows. An Eulerian approach is used with velocity and pressure fields as primary unknowns. These unknowns are to be supplemented with the location of the free surface bordering the flow in the case of a blow molding simulation. A free surface problem always implies a fair amount of remeshing, as the mesh representing the preform geometry will undergo large deformation. The deformation is mostly extension, but a part of if could also be shear. The chosen remeshing method is the Thin Shell Method with Lagrangian Master [1], which is well adapted for blow molding simulations. This method is only available in 2D. The objective is to allow the mesh elements to keep an acceptable shape, as close as possible to a rectangular one. A Lagrangian displacement is imposed for the nodes on the boundary involved in the contact problem, Holding and cooling phase Broach temperature
fn = −k(v − vwall )n
(3)
with n the normal to the surface, k the penalty factor, vwall the wall velocity and v the node velocity. In the tangential direction, the method is similar, allowing the prescription of slip to be more accurate (Eq. 4). ft = −Fslip (vs − vs,wall )
(4)
with Fslip the slip coefficient, vs,wall the tangential wall velocity and v the tangential node velocity. Fslip and k have physical dimensions g.s−1 . A typical value used in air calculation for k is 109 g.s−1 . Slip is allowed when Fslip is lower than k.
Opening of the mold
Broach temperature
Blow molding phase
Broach temperature
normal velocity=0 Flow boundary conditions
contact
Pressure Symmetry axis
Mold Temperature
free convection Symmetry axis
Broach temperature Symmetry axis
free border
free convection
free convection
Mold temperature
Broach temperature
Thermal boundary conditions
Symmetry axis
Fig. 5 geometry 2D of the problem and the conditions imposed in each step of the numerical simulation
this boundary is called master surface. The remaining nodes in the thickness are remeshed following a spine rule, which is described in [1]. This method is robust, even in presence of shear. Such a method requires a structured mesh, with the same amount of nodes along the opposite boundaries of the preform. As defined in [1], the contact with the mold is detected by the use of a local procedure performed at each location along the free surface. The use of the Lagrangian method in the remeshing problem facilitates the contact detection, because this method makes the nodes follow the physical deformation. The contact is then more realistic. Once a node slightly penetrates the mold domain, the discrete momentum equation of the node is modified by a penalty term, which can be different in the normal and tangential directions, allowing slipping for example. In the normal direction, the large value of the coefficient is defined such that the velocity of the node will remain close to the wall velocity (Eq. 3).
Int J Mater Form Fig. 6 Meshed geometry
Yang et al. [24] simulate the injection stretch-blow molding of PET bottles under axisymmetric approach using ABAQUS. The membrane approach allows one to limit the number of nodes in one element, but it constrains the temperature to evolve linearly between nodes. This is not consistent with the experimental observation because the temperature and stress histories are strongly non-linear. Thus it is chosen here to use an axisymmetric approach in the numerical model, so that an accurate determination of the temperature field could be done. The computation is run on a workstation using 4 processors for a typical run time of one hour.
Simulation results One can say that this kind of modelisation is too constrained by the geometry of the mold so that the simulation would not be very sensitive to the chosen material and process parameters. It is shown, in a preliminary study, that this is not quite the case here. Effect of the initial temperature field To illustrate the effect of the initial temperature distribution, a simpler geometry is used here. The simulation only concerns the blow molding stage. Two simulations are run, both with same temperature conditions but with different waiting time before blow molding: – – –
The initial temperature of the material is set at 240 ◦ C; A temperature of 105 ◦ C is imposed on the inner boundary of the preform; A temperature of 50 ◦ C is imposed on the outer boundary of the preform;
–
A constant pressure of 2.5bar is imposed along the inner boundary of the preform;
The temperature boundary conditions have been arbitrary chosen to magnify the temperature gradient. Figure 7 exhibits the results of these computations. The first case (a) has no cooling time. It is immediately blow molded so that no temperature gradient can develop. The second case (b) is blow molded after 10s of cooling, so that the initial temperature of 240 ◦ C decreases, due to the temperatures imposed on the inside and outside of the preform. The upper part of the preform is the thinnest one. When the preform core remains at a uniform temperature, the thinnest part is stretched first, because it is less resilient to the deformation. However, when a non-homogeneous temperature field is present, the local thickness is not the only parameter governing the evolution of the preform. The local temperature matters too, through the thermal dependency of the viscosity. The local stiffness is linked to the local thickness and to the local viscosity. The case (b) exhibits a different stretch pattern of the preform, due to the thinnest part being colder and the thicker part being hotter, and therefore, less viscous. This counter intuitive effect is a direct reminder of the specificity of the Injection Blow Molding Single Stage process. Effect of the shear thinning behavior The chosen model does matter in the final thickness repartition of the simulation. Thus, to highlight the effect of the modelisation of the material behavior, the flask problem have been computed using, on the one hand, the shear thinning law presented above and using, on the other hand, a Newtonian behavior, i.e. with no shear thinning effect. The thermal dependance is assumed in both case by the same
Int J Mater Form
No cooling time
Cooling time of 10s
Temperature [°C] 250.0 235.7 221.4 207.1 192.9 178.6 164.3 150.0 135.7 121.4 107.1 92.86 78.57 64.29 50.00
Temperature [°C] 250.0 235.7 221.4 207.1 192.9 178.6 164.3 150.0 135.7 121.4 107.1 92.86 78.57 64.29 50.00
(a-1) (a-2)
(b-1)(b-2)
Fig. 7 effect of the initial temperature - (a-1) state before blow molding - homogenous temperature - (a-2) first contact - homogenous temperature - (b-1) state before blow molding - inhomogeneous temperature - (b-2) first contact - inhomogeneous temperature
Arrhenius law. A constant pressure of 2.5bar is imposed in both cases, with no evolution during time. Figure 8 presents the initial and finals thickness repartitions of both computations. The relative difference between the two calculations presents a mean value of 16 % with a maximum value of 91 % where the stretch is the higher, thus where the material behavior is expressed the most. The differences cannot be denied and illustrate the need for pertinent material behavior characterizations to aim to better predictability. Effect of the imposed pressure The imposed pressure does matter in the final thickness repartition. A simulation is run with an imposed constant
Fig. 8 Thickness distribution Newtonian case and Cross case
2.5
pressure of 2.5bar. This simulation is compared to the one run with a transient pressure, which is closer to the reality of the process. The pressure evolution is governed by the following law, based on the perfect gas law: P = P (t) − P0 =
P0 V0 + ntRT ˙ air − P0 V (t)
(5)
With P(t) the pressure value at the instant t, P0 the initial pressure value, V0 the initial available volume encapsulated by the preform, V(t) the available volume encapsulated by the preform at the instant t, n˙ the molar rate of the air flow, R the gaz constant, Tair the air temperature. The pressure evolution for the first instant of the blow molding stage is represented in Fig. 9.
Body
Top
Bottom
thickness [mm]
2
1.5
1
0.5 initial Cross final Newtonian final
0
0
0.2
0.4
0.6
normalized length
0.8
1
Int J Mater Form
Figure 10 presents the initial and final thickness repartitions of both computations. The relative difference between the computations presents a mean value of 11 % with a maximum value of 61 %. A transient imposed pressure is smoother than a constant pressure: the first time step of the calculation will see a pressure value far below the constant pressure value, thus lesser strain gradients are involved in the beginning of the calculation. The transient pressure is like a succession of constant pressure problem, the final results are obviously different. This highlights the necessity of approaching the true evolution of the air pressure during the process to improve the accuracy of the simulation. For this work, the transient air pressure is imposed but experimental data is missing to set the value of the air mass flow rate in Eq. 5. Effect of the activation energy The parameters used to represent the material behavior through the Cross-Arrhenius law are determinant on the final thickness repartition obtained by the simulation. To illustrate this sensitivity, the previous simulation with the transient imposed pressure is run under the same conditions, but with an activation energy 30 % lower than the one determined through characterization. Figure 11 shows the comparison between the standard simulation and the one with the lower activation energy. Effect of the activation energy It appears that the differences are more noticeable in the area where the preform is the thickest, thus where the gradient of temperature are more in evidence. These gradients of temperature induce gradients of viscosity through the Arrhenius Fig. 9 Pressure evolution during the blow molding stage
law. A lower activation energy leads to higher viscosity for a given temperature. However, it leads to lower gradient of viscosity for a given gradient of temperature. This explain why the sensitivity to the activation energy express itself the more in this area. Results of the flask simulation The geometry studied here is more complicated, due to the curvature changes along the bottle. This one geometry was set to prove the efficiency of the innovative technology developed in the IS2 project. As one can see in Fig. 12, the temperature decreases a lot during the first stage of the process. The injected material is in the preform mold and the difference between the injection temperature of 240 ◦ C and the mold temperature of 135 ◦ C is large, leading to high cooling rates. As it has been said before, we assume that these high cooling rates do not allow crystallization to occur. A hot spot remains in the preform elbow, a second one can be seen close above the first one. After the opening stage, the temperatures have decreased again, but the hot spots are larger. During the blow molding stage, the temperature increases due to the viscous heating occurring under shear rate as shown in Fig. 13, but not in large amounts. Figure 14 shows the deformation of the preform with time to finally fit the mold surface. Although the thinnest part in the height of the bottle is colder than the thickest one, it is the first to be stretched. This means that the difference in temperature is not high enough to make the thickest part stretches first. Its neighborhood is then stretched, thus its viscosity decrease due to the shear thinning behavior and it will be stretched more easily as stretching goes.
3
2.5
Pressure [bar]
2
1.5
1
0.5
0 16
16.02
16.04
16.06
time [s]
16.08
16.1
Int J Mater Form Fig. 10 Thickness distribution constant pressure versus transient pressure
2.5
Body
Top
Bottom
thickness [mm]
2
1.5
1
0.5 initial Constant final Transient final
0
0
0.2
0.4
0.6
0.8
1
normalized length
Once the nodes are in contact with the wall, they remain in contact with no further slipping, due to the chosen penalty coefficient in the local momentum equation. The nodes in its vicinity are then stretched until contact is made, because they are the less viscous under the local shear. Thickness measurement method Different methods of measurements exist: Non-destructive ones and destructive ones. The non destructive ones, except for tomography X, are point-wise measurements. For example, the magnetic probe method uses a magnetic ball on one
Fig. 11 Thickness distribution 100 % Ea versus 70 % Ea
2.5
side of the part to be measured, and a magnetic sensor on the other side. To obtain the thickness distribution, a sweeping method needs to be set up and this method will be unique for every geometry. Moreover, the diameter of the ball limits the different geometry that can be swept. Steep variation of angle along the geometry could not be reached by the ball. Optical measurement methods could also be invoked, but present the same drawback. The method used here only needs a scanner, which resolution will define the precision of the measurement and digital image analysis, available in a number of software suite. The bottles are cut in two halves and each half is scanned using a resolution of 2400 dpi.
Body
Top
Bottom
thickness [mm]
2
1.5
1
0.5
0
initial Ea final 0.7xEa final
0
0.2
0.4
0.6
normalized length
0.8
1
Int J Mater Form Fig. 12 Evolution of the temperature during the simulation of the Injection Blow Molding process
Temperature [°C] 250.0 232.5 225.0 217.5 210.0 202.5 195.0 187.5 180.0 172.5 165.0 157.5 150.0 142.5 135.0
Temperature [°C] 146.5 145.2 143.8 142.5 141.2 139.9 138.5 137.2 135.9 134.5 133.2 131.9 130.5 129.2 127.9
Temperature [°C] 161.7 159.7 157.8 155.9 154.0 152.1 150.2 148.3 146.4 144.5 142.6 140.7 138.8 136.9 135.0
Holding and Cooling stage
0s
Opening of the mold
10s
Figure 15 represents the scanned image of one of the bottle. A cropping method is used so that the bottle is well centered in the picture. This correction ensure that the applied treatment is the same for all the images. The images are then analyzed: a threshold treatment is applied resulting in binary images. The edges of the bottle are then detected and the local thicknesses are measured. The error of the method is estimated at ± 1 pixels, thus 10μm. Ten thickness repartitions are measured, corresponding to five bottles. A median repartition is deduced from the ten measurements. Fig. 13 Temperature repartition along the outer boundary at the beginning and the end of the blow molding stage
Temperature [°C] 147.0 145.7 144.3 143.0 141.6 140.3 138.9 137.5 136.2 134.8 133.5 132.1 130.7 129.4 128.0
Blow Molding stage
16s
17s
Comparison with simulation is done using a normalized length representation. This normalized length is the curvilinear abscissa along the bottle profile divided by its maximal value. The curvilinear abscissa is calculated as follow: lm =
m dxi2 + dyi2
(6)
i=1
where dxi = xi − xi−1 and dyi = yi − yi−1
140 Initial Final
138
o
Temperature [ C]
136
134
132
130
128
Body
Top 126
0
0.2
0.4
Bottom 0.6
Normalized length
0.8
1
Int J Mater Form Fig. 14 Computation results of the blow molding stage - the color map represents the local thicknesses
(a)
(b)
This is equivalent to the sum of each local distance between the two consecutive points (i-1) and i. Agreement between experience and simulation Figure 16 shows the comparison between experimental and numerical thickness repartitions. The neck region is not considered for comparison here (see grayed areas). The presence of the screw threads at this place has led to a geometrical simplification for the numerical simulation. It would not be relevant to take them into account in the comparison. Considering the other regions, there is 73 % of concordance with a tolerance of 100 μm and 45 % with a tolerance of 50 μm .
(c)
(e)
(f)
(g)
The sharp peak at the edge of the Z5 area in the simulation results in Fig. 16 is not always present in the experimental measurements. On some particular measurements, this peak is indeed observed. The reason why it is not always observed might come from the lack of precision in the cutting operation of the bottle (thickness of the cutting tool, cutting angle, cutting location). Overall, the possible source of error in the simulation can be the following: –
–
– Fig. 15 Scanned picture of a bottle cut in half
(d)
(mm) 2.50 2.27 2.04 1.81 1.58 1.35 1.12 0.89 0.66 0.43 0.20
Simplicity of the constitutive equation: a Cross shear thinning law could not be enough to represent the total behavior of the polymer during the blow molding stage. Crystallization might occur, leading to the reinforcement of the local stiffness and modifying the local deformation. The PPR3221 is a random copolymer which optical properties are enhanced. It is also a semi-crystalline polymer but its transparency is evident. The formulation might refrains the growth of the crystalline structures, leading to a material composed of amorphous polymer charged with tiny crystalline particles, which size are lower than the wavelength of the visible light, allowing good transparency. Moreover, the high cooling rate involved in the process does not promote crystallization. The Cross law is certainly a good approach but it has to be ascertained in further work. Simplification of the thermal and mechanical history before the blow molding stage: residual stress can remain after the injection stage, and the temperature repartition might be different than the one obtained from the simplified simulation of the first two stages. The temperature is important to define the local stiffness of the preform, because it results from the local viscosity coupled to the local thickness. Moreover, in the case of models involving crystallization, the determination of the temperature has to be as accurate as possible. The transient pressure imposed is not close enough to the reality. In situ measurement of the air pressure
Int J Mater Form Fig. 16 Comparison between experience and simulation
2 Z6 Z1 Z2
Thickness [mm]
1.5
Z3
Z1
Z4
1
Z3 Z5
Z4 Z5
0.5
Z6 Z2
0
0
0.2
0.4 0.6 Normalized length
would be needed to increase the accuracy of the law governing the air evolution. In regard of all the sources of error, no definitive conclusion can be made. The experimental measurements presents too many source of dispersion to allow an accurate comparison with the numerical simulation. Nonetheless, the general pattern of the thickness repartition is quite the same. Further works to improve both the experimental measurements and ascertain the numerical simulation assumptions have to be done to validate the model. Currently, we can conclude that a broad agreement have been found between experimental measurement and simulation.
Initial thicknesses
Modification of thicknesses
Experiment
Simulation Experience
0.8
Simulation
1
Optimization Principle of the optimization The optimization of the preform geometry follows the loop described in Fig. 17. The initial geometry is controlled by the use of design points. The simulation of the process, as described above, is then run and the resulting geometry of the blow molded product is tested. If the final geometry does not satisfy the fixed goal, then the initial geometry is modified and the simulation is run again. Once the final geometry is in good agreement with the fixed goal, the optimization is completed. The goal here is to obtain an uniform thickness repartition in the final geometry. Only the thickness in the height of the bottle is optimized. The bottom and the neck are not optimized. The initial geometry is controlled in various
Simulation IBM
node 1 Final thicknesses
node 1 trajectories NO
YES TRIAL
node 2
END
node 2
Optimization Loop
Fig. 17 Principles of the optimization of the process - short dash lined boxes are the parametrized inputs of the loop - straight lined boxes are the running of the calculation - width dash lines boxes are the tested outputs of the loop
Blow molding enini
Fig. 18 Illustration of the local problem
enini
Int J Mater Form
Algorithm
Parametrization
Gain [%]
Number of runs
Elapsed time
Predictor/Corrector Predictor/Corrector
14 points 18 points
97 99
15 15
7h30min 7h30min
locations along the preform. The inner boundary, representing the outer surface of the broach, is fixed. Only the local distances between the broach and the outer side of the preform are parametrized. The parametrized locations are design points which have to be tracked during the blow molding simulation. Knowing the trajectory of the materials points allows to link the final thickness to the concerned initial thickness. The appropriate correction for the next iteration can then be calculated. To do so, the nodes on the outer boundary can be used, as they follow the material deformation due to the Lagrangian remeshing method defined on this side. The objective function to minimize is defined for the n-th iteration as: n = Fobj
N
A Predictor/Corrector algorithm is used to solve the optimization problem. Predictor/Corrector method The Predictor/Corrector consists of considering each design points as an individual problem. The final thicknesses are compared to the targeted thickness and the initial thicknesses are corrected accordingly. The principle of this method is illustrated in Fig. 18. The thickness of the next iteration is then defined: en obj n+1 initial n n einitial = einitial + αrelax efinal − efinal (8) n efinal en
obj
obj
n − ei,final ei,final
2
(7)
i=1
with N the number of points considered in the optimizan the final computed thickness of the tion problem, ei,final obj
i-th point along the preform and ei,final the final objective thickness of the i-th point. For now, each objective thickness is set at the same value so that the optimization aims for an uniform thickness obj obj repartition: ei,final = efinal for i = 1..N Fig. 19 Evolution of the objective function (left axis semi logy scale) and evolution of the coefficient C ruling the relaxation law (right axis)
the local stretch and with efinal the objective thickness, einitial n final α a relaxation parameter following the rule: αrelax = C relax n+1 1 n then . C is initially equal to 0. When Fobj > Fobj 2 C = C + 1. The next Prediction/Correction is done from the n-th solution with the new relaxation parameter. C is incren+1 n . Once F n+1 < F n , its > Fobj mented each time Fobj obj obj value is set to its initial value 0. The results of this method are available in Table 3. This method shows fast convergence and, more importantly, the convergence rate does not seem to be impacted by the increase of the number of design points.
0.12
6 Objective function Coefficient C
Log(Objective function)
0.1
5
0.08
4
0.06
3
0.04
2
0.02
1
0
0 2
4
6
8
Number of iteration
10
12
14
Coefficient C
Table 3 Results of the Predictor/Corrector Algorithm
Int J Mater Form Fig. 20 Optimized profil (blue) versus raw profil (red)
4
raw initial raw final optimized ini raw ini objective thickness
3.5
Top
Body
Bottom
Thickness [mm]
3
2.5
2
1.5
1
0.5
0
0
50
100
150
200
250
300
Mesh Index
Results of optimization The Fig. 19 represents the convergence evolution of the Predictor/Corrector method. The efficiency of this method is well observed: large gain in the firsts runs are then followed by smaller gain in the last runs. At the 8th iteration, the objective function increases, leading to the increase of the exponent C in the relaxation law. When the following iteration leads to an inferior objective function then the coefficient C returns to its initial value. When the following iterations do not give a better objective function than the 9th, the coefficient C keeps increasing. The 9th solution, giving the smallest objective function, is then retained as the optimal one.
Initial thicknesses
Modification of thicknesses
Final thicknesses
Modification of thicknesses
Fig. 21 Complete Optimization principle - short dash lined boxes are the parametrized inputs of the loop - straight lined boxes are the running of the calculation - width dash lines boxes are the tested outputs of the loop
The Fig. 20 compares the initial and final thickness repartition of the raw geometry and the optimized one. One can see that the overall thickness repartition is more uniform after the optimization. In the vicinity of the bottom of the preform, from index 50 to index 80, the optimization of the thickness is not that successful. This is a border effect, these regions corresponding to the limits of the optimization problem where no design point is defined. The same effect can be seen in the neck vicinity, from index 230 to index 250. This is due to the fact that the preform, in the neck and the bottom of the bottle, is very near to the mold. The bottom thickness has to be thick to assure the steadiness of the bottle, and the neck region in the simulation is a simplification of the screw flight. This explains why these regions are not
Mechanical Computation
Mechanical properties
NO
1st loop
Final thicknesses
NO
YES TRIAL
Simulation IBM
Final thicknesses
YES TRIAL
2nd loop
END
Int J Mater Form
involved in the optimization. The fact that the thicknesses at these locations are fixed leads to the border effects, because the design points at the border of the optimized region are in minority, next to a fixed region, thus limiting the effect of their variations. The high variation of thickness in the vicinity of the maximum (index 75) highlights the difficulty to obtain a uniform thickness repartition in the reality of the process. This peak is not a pure numerical effect but is a consequence of the geometry. It appears where the preform is the most stretched before being in contact with the mold. Machining a mold with this huge variation would not be pertinent, and assuming it is, the molding and blow molding of the preform would suffer from this high thickness gradient. The optimization warns that, by imposing the bottom of the preform to be always in contact with the mold, the process is too much constrained and, as things stand, the target of uniform thickness along the bottle cannot be reached. In further study, it would be of interest to define the optimal final thickness repartition in a primary optimization loop, as illustrated in Fig. 21. The objective would be to minimize the weight of the part and keeping it stiff enough to resist to a given load. This optimization loop could be run using the mechanical code of ANSYS Workbench. The Downhill Simplex algorithm, a Nelder-Mead like algorithm presented by [16], could be a good a solution for this study as the the directions of correction are not as evident as in the blow molding process. The results of this first loop would be the goal of the second one.
Conclusion and perspectives The single stage injection blow molding process is an innovative process which purpose is to produce short series of specific hollow part on standard injection machine. Being a new process, the lack of data and background make it difficult to fully understand the parameters range involved. The polymer preform is blown right after being injected. Thus the initial temperature before this blow molding stage cannot completely be controlled, because it is a direct consequence of the injection molding of the preform occurring just a while ago. Thus this limits the available leverage on the process behavior. The optimization of such a process can then only be done on the initial geometry of the preform. A simple approach of the single stage injection blow molding process is proposed. The first step of the simulation, the cooling and injection molding one, should be replaced by a complete fill/cool/pack simulation of the injection stage. A more accurate temperature field could be then exported and used as input in the second stage of the presented simulation.
The use of a simple viscous fluid law coupled to a heat conduction model has proved to be relevant enough to reveal the sensitivity of the process to the material mechanical behavior, to the initial temperature field and to the imposed pressure. Comparison with experimental measurements of the thickness repartition has been made. The experimental measures are performed through the image analysis of the scanned contour of the blown parts. The model seems promising and is able to reproduce the tendencies observed through experimental measurement. With regard to the experimental measurements, the cutting operation of the bottle will have to be more accurate, to lessen the dispersion in the measurement of the thicknesses. Moreover, to have a more accurate median profile, it is necessary to increase the number of measurements. Although this method has its drawbacks, it is still much more reliable than point-wise caliper measurement and much more convenient than full tomography X. Once the experimental method is more robust, the comparison will be more pertinent. The current experimental measurement method is a destructive one due to the part having to be cut. A nondestructive method would be best. Tomography X could be a solution. Another one would be the use of the visible light. The analysis of the primary reflection on the outer side of the part, and the secondary reflection on the inner part of the local thickness would allow one to deduce the given thickness. This method is used by [3] to measure the thickness of the parison in the blow extrusion process. An optimization loop, using a Predictor/Corrector algorithm, shows good efficiency. The optimization script has proven to be robust. The optimized solution highlights constraints imposed by the process on the final thickness repartition of the blown product and that the optimized preform can mitigate them only partially. In further study, a more complete optimization will be run: the determination of the theoretical final thickness repartition in regard to the mechanical properties wanted for the final product would be determined in a first loop using the Downhill Simplex method. Once this final profile is determined, it would be set as the target of the second loop, consisting of an injection molding simulation for the injection stage and a blow molding simulation for the opening stage and the blow molding stage, and determine the optimal initial thickness repartition of the preform leading to the optimal thickness reparation of the final product using the Predictor/Corrector method, in regard to the realistic constraints such as tool tolerance for crafting the mold and manufacturing constraints. DSC trials at different cooling rates have been carried out. They showed that the faster the cooling the lower
Int J Mater Form
the crystallization temperature. These results are consistent with the process observations with the need for the material to remain hot enough to avoid the onset of crystallization which will lead to part defects such as cold drawing. In further work, a more complete rheological law would be used. A model coupling the molten behavior of the polymer and the effect of the onset of the crystallization would be more accurate. This would be compared to the Cross model in process conditions and would ascertain or refute the legitimacy of the use of the Cross law in the numerical model. Acknowledgments This work is a part of the IS2 project supported ´ a french agency for research and industrial development by OSEO (FUI F1111066V). It involves three companies : JP Grosfilley, Euro´ plastiques, Runner Plast together with INSA-Lyon. Funding by OSEO is gratefully acknowledged.
References 1. ANSYS (2010) Ansys polyflow user’s guide. release 13·0 2. Bellet M, Monasse B, Agassant JF (2002) Simulation num´erique des proc´ed´es de soufflage. Techniques de l’ing´enieur (in french) 3. B´ereaux Y, Charmeau JY, Balcaen J (2011) Optical measurement and modelling of parison sag and swell in blow moulding. Int J Mater Form:5 4. Bordival M, Schmidt FM, Maoult YL, Velay V (2009) Optimization of preform temperature distribution for the stretch-blow molding of pet bottles: Infrared heating and blowing modeling. Polym Eng Sci:49 5. Devaux N (2003) Influence d’un cisaillement sur les premiers stades de la cristallisation du polypropyl`ene. PhD thesis, Th`ese de ´ l’Ecole des Mines de Paris (in french) 6. Duffo P, Monasse B, Haudin JM, GSell C, Dahoun A (1995) Rheology of polypropylene in the solid state. J Mater Sci:30 7. Duval C (2004) Polypropyl`enes (pp). Techniques de l’ing´enieur (in french) 8. Fatt MSH, Ouyang X (2008) Three-dimensional constitutive equations for styrene butadiene rubber at high strain rates. Mech Mater:40
9. Groot JAWM, Giannopapa CG, Mattheij RMM (2010) A computer simulation model for the stretch blow moulding process of polymer containers. ASME:2010 10. Hopmann C, Michaeli W, Rasche S (2011) Fe-analysis of stretchblow moulded bottles using an integrative process simulation. The 14th International ESAFORM Conference on Material Forming 11. Huang HX, Yin ZS, Liu JH (2007) Visualization study and analysis on perform growth in polyethylene terephthalate stretch blow molding. Journal of Polymer Applied Science:103 12. Lee DK, Soh SK (1996) Prediction of optimal preform thickness distribution in blow molding. Polym Eng Sci:36 13. Luersen MA, Riche RL (2004) Globalized nelder-mead method for engineering optimization. Computer and Structures:82 14. Menary GH, Tan CW, Armstrong CG, Salomeia Y, Harkin-Jones EMA (2010) Validating injection stretch-blow molding simulation through free blow trials. Polymer Engineering and Science:50 15. Pham XT, Thibault F, Lim LT (2004) Modeling and simulation of stretch blow molding of polyethylene terephthalate. Polym Eng Sci:44 16. Press WH, Teukolsky SA, Vetterling WT, Flannery BP (1994) Numerical recipes in Fortran ? The Art of Scientific Computing, 2nd Edn. Cambrige University Press ISBN 0-51-43064-X 17. Schmidt FM, Agassant JF, Bellet M (1998) Experimental study and numerical simulation of the injection stretch/blow molding process. Polym Eng Sci:38 18. Thevenon A, Fulchiron R (2013) Elongational behavior of amorphous polymers in the vicinity and above the glass transition temperature. Polym Test:32 19. Thevenon A, Fulchiron R (2014) Evolution of poly(propylene) morphology in the rubbery state under uniaxial strain. Macromol Mater Eng:299 20. Thevenon A, Fulchiron R (2014) A thermomechanical modeling approach for the structural changes in semi-crystalline polymers under elongational strain. J Mater Sci:49 21. Thibault F, Malo A, Lanctot B, Diraddo R (2007) Preform shape and operating condition optimization for the stretch blow molding process. Polym Eng Sci:47 22. Wang S, Makinouchi A, Okamoto M, Kotaka T, Maeshima M, Ibe N, Nakagawa T (2000) Viscoplastic material modeling for the stretch blow molding simulation. Int Polym Process XV 23. Yang LM, Shim VPW, Lim CT (2000) A visco-hyperelastic approach to modeling the constitutive behavior of rubber. International Journal of Impact Engineering:24 24. Yang ZJ, Harkin-Jones E, Menary GH, Armstrong CG (2004) A non-isothermal finite element model for injection stretch-blow molding of pet bottles with parametric studies. Polym Eng Sci:44
