Nonlinear Material Modeling for Mechanical Characterization of 3-D Printed PLA Polymer With Different Infill Densities

In additive manufacturing, also called 3-D printing, one of widely used materials is polylactide thermoplastic polymer (PLA) by means of the fused deposition modeling. For weight reduction purposes, infill density is an often used feature in slicing for 3-D printing. We aim at investigating the effect of infill density on the mechanical properties of structures. Therefore, we demonstrate how to prepare tensile specimens and test them by a universal testing machine. Results are collected by a so-called digital image correlation method. As infill density increases, from 10% to 100%, the nominal strain at break decreases from about 2.1% to 1.2%, respectively. In other words, the material becomes more ductile by decreasing the infill density of PLA material, which is possible to justify with an effect of the microstructure created by the infill density. Furthermore, we discuss a possible material model fitting all the presented results and report that a hyperelastic material model is needed for the PLA. We utilize Neo-Hookean, Mooney–Rivlin, and Yeoh models, all for different infill densities. All three models show a fairly good agreement to the experimental data. Neo-Hookean model has an advantage of only one parameter, which increases monotonously with infill density.


Introduction
Additive manufacturing (AM) increases the design freedom to the level not possible by conventional manufacturing methods. Engineers may create or invent more sophisticated geometries leading to tailor-made structures to distinct applications in the future. In AM, also called 3-D printing, the approach is a layer-by-layer production. Hence, the CAD model is prepared by using a "slicer" software for each layer. In this way, not only standard topologies but also internally void structures are possible. Among others, we refer to novel designs by [1][2][3][4][5][6][7] by using advanced fabrication techniques [8][9][10][11][12], for historical remarks, see [13].
The AM gives the possibility to manufacture a material with a microstructure. Hence, this multiscale material opens the possibility to tailor macroscopic behavior. Such material is called metamaterial, for a review, we refer to [14]. Pantographic metamaterials are one of many examples of such mathematically-driven design [15]. With increased design freedom that AM provides, it is possible to experimentally test the theories and models regarding metamaterials. But 3-D printing depends on the processing parameters, we refer to [16] for PLA and [17] for PETG, in the latter a simple formalism is suggested for measuring the effective mechanical properties.
Layer-by-layer fabrication makes internal hollow structures possible with an ease. This feature of the 3-D printing is greatly exploited by a feature is "infill density" or "infill ratio" available in all slicer software solutions out-of-the-box. As the CAD model is imported as a full structure, an outer shell is printed by a chosen wall thickness, whereas inside of this shell is filled by a chosen repeating microstructure. We demonstrate this feature in Fig. 1 by using a (rectangular) grid as the microstructure. This feature is very well accepted by the community since its generation is automatized in slicer programs such as Cura, Slic3r, and similar. Even more importantly, this feature reduces the weight significantly. We know that the mechanical response changes as well; however, this pitfall is often discarded. One possible reason is that there are no predictions how to deal with the mechanical response change in the literature. A common misconception is that the stiffness is simply reduced linearly by the porosity. We clearly stress that experimental observation declares that the relation is more challenging [18]. By incorporating a micrustructure at the same length scale as the geometry itself, infill density may generate a multiscale material with an isotropic behavior [19][20][21] that is under heavy investigation in the literature, among other, we refer to [22][23][24].
The structure of this paper is as follows: First, the concept of AM with specific focus on FDM technique is briefly presented. In Sect. 2, the experimental procedure is explained with PLA specimens for obtaining the stress-strain curves and the derivation of a hyperelastic model. Then in Sect. 4, we present the experimental results in terms of ultimate tensile strength and material parameters of a hyperelastic material model for PLA material. Finally, the main outcome and its connection to microstructure related mechanical response is discussed.

