Nanoparticle Shape Influence over Poly(lactic acid) Barrier Properties by Molecular Dynamics Simulations.

Climate change is leading us to search for new materials that allow a more sustainable environmental situation in the long term. Poly(lactic acid) (PLA) has been proposed as a substitute for traditional plastics due to its high biodegradability. Various components have been added to improve their mechanical, thermal, and barrier properties. The modification of the PLA barrier properties by introducing nanoparticles with different shapes is an important aspect to control the molecular diffusion of oxygen and other gas compounds. In this work, we have described changes in oxygen diffusion by introducing nanoparticles of different shapes through molecular dynamics simulations. Our model illustrates that the existence of curved surfaces and the deposition of PLA around them by short chains generate small holes where oxygen accumulates, forming clusters and reducing their mobility. From the several considered shapes, the sphere is the most suitable structure to improve the barrier properties of the PLA.

Introduction

Climate change is one of the major concerns in recent years; the need to abandon the oil industry and materials from its spread in all sectors is not a simple change. The development of new biodegradable polymers and the improvement of available materials remain a fundamental task, contributing to natural degradation with a lower impact on the environment. In this area, poly(lactic acid) (PLA) has been described as one of the future materials that can replace plastics in the packaging area and medical prostheses.[1−7] Additionally, this material is employed in the chemistry industry for the separation of gases with high energy efficiency[8−10] or in the distribution of medicines for certain treatments.[11]

PLA is a thermoplastic obtained from crops like corn, cassava, sugar beets, and potato starch. It is a promising material for environmental applications, considering that its mechanical and thermal properties are close to petroleum-based products such as poly(ethylene terephthalate) (PET), and it is highly biodegradable and compostable.[12−16] However, its barrier properties need to be improved to allow the extensive commercial application of PLA, along other properties.[14−17] The general improvement of these important properties is still expensive.[1,3,18] Recent studies reported the effect of modifying the PLA structure by the introduction of different materials. Cellulose fibers,[15,19−22] silica nanoparticles,[1,4,9,10,14,20,23−30] and other organic compounds[16,31−39] have remarkably increased the mechanical, thermal, and barrier properties, reaching values close to current plastics.[16,40] Also, some studies analyzed the PLA crystallinity and how its variation produced an improvement in mechanical and barrier properties.[15,31,41,42] From the theoretical point of view, the analyses of the barrier and mechanical properties of PLA and flexible polymers have been performed as a function of content, strength, geometry, and interactions, and how the increased aspect ratio of nanocrystalline cellulose modifies the water permeability.[43]

The potential impact of nanoparticles on the plastics’ average properties is very high due to their extended surface area and the high number of reactive sites with the polymer matrix. Composite materials of polymers and nanoparticles have been studied experimentally, focusing on the changes produced in the thermal and mechanical properties. Only a few of them have measured the overall change in the barrier properties. These studies include organic nanoparticles,[32,33,16] metallic nanoparticles,[14,35] carbon-based materials,[31,44] cellulose nanocomposites,[15,38] and silica nanoparticles.[1,4,9,10,14,20,23−30,34,36] These articles highlighted the importance of the type, quantity, and dimensions of the grafted polymer and studied their effect on the polymer’s different properties. Nevertheless, most of the reported studies are only focused on the experimental characterization of the properties of polymer matrices modified with nanoparticles, not considering the structural influences on the final features of the material. Up to date, only a few works perform systematic analysis of the effect of the shape of the nanoparticles on the final properties of the system. In detail, Knauert et al.(44) analyzed the influence of the shape of carbon materials (nanotubes, fullerene, and graphene sheets) in the viscosity and the ultimate tensile strength of the polymer. Lin et al.(35) observed the nanoparticle’s aggregation as a function of the dimensions, quantity of the grafted polymer, and the shape of the nanoparticles (rod, spheres, and nanoplatelets). This characteristic is remarkable considering that new studies are capable of tailoring the shape of the nanoparticles varying the properties of the full material, as shown for silver by Khodashenas et al.(45)

