The influence of stiffener geometry on flexural properties of 3D printed polylactic acid (PLA) beams

We used finite element analyses (FEA) on Abaqus to study flexural properties of additive manufactured beams using polylactic acid (PLA) polymer. Experimental stress–strain data from flexural testing are used to define elastic–plastic properties of the material in the computation software. The flexural experiments are used to validate the FEA approach suggested. The method provides good results of deflection and stress with errors well below 10% in most of the cases. Therefore, by using the proposed approach, costs related to repeated experimental works can be avoided. In addition, the flexural rigidities of the additive manufactured beams are studied. Five different beam stiffener designs (diamond, honeycomb, square, triangular and wiggle) are studied based on beam bending theory. The force–deflection data from the flexural tests are used to determine the area moments of inertia of the beams. The honeycomb stiffener showed the highest force–deflection behaviour that led to the highest calculated area moment of inertia. However, with the lowest force–deflection behaviour, the square stiffener had the lowest calculated area moment of inertia.


Introduction
The additive manufacturing (AM) is the state-of-the-art technology that changed the conventional approach to manufacturing systems. It is defined as fabricating a 3D model directly from computer aided design (CAD) software without process planning [1]. It also simplified the complexity of manufacturing difficult component geometries. The ease with converting virtual models to products does not always come with the best combinations of strength properties. Particularly, the additive manufactured products are susceptible to process related flaws including pores, surface roughness and geometric deviations from nominal dimensions [2]. Several papers discuss influences of 3D printing parameters on structural rigidity of additive manufactured components.
Porter et al. [3] studied the flexural rigidity of 3D printed slender beams by varying infill densities for a grid pattern. According to the authors, the optimal infill percentage that maximizes specific flexural rigidity is between 10 and 20%. Zuo et al. [4] conducted stress analysis on bidirectional adaptive infill pattern using numerical methods. The method is validated by tensile experiments on additive manufactured PLA samples. Rectilinear, triangular and honeycomb patterns were manufactured using fused deposition modelling (FDM) and tested in tension. Their study concludes that the triangular patters showed the highest tensile strength. Their finite element (FE) model is made by extruding the planar sketch generated by the G-code commands of a single layer. Flexural and tensile strength of acrylonitrile butadiene styrene (ABS) composite is studied using a varying fill speed and layer thickness [5]. Based on the experimental study in [5], low printing speed and low layer thicknesses resulted in maximum tensile and flexural strengths as compared to other process parameters. Fernandes et al. [6] studied the influences of infill density, extrusion temperature and raster angle on mechanical strength. Parameters optimization of FDM process for blended PLA is studied using the Taguchi method [7]. Othman et al. [8] studied the influence of layer thickness on the impact strength of PLA.
A variety of stiffener designs are used to reinforce beams. This paper investigates the flexural properties of 3D printed PLA beams. PLA is the most widely used thermoplastic polymer in FDM processes. It showed a promising beginning in biomedical applications including tissue engineering, scaffolds, drug delivery and implantable devices [9]. In addition, its light weight and relatively higher stiffness are suitable for light weight designs in numerous emerging applications including cellular parts, sacrificial cladding [10,11]. However, the high brittle and low toughness properties have limited its applications. With the advancement in additive manufacturing processes and polymer blending technologies, these mechanical properties will be enhanced to increase its applicability [12].
The objective of this paper is to study the force-deflection behaviours of five different beam stiffener designs under flexural loading. Finite element analyses are made using Abaqus CAE [13] and COMSOL Multiphysics [14] to study these flexural behaviours. The FE calculations are validated by 3-point bending experiments. The standard 3D FEA on Abaqus are also validated using Euler-Bernoulli beam formulation on COMSOL. The experimental stress-strain data are directly used in the FE analyses to configure the elastic-plastic property of the material. The study aims at finding similar flexural behaviour between the FE calculations and experimental investigations.
The different beam stiffeners have different geometrical properties including the area moment of inertia. Therefore, by studying the force-deflection behaviour in linear elastic region, it is possible to examine the area moments of inertia that directly determine the beams' rigidity under bending stresses. In relation to this geometrical property, the force-deflection behaviours are compared to identify the beam stiffener that has the better flexural properties.

