Numerical Modelling of the Thermoforming Behavior of Thermoplastic Honeycomb Composite Sandwich Laminates

1. Introduction

Thermoplastic composite sandwich structure made of honeycomb core and laminated with UD-tape cross-ply cover layers on the top and bottom sides offers reduced mass and high bending stiffness. The honeycomb structure comprises an array of thin-walled hollow hexagonal cells arranged together, resulting in minimal material per unit area. In addition to that, the honeycomb structure provides relatively high out-of-plane compression and out-of-plane shear properties. Meanwhile, cross-ply as cover layers offer a high tensile strength in the fibre direction, resulting in-plane orthotropic properties.

The production of hexagonal honeycombs at an industrial scale is carried out through tube extrusion, expansion, and corrugated processes [1,2,3]. These processes are operated at multiple stages and require long production times, which in total results in high production costs and therefore, limits them to a low volume production [4]. To solve this, ThermHex Waben GmbH has introduced a continuous production process, where a thin thermoplastic polymer film is extruded (1^st Step) and vacuum formed to form a half hexagon (2^nd Step), followed by folding to form a honeycomb (3^rd Step), at the end lamination with cover layers (4^th Step) resulting in composite sandwich structure as shown in Figure 1 (a). The so-produced semi-finished laminate is called as Organosandwich, which is aimed for mass production with maximum production speeds of up to 20 m/min with a panel width of 1.2 m [5].

With the added advantage of thermoplastic material behavior and the suitability for high volume production, the flat semi-finished panels can be further processed into complex shapes using the thermoforming technique so that they can be used in different applications, especially in automotive sector. To form and include functionalizing elements into sandwich components, a multi-stage thermoforming technique known as thermoplastic Sandwich Molding (TSM) is developed at Fraunhofer IMWS, with manufacturing cycle time of 4 s [6]. Simultaneously, with the process technique approach, Fraunhofer IMWS develops the numerical models and perform simulations to adapt the novel technology to the digital process chain.

2. Thermoplastic Sandwich Forming

The TSM technology allows processing of semi-finished thermoplastic sandwich panels into ready-to-use lightweight parts. Focus of this research is to develop a manufacturing process that allows the honeycomb core sandwich laminate to confront the mold shape without wrinkling and core collapsing or buckling. The novel TS-molding manufacturing process of FRTP sandwich components comprises of three main process steps as shown the Figure 1 (b), starting by IR-heating of the semi-finished sandwich laminate, followed by forming of the sandwich laminate into pre-defined shape and then injection molding around the edges [6].

Figure 1. Organosandwich (a) continuous production of honeycomb core and in-line lamination of coverlayers, (b) schematic illustration of the novel manufacturing process of TS-Molding [3,4,23,24].

Figure 1. Organosandwich (a) continuous production of honeycomb core and in-line lamination of coverlayers, (b) schematic illustration of the novel manufacturing process of TS-Molding [3,4,23,24].

Specifically, thermoforming is carried out through a press by closing the two-sided tool in three successive sub-process steps. In the first step of the closing movement, the sandwich structure is formed in areas of the later part in which the core must remain intact. After this forming step, the press stops the closing movement for a short time allowing cooling by tool contact and thereby solidification of the face sheet matrix. After this short break, the tool moves further and fully closes by pressing the sandwich edges into a compact laminate. Pressing is particularly important at the component edge for sealing of the sandwich to protect the core and later for integration of functional elements. The pressed edge also acts as a transition zone between the formed sandwich structure and over-molded edges [4].

3. Materials and Finite Element Modelling

Using mathematical models, physics laws, and computer science, it is essential to reduce the physical work and save the extensive time that is required to conduct numerous experiments through a trial-and-error approach. The finite element analysis (FEA) approach starts with material description and material characterization, followed by FE-modelling with necessary boundary conditions and loading cases. Simulation results are used to identify undesirable deformations in part at the early stages of thermoforming process development, as shown in Figure 2.

3.1. Sandwich Material Description

The Organosandwich used for this study comprised a PP honeycomb core and PP/GF UD-Tape-based Cross-Ply cover layers. During the production of such laminates, a thin PP film is formed to a honeycomb core and laminated with the PP/GF cross-ply using heating and double band press. This lamination allows a strong weld joint between the core and cover layers as they share the same thermoplastic material. There is no requirement for additional adhesive layers, unlike in Nomex or another type of sandwich laminate, making the production more sustainable and economical.

