Finite element simulation of stretch forming of aluminium-polymer laminate foils used for pharmaceutical packaging

Pharmaceutical high barrier blister packages are manufactured from aluminium-polymer laminate foils (e.g. consisting of PA-Al-PVC layers). By a cold stretch forming process cavities are formed. The aim of this work is to determine a homogenized elastic-plastic description of the laminate by micromechanics. Therefore, a microstructural model is developed where the layers are mapped in a representative volume element. The obtained homogenized material model is applied to simulate the stretch forming to gain more insight into the forming process.


Introduction
Blister packs are often used for the packaging of pharmaceutical tablets.Currently, there are two main types of pharmaceutical blister packages on the market.According to the process of manufacturing they can be divided into thermoformed and cold-formed packages.The thermoformed blister packs are produced by thermoforming of a thermoplastic polymer foil into cavity moulds.The other type of package is manufactured by cold stretching of aluminium-polymer laminate foils.The aluminium in this composite material serves as barrier layer.Consequently, the packed tablet is protected against moisture, light, oxygen and other gases.However, the composite foil possesses less formability which leads to bigger cavities than needed for the product being packed.Figure 1 shows the difference in size of thermoformed and cold-formed blisters.
The companies developing forming machines for cold forming blister packs often lack knowledge about the exact specification of foils their costumers will use.Since they have limited information about the deformation behavior of aluminium-polymer laminate foils, they have to design the punches with high safety factors.Because of the endless variety of tablet forms the definition of the optimal punch geometry is an ongoing problem and one of the critical development stages for the developer of cold-forming blister packaging machines.
A direct approach for evaluating the suitability of a punch geometry for forming aluminium-polymer laminate foil is to test the punch on the packaging machine.Actually, this is the current procedure.However, this approach is only possible when the material the customer will use is a e-mail: simon.mueller@hs-karlsruhe.deb e-mail: sabine.weygand@hs-karlsruhe.deknown and available.In addition, this approach requires a lot of time.
A numerical approach through the finite element (FE) method is a possible alternative way to optimize the punch geometry.Such a method allows to design the punches in an early stage of development and no machine time, for setting up the tools and testing, is wasted.One challenge in such a FE-modeling is the description of the mechanical behavior of aluminium-polymer laminate foils, since due to their layered structure their forming behavior and failure behavior is not yet well understood [1].Furthermore, representing all different layers and their interfaces between them in the FE model leads to high computational time when simulating the stretch forming process.Homogenization techniques that help in reducing the computer time would be useful [2].
This work is a first attempt to describe the cold stretch forming behavior of aluminium-polymer laminate foil through the finite element method.As reference material an aluminium-polymer laminate foil, consisting of an outer layer of orientated polyamide (PA), a core layer of aluminium (Al), an inner layer of polyvinyl chloride (PVC) and adhesive and primer layers between is used.This combination of materials is common and finds great use in the pharmaceutical packaging industry.For the numerical investigations a micromechanical model is developed where the layers are mapped with simple material properties.Based on calculations with this model two homogenized material models are presented for the foil.The stretch forming process is simulated with a geometry close to an actually used one.There, three cases are studied.In the first one the foil is subdivided into a layer-by-layer structure.Each layer is assigned with simple material properties like in the micromechanical model.In the second and the third case the material behavior of the foil is described with homogenized behavior.The three cases with different material descriptions are compared by hypotheses regarding the reasons for failure.