Deflection of a beam
The curvature of a Euler-Bernoulli beam subjected to a pure bending moment M due to a normal load P is defined as [15,16] For small deflections, a closed form of solution can be found by integrating Eq. (1) and applying proper boundary conditions [17]. A simply supported beam under a bending load is presented in Fig. 1. For the beam indicated, I is the area moment of inertia with respect to z axis, h is the beam height and L beam length.
The maximum deflection occurring at the location of the bending force (mid-section) is given as [16] (1) Equation (2) is the Euler-Bernoulli beam equation that is used to determine the deflection of the beam under concentrated force acting at the middle. For known values of force and deflection, Eq. (2) can be rearranged to determine the area moment of inertia. Based on the stress-strain relation for a uniaxial loading, the flexural strain of the beam in the Fig. 1 is defined as a function of the maximum deflection, beam length, and thickness [18]. That is, The bending stress at the location of the applied force that is a distance c away from the beam's neutral axis is given as [16] From the theories of the linear elasticity, the equilibrium equation for 3D objects can be presented in terms of displacements [19] In Eq. (5) u and P are displacement and body force vectors in Cartesian coordinates. and are shear modulus and Lamé's constant and are further presented using elastic modulus E and Poisson's ratio v as If the body force is not considered, the equilibrium equation in y direction reduces to a form, In addition to the equilibrium condition, the displacement formulation should fulfil compatibility. Several displacement potential functions can be used to solve the displacement formulations presented in Eq. (7) [16].
In COMSOL, the time dependent Euler-Bernoulli beam bending equation in y direction is given as [14] (2) = In Eq. (8) is density, A cross-section, M z moment about z axis, and f y force per unit length in y.

Finite element modelling and analyses
Abaqus CAE and COMSOL Multiphysics are used for the structural analyses. The five beam models have geometrically different stiffeners. These include triangular, wiggle shaped, square, diamond, and honeycomb geometries. The beams' gross dimensions are; length L = 210mm , width w = 20mm and total height h = 14mm . Each beam has a stiffener height of 12mm that is sandwiched between thin solid slab of 1mm thick at the bottom and top. The wall thicknesses of the stiffeners vary between 0.44and 1.18mm . The average stiffener length represents the side length of the stiffener.
Due to the differences in geometries, the number of stiffener instances vary in each beam. Table 1 presents the details of dimensions. And, Fig. 2 shows the five different structural stiffeners modelled inside the beams.
The boundary load and constraints are defined in a similar manner to the experimental test setup. The stress-strain relations from the experimental investigations are imported to Abaqus to create the material property with a Poisson's ratio v = 0.36 [20]. For the finite element analyses, 10-node 3D quadratic (second-order) elements are used. Tetrahedral mesh, C3D10 is assigned for the element shape. Global element sizes of 0.2 mm on the planar surfaces. However, the local element sizes are edited to 0.01 mm around the stiffener geometries. These mesh refinements are approached by continually increasing the mesh densities. The iterative mesh convergence study indicated that, at these element sizes, the results of the FEA analyses have shown convergence. In addition to the solid (continuum) model analyses in Abaqus, separate FEA are made on COMSOL structural mechanics using a beam physics. This is due to the need to validate the 3D standard analyses. For the analyses, the Euler-Bernoulli beam is selected under beam formulation tool. The beams' polar moment, section moment and area moments about the local axes are calculated and defined as section properties. The experimental data that characterise the stress-strain behaviour of the beams are defined using interpolation function. The convergence of the results is studied by continually decreasing the element sizes from 0.1 mm up to 0.01 mm, see Appendix.
In both analyses, the boundary load is defined using 5 mm/min velocity in the negative y direction. The velocity boundary is defined at the middle of top surfaces of the beams. At this load boundary, the displacements in the remaining axes are constrained. At the bottom surface of the beams, fixed displacement constraints in z axis are  . 2 The five different stiffener patters inside the beam structure placed at two locations that are symmetrically offset from the centre. The displacement constraints represent the support and are 154 mm apart. Even if the models on both analyses are not structurally representative of micropores resulting from additive manufacturing processes, the structural differences are compensated using the experimental stress-strain data in modelling the material properties. Figure 3 shows the 3D printed PLA beams and mounting of the samples for flexural testing. The flexural tests were performed on Testometric M350-5CT machine. While 3D printing the PLA beams, recommended temperature parameters from the manufacturer's technical datasheet were used [20]. The beam models are designed on SolidWorks and 3D printed using MakerBot replicator printers. The infill densities on the stiffeners were adjusted to 100% to print the profiles with solid infill. As the beams were relatively long, they tend to detach from the printing bed at their edges. This resulted in warping that kept the beams out of dimensional tolerances. In order to avoid the warping, the print bed height was adjusted and calibrated first. It was then cleaned with acetone, dried and kept at uniform temperature of 65 ℃ during the printing processes. After setting the temperature, double sided tape was applied on the print bed to make sure the prints adhere to the plate. Additionally, to increase the adherence, the beams were printed on raft. After the prints were completed, the beams' cross-section across their lengths were measured using digital micrometre. The flexural tests of the beams are made according to ISO178:2019 Plastics-Determination of Flexural Properties [18]. Three beams from each beam models are tested and the arithmetic mean of force-deflection data are used to determine the stresses and strains. The span length between the supports is L = 154mm.The flexural tests were performed using a recommended test speed, u = 5mm∕min . The force-deflection data from the flexural experiments are used to estimate the average area moment of inertia of beams in the linear elastic region using Eq. (2). Afterwards, these force-deflection data are converted to stress-strain relations using Eqs. (3) and (4).