Materials and Methods
The AM is based on constructing the structure layer-by-layer by adding material. As opposed to the subtraction based methods such as machining parts by removing material, it uses just the amount of the material, thus, reduces waste in production. Despite technical challenges, AM is beneficial for crafting sophisticated or even hollow designs. Various methods exist for polymers and metals. Among the methods used for polymers, herein, we concentrate on the method of fused deposition modeling (FDM). The process is to fuse layers of materials together according to a defined blueprint to create an object. For this purpose, a commercially available FDM 3-D printer has been fine-tuned and then used for all experiments. FDM deposits molten polymer. Therefore, only thermoplastic polymers are possible to 3-D print with FDM. As the raw material, PLA (polylactic acid) is a type of polyester synthesized from fermented plant starch (sugarcane or corn or similar). Therefore, during melting and printing, no toxic gases are released. Made of renewable sources as well as being non-toxic makes PLA one of the first choices in the community. In addition to its chemistry, its processing is simple with a small interval of temperature range to be printed and it is fully recyclable (if separated correctly from the common trash). PLA and polycaprolactone (PCL) are biomaterials, especially useful as biodegradable polymer-based composites reinforced by hydrogels [25].
We shortly enlist steps to be carried out in order to arrive at the final part: • model the object, • choose the parameters, • prepare the print code, • and finally manufacture the structure. Especially optimization of process parameters is key for a successful 3-D printed structure. Previous research on FDM based manufactured materials, specifically for PLA, demonstrate that layer height, infill density and pattern alter mechanical properties significantly [26].
A review on the most employed AM technologies for polymers, the commonly-used mechanical test standards and a summary of an exhaustive amount of literature regarding the mechanical properties of 3D-printed parts under different loading types is given in [27]. As an alternative to FDM method, Digital Light Synthesis (DLS) technology hugely accelerates the additive manufacturing of soft polymers. A DLS-inspired 3D printer uses a continuous building technique instead of a layer-by-layer approach, where the curing process is activated by an ultra-violet (UV) light. Using DLS approach, mechanical properties of silicone polymer [28] and elastomeric polyurethane polymer [29] is characterised. In [30] a methodology for the prediction of the mechanical properties and mesostructure of FDM polymers is provided. They proposed a computational framework for the simulation of the printing process taking as input data specific manufacturing parameters and filament properties is proposed.
In this study, PLA is used as material for printing and the ISO 527-2 standard is chosen for the geometry of 3-D printed specimens, for a discussion of this standard to be used in 3-D printing [17]. A technical drawing of the specimen is depicted in Fig. 2, showing the chosen dimensions (in mm) according to ISO 527-2 standard.
As the slicer software, we have employed open-source Ultimaker Cura 4.13.1. We stress repeatedly the importance of process parameters, hence, for securing the repeatibility of this study, we compile them in Table 1to-be-used for creating the G-code form the CAD file exported as STL.
As an FDM 3-D printer, we have used Ender 3 v 2 equipped with a bed-leveler and hotend. The firmware has been replaced to Klipper running on Raspberry Pi 3 with an interface of Octoprint for interpreting G-code. We use Klipper libraries to control pressure for  increased accuracy on curvatures. We also run a bed-leveling for each change of the printtable and renew the configuration file. For a better adhesion, we use a standard sticky-pad replaced for each run. All specimens were printed with the same parameters, all contours are avoided in 3-D printing in order to secure layers along the unidirectional configuration of the structure. All layers were unidirectional with 0° orientation, as visible in Fig. 3.
For printing the structure with a unidirectional manner, we used line patterns that steers the layers and thus generates the fibers. Fiber is a specific jargon among the 3-D printing community for describing line pattern. Parameter for connecting infill lines was switched off such that these fibers were not connected to the endpoints.
Prior to mechanical tests, all specimens were preserved in a dehydrator (at 40 °C blowing air) in order to minimize moisture uptake. Uniaxial tensile tests were performed with a Shimadzu AGS-X universal (electromechanical) testing machine. The tensile equipment and the test set-up is visible in Fig. 4. All experiments were conducted by steering the displacement, after an initial test, we observed that 1 mm/min displacement rate was slow enough to neglect viscous effects. In other words, the chosen displacement rate ensures that the material is purely elastic. Also we check if inertia is insignificant in this speed. The procedure is rather simple that we use the ramp at that speed, load, unload, load again. As the force displacement curve is following the same path, the elasticity is guaranteed. Then we double the ramp speed and observe that the same path is resulted such that the inertia may be discarded. So the quasi-static experiment leads to elastic results by the chosen loading rate. Experiments are captured with one camera digital image correlation (DIC) method (also called digital speckle photography). In this way, we have generated strain fields after the experiments have been conducted.
The DIC method provides non-contact measurement of deformation on the surface [31][32][33]. As the thickness is is smaller than other dimensions, without buckling, the surface deformation is representative for the whole body. Hence, digital speckle photography is adequate for measurements on such objects [34]. An open source 2D DIC Matlab code "ncorr" [35] has been used for image analysis in this study.
The basic principle of DIC is recording subsequent images during deformation, and then computing the displacement field by a posteriori image analysis. The strain field is then obtained from the displacement field by taking the derivative numerically. A DIC implementation consists of: • surface preparation, • image recording during measurement, and • a posteriori image processing.
Since the uniaxial test is designed to build-up a constant and planar deformation (for isotropic materials), we used one camera setting leading to a 2D DIC analysis. For a better image analysis, the sample surface is enriched with a random texture before the measurement. This layer has no specific orientation and is non-periodic [36]. Such a layer, thin enough not to add a significant stiffness to the structure, is mainly achieved by spray painting of loose (not connecting) dots. After that, by spraying a black color on a white background, black and white speckles are generated in a stochastic manner. In this way, the software determines the speckles and calculate the displacement based on their motion between subsequent images. For taking photos, a digital camera was used with a constant focus and manual settings. This choice is indeed adequate in the laboratory providing a diffuse, white light in the room. The image processing starts by dividing the image into subsets. The displacement is then calculated separately for each subset by using cross correlation to identify the same subset in two subsequent images recorded during tests. The correlation function is a function of 2D displacement. The location of the peak of the correlation function gives the position of the deformed subset. An iterative cross-correlation algorithm is needed for more complex deformation fields. Deformation of the subset is accounted for by introducing a shape function that translates the pixel coordinates in the reference subset into coordinates in the image after deformation. The grey values are interpolated between the pixels in the image to achieve sub-pixel accuracy [36].