Table 1. Sandwich material description of the used configuration 12THPP120CP820.

Sandwich Entity	Dimension
Sandwich thickness	12 mm
Cover layer thickness	0.5 mm
Honeycomb core thickness	11 mm
Cell size	5 mm
Standard panel dimension	1200 mm x 2500 mm (L x W)
Weight per unit density	3120 – 3240 g/m²

ThermHex Waben GmbH [23].

Table 1. Sandwich material description of the used configuration 12THPP120CP820.

Sandwich Entity	Dimension
Sandwich thickness	12 mm
Cover layer thickness	0.5 mm
Honeycomb core thickness	11 mm
Cell size	5 mm
Standard panel dimension	1200 mm x 2500 mm (L x W)
Weight per unit density	3120 – 3240 g/m²

ThermHex Waben GmbH [23].

Table 2. Mechanical properties at room temperature of the PP core and PP/GF UD-Tape.

Material	Density (kg/m³)	Young’s Modulus (GPa)	Yield stress (GPa)	Poisson’s Ratio	Specific heat capacity (J/kg K)	Thermal Conductivity (W/m K)
PP	980	1.55	0.2	0.41	1794	0.161
PP/GF 60	1540	29.39	0.4	0.40	1160	0.308

Table 2. Mechanical properties at room temperature of the PP core and PP/GF UD-Tape.

Material	Density (kg/m³)	Young’s Modulus (GPa)	Yield stress (GPa)	Poisson’s Ratio	Specific heat capacity (J/kg K)	Thermal Conductivity (W/m K)
PP	980	1.55	0.2	0.41	1794	0.161
PP/GF 60	1540	29.39	0.4	0.40	1160	0.308

3.2. Coverlayers Wrinkle Analysis

The thermal and mechanical loads on sandwich components result in different deformation modes, on the one hand, cover layers exhibit following behavior,

• ply-elongation
• inter-ply slip
• ply shear and
• wrinkle or in-plane buckling [7,8].

On the other hand, the core component exhibits in-plane tension and compression, and out-of-plane buckling due to compression [2]. In thermoforming process development, the focus is to produce a part without any defects. Sometimes, CFRTPs such as woven laminate or non-woven laminate (cross-ply) during thermoforming show surface defects, which usually develop due to buckling or wrinkling. The wrinkle formation in cover layers in non-homogeneous sandwich materials is treated as a local phenomenon. This kind of deformation can be initially estimated through a model of axial loaded flat sandwich structure (see Figure 3), the prediction of onset wrinkle in cover layers is done with calculation of the maximum compressive stresses in the cover layer.

There exist different mathematic models to estimate the maximum compressive stress for initiation of wrinkle formation in the cover layers, based on the homogeneous and non-homogeneous core materials. For non-homogeneous core, isotropic cover layers, and where height of the core is relatively high compared to that of thickness of face sheet (see Figure 1 (a)), Williams [9] made a suitable analysis. This analysis is based on the initial assumption that the axial stress in the direction of load on the core is zero, known as anti-plane stress assumption [8]. The core here is infinitely thick therefore, wrinkles are exponentially decayed on the cover layers away from the core axis. The mathematical expression of wrinkle formation stress for cover layers is derived as,

$${\sigma }_{wr}=\left[\frac{0.825}{{\left(1-{{\nu }_{cl}}^2\right)}^{1/3}}\right]{\left(E_{cl}{E}_c{G}_c\right)}^{1/3},$$

(1)

where ${\sigma }_{wr}$ is the maximum compressive stress for wrinkle formation, ${\nu }_{cl}$ is the Poisson’s ration of the cover layer, $E_{cl}$, $E_c$ are the Young’s modulus of cover layer and core respectively, $G_c$ is the in plane shear modulus of the core. An additional assumption is made in the current study, based on the thermoplastic sandwich lamination, where additional molten polymer is deposited in the core and cover layer bonding regions. Therefore, orthotropic cover layers are considered as locally isotropic in the contact region between core and cover layer, calculating the Equation 1 for Organosandwich at processing temperature, the maximum in-plane compressive stress which result in wrinkle formation is estimated as below,

