Molecular dynamics study on mechanical properties of cellulose with air/nitrogen diffusion behavior

Wang, W., Wang, Y., and Li, X. (2018). "Molecular dynamics study on mechanical properties of cellulose with air/nitrogen diffusion behavior," BioRes. 13(4), 7900-7910.

Abstract

Two hybrid models of cellulose and heat transfer media were designed to simulate changes of cellulose and gas during the heat treatment of wood. The chosen heat transfer media were air and nitrogen, and the range of selected simulation temperature was from 450 to 510 K. The diffusion coefficients of air and nitrogen both increased as the temperature increased, and the diffusion coefficient and amplification of air were larger than that of nitrogen. In addition, the calculation of free volume (fraction) confirmed this conclusion. The mechanical parameters of cellulose were analyzed in terms of the Young’s modulus, Poisson’s ratio, and the ratio of bulk modulus to shear modulus (K/G), which were affected by the temperature and heat transfer medium. Compared with the nitrogen model, the cellulose in the air model decreased in the Young’s modulus and K/G as the temperature increased, while the Poisson’s ratio increased. Young’s modulus and K/G were larger when the heat transfer medium was nitrogen, and the Poisson’s ratio was superior to nitrogen in the air environment. The rate of change of mechanical parameters of cellulose in the air system was greater than that of nitrogen, which was related to the larger diffusion coefficient of air.

Download PDF

Full Article

Molecular Dynamics Study on Mechanical Properties of Cellulose with Air/Nitrogen Diffusion Behavior

Wei Wang,* Yunting Wang, and Xinning Li

Keywords: Wood heat treatment; Molecular dynamics; Diffusion coefficient; Mechanical properties

Contact information: College of Engineering and Technology, Northeast Forestry University, Harbin, 150040, China; * Corresponding author: cedar3@163.com

INTRODUCTION

Wood heat treatment uses air, inert gas, or steam as a heat transfer agent for a short-term pyrolysis treatment at high temperatures between 160 to 250 °C. This process plays an important role in improving wood utilization and quality. Furthermore, it improves the dimensional stability and durability of the wood in wet conditions (Cao et al. 2012; Rowell 2016). The heat treatment process reduces the water absorption capacity by modifying the internal structure and partial components of the wood (Raimo et al. 1996; Gailliot 1998; Sun and Li 2010). Heat-treated wood is better in terms of dimensional stability, shrinkage, and surface quality compared with untreated wood (Kim 2016). In addition, the mechanical properties of cellulose are affected by different heat transfer medium; the properties have been found to decrease less in vacuum (Sibel et al. 2006; Shukla and Sharma 2014). The selection of different heat transfer media during heat treatments has a major influence on the variation of properties. A heat treatment using steam improves the dimensional stability and anti-shrink efficiency (Cao et al. 2012). The mechanical properties of wood treated with inert gas, similar to nitrogen, are lower graded than those treated under a vacuum (Candelier et al. 2013).

Cellulose is the main component of wood and determines its structure and properties. It is a polymeric compound formed by the polymerization of D-glucopyranose, which is composed of crystalline and amorphous regions. Because of its regular arrangement structure, crystalline regions have good mechanical properties even at high temperatures. Small molecules can hardly break the bond between the cellulose chains and do not interact well with crystalline regions (Khazraji and Robert 2013). However, the structure of the amorphous region of cellulose is irregular and easily affected by temperature and small gas molecules, which further affects the mechanical properties of cellulose. In this paper, the amorphous region is chosen as the object to study the wood variations during heat treatment.

The molecular simulation method can be used to construct a model of the amorphous region and to study the diffusion of small gas molecules at different conditions (Miyamoto et al. 2016; Miyamoto et al. 2017; Song et al. 2017). The phenomenon of substance transportation in the system can be called diffusion behavior. It is a spontaneous movement from a high to low concentration, which determines the performance changes of some materials (Najafi and Khademi-Eslam 2011). It is important and difficult to measure the diffusion coefficient of small molecules accurately through experiments. However, through the molecular dynamics simulation, not only can the accurate diffusion coefficient be obtained, but the diffusion behavior of the molecule can be observed from the microscopic level as well. Therefore, the molecular dynamic simulation is considered suitable for the study of diffusing phenomena. Previous studies have established that the behavior of small molecules plays an important role in the performance of the system. The two main diffusion mechanisms of gas in the polymer are the dissolution in polymer phase, and the diffusion in the pores. Diffusion is affected by the mechanism of action, molecular size, volume and other factors (Paunonen 2013). The variation of self-diffusion and inter-diffusion coefficients of hydrogen, carbon monoxide, and water in n-alkanes at different temperatures and pressure has been studied, and the simulation results were in good agreement with the numerical results (Makrodimitri et al. 2011). The size of small molecules, fractional free volume of the system, and the diffusion mechanism at a microscopic level can be better explained by calculating the interaction energy between small molecules and the polymer, which is typically determined by these factors (Hou et al. 2016). It is a general research method to analyze the diffusion of small molecules at the micro level, but there have been few studies on the diffusion behavior of heat transfer media during a wood heat treatment.