Parameter Determination in Hyperelastic Material Models
Stress and strains are linearly related in Hooke's law that is accurate in many engineering metals in the elastic regime (roughly less than 0.1% of strain). Polymers demonstrate larger strains and often a linear relation is not an accurate representation. Nonlinear material models, specifically hyperelasticity in polymers, refers to a nonlinear relation between stress and strain. Determination of hyperelastic material model's parameters is often utilized by using a deformation energy. Typical materials are elastomers such as vulcanized rubber or synthetic polymers, along with some biological materials [37]. One of the most important and interesting aspects of hyperelasticity is that it is the simplest model representation possible to motivate by micromechanisms driving the deformation behavior of polymers [38].
In order to capture the elastic and (approximately) incompressible mechanical behavior of such rubber-like materials, numerous phenomenological and micro-mechanically motivated models have been proposed in the literature. For a review of fourteen selected representatives of these models, we refer to [39]. With regards to application of polymeric material in producing actuators, the mechanical response of VHB 4910 is investigated by [40], where a modified version of the micro-mechanically motivated Bergström-Boyce viscoelastic model has been used along with a finite linear evolution law.
We follow [41] and briefly explain the methodology used herein. Consider the right Cauchy-Green deformation tensor, C = F ⊤ F , as well as the left Cauchy-Green deformation tensor, B = FF ⊤ , where F is the displacement gradient with respect to the reference (initial) frame, where X denotes the position of massive particles in their load-free configuration, According to the Hamilton-Cayley theorem, a symmetric tensor of rank 2 possesses three invariants, I 1 , I 2 , I 3 , in three-dimensional space. We stress that left and right deformation tensors have the same invariants, hence the choice depends on the experimentalist. We need a specially designed experiment, for example in a uni-, bi-, tri-axial tensile test, stretch ratio in one, two, three axes, i (i = 1, 2, 3) are used to connect the invariants to the measurable stretch as follows: As det(C) measures the volumetric change during deformation, we use I 3 = 1 meaning a volume preserving deformation assumption-this assumption is also called incompressibility or also an isochoric deformation. In the case of a uniaxial tensile test shown in Fig. 4, we realize Hence, the deformation gradient and deformation tensor read The latter yields to the following invariants: We use a deformation energy density (also called strain-energy or stored-energy) by means of the invariants of the deformation tensors [42]. In an experiment, this energy density is given in terms of the principal extension ratios (stretches), (2) (3) We take three different models into account. One of the widely employed models is the Mooney-Rivlin (MR) 2 parameters model, which reduces to neo-Hookean model (Neo)-this model is available in all software solutions computing displacement such as Abaqus, Ansys, Comsol. Here, we also use Yeoh model (Yeoh) [43]. Historically, the neo-Hookean model, has been motivated by a network of randomly oriented molecular chains [44][45][46]. This model is obtained by truncating the Mooney-Rivlin series to one order. Hence, a generalization is possible by the general Mooney-Rivlin model: where C ij are the unknown material parameters. They need to be determined by experiments, we emphasize that C 00 = 0 . Such models are used to predict the nonlinear behavior of rubber-like materials [47,48]. Truncating the series gives the widely used forms of Mooney-Rivlin models with 2 (MR2) and 5 parameters (MR5), Yeoh model [43] is based on the first invariant only. It is performing well in cases where triaxility dominates over the shear deformation. We use the Yeoh model with three terms, We examine each model by determining parameters by means of a standard inverse analysis technique. First, we generate from force-displacement data, energy-displacement data by simply numerically integrating the force over displacement. We use the initial volume and obtain a deformation energy density. Second, we fit the energy density as in the aforementioned models and obtain the material parameters. As an indicator of the accuracy of the fit curve, we use R-squared (R2). It is a statistical measure defined as follows: using the total sum of squares, S tot , and the residual sum of squares, S res , i.e., where f i symbolizes the fitted data, for a data set with n values marked by y 1 ...y i ...y n and the mean value reads