$${\sigma }_{wr}=\left[\frac{0.825}{{\left(1-{0.4}^2\right)}^{1/3}}\right]{\left(16.27\ast 130\ast 0.09\right)}^{1/3}=5.03N/mm².$$

(2)

3.3. Spatial Discretization of Organosandwich

TS-molding is a combined process of thermoforming and over-molding, however only thermoforming is discussed in this paper. The demonstrator part used for this study comprises different zones like forming and pressing regions. Therefore, the sandwich is divided into three parts: a honeycomb core, top and bottom face sheets, where each part is sub-modelled using shell element theory. The honeycomb structure is generally modelled using three methods: (1) meso model, (2) volume model, and (3) shell model. Bounding to the manufacturing process, where the honeycomb structure significantly influences the forming locally and globally, the meso model approach is chosen in this study. Using mesh discretization, a complex geometry is divided into smaller entities and simpler subdomains. Then, their respective partial differential equations are solved by applying the principles of boundary value problem. During the selection of element types and their related interpolation functions of the physical quantity for e.g., displacement, temperature, etc., is considered.

In finite element method approach the strong form of the initial boundary value problem for small displacements is written as follows,

$$\nabla ·\mathsf{\sigma }+\mathrm{b}=\rho \ddot{u} \mathrm{in} \Omega ,$$

(3)

where $\rho$ is the density, $\ddot{u}$ is the time derivative of displacement, with the auxiliary conditions like, essential boundary conditions $u_i={\overline{u}}_{i}on{\mathsf{\Gamma }}_{u_i}$ and natural boundary conditions ${\sigma }_{ij}{n}_j={\overline{t}}_{i}on{\mathsf{\Gamma }}_{f_i}$. The weak form is then written as

$$\underset{\Omega }{\overset{}{\int }}\rho \delta u·\ddot{u}d\Omega +\underset{\Omega }{\overset{}{\int }}\delta \epsilon :\mathsf{\sigma }d\Omega =\sum _{i=1}^n\underset{{\Gamma }_f}{\overset{}{\int }}{\delta u}_i{\overline{f}}_{i}d\Gamma +\underset{\Omega }{\overset{}{\int }}{\delta u}_{}·bd\Omega ,$$

(4)

which is used to develop the discrete equations in the Galerkin procedure. The displacement approximation of standard FEM method is defined through shape functions, for a displacement $u_i$ of an isoparametric element, the interpolation function is (Hughes 1987),

𝑢𝑖ℎX,t=𝐼=1𝑛𝑛𝑒𝑁𝐼𝜉𝑋𝑢𝑖𝑗𝑡, X∈Ω,

(5)

where $n_e$ is the number of nodes of elements, for a standard four node element the shape functions $N_I$ are written as,

$$N_I\left(\xi \right)=\frac{1}{4}\left(1+{\xi }_I\xi \right)\left(1+{\eta }_I\eta \right),$$

(6)

here $\xi \left(X\right)$, $\eta \left(X\right)$ are the expressions of mapping of physical domain to its parent domain. [10]

3.4. Thermo-Mechanical Modelling

Thermoforming is a combined problem of mechanical analysis and thermal analysis, where stresses generation due to mechanical deformation as well as temperature variation are vital. In thermo-mechanical coupled finite element model, accurate evaluation of mechanical stresses and thermal stresses at each point of the physical subdomain must be done under the influence temperature. The heat conduction in the structure is governed by Fourier’s law of conduction [18], straining in the structure generates heat in transient and dynamic processes [11]. Thermo-mechanical material model assumed here a three-dimensional non-linear, transient heat conduction problem for a domain $\mathrm{\Omega }$ bounded by $\Gamma$, the governing differential equation written as,

$$\frac{\partial }{\partial x}\left(k_x\frac{\partial T}{\partial x}\right)+\frac{\partial }{\partial y}\left(k_y\frac{\partial T}{\partial y}\right)+\frac{\partial }{\partial z}\left(k_z\frac{\partial T}{\partial z}\right)+\mathrm{Q}\left(x,y,z,t\right)=\rho c\frac{\partial T\left(x,y,z,t\right)}{\partial t},$$

(7)