In the past decades, several works have applied simulations to understand the changes in polymers when introducing nanoparticles.[3,18,46,47] Computational approximations appear as a key instrument to verify and characterize the chemical interactions at the atomic scale, delivering relevant energetic and structural information to complement the experimental analysis. Multiple studies have assessed, at different simulation scales of dimension and time, the concentration effect of a modifier on the interface between the polymer and the nanoparticles.[29,48−53]

These studies revealed the importance of a higher or lower density interface,[48−51] the effects of the size, agglomeration, and spherical nanoparticles’ curvature.[29,47,48,52−54] These theoretical reports, as most of the experimental ones, are focused on evaluating the changes in thermal and mechanical properties of the final nanocomposite and, to a lesser extent, studying the barrier properties. To the best of our knowledge, a key factor has not been addressed yet: the effect of the nanoparticles’ shapes.

In this work, we focused on studying the barrier properties of PLA and how the introduction of a modifier agent, such as a nanoparticle, affects the final properties of the polymer matrix. The shape of the modifier was systematically analyzed to evaluate its influence on the chemical configuration and properties of PLA. The diffusion of molecular oxygen (O2) in PLA and modified PLA was studied based on atomistic models described by molecular dynamics. These models highlighted the influence of the shape of the nanoparticle on the barrier properties of PLA. We observe the curvature effect of the different morphologies, identifying the sphere like-shape as the best solution to improve the barrier properties of the polymer.

Methods

Modeling of Neat PLA