At present, wood heat treatment research is mainly focused on the experimental level, and the interaction mechanism between the heat transfer medium and wood has not been precisely explained. With the development of heat treatment technology, it is essential to study the relationship between a heat transfer medium and temperature to improve the quality and availability of wood. In this paper, the effects of air and nitrogen on the mechanical properties of cellulose during heat treatment using molecular dynamics simulation were studied. The change in macroscopic properties was explained from the microscopic level, which provided theoretical support for the wood heat treatment process.

EXPERIMENTAL

Methods

Computational details

Simulations were performed using Materials Studio 8.0 software (Accelrys, USA) to simulate the change of gas molecules and the amorphous polymer. The amorphous region of cellulose was established by constructing the amorphous polymer (Theodorou and Suter 1985). A cellulose chain was created with DP = 20 (Mazeau and Heux 2008). The amorphous cell (AC) tools module was used to simulate the change of cellulose in nitrogen and air environment. The hybrid models of cellulose and heat transfer medium involved polymers and covalent organic small molecules, and the polymer consistent force field (PCFF) had been shown to be very suitable for the calculation of organic compounds and carbohydrate, so PCFF was chosen as the force field. (Sun 1995; Zhang et al. 2014; Wang et al. 2017; Sun et al. 2018). The Ewald method was used as the charge summation method, and the Van der Waals sum method was based on the atomic method. Both the nitrogen model and the air model were established by setting the target density as 1.5 g/cm³ (Brandrup and Immergut 1999), which was similar to cellulose. In the nitrogen model, 20 nitrogen molecules and a cellulose chain with DP = 20 were added. In the air model, 16 nitrogen molecules, 4 oxygen molecules, and a chain of cellulose with DP = 20 were added. The size of this box was 19.6 × 19.6 × 19.6 (A³).

Fig. 1. Cellulose amorphous region model

Dynamic simulation

Followed by the model construction, the initial geometric structure was optimized, as shown in Fig. 1. This was created through calculations in the Forcite module, which is a classical molecular mechanics tool to make geometric optimization, dynamics calculation, mechanical properties calculation, etc. This module supports universal, Dreiding, PCFF, and other force fields. The PCFF was chosen to study the dynamic process as well, and the smart algorithm was selected to calculate 5000 steps during the optimization. The free motion of molecules in the mixed system reached the equilibrium state and minimized its energy. Then, the initial dynamic relaxation of the model was carried out at 300 K in the canonical ensemble (NVT). A random velocity was selected for the initial velocity, and integrated the equation of movement with 1 fs step which resulted in 1 ns of the total simulation length. The output frame was exported every 5000 steps, and the temperature control method was Nose (Shuichi 1991). After relaxation, the initial internal stress in the system decreased and the system became gradually stable, which could be used as the initial structure of the dynamic simulation (Paajanen and Vaari 2017). The temperature of the wood heat treatment was taken in 10 K steps ranging from 450 K to 510 K to observe the effects of the temperature on the mixed system.

The simulation was carried out in the isothermal isobaric ensemble (NPT) with a time of 1 ns because it is enough to equilibrate both systems. The pressure was controlled by the Berendsen method (Berendsen et al. 1984).

RESULTS AND DISCUSSION

Equilibrium of System

In a molecular dynamics simulation, temperature and energy can be used as the basis for the system’s equilibrium (Allen and Tidesley 1989; Sangiovanni et al. 2016.). The system tends to be stable when their values are in the range of 5 to 10% fluctuation. Through the molecular dynamics simulation of 1 ns at various temperatures, the energy fluctuation of the last 200 ps was stable and the fluctuation value was about 9.1%. Therefore, the total energy of this mixed system was considered to be in an equilibrium state after the 1 ns molecular dynamics simulation (Fig. 2).