where $T\left(x,y,z\right)$ is the temperature at point and for time $t$, $\mathrm{Q}\left(x,y,z,t\right)$ is the internal strength of heat source, $\rho$ is the density, $c$ is the specific heat and $k_x,k_y,k_z$ are thermal conductivities in the $x,y,z$ directions respectively, at initial conditions $t=0$, temperature is $T_0$ [12].

Standard mesh FE model has an inherent, overly stiff phenomenon that causes poor accuracy in solving complex problems, to address this Liu [19] applied the smoothing gradient technique and established the smoothed Galerkin weak form. Further, a node based smoothing finite element method was introduced to construct the smoothing domain of transient thermal mechanical analysis to obtain better results than standard FEM [12]. The balance equations for the coupled problem are composed of force field and thermal field balance equation [13]. By substituting the thermal shape functions of element nodes in the Galerkin weak form for transient heat conduction, the discretized system of equilibrium equations is expressed as following matrix form,

$${\mathit{N}}_T\left\{\frac{\partial \mathit{T}}{\partial t}\right\}+\left[{\mathit{K}}_T+{\mathit{K}}_T^h\right]\left\{\mathit{T}\right\}=\left\{\mathit{P}\right\},$$

(8)

where $N_T$ is the thermal capacitance matrix, $K_T$ is the conduction matrix, $P$ is the thermal load vector of the predefined thermal shape functions. After defining the initial temperature field, the respective values at an arbitrary time $t$ can be approximated by applying the backward difference technique to the Equation 8.

The governing equations for thermo-mechanical problem over the specific domain $\mathrm{\Omega }$ bounded by $\Gamma$ are written as,

$$\mathrm{Equilibrium}\mathrm{equation}:{\sigma }_{ij,j}+b_i=0,$$

(9)

$$\mathrm{Displacement}\mathrm{boundary}:u_i=u_{\Gamma }\mathrm{on}{\Gamma }_u,$$

(10)

$$\mathrm{Stress}\mathrm{boundary}:{\sigma }_{ij}{n}_i=f_{\Gamma }\mathrm{on}{\Gamma }_f,$$

(11)

$$\mathrm{Strain}\mathrm{displacement}\mathrm{relations}:{\epsilon }_{ij}=\frac{1}{2}\left(u_{i,j}+u_{j,i}\right),$$

(12)

where ${\sigma }_{}$, $\epsilon$ and $u$ are the stress, strain, and displacement, $b$ is the body force, $u_{\Gamma }$ and $f_{\Gamma }$ are the displacement and traction on the essential and natural boundaries. The above stress and strain are an expression of thermal expansion coefficient,

$${\sigma }_{ij}={\delta }_{ij}\lambda {\epsilon }_{kk}+2\mu {\epsilon }_{ij}-{\delta }_{ij}\left(3\lambda +2\mu \right)\alpha ∆T\left(t\right),$$

(13)

in the above equation $\lambda$ and $\mu$ are Lame’s constants, they can be derived from the Poisson’s ratio and Young’s modulus, ${\delta }_{ij}$ is the Kronecker delta, $\alpha$ is the thermal expansion coefficient, and $∆T\left(t\right)$ is the change in temperature with the time $t$. In the Equation 10, ${\epsilon }_{ij}$ is the total strain, which is sum of elastic, plastic, and thermal strain,

$${\epsilon }_{ij}={\epsilon }_{ij}^e+{\epsilon }_{ij}^p+{\epsilon }_{ij}^T.$$

(14)

With the temperature field obtained, the thermal stress and strain analysis can be performed using the standard Galerkin weak form,

$$\underset{\Omega }{\overset{}{\int }}{\left\{\delta \epsilon \left(u\right)\right\}}^T\left\{\mathsf{\sigma }\right\}d\Omega -\underset{\Omega }{\overset{}{\int }}{{\left\{\delta u\right\}}^T\left\{b\right\}}_{}d\Omega -\underset{{\Gamma }_f}{\overset{}{\int }}{\left\{\delta u\right\}}^T\left\{f_{\Gamma }\right\}d\Gamma =0,$$

(15)

where $u$ is the displacement field, which is interpolated using linear shape functions as below,

$$u=\sum _{i=1}^n{\mathsf{\Phi }}_i{d}_i,$$

(16)

$d_i$ is the nodal displacement of the node $i$, and ${\mathsf{\Phi }}_i$ is the linear shape function. By solving the Equation 15 using Equation 16, the discretized system of equations for a thermo-mechanical problem can be expressed in the form of matrix as follow,