Material models for the individual layers
In our numerical study each layer is described as a linear elastic, rate independent plastic material with isotropic hardening.To identify the material parameters of each layer uniaxial tensile tests were carried out on the individual layers.Thereby, nominal stress-strain curves were measured for PA, Al and PVC and converted into true stressstrain curves.Because this translation is only possible up to the tensile strength the curves are tangentially linearly extended above this point up to a true strain ε tr = 1.The obtained flow curves are shown in Figure 2 and they provide the starting point for the FE simulations.It has to be noted that PVC shows strain softening right after yield initiation.That means that after the elastic deformation the nominal stress-strain curve includes a stress overshoot [3].In this work, softening is neglected.This means that the stress overshoot of the measured nominal stress-strain curve is cut off.Because of that, we will denote the material model for the modified PVC with an asterisk.Additionally, the values for the Young's modulus E for each layer were determined from the tensile tests.For the Poisson's ratio ν of aluminium the well-known value is used [4].For the Poisson's ratio of PVC and PA assumptions are made.Table 1 summarizes the applied elastic constants of the three materials.  .True stress-strain curves of the layer materials of the aluminium-polymer laminate.These three curves were derived from nominal stress-strain curves and represent the base for the numerical calculations.Because of simplifications regarding the mechanical behavior of PVC it is denoted with an asterisk.

Numerical modeling
To obtain the overall behavior of the aluminium-polymer laminate foil a micromechanical model has been developed.In this so called representative volume element (RVE) all three single layers are mapped individually, and the previously described material models are assigned to them.This model allows to study the influence of the layered structure on the elastic and plastic behavior of the foil.Furthermore, it enables to apply load cases, which are not accessible in experiments, due to the fact that the composite material is a thin foil.
In addition, a FE stretch forming model close to real process was designed.In the first simulation (called reference case) the single layers are again represented with their individual material behavior.Then the simulations are repeated with homogenized material models obtained from the studies on the representative volume element.For the modeling and the following simulations the commercial FE software Abaqus is used.

Representative volume element
The developed microstructural model is based on the representative volume element technique.The RVE employed in this work is shown in Figure 3.It was build up as a cube according to the geometry given in the figure.As can be seen, the model represents the laminar structure with its different layers (PA-Al-PVC).The material models presented in the previous section were assigned according to the labels.It is assumed that the individual layers are perfectly bounded and that the adhesives and primer layer thicknesses are negligible.Therefore, they are not included into the model.The foil can be described as a periodic arrangement of a set of RVEs in the (x, y)-plane.To ensure this behavior under load, periodic boundary conditions have to be applied via the following constraints.The difference of the displacement vectors in two corresponding points situated at opposite volume element boundaries is a constant displacement vector corresponding pairs of points.If one looks at Figure 3, these constraints can be expressed by the relations u(x, y, z) is the displacement vector of a point with the coordinates x, y and z.The vectors k 1 , k 2 and k 3 are constant displacement vectors which may be predefined or follow from the global equilibrium conditions of the RVE.
Instead of using displacement boundary conditions for hindering the rigid body motion, spring elements with a very high stiffness are employed.Here, the benefit lies in a more homogeneous stress distribution in the model than using displacement boundary conditions.The RVE was meshed with first-order, full-integrated, hexahedral elements.The resulting model has a mesh consisting of 21952 elements with the same dimensions.

Stretch forming model
The stretch forming model is designed according to a tool for industrial use.It consists of one punch, a die plate and the aluminium-polymer laminate foil.The foil has a length of 32 mm, a width of 26 mm and again a thickness of 0.14 mm.The dimensions of the single layers are the same and will be included into the model, too.The adhesive and primer layers are neglected.The dimensions of the other employed parts are shown in Figure 4. Due to geometrical symmetry only a quarter of the tool and the foil have to be realized.Furthermore, the punch and the die plate are assumed to be rigid and are meshed with threedimensional bilinear rigid quadrilateral elements to reduce the computation time.
The foil is meshed with first-order, full-integrated, hexahedral elements.Thereby, it is ensured that over the thickness of each layer five elements are generated.The foil is placed with the PA layer facing the die plate as can be seen in Figure 4. To establish the contact between the foil and the other parts, a hard surface to surface contact is used.Frictionless behavior between the different parts is assumed.
The applied boundary conditions are such that all translational and rotational degrees of freedom of the die plate are blocked.The punch is only allowed to do translational movements in z-direction.Thus, the loading is applied there.Because only a quarter model of the tool was created, symmetric boundary conditions have to be defined on the foil.In addition, the edges of the laminate are fixed in all directions.