Fig. 2. Energy – time variation chart

The accuracy of this simulation process can be expressed as an energy convergence parameter (δ E), which is defined in Eq. 1,

(1)

where N is the total step number of the simulation, E₀ represents the initial energy of the simulation system, and E_i denotes the energy of system after the simulation steps (i). When δE ≤ 0.001 ~ 0.003, the simulation results were reliable. Through the calculation, δE in this system was 0.0021, which proved that the energy convergence was reliable.

Figure 3 shows the temperature-time variation of the cellulose-air system at 450 K. As shown, the temperature varied slightly with time and was controlled in the range of ±25 K, which indicated that the system achieved equilibrium during relaxation. The temperature-time curves of this system and nitrogen system at other temperatures were similar to this and all of the systems were in equilibrium.

Fig. 3. Temperature – time variation chart

Coefficient of Diffusion

The mean square displacement (MSD) refers to the mean square sum of the position vector increments of gas molecules after a period of simulation time in the system. The MSD can effectively characterize the migration path and behavior of molecules. It plays an important role in understanding the diffusion changes of gas molecules, represented by Eq. 2,

(2)

where < > denotes the average value, denotes the centroid displacement of the i^th particle and the particle at initial time after a certain simulation time, n is the number of diffusion particles, r is the position of gas molecular, and t is the elapsed time.

The molecular diffusion coefficient D is often used to describe the diffusion behavior of small molecules of gas. Its expression is shown in Eq. 3,

(3)

which can be further simplified in Eq. 4,

(4)

where m is the slope obtained by fitting the MSD curve with the least square method

Based on the analysis of the two models in which either air or nitrogen is the heat transfer medium, the MSD of the two gases at different temperatures were obtained (Fig. 4). The MSD increased linearly with time in the simulation time of 0 to 800 ps, which indicated that the diffusion coefficient of air and nitrogen was constant at the temperature of 450 to 510 K. The slope can be obtained by the origin fitting and then calculated by Eq. 4 to the get diffusion coefficient D at each temperature. We used Origin 8 software to fit the obtained data, and then calculated the goodness of fit by it. (as shown in Table 1). The curves in the range of 900 to 1000 ps presented different degrees of deviations because statistical errors were included in the calculations.

Fig. 4. MSD of gas molecules at different temperatures

Table 1. Diffusion Coefficients of Two Models at Different Temperatures

As shown in Table 1, the goodness of fit was above 0.9, which was reliable. Therefore, the obtained diffusion coefficient can be further analyzed. In the temperature range of 450 to 510 K, the diffusion coefficient of air was greater than that of nitrogen, and both kinds of gas diffusion coefficients increased with a rise of temperature, while the diffusion coefficient of air had a larger increase than that of nitrogen. This is because the molecular dynamics radius of oxygen in the air was smaller than that of nitrogen and the required activation energy was lower. The thermal movement of oxygen was more intense with an increase of temperature than that of the nitrogen, thus affecting the change of the diffusion coefficient.

Free Volume

According to the free volume theory, lattice volume consists of polymer occupied volume and free volume in the form of a “cavity”. The fractional free volume (FFV) reflects the distribution of the free volume. The diffusion coefficient of air and nitrogen in the simulation system is closely related to the free volume, which provides necessary space for the diffusion of molecules. Therefore, analyzing the size and distribution of free volume at the microscopic level is beneficial to understand the difference between the air and nitrogen diffusion coefficients. In this paper, the hard sphere probe method (Swol and Petsev 2014) was used to calculate the free volume and their fractions of both systems with cellulose and the heat transfer medium as air and nitrogen (Fig. 5). The Van der Waals force radius (R_p) was set to 1.52 nm for oxygen atoms and 1.57 nm for nitrogen atoms, respectively.

Fig. 5. Cellulose–air / nitrogen system free volume

The grey part of Figure 5(a) and (b) is the free volume. As shown, the free volume of the air system is larger than for nitrogen. As Fig. 5(a) shows, the free volume of cellulose-air system was 670.53 Å³ and its FFV was 8.9%. In the cellulose-nitrogen system (b), the value is 427.60 Å³ and 5.7%, respectively. It is indicated that the R_p of different gases had an effect on the calculation of free volume, and an increase of the gas’s atomic radius led to the decrease of free volume (fraction), which made the diffusion space of the gas molecules smaller. The phenomenon that the diffusion coefficient of air in cellulose system is larger than that of nitrogen is explained reasonably from the microcosmic point of view.