$${\mathit{K}}_d\left\{\mathit{d}\right\}=\left\{\mathit{F}\right\},$$

(17)

where ${\mathit{K}}_d$ is the stiffness matrix comprised of product of elastic coefficient matrix and strain displacement matrix, $\mathit{F}$ is the load vector. For time-dependent non-linearity of materials the Equation 11 is solved by using Newton-Raphson iteration method at each time step.

3.5. Thermoforming Model Development

The thermoforming model developed here is aimed to describe all three, geometric, material, and boundary non-linearities of thermoforming process of sandwich laminate. For this a simplified geometry model is developed considering the actual geometry. A unit cell representation of the sandwich and honeycomb core is shown in the Figure 5 (a). thereby, L is the length orientation, W is the width orientation, and H is the sandwich thickness and out-of-plane direction. The shell element formulation for cell walls as well as face sheet is done using Reissner-Mindlin kinematic assumption. The thermo-mechanical properties like temperature dependent E-modulus from DMA tests and temperature dependent Stress-Strain behavior from tensile tests (DIN EN ISO 527-1) are presented in Figure 4 (a) and (b) respectively. The thermo-physical properties like temperature specific heat capacitance and thermal conductivity are obtained through respective testing methods. Further, validated test results are presented by Minupala [4]. These temperature dependent material properties are taken as the input for the respective material models.

As there exists two different materials, face sheet is modelled using MAT_Reinforced_Thermoplastic, and core is modelled using MAT_Elastic-Viscoplastic_Thermal, both the materials are combined with the MAT_Thermal_Isotropic to form a thermoforming model in LS-DYNA interface.

Figure 5. (a) Representative idealized sandwich honeycomb cell (b) Initial thermal state of the sandwich after heating step in thermoforming process.

Figure 5. (a) Representative idealized sandwich honeycomb cell (b) Initial thermal state of the sandwich after heating step in thermoforming process.

The thermal state of the sandwich laminate is achieved through IR-heating process, where the sandwich laminate is heated from top and bottom sides so there is a symmetrical thermal state through the thickness. Numerous experiments were conducted by controlling intensity of the heat and heating time, Glaesser explained this technique in detail in his thesis work [6]. The finalized temperature profile of the sandwich, where face sheet is brought to matrix melt temperature T_m and core is at varying temperatures from T_g and below as shown in Figure 5 (b).

4. Results and Discussion

4.1. Thermoforming

Thermoforming results of the PP and PP/GF based Organosandwich components are shown in the Figure 6 (a), which is a simplified sandwich model formed between positive and negative tool forms. The form structure is a demonstrator part, which is developed to analyze different forming mechanisms and forming parameters.

Here, we show the optimized thermoforming results with controlled laminate temperature and forming velocity. Owing to sandwich bending principle, the top face sheet experiences in-plane tension, on the other side bottom face sheet experiences in-plane compression [8]. With tension and compression on top and bottom sides, the core exhibits pure out-of-plane shear deformation under bending [14]. The cross-ply as face sheet which has fiber orientations in 0° and 90°, in Y and X directions respectively, show in-plane shear stresses in 3D forming area as shown in Figure 6 (b). In the 3D forming region as shown in Figure 6 (b), it is observed that shear stresses are highly concentrated and larger in magnitude compared to 2D forming region. The in-plane shear stresses will develop axial compression in the transition between shear compression and shear tension. Therefore, unwanted deformations like wrinkle formation as in Figure 7 occur only in distinctive areas of the cross-ply laminate depending on the fiber orientation [7].

4.2. Critical Stress Analysis

Sandwich bending for a global tension on the top side exerts global compression on the bottom side. Examining a cross-ply under forming or draping, whilst fibres in the principal axis are subjected to tension, the transverse fibres show compression behavior. Combining these two effects, the in-plane instabilities likely to happen on the bottom face sheet are wrinkling or buckling. For a chosen sandwich construction and selected 3D part geometry, the thermoforming results show a surface waviness on the bottom cover layer in 3D forming region as shown in Figure 7. With the detailed study of the results and stress plots in the specific region, the outer ply oriented at 90° on the bottom face sheet is subjected to local tension in X-direction, and the inner ply oriented at 0° on the bottom face sheet is subjected to local compression in Y-direction.