Finite element analyses and results
To obtain the overall mechanical behavior of the aluminium-polymer laminate foil and to gain insight into the internal stress variation, the RVE was applied to simulate various different load cases.This was done by defining force vectors F i on reference nodes, which are coupled to the displacement vectors k i .
In the first set of calculations the RVE was applied to simulate uniaxial loading in three different directions, namely in x-, yand z-direction.As can be seen in Figure 3, the x and y-directions lie in the foil plane while the z-direction corresponds to the transverse direction.The resulting three nominal stress-strain curves in tension are shown in Figure 5.It can be clearly seen, that loads in xand y-directions lead to identical curves.In contrast, a load in z-direction shows a slightly different behavior.Here, the low value for the Young's modulus is noticeable.Another point is the elastic limit, which is approximately 21 MPa when loaded in xand y-directions and 38 MPa when loaded in z-direction.After a strain of 0.02, the nominal stress-strain behavior is almost similar.The hardening only differs in a few units between the load cases.In the second set of calculations the RVE was applied to simple shear.The results are shown in Figure 6.Shearing on the x-plane in z-direction and shearing on the y-plane in z-direction shows an identical mechanical behavior.In contrast, shear load on the x-plane in y-direction shows a different trend.In this load case the shear modulus is higher than in the other two cases.The work hardening behavior differs, too.
The computed stress-strain curves together with the geometrical changes transverse to the loading direction are used to determine the elastic moduli of the composite foil.The whole 12 material constants are summarized in Table 2. E i stands for the Young's modulus along axis i. G i j is the shear modulus in direction j on the plane whose normal is in direction i. ν i j is the Poisson's ratio that corresponds to a contraction in direction j when an extension is applied in direction i.It is noticeable that a lot of these constants are similar.As the xand y-directions lie in the foil plane, their corresponding Young's moduli E 1 and E 2 are identical.Its value roughly corresponds to the arithmetical average of the individual elastic moduli of the layers (upper bound).Is the foil loaded transverse to its plane (i.e. in z-direction), the Young's modulus E 3 is smaller as the layers are loaded in row (lower bound).
As our investigated foil is symmetric with respect to a rotation about the z-axis, its overall elastic behavior is transversely isotropic.Therefore, the number of indepen- dent constants in the elasticity tensor is reduced to 5.These will be used in the following stretch forming simulations.The simulations of the stretch forming process have been grouped in three different cases.These cases differ in the material description of the foil.In the first case the laminar structure of the foil was maintained.Each layer was assigned to its individual material properties given in Figure 2 as done in the calculations with the RVE.As in this case the layered structure is represented, it is closest to reality.That's why it is considered as reference case.Two more cases are investigated, where the laminar structure is not taken into consideration.Here the material behavior of the foil is described with homogenized elastic-plastic material models.
Table 2. Computed elastic moduli of the foil.The indexes 1, 2 and 3 correspond to the x-, y-, z-direction defined in Figure 3.In the second case the material behavior of the foil is described with linear elasticity with isotropic hardening, rate independent Mises plasticity.To derived the true stress-strain behavior the nominal stress-strain curve from the uniaxial loading of the RVE in x-direction is taken and transferred to the true stress-plastic strain curve (hardening curve) shown in Figure 7.The plastic strain ε pl was calculated subtracting the elastic strain by using the Young's modulus E 1 .To describe the transverse contraction the Poisson's ratio ν 12 is used.
In the third case the elastic material behavior of the foil is defined by the 5 parameters of the transversely isotropic elasticity, which had been determined from the calculations on the RVE.The plastic material behavior is specified by the true stress-plastic strain curve in Figure 7.
All stretch forming simulations started with the punch touching the foil.Then the load was applied through a translational punch displacement of 7.14 mm in negative zdirection.The detailed analysis of the stress and strain distribution is performed at maximum loading (i.e. at punch displacement of 7.14 mm).At the moment, the crucial failure mechanism in stretch forming aluminium-polymer laminate foils is not yet known.To obtain a better understanding of possibly involved mechanisms some hypotheses regarding the reasons for failure are made.These hypotheses will be evaluated with the results from experimental investigations in a future work.
The first hypothesis is that a crack occurs at the location of the maximal principal stress.The interest does not only lie in the maximum value, but also in its location and whether the more simple material descriptions in case 2 and case 3 predict the same location.The outcomes of the simulations are shown in Figure 8.The two cases with the homogenized material behavior show nearly the same maximal principal stress at the same location.This was expected since in both cases the hardening is described by the same curve.A closer examination shows that the maximal value occurs in the area where the punch has an edge.In the second case the maximal principal stress is approximately 241 MPa.In the case with the transversely isotropic elasticity the value is greater by 270 MPa.Both stresses appear on the surface of the foil.In the first case, where the foil is described by a multilayer model, the maximal principal stress occurred with 325 MPa in the area where the punch has an edge.However, a detailed investigation into the laminate structure shows that this stress doesn't occur on the surface of the foil.Here the maximum stress value appears in the aluminium core layer.The second hypothesis is the appearance of failure due to a critical thickness.Is the Al layer too thin, permeation can occur which is equal to failure.Therefore, the strain in thickness direction is analyzed.The results of the simulations are visible in Figure 9.The given strain component in the legend refers to the material orientation.It was chosen such that the shown component is equal to the strain in thickness direction.Visible is the fact that in all three cases high negative strains occur.Furthermore, they appear at nearly the same location.In the first case, which is the model with the laminate structure, a minimal strain value of -0.911 is reached.The location of appearance is at the tran- sition from the flat area at the bottom of the cavity to the curvature.A detailed analysis of the different layers in this area shows, that nearly the same strain occurs in every material layer.This was expected since the foil has to deform with the displacement movement of the punch.By contrast, the stresses are different in every layer.In the second case, where linear elasticity and rate independent plastic behavior with isotropic hardening was defined, a minimal strain of -0.927 is reached.In the case with the transversely isotropic elasticity the lowest strains in thickness direction occur with a value of -0.944.
The third hypothesis refers to the forming limit of the foil.Therefore, the maximal principal strains are evaluated.The outcome of the simulation is shown in Figure 10.Here, the location of the highest value is not that clearly as for the other hypotheses.In the first case the maximal value of strain reaches 0.536 and is located at the area, where the punch has an edge.But on the bottom of the cavity and near the edge of the die plate high strains occur, too.The second case has a maximal principal strain of 0.512.The location is the same as in the first case.Again high strains on the bottom of the cavity and near the edge of the die plate occur.In the third case the maximal principal strain of 0.515 appears on the bottom of the cavity near the radius.At the edge of the punch and at the edge of the die plate are other high values, too.