Mechanical Property

Mechanical parameters characterize the mechanical properties of wood. The cellulose kinetic trajectories of cellulose at various temperatures were calculated to understand the relationship on the macro and micro level. Because the amorphous region of the cellulose could be considered as isotropic, its mechanical properties can be calculated by C11 and C12 in the elastic coefficient matrix of cellulose, as shown in Eq. 5,

(5)

where λ = C12, μ = C11-C12, and λ and μ are called Lame constants.

Lame constants can be used to calculate mechanical properties such as Young’s modulus (E), shear modulus (G), bulk modulus (K), and Poisson’s ratio (γ). The formulas are as follows:

(6)

(7)

(8)

(9)

The parameters of mechanical properties of cellulose-air and nitrogen models at different temperatures were calculated by MD, and the mechanical properties were obtained by the above formulas (Table 2). Poisson’s ratio is the transverse modulus of elasticity, which characterizes the plasticity of materials. A greater value resulted in a stronger plasticity. The Poisson’s ratio of air system at each temperature was higher than that of nitrogen. It increased when the temperature rose, and the rate of increase in the air and nitrogen was 41.4% and 37.8%, respectively. Young’s modulus is the longitudinal elastic modulus, which represents the rigidity of the material. A higher value resulted in a stronger resistance to deformation. The Young’s modulus in the air system was smaller than that in the nitrogen system. It decreased when the temperature increased in both the air and nitrogen environment. The mechanical strength of the macroscopic wood decreased gradually by 20.2% and 16.6%, respectively. The ratio of bulk modulus to shear modulus (K/G) reflects the toughness of materials. A higher value resulted in better toughness. Obtained data in the table shows that the K/G in the nitrogen environment was higher than in the air environment. The K/G value revealed a downward trend with increasing temperature, and the rate of change was not substantial as 2.1% and 1.1%, respectively.

Table 2. Mechanical Parameter of Two Models at Different Temperatures

The value of Young’s modulus and K/G of cellulose was larger in the nitrogen environment while Poisson’s ratio was larger in the air environment in the range of heat treatment temperature. Therefore, if the plasticity of wood is needed, air could be selected as the heat transfer medium. The change of mechanical properties of cellulose in the air system was larger than that of nitrogen; this was because air molecules have a higher diffusion coefficient, which was influential on the change of mechanical properties. In general, the mechanical properties of cellulose are well protected in a nitrogen environment during heat treatment.

CONCLUSIONS

Models of cellulose-air and cellulose-nitrogen were stabled through the MD simulations to explore the property evolution during wood heat treatments.
The diffusion coefficient of air and nitrogen in cellulose increased with an increase of temperature, and the diffusion coefficient of air was larger than that of nitrogen. This was due to the fact that the Van der Waals force radius of oxygen atom in air is larger than that of nitrogen atom, and the activation energy required is lower. High temperatures can make oxygen atoms more active in diffusion. The free volume and fraction of the cellulose system were affected by different Van der Waals force radius of oxygen atom and nitrogen atom in air. Because the kinetic radius of oxygen atom was smaller than that of nitrogen atom and the diffusion space was large, the air diffusion coefficient was larger than that of nitrogen.
The diffusion behavior of the heat transfer medium was shown to be related to the mechanical properties of the cellulose, and the performance of temperature and type of medium on the properties were observed. With the increase of temperature, the Young’s modulus and ratio of bulk modulus to shear modulus (K/G) value of cellulose presented an overall downward trend and Poisson’s ratio showed an entirely rising tendency. The range of mechanical properties varied in different heat transfer media, which was caused by their different diffusion behavior. When the heat transfer medium was air, the rate of change of mechanical properties was larger than that of nitrogen, and the large diffusion coefficient of air had a certain influence on it. When the heat transfer medium was nitrogen, the mechanical properties of cellulose were better protected.

ACKNOWLEDGMENTS

The authors are grateful for the support of the research funds for central universities.

REFERENCES CITED

Allen, M., and Tidesley, D. (1989). Computer Simulation of Liquids, Oxford University Press, Oxford, UK.