In Figure 8, the inner UD-ply oriented in 0° to Y-axis shows maximum tensile stress 17 MPa and maximum compressive stress 5.2 MPa. The axial compressive stresses in the respective ply are above the maximum critical compressive stress for wrinkle formation, as mentioned in the Equation 2, ${\sigma }_{wr}$ is 5.03 MPa. As the maximum compressive stresses in the ply exceeds the limit the wrinkle formation on the coverlayer is inevitable. The prediction of mathematical model and simulation results show critical regions in the thermoforming of Organosandwich, which can be used for the design or process optimization.

4.3. Considerations and Product Optimization

From thermoforming simulations, the identified critical regions can be rectified with different approaches like optimization of process parameters, boundary conditions and design modification. Although process parameters like laminate temperature and tool velocity influence both in-plane and out-of-plane buckling, the wrinkle effect in Figure 7 is not their influence. From numerous experiments, the thermal state of the laminate and velocity ramp of the tool are optimally finalized; further modifications in these parameters lead to adverse effects like buckling. The second approach is modification of boundary conditions, which means the holding position of the laminate inside tool with grippers and clamps. There exists already a clamping system in the tool around U-profile (acts like clamping ring in classical forming technologies [15]), which minimizes the bending stresses, its influence was discussed in the previous work [16]. In composite laminate forming, this type of wrinkles are rectified by using a spring-loaded frame gripper system; the springs act as tensioners and generate in-plane tension during forming in the respective direction [17]. The pre-consolidated laminate is placed inside a frame with specific springs and transported to the heating station and forming machine, this kind of gripper system inside the process slows down the manufacturing cycle as hitching and removing the laminate will add additional steps. With time constraints, the idea of including additional steps to the process is not taken into consideration. The existing vacuum gripper system inside the tool middle of the form cavity stays as the standard.

The direct basis of optimization is to get rid of critical axial compressive stresses, when existing process parameters does not allow further modifications, the remaining option is designing modification. To make any changes to the existing part design, one needs to investigate the forming mechanisms at first and their influence on the structural performance subsequently. As the stress concentration is limited to 3D region, a minor change is made to the existing shell structure as shown in Figure 9.

The existing design as shown Figure 9 (a) of uniform thickness through the structure is modified to a reduced thickness of approximately 10% in the three-dimensional bending region, as shown in Figure 9 (b). This reduced thickness generates out-of-plane compression inside the sandwich, and as a reaction, an in-plane local tension is developed, which counters the previously existing critical wrinkle stress. The results of the modified part design are shown in Figure 10. The maximum compressive stresses in Y-direction for UD-ply orientation in 0° are 3.5 MPa, whereas average maximum compressive stresses are around 1 MPa, see Figure 10 (a), which are way below the limit of ${\sigma }_{wr}$ 5.20 MPa. With this design modification, there is a drop in maximum axial compressive stresses by 33% during forming.

Cross-ply of 0°/90° orientation is good at tensile loads on the sandwich structure. Therefore, the top cover layer on the demonstrator part shows no defects in both designs. Thermoforming trials, with the Figure 9 (a) design, have always resulted in wrinkles on the bottom cover layer, with the modified design as Figure 9 (b) does not show any wrinkles, see Figure 10 (b).

5. Conclusions

The PP-based sandwich laminates tested for thermoforming bahavior under different conditions described in this paper. All the thermoforming parameters tested through numerical modelling validated with experimentally produced results. At the end, thermoplastic honeycomb sandwich laminate’s limitations for thermoformability are described. A brief introduction to the thermoforming of thermoplastic sandwich structure was given. With advanced numerical tools, a finite element model is developed to describe the complex thermoform process, making the TSM technology suitable for industrial applications. As the novel thermoforming process is aimed for a high-volume production process like automotive applications, critical deformations can be analyzed with this numerical method at early stages of part design and optimizations are provided by an economical way. Because the wrinkle deformations in the CFRTP could influence the structural properties of a part and decrease the tensile and compressive strength up to 20 – 30 % [25]. Therefore, it is necessary to identify the critical regions in the part through its manufacturing process and eradicate them with optimal solutions. This research work provides the analysis of such deformations in the sandwich structure and helps in optimizing the thermoforming process and further contributes to structural performance of a sandwich component.