Conclusion
For numerical investigations of the deformation behavior of aluminium-polymer laminate foils a micromechanical FE model was developed.This model is based on the repre-sentative volume element technique.It consists of a layered structure with periodic boundary conditions.Different linear elastic, rate independent plastic material behaviors with isotropic hardening have been assigned to the single layers of the RVE.
By simulations with the RVE nominal stress-strain curves and shear stress-shear strain curves were gained.These curves give an insight into the anisotropy that results only through the laminate structure.With the help of these curves the 12 parameters of the orthotropic linear elasticity have been determined.The values of these parameters show that the elastic behavior of the RVE can be described with the special case of transversely isotropic elasticity.Furthermore, the RVE was used to develop two simple homogenized material models on the basis of the composite behavior under tension load and the anisotropic linear elasticity.
In addition a stretch forming model close to the real cold stretch forming process of blister packaging was designed.That model offers the possibility to define a foil with single layers with different material models.Thereby, three simulations of the forming process were done.In the first case the three layer of the foil were assigned with the same material properties used in the simulations with the RVE.For the second case a homogenized material behavior was used which corresponds to a linear elastic, rate independent plastic material behavior with isotropic hardening.In the third case the linear elasticity of the homogenized material model was changed to transversely isotropic elasticity with determined parameters from the calculations with the RVE.
Three hypothesis regarding the reason for failure while stretch forming the foil are made.The first hypothesis is that a crack may nucleate at the location of the maximal principal stress.All simulations showed that the location of the maximal principal stress is where the punch has an edge.In the second case and in third case these stresses appear at the surface of the foil.In the first case it occurs in the aluminium core layer of the laminate foil.The stress values differ in every case.The second hypothesis is that a critical thickness reduction is responsible for failure.Regarding this the strain component in thickness direction has been evaluated.The results of the simulations show that high negative strains occur between the bottom of the cavity and the curvature.It is worth noting that in the first case in every point of the foil the different layers underlie the same strains.The third hypothesis refers to failure of the foil at the area of the maximal principal strain.The outcome of the simulations shows that the maximal values are similar.But the locations differ between the cases.To mention is that in every calculation three areas with high maximal principal strains occur.One is the transition from the bottom of the cavity to the curvature.The other are at the edge of the punch and at the edge of the die plate.