Berendsen, H. J. C., Postma, J. P. M., Gunsteren, W. F. V., Dinola, A., and Haak, J. R. (1984). “Molecular dynamics with coupling to an external bath,” Journal of Chemical Physics 81(8), 3684-3690. DOI: 10.1063/1.448118

Brandrup, J., and Immergut, E. H. (1999). Polymer Handbook, Wiley, Hoboken, NJ, USA.

Candelier, K., Dumarçay, S., Pétrissans, A., Gérardin, P., and Pétrissans, M. (2013). “Comparison of mechanical properties of heat treated beech wood cured under nitrogen or vacuum,” Polymer Degradation & Stability 98(9), 1762-1765. DOI: 10.1016/j.polymdegradstab.2013.05.026

Cao, Y., Lu, J., Huang, R., Zhao, X., and Jiang, J. (2012). “Effect of steam-heat treatment on mechanical properties of Chinese fir,” BioResources 7(1), 1123-1133. DOI: 10.15376/biores.7.1.1123-1133

Gailliot, F. P. (1998). “Initial extraction and product capture,” in: Natural Products Isolation, R. J. P. Cannell (ed.), Springer, New York, pp. 53-89.

Hou, Y., Chi, J., Min, L., Liu, Y., Peng, L. I., and Niu, Q. (2016). “Diffusion behavior of the model gasoline components in different polymer membranes studied by molecular dynamic simulation,” Acta Petrolei Sinica 32(6), 1121-1127. DOI: 10.3969/j.issn.1001-8719.2016.06.006

Huang, S., Sun, Y. Z., Xu, Y., and Meng, Z. “Studies on influence of ammonia on properties of cellulose iβ based on molecular dynamics simulation,” Acta Polymerica Sinica 014.2(2014), 188-193. DOI: 10.3724/SP.J.1105.2014.13185

Khazraji, A. C., and Robert, S. (2013). “Interaction effects between cellulose and water in nanocrystalline and amorphous regions: A novel approach using molecular modeling,” Journal of Nanomaterials 2013, Article ID 409676. DOI: 10.1155/2013/409676

Kim, Y. S. (2016). “Research trend of the heat-treatment of wood for improvement of dimensional stability and resistance to biological degradation,” Journal of the Korean Wood Science & Technology 44(3), 457-476. DOI: 10.5658/WOOD.2016.44.3.457

Makrodimitri, Z. A., Unruh, D. J. M., and Economou, I G. (2011). “Molecular simulation of diffusion of hydrogen, carbon monoxide, and water in heavy n-alkanes,” Journal of Physical Chemistry B 115(6), 1429-1439. DOI: 10.1021/jp1063269

Mazeau, K., and Heux, L. (2008). “Molecular dynamics simulations of bulk native crystalline and amorphous structures of cellulose,” Journal of Physical Chemistry B 107(10), 2394-2403. DOI: 10.1021/jp0219395

Miyamoto, H., Schnupf, U., Crowley, M. F., and Brady, J. W. (2016). “Comparison of the simulations of cellulosic crystals with three carbohydrate force fields,” Carbohydrate Research422, 17-23. DOI: 10.1016/j.carres.2016.01.001

Miyamoto, H., Rein, D. M., Ueda, K., Yamane, C., and Cohen, Y. (2017). “Molecular dynamics simulation of cellulose-coated oil-in-water emulsions,” Cellulose 24(7), 2699-2711. DOI: 10.1007/s10570-017-1290-1

Mohammad, A. S., Biernacki, J. J., Northrup, S., and Adenson, M (2018). “Diffusion of CO₂, and fractional free volume in crystalline and amorphous cellulose,” Journal of Analytical & Applied Pyrolysis 134, 43-51.

Najafi, A., and Khademi-Eslam, H. (2011). “Lignocellulosic filler/recycled HDPE composites: Effect of filler type on physical and flexural properties,” BioResources 6(3), 2411-2424. DOI: 10.15376/biores.6.3.2411-2424

Paajanen, A., and Vaari, J. (2017). “High-temperature decomposition of the cellulose molecule: A stochastic molecular dynamics study,” Cellulose 24(7), 2713-2725. DOI: 10.1007/s10570-017-1325-7