The study of the barrier properties of PLA with nanoparticles of different shapes was carried out using an atomistic model based on molecular dynamics (MD). The PLA samples were built using the SCIENOMICS MAPS code [https://www.scienomics.com] and a Dreading-type interaction.[55] The MD simulations were performed using the LAMMPS code[56,57] and analyzed by the OVITO program.[58] The CHARMM potential represents the interaction of PLA atoms and molecular oxygen.[59,60] In all the simulations, the time step was 2 fs and an rRESPA-type integrator with two additional time steps, the first half of the original and the second a quarter, likewise the interaction of pairs with one step is imposed of time equal to the original and a k-space with half the initial time step.

Following the criteria of minimum energy and non-size effects described by McAliley and Bruce[59] and Zhang et al.,[28,61] we elaborated PLA samples ([C3H4O2]) of 32 chains with 50 monomers each. The samples were arranged in the MAPS code using the Dreading potential. These samples were thermalized by a simulation process in LAMMPS using the CHARMM potential with the method described by Wensink et al.(62) for the relaxation of polymers at room temperature and pressure and used to obtain a good density approximation. In this method, the sample is under successive simulations with a constant number of atoms, volume, and temperature (NVT). Alternatively, the model is treated with a constant number of atoms, pressure, and temperature simulation (NPT). Table describes in detail the performed simulations of the sample in the box relaxation.

Table 1

Summary of the Conditions Used in the MD Simulations to Build the Initial Box

ensemble	pressure (bar)	temperature (K)	run time (ns)
NVT		500	0.05
NVT		300	0.05
NPT	1000	300	0.05
NVT		500	0.1
NVT		300	0.05
NPT	5000	300	0.05
NVT		500	0.1
NVT		300	0.05
NPT	30,000	300	0.05
NVT		500	0.1
NVT		300	0.05
NPT	1	300	1

The relaxed PLA samples were analyzed in terms of the density, the glass transition temperature (Tg), the melting temperature (Tm), and the self-diffusion coefficient. Tg and Tm were calculated by the method of McAliley and Bruce,[59] which describes a variation in the glass transition temperature with the cooling speed when using high speed rates in MD simulations. In this case, the temperature was raised to 750 K with an NPT ensemble at 1 bar during 100 ps. Next, the temperature was reduced at different speeds (500, 300, 100, and 50 K·ns–1), using the same pressure. The variation in density with respect to the temperature for each case is shown in Figure , where two different changes of the slope are observed: a constant in all cases corresponding to the Tm and a variable for each case corresponding to the Tg. Later, the glass transition temperature variation with the cooling rate was represented, and the experimental value of 10 K·s–1 was extrapolated.

Figure 1

(a) Density as a function of the temperature for molecular dynamics simulation of PLA cooling at different rates from 750 to 300 K. (b) Zoomed in view of the glass temperature (Tg) transition section (red). (c) Zoomed in view of the melting temperature (Tm) transition section (blue).

The study of the self-diffusion coefficient was calculated by the analysis of the mean-square displacement (MSD) of the PLA atoms as a function of time, according to the Einstein equation:[4,60,63−67]where D is the diffusion coefficient, N is the number of molecules, and r⃗(t) and r⃗(0) are the final and initial positions of the PLA chain in a time interval. Thus, an NVT simulation at 300 K for 4 ns was performed. The MSD as a function of time was plotted for the last nanosecond of each simulation. Additionally, another five simulations were run for each case, under the same conditions with different initial velocity distributions, to calculate the error of the diffusion coefficient. Finally, some additional restrictions were imposed on the PLA system to mimic the real behavior of PLA. The central atoms of the PLA (*–COC–*) were constrained to reduce the self-diffusion, following the models described by Suarez et al.,[68] Takahashi et al.,[69] and Maskey et al.(70) Methyl groups, double-bond oxygen, and hydrogens were not constrained. A test without constrain showed that the MSD function has many fluctuations, and the PLA self-diffusion coefficient is far from the real value.

Modeling of PLA with Nanoparticles

The simulations with nanoparticles were carried out defining regions in space with four different shapes (sphere, cylinder, cone, and flake), enclosed by repulsive surfaces to the atoms of the system, that is, to the atoms forming the PLA and to the oxygen molecules. Therefore, these regions are prohibited areas to occupy by an atom. We did it this way to focus exclusively on the shapes of nanoparticles rather than their types. The flake shape was modeled as a cylinder with a high radius and a low height, while the rest of the shapes were modeled following their geometry. The LAMMPS fix wall command was set, which imposes a potential in a region centered on the box with the desired shape, in our case, Lennard-Jones with ε = 1.0 kcal·mol–1 and σ = 1.0 Å. The samples with the gap were obtained with this command starting from a region 1000 times smaller (an equivalent volume of 10–9 Å3) and increasing its dimensions one time every thousand steps. During each interval of 1000 steps, an NVT simulation was run at 300 K and at the end of the growth for one million steps. These simulations permit atom evolution and relaxation for the non-equilibrium process.

Diffusion and Barrier Analysis

The barrier properties were studied by the diffusion of oxygen through each sample. Traditionally, experimental permeation tests are performed in an atmosphere saturated with oxygen. Therefore, in the simulations, the change of volume was analyzed while increasing the number of oxygen molecules inserted randomly in the voids of PLA and outside of the nanoparticle region, seeking the highest oxygen concentration without producing volume expansion without control. By this strategy, a maximum of 15% in mass of oxygen (652 molecules) was added randomly to the PLA system, with a volume variation below 10%. Percentages higher than this value produced the destabilization of the original sample. Simulations with oxygen molecules were run in an NVT ensemble at 300 K for 4 ns. As in the PLA sample, five other systems with the same conditions were run with different starting velocity conditions. The MSD of oxygen during the last nanosecond of each simulation was used to calculate the diffusion coefficient by applying eq . Simulations of 8 ns were performed to observe that the diffusion coefficient takes a similar value and a constant order of magnitude compared with the case of 4 ns. The comparison between the different nanoparticle shapes was made with an equivalent volume fraction.

Finally, a detailed analysis of the oxygen interaction with the PLA sample was performed to observe which atoms are responsible for the molecular oxygen interaction reducing the diffusion. We studied the oxygen interactions with other atoms of the sample, taking special attention to the duration of each interaction. The oxygen trajectory was followed by observing the nearest atom in each case. A possible interaction was considered when the oxygen atom was found closer than 2 Å from another atom X. Then, we observed that the most common possible interactions were with hydrogen and other oxygen molecules. Therefore, the oxygen molecules were labeled as interacted if the distance was lower than 2 Å and the angle was 104 ± 5°, applying only the last condition for the interactions with hydrogen atoms.

Results

Analysis of Neat PLA

PLA samples simulated by molecular dynamics have a density of 1.09 ± 0.01 gr·cm–3 (box dimensions of 5.65 × 5.65 × 5.65 nm3), an approximate glass transition temperature of 364 K, a melting temperature of 600 K (as observed in Figure ), and a self-diffusion coefficient of around 10–7 cm2·s–1. These values are coincident with those described for the CHARMM potential in a pure PLA study[59] and can be considered a good approximation to the experimental values of density = 1.25 gr·cm–3, Tg = 300 K, and Tm = 500 K.[4] However, the value of the self-diffusion coefficient is higher than the value described in the experiments ∼10–8–10–9 cm2·s–1. To mimic this value, a series of restrictions on the central atoms of the PLA chain (*–C–O–C–*) were imposed, which reduced the mobility of the chains and allowed the representation of the movement of micrometric chains that occur in reality. Figure shows how the average square displacement of the PLA atoms varies as a function of time with and without restrictions on the chain’s central atoms. The behavior of the PLA without any restriction is quite similar to that of a liquid material. When imposing the restriction, the system tends to a stable value. Additionally, the fictitious self-diffusion coefficient of our fixed system would have a value around 10–9 cm2·s–1, which is in the order of the experimentally expected magnitude. Based on these results, the simulated PLA sample has density, melting and glass transition temperatures, and a self-diffusion coefficient close to the real value, which allowed us to analyze the barrier properties of the system.

Figure 2

Mean-square displacement of the PLA atoms in a simulation with the force of the central atoms established as zero (black line) or without any restriction (red line). The arrows indicate to which scale each curve corresponds.

Diffusion Analysis

The barrier properties of packaging materials mainly focus on the permeation of oxygen, water vapor, and carbon dioxide. The permeation of a gas through a material is described as the product of the diffusion coefficient by the solubility coefficient. This work was focused on studying the diffusion of oxygen through PLA, observing how it varies when introducing structures of different shapes in its interior. The diffusion coefficient variation is directly related to the permeation coefficient since the sample of PLA was immersed in a supersaturated oxygen ambient to reach the maximum solubility inside it. This super-saturation will increase the internal pressure in the sample, which will produce an overestimation of the expected diffusion coefficient. Still, it will not affect the study of the shape influence on the barrier properties of the nanoparticles.

First, we studied the oxygen diffusion coefficient in pure PLA without the insertion of nanoparticles (D0), obtaining an approximate value of 1.75 × 10–5 ± 1.34 × 10–12 cm2·s–1. Figure shows the reduction in the diffusion coefficient () for the different inserted nanoparticles (sphere, cylinder, cone, and flake) with a volume equivalent to a sphere of 1 nm radius (black) and a sphere of 1.5 nm radius (red). As can be seen, the nanoparticles with a smaller volume have a lower impact on decreasing the diffusion coefficient. In contrast, for larger volumes, the effect of the shape on the diffusion coefficient is similar for all the evaluated nanoparticle shapes. This difference can be related to the size of the nanoparticle and the total system. Once nanoparticle volume is increased and the shape approximates to the limits of the simulation box, unrealistic effects could be generated that change the structure of the sample.

Figure 3

Percentage of the decrease of the diffusion coefficient for the different nanoparticles introduced in the PLA. Black bars correspond to an equivalent volume of 4.19 nm3, and the red to 14.14 nm3.

On the other hand, in the case of nanoparticles with smaller volumes (V = 4.19 nm3), the sphere is the one that offers a greater reduction of the diffusion coefficient, followed closely by the cone. Both structures would then produce a higher reduction in the permeation coefficient and improve the PLA for packaging. Finally, the cylinder and the flake cases are the least favorable, and by their design, we are comparing a higher cylinder (cylinder) with a wider one (flake). Therefore, a higher cylinder would always be more favorable, a shape with a larger curved surface. From this analysis, it can be observed that the sphere and the cone are suitable forms to reduce the oxygen permeation in the PLA. Probably this effect is related to the greater amount of curved surface on both structures. If we focus on the curved surface, the sphere has the highest value (Sc = 12.57 nm2, data for structures of smaller volume) followed by the cone (Sc = 11.01 nm2), the cylinder (Sc = 7.60 nm2), and the flake (Sc = 6.48 nm2), coinciding with the percentage decrease in the diffusion coefficient.

Barrier Analysis

Detailed analysis of the oxygen occupation was performed inside the PLA sample averaged over the last nanosecond of simulation and considered to calculate the mean-square displacement and diffusion coefficient. The degree of occupation of the sample allows us to describe the influence of each of the nanoparticles and observe if there are areas of greater concentration. Figure shows a color map of the oxygen density averaged for the last nanosecond in the PLA samples with and without nanoparticles. It is observed that, indeed, for the sphere and the cone, the concentration is much greater around the nanoparticle, i.e., they produce an oxygen capture region that reduces the diffusion of oxygen molecules. In the case of the cylinder and the flake, an increase is also generated in the capture regions, but they present a dispersion very similar to the one seen for pure PLA.

Figure 4

Color map of the oxygen occupation (average of the last nanosecond of the molecular simulation) in pure PLA and the different forms of the introduced nanoparticles.

Additionally, when the diffusion coefficient is greatly reduced, it is found that oxygen atoms are concentrated in fewer regions than pure PLA, that is, regions with a higher concentration of oxygen. In summary, we observe that oxygen is deposited preferentially around the different nanoparticles and concentrates more than what is seen for pure PLA, which decreases oxygen mobility. This concentration is uniform throughout the nanoparticle, although higher values are observed in the curved regions.

Finally, interaction analysis of the oxygen with the different PLA atoms was carried on. We observe that oxygen is mainly interacting with hydrogen in the PLA chain and other oxygen molecules, as shown in Figure , where the oxygen interactions are shown for a simulation of 150 fs, every 50 fs. The oxygen interactions with hydrogen are normally produced with the methyl group (−CH3), with more than 45,000 interactions in 20 ps, followed by the CH group, with more than 9000 interactions in 20 ps, and the hydroxyl group (−OH), with less than a thousand interactions in 20 ps, as expected from the number of hydrogen atoms per group. The contact duration is similar in the three cases (around 30 fs), producing an effect of reducing the oxygen diffusion by accumulative short-term interactions with the polymer matrix. The interactions between oxygen molecules are very common through the simulation (more than 2 million in 20 ps), with a long duration (around 260 fs), much higher than the interaction with the PLA groups. Therefore, the main cause of oxygen retention is the formation of oxygen clusters that cannot pass through PLA with the same mobility as individual oxygen molecules. In summary, more porosity in the PLA samples could be a favorable method to reduce the diffusion coefficient and increase the barrier properties of PLA.

Figure 5

Oxygen evolution in pure PLA simulation every 50 fs. Color code: hydrogen is represented in white, carbon in gray, and oxygen from PLA in red. Oxygen molecules are represented in blue, and the one that is followed in cyan. Hydrogen interaction and oxygen molecule interactions are represented in black with dashed lines. In the top series, we observe an oxygen interaction and duration with different hydrogen atoms, and in the bottom series, we observe the formation and duration of an oxygen cluster.

Discussion

In the literature, it has been detailed through experiments and simulations how the barrier properties of PLA are improved in the presence of nanocomposites of different materials. Petersson and Oksman[20] describe that the presence of bentonite sheets causes a reduction in oxygen permeation, Sanchez-Garcia and Lagaron[21] observe the same behavior in cellulose nanothreads for water steam permeation, Sepulveda et al.(30) explain that the reduction in oxygen permeation is due to the insertion of silica nanoparticles with cinnamaldehyde, Mulla et al.(14) observe a great reduction in water permeability when ZnO nanoparticles are included in PLA, and Mohammadalinejhade et al.(38) describe that PLA with nanohybrids of silver and lignocellulose reduces the water vapor permeability, due to the highly hydrophobic lignocellulose and the lower mobility of PLA chains due to the nanohybrids’ occupation of the cavities between polymer chains. All of them indicate that the presence of a certain material prevents the passage of oxygen or water through it and generates a tortuous path for the gas. This type of behavior can describe our simulations of greater volume, in which independent of the shape, a reduction in the value of the permeation occurs in a similar percentage. However, the amount of material introduced in this case can be counterproductive and generate an undesired increase in biodegradability, as described by Fukushima et al.(12) and Mulla et al.(14) In smaller nanocomposites such as spherical nanoparticles, Ndoro et al.(53) indicate that the curvature influences the arrangement of the polymers. A surface with smaller curvature (flatter) prevents polymers from distributing homogeneously on their surface, interfering with each other, producing lower polymer density and greater free volume. Merkel et al.(9) and Zhou et al.(29) add that gas permeation is reduced when the spherical nanoparticles have a smaller size and therefore a smaller curvature. Both clarify that the polymer’s density on the surface has increased, reducing the free-volume fraction and reducing the permeation. Finally, Wen et al.(26,27) and Zhou et al.(29) support this theory by studying the effect of nanoparticle concentration on permeation. They observe that more nanoparticles reduce the permeation until the agglomeration is so high that it causes the opposite effect. The agglomeration produces porosity where the polymer does not enter, increasing the fraction of free volume through which the gases move. Also, Singh,[15] in his PhD thesis, describes that chitin nanocrystals and cellulose nanocrystals produce an increment in PLA crystallinity, forming small spherulites that increase the barrier properties more than the bigger ones. Huang et al.(31) observe a similar variation for graphene nanosheets due to the tortuous path produced. In our model, we only study the variation in amorphous PLA, to reduce the number of possible variables that affect the result; however, we expect a major reduction when crystallinity is considered. The study presented in this work on how the shape of the nanoparticles introduced in a PLA polymer matrix affects the barrier properties agrees with what was observed previously. We have observed that the sphere improved, to a larger extent, the barrier properties of PLA compared to the other shapes. It has a lower curvature than the rest of the forms. Therefore, we can assume that this effect is purely caused by the shape of the nanoparticles, and it is not related with their nature, in which our model construction is based.

Additionally, the curvature effect is more evident when studying two similar structures, such as the cylinder and the flake, since the former presents a greater amount of surface with less curvature than the latter with larger flat surfaces. In Figure , we have observed that indeed curved sections are points of agglomeration of oxygen, that is, having a lower fraction of free volume in this region, oxygen is trapped in these areas, accumulating and moving with difficulty. Also, we observe that the interaction with other oxygen molecules is the main factor of oxygen retention. Hence, an accumulation zone with an intermediate porosity is the best way to reduce oxygen diffusion.

Conclusions

In this article, we present a model based on molecular dynamics to study the nanoparticle shape influence on the barrier properties of a PLA matrix, isolating this effect from the nature of the nanoparticles. First, we validate the CHARMM potential for reproducing oxygen diffusion in a PLA matrix. Concerning the shape effect, we have observed that the most significant reductions in the permeation values are obtained by establishing a tortuous path through large blocks of materials, always limited by the biodegradation time of the same, and using nanoparticles with a larger curvature like the spheres. We expect that these findings could be useful in designing new nanocomposite materials with improved barrier properties not only for oxygen molecules but also for water or other gases and as a piece of extra information to describe the difference observed in the experimental results.