Results and Discussions
Experiments for several infill ratios have been conducted. PLA is used as material with infill densities from 10% to 100%. The grid infill pattern was selected as presented in Fig. 1. The orientation has been directed along the axis of the sample, thus, along the loading. We demonstrate exemplary results for 4 infill ratios in Fig. 5. All specimen have been printed by using the same process parameters and the same batch of material. Force, F, is directly measured as well as displacement in DIC. From the displacement gradient, we obtain stretch, , as a mean value within the region of constant cross-sectional area of length 50 mm in the specimen shown in Fig. 2. The outcome is force versus stretch-force and stretch are measurable quantities. As obvious in Fig. 5, the repeatability is adequate. We stress that this quality is owing to the used firmware in the system. We have tested the original firmware and observed (not shown here) that Klipper is performing significantly better in terms of repeatability.
For determining parameters, we acquire the energy density in Eq. (6). For this purpose, we integrate numerically force, F, per (effective) cross-sectional area, A, over stretch, , as follows: for calculating the energy density (per volume). Irrespective of the chosen infill ratio, we use the same cross-sectional area defined by the geometry in Fig. 2, namely 4 × 10 mm. The energy density per stretch is then a derived variable from measured variables and used for determining the material parameters in the inverse analysis (curve fitting procedure).
We emphasize that infill ratio is often employed as weight reduction in 3-D printing. The consensus is that the material characteristics fail to change, in other words, a linear elastic material performs the same way for all infill ratios. We show the directly measurable and derived quantities in Fig. 6.
As infill density decreases, from 100% to 10%, the response is getting weaker or material becomes softer, as expected. We comprehend that this materials response change is not only for PLA material, but more because of the structure itself. In other words, we predict that this response change is directly related to the microstructural change since the introduced microstructure is at the same length scale of the geometry. We understand that the stiffness increase and denser material-increasing the infill density-carries more load. We emphasize that the load bearing (cross-sectional) area changes as well. We do not account it here and the value of normalized stress is left undetermined, we refer to [26] for a similar behavior in PLA, where fracture behavior has been investigated. Infill ratio alters the strength; however, it is misleading to consider this property as a porosity of a bulk material. A porosity is rather simple to be involved in a material equation by a mixture rule from classical laminate theory. In this way, a modulus value of 100% would be reduced to 50% linearly. If we compare 70% infill ratio and 50% infill ratio, their Ultimate Tensile Strength (UTS) values are roughly equivalent. Thus, we circumvent ourselves from introducing an infill ratio or porosity related material properties variation. Indeed, 50% has less material and performs softer, but the total force bearing capacity is the same. We stress that the force is compared such that we consider the same cross-sectional area for both cases. Below 50% infill ratio, such an effective cross-sectional area is questionable and in many cases other effects like local rupture under tension or buckling under compression arise. In general, we suggest not even to consider to use less than 50% infill ratio for design.
Interestingly, the material characteristics start changing. For 100% infill ratio, the material behavior is linear elastic (and brittle) that is also the bulk material response of the PLA. We decrease the infill ratio by introducing a grid structure and observe that 50% infill ratio becomes a nonlinear (hyperelastic) material behavior. We understand that this change is purely microstructure related since the same material is used in all specimen. For supporting this observation, we demonstrate an extreme case of 10% infill ratio that this hyperelastic behavior remains as a system response. By excluding less than 50% infill ratio - Fig. 1 indicates a rule-of-thumb that the microstructure turn into the same length-scale as the bulk -we may compare 70% and 50% and understand that the load bearing capacity remains the same owing to this response change. A possible understanding is that the microstructure starts dominating the bulk properties. The simplest model is a hyperelastic material model for covering all infill ratios.
For modeling this hyperelastic material behavior, we utilize three models. Namely, neo-Hookean, Mooney-Rivlin, and Yeoh models are modeling different infill densities as shown in Fig. 7. The corresponding material parameters, C ij , and the accuracy of fitted curves, R 2 , are summarized in Table 2.