Figure 1 .
Figure 1.A thermoformed blister is shown in a) and a cold stretch formed blister is displayed in b).The difference of the dimensions of the cavities is clearly visible.

Figure 2
Figure2.True stress-strain curves of the layer materials of the aluminium-polymer laminate.These three curves were derived from nominal stress-strain curves and represent the base for the numerical calculations.Because of simplifications regarding the mechanical behavior of PVC it is denoted with an asterisk.

Figure 3 .
Figure 3. Model of the representative volume element (RVE).Depicted are the three layers of the aluminium-polymer laminate with their dimensions given in mm.The spring elements hinder the rigid body motion.

Figure 4 .
Figure 4. Side view and top view of the stretch forming model.The die plate, the punch and the laminate foil are shown.Their boundary conditions are indicated in the side view.The dimensions of the die plate and the punch are given in mm.

Figure 5 .
Figure 5. Nominal stress-strain curves of tensile loads on the RVE.The diagram shows that the mechanical behavior under load in xand ydirection shows similar behavior.The small diagram shows the stress-strain behavior at small strains.

Figure 6 .
Figure 6.Nominal shear stress-strain of shear loading on the RVE.The diagram shows that shearing on the x-plane in zdirection and shearing on the y-plane in z-direction leads to the same mechanical behavior.

Figure 7 .
Figure7.Hardening curve determined from the nominal stressstrain curve which was calculated by tensile loading of the RVE in x-direction.For translation the Young's modulus E 1 was used.This curve describes the plastic material behavior in the second and the third case of the stretch forming simulations.

Figure 8 .
Figure8.Computed maximal principal stress in the simulation of the stretch forming process at highest load.The first case, where the foil is described by a laminate structure, is shown in a).The second case, where the foil was defined with linear elasticity and rate independent plastic material behavior with isotropic hardening, is visible in b).The third case is shown in c).Here, the linear elasticity was replaced with transversely isotropic elasticity.The marks show the location of the maximal principal stresses.

Figure 9 .
Figure 9. Computed minimal strain in thickness direction in the stretch forming simulations.The first case is shown in a), the second in b) and the third case in c).The marks show the location of the minimal strain in thickness direction.

Figure 10 .
Figure 10.Computed maximal principal strain in the simulation of the stretch forming process at highest load.The first case is shown in a), the second in b) and the third case in c).The numbers with the boxes around indicate the highest values.

Table 1 .
Linear elastic properties of the three material models