On the Springback and Load in Three-Point Air Bending of the AW-2024 Aluminium Alloy Sheet with AW-1050A Aluminium Cladding

This article presents the results of an analysis of the bending load characteristics and the springback phenomenon occurring during three-point bending of 1.0 and 2.0 mm thick AW-2024 aluminium alloy sheets with rolled AW-1050A cladding. A new proprietary equation was proposed for determining the bending angle as a function of deflection, which takes into account the influence of the tool radius and the sheet thickness. The experimentally determined springback and bending load characteristics were compared with the results of numerical modelling using different models: Model I, a 2D model for a plane deformation state, disregarding the material properties of the clad layers; Model II, a 2D model for a plane deformation state, taking into account the material properties of the cladding layers; Model III, a 3D shell model with the Huber–von Mises isotropic plasticity condition; Model IV, a 3D shell model with the Hill anisotropic plasticity condition; and Model V, a 3D shell model with the Barlat anisotropic plasticity condition. The effectiveness of these five tested FEM models in predicting the bending load and springback characteristics was demonstrated. Model II was the most effective in predicting bending load, while Model III was the most effective in predicting the amount of springback after bending.


Introduction
In many industries, especially in aviation and automotive production, the selection of material is based on the ratio of its density to strength [1]. In the aviation industry, the materials used for aircraft coatings are high-strength aluminium alloys; however, the disadvantage of these materials is their low corrosion resistance. Corrosion damage occurring during operation is a serious threat, especially for aging aluminium aircraft structures [2,3]. A way to improve the corrosion resistance of the materials is to apply protective coatings made of corrosion-resistant materials. In the case of high-strength aluminium alloy sheets, technically pure aluminium is used as a material for coatings, and applied on both sides of the parent material in the rolling process [4].
Bending is one of the processes of shaping sheet metal elements with unfolded surfaces. Regardless of the bending method and operation used, three basic stages can be distinguished in this process, i.e., elastic bending, plastic bending and springback. One of the features determining the quality of bent elements is their dimensional and shape accuracy. Due to this reason, the main problem during the analysis and design of this process is the determination of the amount of springback of the material after bending and the bending load.
An analysis of the bending process can be carried out with the use of analytical, experimental and numerical methods. The springback value is most often determined by the dimensionless springback factor, K. This factor is calculated as the quotient of the bending radius under load, ρ g , and the bending radius after unloading, ρ s , or as the quotient of the bending angle after unloading, γ s , and the bending angle under load, γ g ( Figure 1). where: γ g and γ s -band deflection angles under load and after springback, respectively, ρ g and ρ s -band radii under load and after springback, respectively. Knowing the value of this factor for specific bending cases enables calculating the radius or angle that should be used during bending in order to obtain a product of the desired shape after bending. In other words, it makes it possible to eliminate or at least minimize the adverse effects of the springback phenomenon after bending in the process and technological equipment design stage. On the other hand, the knowledge of the course and values of the bending load of components is necessary when designing dies and tools, as well as in the selection of bending presses.
Depending on the accuracy of the results of the analysis of the bending process in the design stage, in engineering practice, the correction of the shape of prototype bending tools specified in the design is often performed by trial and error, with the number of trials decreasing as the accuracy of the calculation result increases. From this aspect, it is extremely important not only to know the computational methods used in this area, but also the potential of these methods and their effectiveness depending on the models used and the degree of inclusion in these models of factors influencing the analysed parameters of the analysed process.
In the scientific literature, one can find many works devoted to analyses of the sheet bending process, including the process force parameters and sheet springback for various materials, taking into account various models of strain hardening, various plasticity conditions, changes in the value of Young's modulus, the Bauschinger effect and others. Various approaches have been used to solve the issues under consideration, e.g., analytical [5,6] and numerical [7,8], as well as using neural networks [9]. In most numerical calculations related to the bending process, both implicit and explicit methods are used using commercial software. Much less information is available on the effect of cladding layers in predicting the springback or bending load of clad plates. One can find only the results of research on bending of laminated (double layer) sheets, conducted mainly in terms of determining the impact of the thickness of individual layers on the value of the bending load and springback [10] or determining the bendability limit of such materials [11,12]. However, there are no studies available on clad (three-layer) sheets, which is the subject of this study. In a study on the analysis of the bending process of double clad (three-layer) Cu/Al/Cu sheets [13], the focus was mainly on determining the cohesion of individual layers after deformation. This paper presents the results of experimental research and numerical analyses of the three-point air bending process of AW-2024 aluminium alloy with rolled AW-1050A cladding for anti-corrosion protection which, which, among others, is used to make the outer plating of aircraft structures. The bending load and springback characteristics determined experimentally were compared with those obtained in the process of numerical modelling with the use of various plasticity conditions.