Experimental results
A total of 15 flexural tests were made, three for each beam. The results of the flexural experiment indicate that the beams with honeycomb stiffener showed the highest force-deflection relation compared to the others. The beams with square stiffener had the lowest one. The beams with triangular stiffener exhibit brittle behaviour snapping after a deflection of 2.53mm . To the contrary, the beams with diamond stiffener showed a good toughness property sustain a deflection of 6.97mm . The wiggle stiffener exhibited a better force-deflection than the square one. However, the stiffener instances failed prematurely and their cross supporting after failure influenced the force-deflection data. Figure 4 shows the experimental force-deflection results of the 5 beam models.
The force-deflection results in the linear regions are used to estimate the average area moments of inertia. These geometrical properties are related to the force-deflection behaviours of the beams. The beam with the larger area moment of inertia has higher force-deflection behaviour and the beam with the smaller area moment of inertia has lower force-deflection behaviour. Figure 4 and Table 2 present these relationships.
The force-deflection data are further converted to flexural stress-strain relations for modelling material property in Abaqus CAE and COMSOL Multiphysics. These stress-strain relationships are the material's inherent behaviour and show similarity among the different beam models. Figure 5 shows the experimental stress-strain relations.

Finite element analyses results
Abaqus standard for the 3D models and COMSOL structural mechanics for the beams are used in the FE analyses. Converged solution of displacement, flexural stresses and flexural strains are obtained from the two analyses software. The solutions of the solid 3D models and beam structures showed good agreement. This indicated the validation of the results from the Abaqus analyses. Figures 6,  7 and 8 present the contour plots from Abaqus for the displacement, bending stress and strain respectively. The results of the Euler-Bernoulli beam analyses on COMSOL are also presented, see the Appendix.
Elements that are located at the mid-sections of the beams are selected to get the stress and strain data from the FE calculations, see Figs. 7 and 8. The average stresses and strains on the selected tetrahedral elements are used to plot the stress-strain curves. Figure 9 shows the FE calculated stress-strain relations of the beam models.
The deflections and stress-strain results from the experiment and finite element analyses on Abaqus are compared. The comparisons are made by selecting deflection and stress data at selected time frames on the finite element study. The average errors in deflection and stress are calculated as, The average percent errors of the finite element calculations as compared to the 3-point bending experiment are calculated using Eq. (9) and presented in Table 3.
The deviations in deflections and stresses through the progress in time are also clearly presented on Figs. 10 and 11.