Paunonen, S. (2013). “Strength and barrier enhancements of cellophane and cellulose derivative films: A review,” BioResources 8(2), 3098-3121. DOI: 10.15376/biores.8.2.3098-3121

Rowell, R. (2016). “Dimensional stability and fungal durability of acetylated wood,” Drewno 59(197), 139-150. DOI: 10.12841/wood.1644-3985.C14.04

Sangiovanni, D. G., Hellman, O., Alling, B., and Abrikosov, I. A. (2016). “Efficient and accurate determination of lattice-vacancy diffusion coefficients via non equilibrium ab initio molecular dynamics,” Physical Review B 93(9). DOI: 10.1103/PhysRevB.93.094305

Shuichi, N. (1991). “Constant temperature molecular dynamics methods,” Progress of Theoretical Physics Supplement 103, 1-46. DOI: 10.1143/PTPS.103.1

Shukla, S. R., and Sharma, S. K. (2014). “Effect of high temperature processing under different environments on physical and surface properties of rubberwood (Hevea brasiliensis),” Journal of the Indian Academy of Wood Science 11(2), 182-189. DOI: 10.1007/s13196-014-0135-x

Sibel, Y., Engin, D. G., and Umit, C. Y. (2006). “Mechanical and chemical behavior of spruce wood modified by heat,” Building and Environment 41(12), 1762-1766. DOI: 10.1016/j.buildenv.2005.07.017

Song, Y., Jiang, B., and Li, F. L. (2017). “Molecular dynamic simulations of selective self-diffusion of CH₄/CO₂/H₂O/N₂ in coal,” IOP Conference Series: Materials Science and Engineering 213(2017), Article ID 012014. DOI: 10.1088/1757-899X/213/1/012014

Sun, H. (1995). “Ab initio calculations and force field development for computer simulation of polysilanes,” Macromolecules 28(3), 701-712. DOI: 10.1021/ma00107a006

Sun, W. L., and Li, J. (2010). “Analysis and characterization of dimensional stability and crystallinity of heat-treated Larix spp.,” Scientia Silvae Sinicae 46(12), 114-118. DOI: 10.11707/j.1001-7488.20101219

Sun, Y., Chen, L., Cui, L., Zhang, Y., and Du, X. (2018). “Molecular dynamics simulation of cross-linked epoxy resin and its interaction energy with graphene under two typical force fields,” Computational Materials Science 143, 240-247. DOI: 10.1016/j.commatsci.2017.11.007

Swol, F. V., and Petsev, D. N. (2014). “Molecular dynamics simulation of binary hard sphere colloids near the glass transition,” RSC Advances 4(41), 21631-21637. DOI: 10.1039/c4ra02391a

Tang, C., Zhang, S., Wang, Q., Wang, X., and Hao, J. (2017). “Thermal stability of modified insulation paper cellulose based on molecular dynamics simulation,” Energies 10(3), 397. DOI: 10.3390/en10030397

Theodorou, D. N., and Suter, U. W. (1985). “Detailed molecular structure of a vinyl polymer glass,” Macromolecules 18(7), 1467-1478. DOI: 10.1021/ma00149a018

Wang, X., Tang, C., Wang, Q., Li, X., and Hao, J. (2017). “Selection of optimal polymerization degree and force field in the molecular dynamics simulation of insulating paper cellulose,” Energies 10(9), 1377. DOI: 10.3390/en10091377

Wu, X. W., Moon, R. J., and Martini, A. (2014). “Tensile strength of I β, crystalline cellulose predicted by molecular dynamics simulation,” Cellulose 21(4), 2233-2245. DOI: 10.1007/s10570-014-0325-0

Zhang, S., Li, X., and Tang, C. (2014). “Analysis on the thermal stability of the amorphous region in insulation paper,” International Conference on Mechatronic Sciences, Electric Engineering and Computer IEEE, 2014, pp. 3581-3585. DOI: 10.1109/MEC.2013.6885622

Zhao, Y., Liu, X., Wang, J., and Zhang, S. (2013). “Effects of anionic structure on the dissolution of cellulose in ionic liquids revealed by molecular simulation,” Carbohydrate Polymers94(2), 723-730. DOI: 10.1016/j.carbpol.2013.02.011

Article submitted: June 8, 2018; Peer review completed: August 17, 2018; Revised version received and accepted: August 22, 2018; Published: August 30, 2018.

DOI: 10.15376/biores.13.4.7900-7910