Material and Experimental Procedure
The tests were carried out with the use of AW-2024 aluminium alloy sheets with the thicknesses of 1.0 and 2.0 mm used as-received, i.e., after softening annealing. The chemical composition of the tested sheet material is presented in Table 1. The thickness of the clad layers applied on both sides was 50 µm for a 1.0 mm sheet and 55 µm for a 2.0 mm sheet (Figure 2), which corresponds to 10.0% and 5.5% of the sheet thickness, respectively. The rolled AW-1050A cladding layer closely adheres to the AW-2024 sheet material, with no visible transition zones. From the supplied sheets, samples were taken for uniaxial stretching in the 0 • , 45 • and 90 • directions according to the rolling direction of the sheets. The shape and dimensions of samples for testing in accordance with ISO 6892-1 are shown in Figure 3.   Quasi-static uniaxial tensile tests of individual samples were carried out on a Zwick/ Roell Z030 testing machine by measuring the elongation and change in the sample width using a multi extensometer. Three test samples were taken for each of the rolling directions, resulting in nine samples in total for each of the tested sheets. Tensile (and bending) tests were carried out at the room temperature of 293 K. The tensile speed of the samples was 30 mm/min. On the basis of the results of these tests, strain hardening curves were developed. They were used to determine the material parameters for the four models of the strain hardening curve described later in this work and presented in Section 4. The plastic anisotropy factor, r, was determined on the basis of the relationship between the width strain and thickness strain in the whole range of straining using the method proposed by Welch et al. [15]. For both the 1.0 mm thick sheet and the 2.0 mm thick sheet, the value of the yield strength of the base material was more than twice as high as the yield strength of the cladding material, while the value of the normal anisotropy of the cladding material was much higher than the value of this parameter for the base material ( Table 2). For practical reasons, the curves of strain hardening are presented in the form of constitutive equations of the so-called function of flow stress. Such equations were used, among others, for the analysis and simulation of cold forming processes at relatively low strain rates, when their influence on the flow stress can be neglected. In this paper, four models of flow stress of various complexity levels were selected to describe the course of the material strain hardening of the tested sheet: -Hollomon [16]: -Swift [17]: -Voce [18]: -El-Magd [19]: where σ p -flow stress, ε p -equivalent plastic strain, K 1 -K 4 , A 3 , A 4 , B 4 , ε 0 and n 1 -n 4 -material constants determined experimentally. The Hollomon model, which is the simplest and most often used in engineering, is the strain hardening model, and provides a good description of the hardening curve in a wide range of deformations, which is why it is commonly used in modelling plastic forming processes, especially those with large deformations, such as forging, extrusion, punching, etc. The Swift model, similar to the Hollomon model, due to its versatility but also greater accuracy in the description of the initial course of the strain hardening curve, is very often used in the numerical modelling of a wide range of plastic forming processes in the field of small and large deformations. The Voce model is also often used to describe the course of the strain hardening curve, which, similar to the Swift model, requires knowledge of three material constants. The most complex of the selected models is the Voce model, extended with an additional linear component, and its use requires determination of as many as four material constants. In the literature, the El-Magd model is also referred to as the extended Voce model [20].
Bending specimens that were 40 mm wide and 100 mm long were taken parallel to the rolling direction. Experimental bending tests were carried out on a Zwick/Roell Z030 testing machine in a three-point tool system with the same radii (r = 5 mm) ( Figure 4). For each of the tested sheets, the curves of the bending load as a function of the deflection arrow of the samples were determined in a one load cycle in the range of the deflection arrow (0-20 mm) and in eight load cycles for the deflection, which were successively under load (2, 4, 6, 8, 11, 14, 17 and 20 mm). Then, based on the experimental data, the numerical values of the deflection arrow of the sample after unloading in subsequent bending cycles were determined. The values of the bending angle and the springback angle in subsequent cycles were calculated, together with numerical values of the springback factor corresponding to these cycles. To calculate the bending angles, γ g , and springback, γ s , in individual bending cycles, Equations (6) and (7), respectively, are most often used. However, these formulas do not take into account the radii of the bending tools and the thickness of the bent strand [9].
where • f g and f s -band deflection arrows under load and after unloading, • w-spacing of external bending rollers.
Simultaneously, it can be observed that as the deflection arrow grows, the location of the contact point of the sample with the surface of the shaping tools changes ( Figure 5), which has a decisive impact on the value of the calculated bending angle (Figure 5a) and springback angle (Figure 5b). Therefore, original formulas were proposed for calculating the angle of bending (Equation (8)) and springback angle (Equation (9)), which take into account the observed phenomenon and also consider the influence of the thickness of the band: where A-dimensionless coefficient depending on the thinness coefficient, t/r, and the size of the spacing of the external bending rollers, w. The value of parameter A was determined on the basis of geometrical relationships using the least squares method with the use of the Logger Pro program. The data used to determine the value of coefficient A for w = 50 mm and w = 100 mm are shown in Tables 3 and 4, respectively. The values of the angles γ g3 or γ s3 for individual deflection arrows f g or f s , respectively, were determined using CAD software. . γ g1 and γ s1 -bending angles without taking into account tool radii and sheet thickness, calculated from Equations (6) and (7), γ g2 and γ s2 -bending angles taking into account tool radii without taking into account the sheet thickness, γ g3 and γ s3 -bending angles taking into account tool radii and sheet thickness, calculated from Equations (8) and (9).  The dependence of the value of coefficient A on the value of the ratio of the sheet thickness to the radius of the bending rollers, t/r, is more greatly emphasised with the smaller spacing of the external bending rollers ( Figure 6).