Conclusion
Five different beam stiffener designs that are 3D printed from PLA were analysed. Finite element analyses and experimental methods were used to study the beams' flexural rigidity and force-deflection behaviours. The experimental force-deflection results were converted to stress-strain relations and used to build the PLA elastic-plastic material property in Abaqus and COMSOL. We used two FEA methods and analysed the beams as 3D structure on the Abaqus and Euler-Bernoulli beam model on COMSOL. The FEA results showed good agreement among results produced.
Comparing the FE calculations with the experimental results, the errors were well below 10% in most of the Fig. 9 Flexural stress-strain relations of the beam models by finite element calculations by Abaqus CAE  cases. In deflection, the beam with square stiffener showed a maximum of 4.29% error. However, the stress comparison showed that the beam with wiggle stiffener had a maximum of 9.5% error. The experimental force-deflection data of the wiggle beam near fracture point were influenced by the effect of failed stiffeners. After break, the stiffener instances leaned on each other and reinforce. Hence, the force-deflection data past the linear region were not purely from beam's inherent flexural rigidity and are found being the sources of the maximum error observed. However, in most cases the deviations are close to 1% in deflection and below 5% in stress. These are tolerable error margins in structural analyses and can be used as reasonable approximation of the stress-strain behaviour of the polymer. The property definition using the experimental data mostly compensates the geometrical differences arising from the solid models in finite element and porous models in 3D. However, the reason for the deviations of results presented by Table 3 can still be attributed to the geometrical differences between the 3D printed beams and the FE modelling. While 3D printing the beams structures, significant efforts were made to optimize printing parameters including glass and melt temperatures, layer thickness, print speed and so forth. To maintain process uniformity, the three samples from each geometry were printed at once. After the beams were produced, measurements were made to inspect their dimensions. However, with all these efforts made, the prints had process related flaws. The wiggle stiffener beam with the least number of instances showed higher warping at edges. Process related flaws like these are additional sources of deviations between the FEA models and 3D printed beam structures.
In this paper, the influence of the beams' geometry on their force-deflection behaviour is studied. Having the highest area moment of inertia, the beam with honeycomb stiffener showed a good force-deflection relation. The square stiffener with the smallest area moment of inertia showed a poor force-deflection relation. These relationships agree with the classical beam theory. The beam with diamond stiffener, even if with lower force-deflection slope than that of the honeycomb and triangular stiffeners, it had the highest deflection before fracture. The fracture occurred on the stiffener instances and the beam deflected higher before the stiffeners are broken all across. The wiggle stiffener, even if better than the square one, the linearity of its force-deflection curve is maintained up to 4.5 mm in deflection. The force-deflection reading beyond the 4.5 mm deflection was influenced by the cross bracing of failed stiffeners.
Porter et al. [3] investigated the influences of infill densities for 3D printed beams, however they only considered a single grid pattern. Geometrical influences were also studied by Zoe et al. [4] using FEA, however their method lacks real material behaviour and characterizes only tensile properties. The impact strength of PLA was studied by Othman et al. [8], however their experimental study does not address the influences of infill geometries. Our work incorporates flexural rigidity analyses of wide variety of stiffener models that are validated using experiments. This work showed the influences of different stiffener geometries on flexural rigidity of PLA beam structures using FE approaches. By using the developed FE methods, costs related to repeated experimental works can be avoided. Following the path, our future researches will focus on studying the influences of the geometrical variation caused by additive manufacturing processes at microscale. Sciences for providing accesses to materials, equipment and software used in this research.

Funding Open access funding provided by Arcada University of Applied Sciences (UAS).
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creat iveco mmons .org/licen ses/by/4.0/.

Appendix
Euler-Bernoulli models model analyses results of deflection and bending stress using Comsol Multiphysics.
See Appendix Figs. 12 and 13