(c)
All three models demonstrate an adequate agreement to the experimental data as evident from Fig. 7 as well as from the given value of R 2 in Table 2. In comparison, the lower number of parameters of Mooney-Rivlin with 2 parameters justifies its advantage to the Yeoh model with three terms. Furthermore, neo-Hookean has an advantage of only one parameter, which increases monotonously with infill density. We chose neo-Hookean model to fit its parameter, C 10 , over infill densities. The fitted curve is shown in Fig. 8.
By using a simple second order polynomial, the fit function for C 10 parameter regarding infill ratio, , reads This function may be used in order to obtain a functionally gradient infill ratio hyperelastic material. In simpler words, infill ratio may change in space and the material response is simulated by this fit function in Eq. (15) smopthly altering C 10 in space.

Conclusions
PLA is dominating the 3-D printing raw materials market because of its non-toxic and relatively low temperature manufacturing capability. One key feature in 3-D printing is an increased design freedom incorporating hollow structures. Hence, infill density is in standard use that is based on introducing a repetitive pattern as a microstructure and reducing amount of material inside of the structure. By experimentally investigating the change of material behavior by varying infill density, we have found out • a transition of material behavior from a linear elastic to a nonlinear elastic relation, • a possibility to optimize the infill ratio between 70% and 50% since their load bearing capacity are equal but deformation patterns differ.
Specifically, a hyperelastic material model of the PLA is used to fit the nonlinear material relation. Namely neo-Hookean, Mooney-Rivlin, and Yeoh models have been used for different infill densities. All three models are adequate, yet we propose to use neo-Hookean model because of advantage of only one parameter and monotonous change in material parameter with infill density. We emphasize that we have used a uniaxial tensile test and assume an isochoric behavior. This response is very well known in elastomers, but it is more an assumption to-beused in thermoplastics. Therefore, we stress that the presented outcome herein needs to be generalized by using biaxial tests or other type of deformations. In addition, we propose to perform numerical simulations of a component made of this material under different types of loading and the solution's validation by experiments. We see as a possibility to introduce functionally graded materials where the infill ratio varies in space. Commercial software solutions use such features in topology optimization yet their models use infill ratio as a porosity value and alter material parameters linearly. We demonstrate herein that the material response changes from linear to nonlinear as well.
Funding Open access funding provided by Uppsala University.

Data Availability
The datasets generated during and/or analysed during the current study are available from the corresponding author on reasonable request.
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/.