Numerical Modelling
An analysis of the bending process was performed with the use of a nonlinear FEM and the commercial MSC.MARC/Mentat 2020 software. The aim of this analysis was to determine the influence of five various variants of FEM models on the effectiveness of forecasting the springback value and the bending load of clad sheets. The numerical models were developed based on the experiments, taking into account the dimensions of samples and tools as well as the pattern of loading the samples during bending. In all models, the deformed strand of material was defined as a deformable body and the tools of the bending device were defined as perfectly rigid bodies. The friction model was described by Coulomb's law. The friction coefficient between the sample material and the tools was 0.23 [21]. The El-Magd model [19] was used to describe the course of the strain hardening curves in all numerical models. The evolution of the surface of plasticity as a result of the phenomenon of strain hardening was described using the isotropic model.
In the calculations, the associated law of plastic flow of Prandtl-Reuss and the implicit scheme of integration of differential equations over time with the Newton-Raphson method were used.
The simulations of the bending process were carried out with the use of two 2D numerical models (Figure 7) and three 3D shell models (Figure 8): model I-a 2D model analysed in a plane deformation state disregarding the material properties of the clad layers (Figure 7, left), model II-a 2D model analysed in a plane deformation state taking into account the material properties of the clad layers (Figure 7, right), model III-a 3D shell model with the Huber-von Mises isotropic plasticity condition, model IV-a 3D shell model with the Hill anisotropic plasticity condition, V model-a 3D shell model with the Barlat anisotropic plasticity condition. For model II, rigid connections between the cladding and the base material were assumed. The cladding material properties were defined for the outer layers of the mesh. For models III-V, due to the inability to use shell elements in tasks requiring two-sided contact [22], the different material properties of the cladding layers over the sheet were not taken into account. The process of unloading the sample was modelled by switching off the contact between the sample and the bending punch, which caused springback to occur.  For the discretization of the deformable body in 2D models, quad 4 type 11 elements were used. Element type 11 is a four-node, isoparametric, arbitrary quadrilateral written for plane strain applications. For these elements, a modified interpolation scheme was used to improve the bending characteristics of the elements. This allows the ability to capture pure bending using a single element through the thickness. This substantially improves the accuracy of the solution, although the stiffness assembly computational costs increase [22]. The finite element meshes in model I and model II were the same. Four hundred finite elements were used along the length and nine finite elements were used in the thickness of the bent sample (Figure 7). The isotropic Huber-von Mises plasticity condition was used in 2D models. In contrast to discretization of the deformable body in the 3D shell models, 4-node quad 4 bilinear elements of type 75 were used in the formula, considering the effects of transverse shear and improving the behaviour of the shell elements during bending [20]. For these shell elements, the number of integration points on the thickness was 11. The finite element meshes in all 3D shell models were the same (Figure 8a). The size of the mesh elements was 1 mm in the width and length of the bent sample. Additionally, the finite element mesh along the length of the bent specimen in the areas of contact with the tools was densified and amounted to 0.5 mm. Ultimately, each of the 3D shell models consisted of 5120 finite elements. In order to make the results of numerical calculations independent of the mesh size of bodies modelled as deformable, during the development of the FEM model, the influence of the finite element ship size was analysed by densifying the mesh. After each compaction of the mesh, the bending load course was compared with the course of this force obtained with the use of the mesh before refinement. In this way, the nets were compacted until their further compaction had no greater impact on the compared parameter. The material parameters used to calculate the required coefficients in the Hill [23] and Barlat [24] models are shown in Table 1. Additionally, the value of the exponent m in the Barlat equation was assumed to be m = 8 [24].

Results and Discussion
In order to describe the material properties of a deformable sheet, an elastic plastic material model with nonlinear strain hardening was adopted. The material parameters of the tested sheets and clad layer were determined experimentally based on a uniaxial tensile test for four models of strain hardening. The material constants in Equations (2)-(5) were determined for the individual strain hardening curves using the least squares method by applying the Logger Pro program. The degree of matching of the individual equations of the strain hardening curve to its course obtained during the tensile test was assessed using the B f index, which was calculated by relating the root mean square error, RMSE, to the mean feature level of σ p : The values of material parameters in models (2)- (5) and the values of the fit error ( Table 2) ware calculated as the average of the three sample orientations with respect to the rolling direction from the formula: where X-the value of selected material parameter, the subscripts denote the orientation of the specimen with respect to the rolling direction of the sheet. A comparison of the values of errors in fitting the strain hardening curves with the use of individual constitutive models shows that the El-Magd model was the best for both the sheet material and the clad layer material (Table 5). For example, for a 1.0 mm thick sheet, the error of fit with the El-Magd model was approximately 5 times smaller than for the Hollomon model and approximately 4 times smaller than for the Swift and Voce models. Even greater differences were found in the case of the 2.0 mm thick sheet metal. The best degree of fit of the El-Magd model to the results of the experiment in the entire range of sample elongation is also documented in Figure 9. Additionally, it should be noted that the degree of matching of individual Equations (2)-(5) was estimated for the range of tensile test deformation was close to the range of deformations observed in the bending process.  In order to assess the validity of introducing the relationships (8) and (9) to determine the sample deflection angles, a comparison of the values of the deflection angles measured for the samples after springback, γ s , with the results of the angles calculated based on the deflection arrow value was conducted. On the basis of the obtained characteristics ( Figure 10), it can be seen that the values of the bending angles calculated with the use of Equation (9) are definitely closer to the values determined by the measurements than those calculated using Equation (7), especially for larger deflection values of the samples. As mentioned earlier, in order to determine the dependence of the bending load on the deflection value, tests were carried out according to two schemes: multi-cycle bending with sample unloading and single-cycle bending. It is worth noting here that during bending, the surface of the samples did not reveal any changes that could be caused by damage to the clad layer. Based on the bending characteristics obtained (Figure 11), it can be seen that for both sheet thicknesses, the curves of single-cycle and multi-cycle bending practically coincide with each other. Such bending characteristics can be useful in fairly common cases when to obtain the desired curvature of the finished product it is necessary to re-bend it to compensate for springback. The numerical values of the band deflection arrows f g and f s both in the experiment and in the modelling were determined on the basis of the measurement points of the bending characteristics for eight load cycles of the sample. The method of reading the band deflection arrows f g and f s for the example of the second load cycle of the sample is shown for the load characteristic in Figure 11. The values of deflection arrows f g and f s determined experimentally and calculated from individual numerical models are presented in Table 6 for a sheet of 1 mm thickness and Table 7 for sheet metal of 2 mm thickness.  On the basis of the conducted research, the influence of the investigated numerical models on the effectiveness of prediction of springback and the course of the bending load characteristics was determined. Graphs (Figures 12 and 13) show the bending load courses determined experimentally and calculated for the five numerical models tested, respectively, for sheets with thicknesses of 1.0 and 2.0 mm. These results indicate that the bending load characteristic determined numerically with the use of model II is the closest to the bending load characteristic determined experimentally. In order to estimate which of the proposed numerical models I-V is most useful for determining the value of the bending load, the relative error value of the maximum bending load was determined based on measurements and numerical calculations according to the relationship: where X EXP and X FEM are the maximum values of the analysed parameter determined experimentally and numerically, respectively, with the use of various constitutive equations.  In the case of both tested sheets, the best agreement of the numerical modelling results with the results of experimental measurements was found for model II, which takes into account the presence of the cladding layers ( Figure 14). The observed higher values of determined error for the sheet with a thickness of 1.0 mm according to all applied numerical models result from the fact that for this sheet the percentage of cladding layer thickness is almost twice as high as for the sheet with a thickness of 2.0 mm. Figure 14. Comparison of the relative calculation error of the maximum bending load for the tested sheets depending on the numerical models used.
On the basis of the measurement of the amount of elastic deformation of the tested samples, the springback characteristics were determined. The springback characteristics are presented as the dependence of the springback factor, K, on the relative deflection, w/f g (w-distance between supports and f g -deflection arrow under load). In the case of both tested sheets, we obtained linear relationships with a high degree of correlation ( Figure 15). The next step in the research was to determine which of the analysed numerical models correlated best with the results of the experimental measurements. To assess the correctness of the models used, the value of the relative error determined according to Equation (12) was used, as in the case of the bending load. In the case of springback characteristics, Model III (3D shell model with the Huber-von Mises isotropic plasticity condition) can be considered as the best numerical model ( Figure 16). The values of the relative error for the entire range of bending curvatures of the 1.0 mm thick sheet metal ranged from 0.96 to 3.24%, while for the 2.0 mm sheet it ranged from 0.87 to 3.44%. The numerical model II (which takes into account the presence of the cladding layers), which gave the best degree of fit to the experimental results in the case of the bending load, showed a slightly lower fit in the case of springback characteristics, but the relative error values were still within the acceptable range. For this model, the values of the relative error for the entire range of bending curvatures of the 1.0 mm thick sheet metal ranged from 1.54 to 5.49%, while for the 2.0 mm sheet it ranged from 1.25 to 4.77%. Figure 16. The values of the relative error in the calculation of the springback factor, K, of a 2.0 mm thick sheet depending on the quotient (w/f g ) for the tested FEM models.

Conclusions
This work presents the results of an analysis of the bending load and the springback phenomenon observed during elastic plastic cold shaping of AW-2024 aluminium alloy with AW-1050A aluminium cladding. The main conclusions of this work can be summarized as follows: -Among the four popular models of yield stress, the El-Magd model proved to be the most suitable for describing the course of the material hardening curve of the tested sheets, for which the average error of the curve fitting was approx. 0.5%, while for the Hollomon model it was about 3%. - Among the five numerical models tested, Model II, taking into account the cladding layer, turned out to be the most effective in forecasting the bending load, for which the relative error value was below 1% regardless of the tested sheet. -For all the tested models, a significant influence of the ratio of the thickness of the cladding layer to the thickness of the sheet on the accuracy of the calculated bending load was found. With the increase in the proportion of the cladding layer in the sheet thickness, the accuracy of bending load prediction decreases. In the tested case, an increase in the thickness of the cladding layer from 5.5% to 10% of the sheet thickness resulted in an approximately two-fold increase in the relative error of the maximum bending load for all models, except for Model V, where the increase was approximately four-fold. - The developed method for calculating the bending angle on the basis of the deflection arrow is much more accurate than the previously known one, as it additionally takes into account the influence of the tool radii and the thickness of the bent strand on the bending angle value. -Out of all the numerical models tested, Model III turned out to be the most effective in predicting the amount of springback after bending.
Author Contributions: Conceptualization, S.K.; Methodology, S.K. and F.S.; Experimental investigation, S.K., F.S. and G.P.; Numerical calculation, S.K. and G.P; Analytical and formal analysis, S.K., F.S. and G.P.; review and editing F.S. and S.K. All authors have read and agreed to the published version of the manuscript.