journal of the mechanical behavior of biomedical materials 42 (2015) 198 –206
Available online at www.sciencedirect.com
www.elsevier.com/locate/jmbbm
Research Paper
Molecular deformation mechanisms of the wood cell wall material Kai Jin, Zhao Qin, Markus J. Buehlern Laboratory for Atomistic and Molecular Mechanics (LAMM), Department of Civil and Environmental Engineering, Massachusetts Institute of Technology, 77 Massachusetts Ave., Cambridge, MA 02139, USA
ar t ic l e in f o
abs tra ct
Article history:
Wood is a biological material with outstanding mechanical properties resulting from its
Received 18 July 2014
hierarchical structure across different scales. Although earlier work has shown that the
Received in revised form
cellular structure of wood is a key factor that renders it excellent mechanical properties at
7 November 2014
light weight, the mechanical properties of the wood cell wall material itself still needs to be
Accepted 13 November 2014
understood comprehensively. The wood cell wall material features a fiber reinforced
Available online 21 November 2014
composite structure, where cellulose fibrils act as stiff fibers, and hemicellulose and lignin
Keywords:
molecules act as soft matrix. The angle between the fiber direction and the loading
Wood cell wall
direction has been found to be the key factor controlling the mechanical properties.
Molecular dynamics
However, how the interactions between theses constitutive molecules contribute to the
Yielding of matrix
overall properties is still unclear, although the shearing between fibers has been proposed
Slip-stick mechanism
as a primary deformation mechanism. Here we report a molecular model of the wood cell wall material with atomistic resolution, used to assess the mechanical behavior under shear loading in order to understand the deformation mechanisms at the molecular level. The model includes an explicit description of cellulose crystals, hemicellulose, as well as lignin molecules arranged in a layered nanocomposite. The results obtained using this model show that the wood cell wall material under shear loading deforms in an elastic and then plastic manner. The plastic regime can be divided into two parts according to the different deformation mechanisms: yielding of the matrix and sliding of matrix along the cellulose surface. Our molecular dynamics study provides insights of the mechanical behavior of wood cell wall material at the molecular level, and paves a way for the multiscale understanding of the mechanical properties of wood. & 2014 Elsevier Ltd. All rights reserved.
Introduction
bearing systems, for example, houses, bridges, ships and even aircrafts, since it also possesses excellent mechanical
Wood, with its broad availability, has served as an important bulk material since ancient times and still plays a major role today. One important application area of wood is in the load n
Corresponding author. Tel.: þ1 617 452 2750; fax: þ1 617 324 4014. E-mail address: [email protected] (M.J. Buehler). URL: http://web.mit.edu/mbuehler/www/ (M.J. Buehler).
http://dx.doi.org/10.1016/j.jmbbm.2014.11.010 1751-6161/& 2014 Elsevier Ltd. All rights reserved.
properties. In the sense of specific mechanical properties (normalized by mass density), wood material parallels and even behaves better than many other engineering materials
journal of the mechanical behavior of biomedical materials 42 (2015) 198 –206
such as steel (Ashby et al., 1995; Gibson and Ashby, 1999; Gibson et al., 1995). In fact, the origin of most woods, the tree, is also a load bearing system, where the smart design of the wood structure renders a tree capable of supporting its own weight and sustaining various climate conditions (Fratzl and Weinkamer, 2007; Qin et al., 2014). The unique combination of composite material and microstructure design in wood material has attracted research interests in comprehensive understanding the mechanical properties of wood and designing inspired materials and structures (Compton and Lewis, 2014; Gibson and Ashby, 1999; Gibson et al., 2010). Similar to other outstanding biological materials, for example, spider silk (Eisoldt et al., 2011; Keten et al., 2010), bone (Launey et al., 2010; Nair et al., 2013; Weiner and Wagner, 1998) and sponge (Aizenberg et al., 2005; Meyers et al., 2013), wood also features a hierarchical structure (Ali and Gibson, 2013) (Fig. 1). At the macro-scale, we can clearly see the annual rings in a cross-section of a wood trunk (Fig. 1a), which is a result of the different growing rates during different seasons (Fratzl and Weinkamer, 2007). Zooming into the micrometer scale, we can clearly see that the wood material is not solid but features a cellular structure (Fig. 1b, c), where hollow space is enclosed by wood cell walls. Studies have shown that the cellular structure is the underlying feature that endows wood as an excellent mechanical material: it makes wood light weight and perform efficiently when it is used as engineering components (Ashby et al., 1995; Gibson and Ashby, 1999; Gibson et al., 1995; Vural and Ravichandran, 2003). The significant hollow space in wood also gives it high deformability. For example, cork is an excellent material for energy absorption and bottle bungs (Gibson and Ashby, 1999).
199
As shown in Fig. 1c-e, the wood cell wall material is a composite, combing stiff and soft polymer components. Studies have shown that the wood cell wall is composed of multiple layers and it is a kind of fiber-reinforced composite with stiff cellulose fibrils embedded in soft matrix consisting of hemicellulose and lignin (Booker and Sell, 1998; Fahlen and Salmen, 2002) (Fig. 1c, d). Among these layers, the middle layer of the secondary cell wall (S2 layer), which occupies more than 80% of the thickness of the cell wall, is the thickest and mainly controls the mechanical properties of the wood cell wall (Booker and Sell, 1998). A unique feature of S2 layer is that the stiff cellulose fibrils lie parallel with each other and form an angle (micro-fibril angle, MFA) with the longitudinal axis of wood cell (Fig. 1d) (Booker and Sell, 1998). As the longitudinal direction is the main load-bearing direction in wood, the mechanical properties of wood mainly relate to the cell wall's properties along this direction. Since MFA measures how close are the stiff cellulose fibrils aligned with this load-bearing direction, it plays a key role in controlling the mechanical properties of wood cell wall: smaller MFA indicates that the stiff cellulose fibrils are more aligned with the load-bearing direction, and the material will behave stiffer. (Booker and Sell, 1998; Navi et al., 1995; Reiterer et al., 2001; Reiterer et al., 1999). Since the stiff cellulose fibrils spirally wind the wood cell (Booker and Sell, 1998), tensile loading along the longitudinal direction will result in the unwinding of the fibrils and the decreasing of the MFA (Keckes et al., 2003). Theoretical and simulation works have proposed that this reorientation of the fibrils generates relative shearing between adjacent fibrils, and that the shear deformation on the matrix is the underlying mechanical response that is responsible for the mechanical
Fig. 1 – The hierarchical structure of wood material. (a) A drawing shows the annual rings structure at the cross-section of a wood trunk. (b) The cellular structure of wood (reproduced from (Vural and Ravichandran, 2003), with permission from Elsevier). (c) Multiple layers in wood cell wall (reproduced from (Booker and Sell, 1998), with permission from Springer). M: middle layer, P: primary layer, S1: the first secondary layer, S2: the second secondary layer, S3: the third secondary layer. (d) The fiber reinforced structure of wood cell wall. The stiff fibrils form an angle θ (micro-fibril angle, MFA) with the longitudinal direction of cell. (e) The molecular model used in the present study. Black arrows indicate the crystal orientations of cellulose.
200
journal of the mechanical behavior of biomedical materials 42 (2015) 198 –206
behavior of wood cell wall (Adler and Buehler, 2013; Altaner and Jarvis, 2008; Burgert, 2006; Keckes et al., 2003). However, the studies of the mechanical behaviors of the constituting molecules are out of current capability of experiments. Some earlier simulation work (Adler and Buehler, 2013) explored the mechanics of wood using a coarse-grained model, but lacked atomistic scale resolution and chemical details of the constituents. This motivates us to apply molecular dynamics (MD) in the present study, as it enables us to build a bottom-up model starting from the chemical structure of the basic constituents of wood. Because studies have revealed the micro-structural design of wood cell wall material (Salmén, 2004) and the force field parameters for polysaccharides (Guvench et al., 2008; Guvench et al., 2009; Hatcher et al., 2009; Raman et al., 2010) and lignin (Petridis and Smith, 2009) have been developed, we are able to build a molecular model of wood cell wall material and perform MD simulations to achieve the goal of this study. A layered molecular model of the wood cell wall is built (Fig. 1e) and shear loading, as suggested earlier, is applied to the model. By analyzing the stress-strain response of the model and simulation trajectories during deformation, this study aims to recognize the deformation mechanisms at the molecular level.
Model A general description As shown in Fig. 1e, the molecular model in the present work includes the three major components of wood cell wall at atomistic resolution: cellulose, hemicellulose and lignin, among which cellulose has the highest mass ratio (Lewin and Goldstein, 1991; Terashima et al., 2009). Although these three components constitute the basic structure of wood cell wall material, water molecules, another factor that it not included here, act as a regulator of the material's mechanical properties. Studies have shown that the bound water within the cell wall can soften the material (Skaar, 1988). Since the bound water in the cell wall can exist at the cellulose surface and within the matrix (Skaar, 1988), the softening effects are expected to be lubricating the cellulose-matrix interface and softening the matrix. In order to comprehensively unravel water's role from these two different aspects, simulations are needed to study water's effects on matrix and fibril-matrix interface separately and it will also involve a study of water's distribution within the cell wall. These studies are out of the scope of our present modelling work. We limit our exploration in the present work on the deformation mechanisms and all the simulations are performed on the dry model as show in Fig. 1e. The cellulose fibril has a crystalline structure that has been revealed by experiments (Nishiyama et al., 2002; Nishiyama et al., 2003). The Iβ phase cellulose is chosen in the present model since it is rich in wood (Horii et al., 1987; Sugiyama et al., 1991). Linear cellulose chains, in which glucose rings are covalently connected, are held together through inter-chain hydrogen bonds (H-bonds) to form cellulose layers and these layers stack in [100] direction (Nishiyama et al., 2002). The crystalline structure and the covalent bonding along the chain
direction render the cellulose fibrils high modulus. Experiments and simulation works have reported the modulus along the chain direction ranging from 150 to 200 GPa (Diddens et al., 2008; Dri et al., 2013; Kulasinski et al., 2014; Tanaka and Iwata, 2006; Tashiro and Kobayashi, 1991). Different morphologies have been proposed for the cellulose fibrils, including square and hexagonal shape (Ding and Himmel, 2006; Fernandes et al., 2011; Terashima et al., 2009). This means that the cellulose fibrils are able to interact with matrix molecules on different crystalline planes. In the present model, the (110) crystalline surfaces (Charlier and Mazeau, 2012; Mazeau, 2011) are chosen to be the fibril surfaces. The size of cellulose fibrils ranges from several nanometers to tens of nanometers (Sinko et al., 2014), which makes it too expensive to include a full cellulose fibril in the modelling. Considering that the deformation of cellulose fibril is much smaller than that of matrix molecules because of its high stiffness, it contributes little to the overall deformation and it is not necessary to model a full cellulose fibril. As shown in Fig. 1e, only two layers (thickness is approximately 1.1 nm) of cellulose in (110) plane are chosen to represent each cellulose fibril. This two-layer setting meets the minimum requirements: it provides one layer to be added loading condition (restricted movement) and the other layer to interact with matrix molecules without external restriction. The thickness, which is beyond the cutoff for non-bond interaction, prevents the matrix molecules interact with vacuum space. Therefore, this simplification significantly reduces the computation consumption while it is still capable of exploring the deformation response. Among the two soft matrix components, hemicellulose directly attaches on the cellulose surface (Lewin and Goldstein, 1991). Although hemicellulose has similar composition as cellulose, it is in an amorphous state (Lewin and Goldstein, 1991; Terashima et al., 2009). Lignin molecules fill the remaining space between hemicellulose molecules. The thickness of matrix between cellulose fibrils is proposed to be approximately 6nm (Terashima et al., 2009). Together with the volume ratio between hemicellulose and lignin, which is proposed to be 2:3, and the density of each component (1.4 g/ cm3 and 1.36 g/cm3 respectively) (Terashima et al., 2009), the amount of both hemicellulose and lignin in the model can be calculated. Due to the limitation of computation consumption, the modelling work does not attempt to model a fibril surrounded by matrix but instead focuses on the behavior of matrix between adjacent fibrils. Periodic boundaries are applied in both directions of the cellulose layer. Thus, the present model represents actually two infinite cellulose slabs separated by an infinite layer of matrix molecules.
Building the molecular model The first step of building the model is to generate the (110) cellulose surfaces (Fig. 2a). The crystalline structure revealed by synchrotron X-ray and neutron fiber diffraction (Nishiyama et al., 2002) is used to build a cellulose block. Cellulose unit cells are copied to form a block of 4 layers of cellulose in the (110) plane. Each layer has 6 chains and each chain consists of 10 covalently bonded glucose rings. Periodic boundaries are applied in all directions and each chain is covalently connected across the boundary in the chain extending direction. Equilibration of the block is performed: 1 ns under NPT ensemble
journal of the mechanical behavior of biomedical materials 42 (2015) 198 –206
201
Fig. 2 – Overview of the model building procedure, consisting of creating (110) cellulose crystal surfaces (a), and then adding hemicellulose (b) and lignin (c) molecules. (d) A dynamic bond-forming procedure of generating linkages between lignin units is applied subsequently. The green dash line indicates the formed bond. (e) After the building of the model, it is loaded by moving the outer-most cellulose layers in opposite directions (indicated by flat green arrows). (f) The lignin units and linkages included in the present model. The reactive atoms that can form bonds are highlighted by circles. (g) The density profiles of the three components across the model's thickness direction.
(1 a.t.m. and 293 K) and another 1 ns under NVT ensemble (293 K). Once the cellulose block is equilibrated, it is cleaved into two parts, each has two (110) layers, and moved apart (Fig. 2a). In this way, two (110) cellulose surfaces facing each other are created and the space between them will be filled by matrix molecules in the following steps. The second step is to add hemicellulose molecules onto the cellulose surfaces (Fig. 2b). Although hemicellulose has various composing units (Lewin and Goldstein, 1991), we choose the Xylan type cellulose in the present study. Significant disparity is not expected for other types of cellulose, since they share similar structure. For each cellulose surface, in total 22 hemicellulose segments are added. These include 12 segments with 5 D-xylopyranosyl residues covalently connected and another 10 segments with additional L-Arabinofuranosyl (Araf) linked to the third the xylosyl residues through α-(1-3) linkage (Charlier and Mazeau, 2012). The Araf residues are used for the further connection with lignin molecules. The segments are added sequentially with random orientation and an equilibration is performed after every adding step (Charlier and Mazeau, 2012; Mazeau and Charlier, 2012). Additional equilibrations under high temperature and following low temperature are performed to obtain the final hemicellulose phase. The final step is to generate the lignin molecules in the model (Fig. 2c-e). For each hemicellulose layer, 3 of the 10 Araf residues are selected to be connected with ferulate moieties that bridging the hemicellulose with lignin (Charlier and Mazeau, 2012). Lignin molecules have various composing units
and linkages between the units (Lewin and Goldstein, 1991; Pu et al., 2008; Vanholme et al., 2010). In the present model, the two major units, G and S unit, and part of the major linkages (shown in Fig. 2f) are included. It is believed that lignin has very complex 3D structure resulting from its various types of units and linkages and possible branching structure. Experiments on the structure of lignin usually provide the ratios of different units and linkages (Pu et al., 2008). In contrast to applying pre-defined sequence of a lignin molecule (Charlier and Mazeau, 2012; Petridis and Smith, 2009), we apply a dynamic bond-forming procedure to build the amorphous lignin structure. This method follows the bio-synthesis pathway of lignin, in which the lignin units are oxidized to form phenolic radicals and the radicals can further be bonded through forming the linkages (Vanholme et al., 2010). In order to reduce the complexity, only the dominant β-O-4 linkages (Pu et al., 2008; Vanholme et al., 2010) are formed in the procedure. Other types of linkages are introduced by adding corresponding dimers (Fig. 2f) into the system. This simplification can also precisely control the ratios of different units and linkages by controlling the number of units and dimers added into the system and the number of formed β-O-4 linkages. Prior to the bond-forming procedure, the units and dimers are set such that the Cβ and O4 atoms (Fig. 2d, f) can form bonds. They are packed into the system sequentially with random positions and orientations (Abbott et al., 2013) and an equilibration with high temperature (800 K) is then performed to thoroughly distribute the units. The tool developed by L. J. Abbott et al. (Abbott et al., 2013) is applied to accomplish this
202
journal of the mechanical behavior of biomedical materials 42 (2015) 198 –206
bond-forming procedure. This tool mimics a polymerization procedure that generates polymer molecules from input monomers. It is noted that when the cellulose block is separated, the space between them is larger than the anticipated 6nm so that the lignin units can have more space to adjust their positions and orientations in order to relax during the bond-forming procedure. An energy minimization is performed after each bond-forming and once every five bonds are formed an equilibration under NVT ensemble with high temperature (800 K) is performed in order to sample the structure. Additional equilibration is performed after the forming of desired number of β-O-4 linkages and then the model is compressed to the desired thickness. In order to study the effects of different lignin structures, a bunch of models with different lignin structures: different units’ ratios (S: G¼ 0: 1, 1: 1, and 1.2: 1), different linkages’ ratios ((β-O-4): (α-O-4): (β-5): (55): (5-O-4) ranges from 50%: 12.5%: 12.5%: 12.5%: 12.5% to 76%: 6%: 6%: 6%: 6%), with and without branching in molecule, different molecule weight (average molecule size ranges from consisting of 6.5 monomers to consisting of 65 monomers), are built and studied. Once the model is built according to the above steps, it is equilibrated in order to reduce the artificial effects introduced during the building procedure. The shear loading is added by moving the outer-most cellulose layers in opposite directions (Fig. 2e). We assume that the two-layer cellulose fibril here represents a fibril of about 12 nm in thickness (Terashima et al., 2009). Thus, the equivalent distance between two adjacent cellulose fibrils is approximately 18 nm (counting half thickness of each cellulose fibril and the thickness of the matrix) and a shear strain can be calculated as dividing the relative displacement of the two moving cellulose fibrils by 18 nm. The stress in the matrix can also be recorded and further we can obtain the stress-strain response on the model. All the simulations are performed by the large-scale atomic/molecular massively parallel simulator (LAMMPS) (Plimpton, 1995). The timestep is chosen as 1 fs and the cutoff for the non-bonded interactions is 1 nm. The particle–particle particle–mesh (PPPM) method is used to compute the long-range columbic interactions and the SHAKE algorithm is applied to constrain highfrequency dynamics from hydrogen-related energy terms.
these two types of molecules. In general, the profiles here are in agreement with that in literature (Charlier and Mazeau, 2012). However, the density profile of lignin here is smoother than that suggested in the literature (Charlier and Mazeau, 2012). This is a result of applying different lignin building methods. The random packing and equilibration before the bond-forming procedure in the present study can provide more uniform distribution of lignin mass than that built by the predefined sequence method (Charlier and Mazeau, 2012). An average value of the lignin density can be calculated to be approximately 1.2 g/cm3, which is close to the experimental value of 1.36 g/cm3 (Terashima et al., 2009).
The elastic-plastic behavior Fig. 3 shows a representative stress-strain curve of the model response. The curve shows elastic behavior at the initial regime and yielding at the shear strain of approximately 0.03. This value can be applied into a theoretical model (Keckes et al., 2003) to predict the tensile yielding strain of wood cell wall materials. Keckes et al.'s study (Keckes et al., 2003) reported the stress-strain curve of a wood cell wall with initial MFA of about 46 degrees and it yielded at the tensile strain of 0.015. Our simulation result, accompanied with the theoretical model, predicts that it happens at the strain of 0.016, which agrees quite well with the experimental result. After yielding, the stress continues to gradually increase with occasional drops. In order to show the plastic deformation, unloading and reloading are performed. As shown in the figure, the strain cannot return back to the original position after unloading, indicating that irreversible deformation after yielding has occurred. During the reloading procedure, the initial elasticity recovers and once the stress reaches the plastic level, the yielding behavior restarts again. At a strain of approximately 0.18, the stress shows a significant drop, and at larger deformation the stress evolves in a “saw-tooth” shape with large amplitude fluctuations. Although this is an extension of the post-yielding behavior, the significant fluctuations indicate that the underlying deformation mechanism might be different from the one that is responsible for the deformation before the strain of 0.18 (it will be discussed in more depth in the next section).
Results and discussion The molecular model The final thickness of the matrix after equilibration is calculated as approximately 5.7 nm, which is in agreement with the desired value of 6 nm (Terashima et al., 2009). As shown in Fig. 2g, the density profiles are calculated across the model's thickness direction for each component. The curves of these components are roughly separated, showing the layered structure of the model: hemicellulose layers lying on cellulose surface and then lignin molecules fill the remaining space. There are two peaks for each cellulose part, which is a result of the crystalline layered structure. The density profile of hemicellulose also shows peaks but they are not intense, showing it is in amorphous state. The overlap regions between the profiles of hemicellulose and lignin indicate the interpenetration of
Fig. 3 – A representative stress-strain response (blue curve) of the model. Three regimes are distinguished: one initial elastic regime and two plastic regimes with different underlying mechanisms. The unloading (red curves) and reloading (green curves) tests show irreversible deformation after yielding.
journal of the mechanical behavior of biomedical materials 42 (2015) 198 –206
Since the stress level before and after the strain of 0.18 does not show significant change, it cannot show clear difference in the experimental stress-stain curves to identify when the mechanism-transition happens. Thus, this mechanismtransition strain cannot be verified by experimental data. Actually, this is also a proof that shows the advantages of molecular simulations over experiments to explore the molecular level deformation mechanisms. Unloading and reloading tests are also performed in this regime. Similarly, the non-zero crossing at the strain axis indicates irreversible deformation and the elasticity recovers during the reloading. Once the stress reaches the yielding stress, the “saw-tooth” behavior is observed again. According to the above description, the stress-strain response can be divided into three regimes: an initial elastic regime and two plastic regimes that are separated by the strain of approximately 0.03 and 0.18 (indicated by the dash lines in Fig. 3).
Deformation mechanisms In order to distinguish the underlying mechanisms for different deformation regimes, snapshots of the model under loading are analyzed. Fig. 4a-c show the snapshots at the initial state, the first plastic regime and the second plastic regime, respectively. The color contours correspond to the initial coordinates in the loading direction. As shown in the figures, the color contour inclines under loading, which indicates shear deformation is generated on the matrix. In contrast, the color contours within the cellulose parts keep vertical, indicating that the adjacent cellulose layers stick with each other firmly during deformation. This result verifies our assumption of the small deformation of cellulose and indicates that the simplification of using two-layer cellulose to
203
represent a thicker cellulose fibril is reasonable. Although the inclined color contours indicate shearing in both plastic regimes, there is a difference of the displacement continuity at the matrix-cellulose interface: while it is continuous during the first plastic regime (Fig. 4b), the displacement shows a jump at the matrix-cellulose interface during the second plastic regime (Fig. 4c). This difference is most likely a result of distinct deformation mechanisms. For the first plastic regime, the irreversible deformation after yielding should be attributed to the deformation within the matrix, since the matrix keeps sticking with the cellulose surface. A close observation of the matrix molecules shows some extent of irreversible reorganizing of matrix molecules. Fig. 4d shows an example in which one molecule squeezed between two other molecules after yielding and cannot return back to the original position when the model is unloaded. Therefore, the mechanism that is responsible for the first plastic regime is the yielding of matrix molecules. In contrast, the discontinuous displacement during the second plastic regime indicates that the matrix molecules slide along the cellulose surface. The movie of the deformation procedure shows that the sliding does not happen in a smoothly continuous way but with a “slip-stick” manner. The slip of matrix happens in a very short time that the relative displacement of the two moving cellulose and further the shear strain nearly keeps the same value. This results in the sudden drops of stress, which are nearly vertical in the stress curve with respect to shear strain (Fig. 3). After the slip event, the matrix will immediately stick with the cellulose surface at a new position and the stress starts to increase again. As mentioned in the model section, hemicellulose has similar units and structure with cellulose. The rich hydroxyl groups in both hemicellulose and the (110) cellulose surface help to
Fig. 4 – (a–c) The snapshots at the initial state, during the first plastic regime and during the second plastic regime, respectively. The color contours correspond to the initial coordinates in the loading direction (indicated by green arrows). (d) Snapshots show the reorganization of molecules within the matrix after yielding. (e) The dynamic H-bonding (the short red line) between hemicellulose and cellulose. (f) The number of H-bonds at the matrix-cellulose interface with respect to the shear strain on the model.
204
journal of the mechanical behavior of biomedical materials 42 (2015) 198 –206
construct a hydrogen bonding network, an important portion of the non-boned interactions at the matrix-cellulose interface. Therefore, the stress-strain curve can be connected with the dynamics of the interfacial H-bonds. In every slip event, the matrix molecules draw back a little bit and the former Hbonds network is broken. This is followed by a reforming of the H-bonds network in a new position once the matrix sticks with the cellulose surface again. The reformed H-bonds help the interface to regain the stiffness, and the stress thus restarts to increase. Fig. 4e shows an example of the dynamic H-bonds at the interface: a hydroxyl group of hemicellulose sequentially forms H-bonds with the acceptor oxygen atoms in cellulose. The number of H-bonds at the interface is recorded with respect to the shear strain (Fig. 4f), where the big amplitudes of “saw-tooth” indicate the dynamic breaking and reforming procedure. Since the matrix draws back, it is unloaded at the slip events and the following stick reloads the matrix again. As an indication, the slope of the “sawtooth” during the second plastic regime is close to the initial elasticity. Meanwhile, the actual deformation of the matrix is limited because of the drawing-back and thus the sliding of matrix along the cellulose surface is the dominant mechanism during the second plastic regime.
In sum, two different mechanisms are recognized for the plastic behavior of the model. Although the yielding of matrix results into the first plastic regime, the “slip-stick” mechanism is dominant during the second plastic regime. In fact, this kind of “self-healing” interface is a common strategy utilized by nature to achieve excellent mechanical properties such as the dynamic H-bonds between collagen and minerals (Qin et al., 2012). This strategy can also be applied in the interfacial functionalization for synthetic high-performance materials (Nair et al., 2012; Naraghi et al., 2013) and potentially many other areas of bio-inspired materials design.
Effect of lignin structure As mentioned in the Model section, models with different lignin structures are built and studied in order to understand its effect on the mechanical response of the model. Fig. 5 shows the stress-strain curves of 9 cases with different lignin structures. As is apparent from the data, all results share the above presented 3-regime behavior; one initial elastic regime and two plastic regimes with different underlying mechanisms, although slight difference exists between each case. In fact, the dominant interaction within the lignin matrix is the squeezing and separation between atoms and the effects of structural level difference are not apparent in the mechanical behavior.
Effect of loading rate
Fig. 5 – Stress-strain responses on models with different lignin structures. All responses reproduce the three-regime behavior. The shadow area indicates the range of the first sliding event for different models.
It is noted that the effect of high-strain rate ( 108 per second in the present model) is unavoidable due to the limitation of computational resources. Two cases with lower strain rate ( 0.2 108) are performed to reveal the effects of loading rate. The stress-strain response of one of them is shown in Fig. 6a, accompanied with the case of the same lignin structure but with higher strain rate. The results show that the 3-regime deformation procedure preserves for the case of lower strain rate. The main difference lies in the initial elasticity and the stress level during the first plastic regime.
Fig. 6 – (a) A comparison of stress-strain responses with different loading rates. (b) The average initial modulus of the nine cases with different lignin structures (the red error bar indicates the standard error) in comparison with the initial modulus of the two cases (indicated by different colors) with lower strain rate.
journal of the mechanical behavior of biomedical materials 42 (2015) 198 –206
Obviously, the matrix behaves stiffer under higher loading rate. Fig. 6b shows the average initial modulus of the 9 cases shown in Fig. 5 and the values of the two cases with lower strain rate. It has an approximate 38% drop from 2.89 GPa. The difference during the second plastic regime is not obvious, probably because of the big fluctuations. Although it is not reliable to extract the stress-strain response quantitatively, the modelings here provide the insights of the molecular-level deformation mechanisms qualitatively.
Conclusion and outlook Through the study with a molecular model of the wood cell wall material, at atomistic resolution, several deformation mechanisms at the molecular level are identified. The yielding of the matrix and sliding of matrix along cellulose surface, which start at the shear strain of approximately 0.03 and 0.18 respectively, sequentially contribute to the post-yielding deformation of wood cell wall material at molecular level. The present MD simulation work provides foundation for both the further study on water's effects and further coarse-graining modelling of wood cell wall material, which could enable us to study the effects of various factors, for example, MFA, fibril length, fibril mass ratio, defects in the material, on the mechanical properties at the larger mesoscale. Further investigation can also include experiments on the inspired wood cell wall structure created using 3D printing techniques.
Acknowledgment We acknowledge the support from the BASF-NORA program, and helpful discussions with Professor Lorna Gibson (MIT).
r e f e r e n c e s
Abbott, L.J., Hart, K.E., Colina, C.M., 2013. Polymatic: a generalized simulated polymerization algorithm for amorphous polymers. Theor. Chem. Acc. 132, 1–19. Adler, D.C., Buehler, M.J., 2013. Mesoscale mechanics of wood cell walls under axial strain. Soft Matter 9, 7138. Aizenberg, J., Weaver, J.C., Thanawala, M.S., Sundar, V.C., Morse, D.E., Fratzl, P., 2005. Skeleton of euplectella sp.: Structural hierarchy from the nanoscale to the macroscale. Science 309, 275–278. Ali, Z.M., Gibson, L.J., 2013. The structure and mechanics of nanofibrillar cellulose foams. Soft Matter 9, 1580. Altaner, C.M., Jarvis, M.C., 2008. Modelling polymer interactions of the molecular velcro type in wood under mechanical stress. J. Theor. Biol. 253, 434–445. Ashby, M.F., Gibson, L.J., Wegst, U., Olive, R., 1995. The mechanical properties of natural materials. I. Material property charts. P. Roy. Soc. A: Mat. 450, 123–140. Booker, R.E., Sell, J., 1998. The nanostructure of the cell wall of softwoods and its functions in a living tree. Holz. Roh. Werkst. 56, 1–8. Burgert, I., 2006. Exploring the micromechanical design of plant cell walls. Am. J. Bot. 93, 1391–1401. Charlier, L., Mazeau, K., 2012. Molecular modeling of the structural and dynamical properties of secondary plant cell walls: influence of lignin chemistry. J. Phys. Chem. B 116, 4163–4174.
205
Compton, B.G., Lewis, J.A., 2014. 3D-Printing of lightweight cellular composites. Adv. Mater. 26, 5930–5935. Diddens, I., Murphy, B., Krisch, M., Muller, M., 2008. Anisotropic elastic properties of cellulose measured using inelastic X-ray scattering. Macromolecules 41, 9755–9759. Ding, S.Y., Himmel, M.E., 2006. The maize primary cell wall microfibril: a new model derived from direct visualization. J. Agr. Food Chem. 54, 597–606. Dri, F.L., Hector, L.G., Moon, R.J., Zavattieri, P.D., 2013. Anisotropy of the elastic properties of crystalline cellulose Iβ from first principles density functional theory with Van der Waals interactions. Cellulose 20, 2703–2718. Eisoldt, L., Smith, A., Scheibel, T., 2011. Decoding the secrets of spider silk. Mater Today 14, 80–86. Fahlen, J., Salmen, L., 2002. On the lamellar structure of the tracheid cell wall. Plant Biol 4, 339–345. Fernandes, A.N., Thomas, L.H., Altaner, C.M., Callow, P., Forsyth, V.T., Apperley, D.C., Kennedy, C.J., Jarvis, M.C., 2011. Nanostructure of cellulose microfibrils in spruce wood. Proc. Natl. Acad. Sci. U.S.A. 108, E1195–E1203. Fratzl, P., Weinkamer, R., 2007. Nature’s hierarchical materials. Prog. Mater. Sci. 52, 1263–1334. Gibson, L.J., Ashby, M.F., 1999. Cellular solids: structure and properties, 2nd ed. Cambridge University Press, Cambridge. Gibson, L.J., Ashby, M.F., Harley, B.A., 2010. Cellular materials in nature and medicine, 1st ed. Cambridge University Press, Cambridge. Gibson, L.J., Ashby, M.F., Karam, G.N., Wegst, U., Shercliff, H.R., 1995. The mechanical properties of natural materials. II. Microstructures for mechanical efficiency. P. Roy. Soc. A: Mat. 450, 141–162. Guvench, O., Greene, S.N., Kamath, G., Brady, J.W., Venable, R.M., Pastor, R.W., Mackerell Jr., A.D., 2008. Additive empirical force field for hexopyranose monosaccharides. J. Comput. Chem. 29, 2543–2564. Guvench, O., Hatcher, E.R., Venable, R.M., Pastor, R.W., Mackerell, A.D., 2009. CHARMM additive all-atom force field for glycosidic linkages between hexopyranoses. J. Chem. Theory Comput. 5, 2353–2370. Hatcher, E., Guvench, O., MacKerell, A.D., 2009. CHARMM Additive all-atom force field for aldopentofuranoses, methylaldopentofuranosides, and fructofuranose. J. Phys. Chem. B 113, 12466–12476. Horii, F., Hirai, A., Kitamaru, R., 1987. Cp/Mas C-13 Nmr-spectra of the crystalline components of native celluloses. Macromolecules 20, 2117–2120. Keckes, J., Burgert, I., Fruhmann, K., Muller, M., Kolln, K., Hamilton, M., Burghammer, M., Roth, S.V., Stanzl-Tschegg, S., Fratzl, P., 2003. Cell-wall recovery after irreversible deformation of wood. Nat. Mater. 2, 810–814. Keten, S., Xu, Z., Ihle, B., Buehler, M.J., 2010. Nanoconfinement controls stiffness, strength and mechanical toughness of beta-sheet crystals in silk. Nat. Mater. 9, 359–367. Kulasinski, K., Keten, S., Churakov, S.V., Derome, D., Carmeliet, J., 2014. A comparative molecular dynamics study of crystalline, paracrystalline and amorphous states of cellulose. Cellulose 21, 1103–1116. Launey, M.E., Buehler, M.J., Ritchie, R.O., 2010. On the mechanistic origins of toughness in Bone. Annu. Rev. Mater. Res. 40, 25–53. Lewin, M., Goldstein, I.S., 1991. Wood structure and composition. Marcel Dekker, New York, Basel, Hong Kong. Mazeau, K., 2011. On the external morphology of native cellulose microfibrils. Carbohyd. Polym. 84, 524–532. Mazeau, K., Charlier, L., 2012. The molecular basis of the adsorption of xylans on cellulose surface. Cellulose 19, 337–349.
206
journal of the mechanical behavior of biomedical materials 42 (2015) 198 –206
Meyers, M.A., McKittrick, J., Chen, P.Y., 2013. Structural biological materials: critical mechanics-materials connections. Science 339, 773–779. Nair, A.K., Gautieri, A., Chang, S.W., Buehler, M.J., 2013. Molecular mechanics of mineralized collagen fibrils in bone. Nat. Commun. 4, 1724. Nair, A.K., Qin, Z., Buehler, M.J., 2012. Cooperative deformation of carboxyl groups in functionalized carbon nanotubes. Int. J. Solids Struct. 49, 2418–2423. Naraghi, M., Bratzel, G.H., Filleter, T., An, Z., Wei, X.D., Nguyen, S. T., Buehler, M.J., Espinosa, H.D., 2013. Atomistic investigation of load transfer between DWNT bundles crosslinked by PMMA oligomers. Adv. Funct. Mater. 23, 1883–1892. Navi, P., Rastogi, P.K., Gresse, V., Tolou, A., 1995. Micromechanics of wood subjected to axial tension. Wood Sci. Technol. 29, 411–429. Nishiyama, Y., Langan, P., Chanzy, H., 2002. Crystal structure and hydrogen-bonding system in cellulose 1 beta from synchrotron X-ray and neutron fiber diffraction. J. Am. Chem. Soc. 124, 9074–9082. Nishiyama, Y., Sugiyama, J., Chanzy, H., Langan, P., 2003. Crystal structure and hydrogen bonding system in cellulose 1(alpha), from synchrotron X-ray and neutron fiber diffraction. J. Am. Chem. Soc. 125, 14300–14306. Petridis, L., Smith, J.C., 2009. A molecular mechanics force field for lignin. J. Comput. Chem. 30, 457–467. Plimpton, S., 1995. Fast parallel algorithms for short-range molecular-dynamics. J. Comput. Phys. 117, 1–19. Pu, Y., Zhang, D., Singh, P.M., Ragauskas, A.J., 2008. The new forestry biofuels sector. Biofuel Bioprod. Bior. 2, 58–73. Qin, Z., Dimas, L., Adler, D., Bratzel, G., Buehler, M.J., 2014. Biological materials by design. J. Phys. Condens. Matter: Inst. Phys. J. 26, 073101. Qin, Z., Gautieri, A., Nair, A.K., Inbar, H., Buehler, M.J., 2012. Thickness of hydroxyapatite nanocrystal controls mechanical properties of the collagen-hydroxyapatite interface. Langmuir: ACS J. Surf. Colloid 28, 1982–1992. Raman, E.P., Guvench, O., MacKerell, A.D., 2010. CHARMM additive all-atom force field for glycosidic linkages in
carbohydrates involving furanoses. J. Phys. Chem. B 114, 12981–12994. Reiterer, A., Lichtenegger, H., Fratzl, P., Stanzl-Tschegg, S.E., 2001. Deformation and energy absorption of wood cell walls with different nanostructure under tensile loading. J. Mater. Sci. 36, 4681–4686. Reiterer, A., Lichtenegger, H., Tschegg, S., Fratzl, P., 1999. Experimental evidence for a mechanical function of the cellulose microfibril angle in wood cell walls. Philos. Mag. A 79, 2173–2184. Salme´n, L., 2004. Micromechanical understanding of the cell-wall structure. C. R. Biol. 327, 873–880. Sinko, R., Mishra, S., Ruiz, L., Brandis, N., Keten, S., 2014. Dimensions of biological cellulose nanocrystals maximize fracture strength. ACS Macromol Lett. 3, 64–69. Skaar, C., 1988. Wood-water relations, 1st ed. Springer-Verlag, Berlin, Heidelberg, New York. Sugiyama, J., Persson, J., Chanzy, H., 1991. Combined infrared and electron-diffraction study of the polymorphism of native celluloses. Macromolecules 24, 2461–2466. Tanaka, F., Iwata, T., 2006. Estimation of the elastic modulus of cellulose crystal by molecular mechanics simulation. Cellulose 13, 509–517. Tashiro, K., Kobayashi, M., 1991. Theoretical evaluation of 3-dimensional elastic-constants of native and regenerated celluloses - role of hydrogen-bonds. Polymer 32, 1516–1530. Terashima, N., Kitano, K., Kojima, M., Yoshida, M., Yamamoto, H., Westermark, U., 2009. Nanostructural assembly of cellulose, hemicellulose, and lignin in the middle layer of secondary wall of ginkgo tracheid. J. Wood. Sci. 55, 409–416. Vanholme, R., Demedts, B., Morreel, K., Ralph, J., Boerjan, W., 2010. Lignin biosynthesis and structure. Plant Physiol. 153, 895–905. Vural, M., Ravichandran, G., 2003. Microstructural aspects and modeling of failure in naturally occurring porous composites. Mech. Mater. 35, 523–536. Weiner, S., Wagner, H.D., 1998. The material bone: structure mechanical function relations. Annu. Rev. Mater. Sci. 28, 271–298.
Available online at www.sciencedirect.com
www.elsevier.com/locate/jmbbm
Research Paper
Molecular deformation mechanisms of the wood cell wall material Kai Jin, Zhao Qin, Markus J. Buehlern Laboratory for Atomistic and Molecular Mechanics (LAMM), Department of Civil and Environmental Engineering, Massachusetts Institute of Technology, 77 Massachusetts Ave., Cambridge, MA 02139, USA
ar t ic l e in f o
abs tra ct
Article history:
Wood is a biological material with outstanding mechanical properties resulting from its
Received 18 July 2014
hierarchical structure across different scales. Although earlier work has shown that the
Received in revised form
cellular structure of wood is a key factor that renders it excellent mechanical properties at
7 November 2014
light weight, the mechanical properties of the wood cell wall material itself still needs to be
Accepted 13 November 2014
understood comprehensively. The wood cell wall material features a fiber reinforced
Available online 21 November 2014
composite structure, where cellulose fibrils act as stiff fibers, and hemicellulose and lignin
Keywords:
molecules act as soft matrix. The angle between the fiber direction and the loading
Wood cell wall
direction has been found to be the key factor controlling the mechanical properties.
Molecular dynamics
However, how the interactions between theses constitutive molecules contribute to the
Yielding of matrix
overall properties is still unclear, although the shearing between fibers has been proposed
Slip-stick mechanism
as a primary deformation mechanism. Here we report a molecular model of the wood cell wall material with atomistic resolution, used to assess the mechanical behavior under shear loading in order to understand the deformation mechanisms at the molecular level. The model includes an explicit description of cellulose crystals, hemicellulose, as well as lignin molecules arranged in a layered nanocomposite. The results obtained using this model show that the wood cell wall material under shear loading deforms in an elastic and then plastic manner. The plastic regime can be divided into two parts according to the different deformation mechanisms: yielding of the matrix and sliding of matrix along the cellulose surface. Our molecular dynamics study provides insights of the mechanical behavior of wood cell wall material at the molecular level, and paves a way for the multiscale understanding of the mechanical properties of wood. & 2014 Elsevier Ltd. All rights reserved.
Introduction
bearing systems, for example, houses, bridges, ships and even aircrafts, since it also possesses excellent mechanical
Wood, with its broad availability, has served as an important bulk material since ancient times and still plays a major role today. One important application area of wood is in the load n
Corresponding author. Tel.: þ1 617 452 2750; fax: þ1 617 324 4014. E-mail address: [email protected] (M.J. Buehler). URL: http://web.mit.edu/mbuehler/www/ (M.J. Buehler).
http://dx.doi.org/10.1016/j.jmbbm.2014.11.010 1751-6161/& 2014 Elsevier Ltd. All rights reserved.
properties. In the sense of specific mechanical properties (normalized by mass density), wood material parallels and even behaves better than many other engineering materials
journal of the mechanical behavior of biomedical materials 42 (2015) 198 –206
such as steel (Ashby et al., 1995; Gibson and Ashby, 1999; Gibson et al., 1995). In fact, the origin of most woods, the tree, is also a load bearing system, where the smart design of the wood structure renders a tree capable of supporting its own weight and sustaining various climate conditions (Fratzl and Weinkamer, 2007; Qin et al., 2014). The unique combination of composite material and microstructure design in wood material has attracted research interests in comprehensive understanding the mechanical properties of wood and designing inspired materials and structures (Compton and Lewis, 2014; Gibson and Ashby, 1999; Gibson et al., 2010). Similar to other outstanding biological materials, for example, spider silk (Eisoldt et al., 2011; Keten et al., 2010), bone (Launey et al., 2010; Nair et al., 2013; Weiner and Wagner, 1998) and sponge (Aizenberg et al., 2005; Meyers et al., 2013), wood also features a hierarchical structure (Ali and Gibson, 2013) (Fig. 1). At the macro-scale, we can clearly see the annual rings in a cross-section of a wood trunk (Fig. 1a), which is a result of the different growing rates during different seasons (Fratzl and Weinkamer, 2007). Zooming into the micrometer scale, we can clearly see that the wood material is not solid but features a cellular structure (Fig. 1b, c), where hollow space is enclosed by wood cell walls. Studies have shown that the cellular structure is the underlying feature that endows wood as an excellent mechanical material: it makes wood light weight and perform efficiently when it is used as engineering components (Ashby et al., 1995; Gibson and Ashby, 1999; Gibson et al., 1995; Vural and Ravichandran, 2003). The significant hollow space in wood also gives it high deformability. For example, cork is an excellent material for energy absorption and bottle bungs (Gibson and Ashby, 1999).
199
As shown in Fig. 1c-e, the wood cell wall material is a composite, combing stiff and soft polymer components. Studies have shown that the wood cell wall is composed of multiple layers and it is a kind of fiber-reinforced composite with stiff cellulose fibrils embedded in soft matrix consisting of hemicellulose and lignin (Booker and Sell, 1998; Fahlen and Salmen, 2002) (Fig. 1c, d). Among these layers, the middle layer of the secondary cell wall (S2 layer), which occupies more than 80% of the thickness of the cell wall, is the thickest and mainly controls the mechanical properties of the wood cell wall (Booker and Sell, 1998). A unique feature of S2 layer is that the stiff cellulose fibrils lie parallel with each other and form an angle (micro-fibril angle, MFA) with the longitudinal axis of wood cell (Fig. 1d) (Booker and Sell, 1998). As the longitudinal direction is the main load-bearing direction in wood, the mechanical properties of wood mainly relate to the cell wall's properties along this direction. Since MFA measures how close are the stiff cellulose fibrils aligned with this load-bearing direction, it plays a key role in controlling the mechanical properties of wood cell wall: smaller MFA indicates that the stiff cellulose fibrils are more aligned with the load-bearing direction, and the material will behave stiffer. (Booker and Sell, 1998; Navi et al., 1995; Reiterer et al., 2001; Reiterer et al., 1999). Since the stiff cellulose fibrils spirally wind the wood cell (Booker and Sell, 1998), tensile loading along the longitudinal direction will result in the unwinding of the fibrils and the decreasing of the MFA (Keckes et al., 2003). Theoretical and simulation works have proposed that this reorientation of the fibrils generates relative shearing between adjacent fibrils, and that the shear deformation on the matrix is the underlying mechanical response that is responsible for the mechanical
Fig. 1 – The hierarchical structure of wood material. (a) A drawing shows the annual rings structure at the cross-section of a wood trunk. (b) The cellular structure of wood (reproduced from (Vural and Ravichandran, 2003), with permission from Elsevier). (c) Multiple layers in wood cell wall (reproduced from (Booker and Sell, 1998), with permission from Springer). M: middle layer, P: primary layer, S1: the first secondary layer, S2: the second secondary layer, S3: the third secondary layer. (d) The fiber reinforced structure of wood cell wall. The stiff fibrils form an angle θ (micro-fibril angle, MFA) with the longitudinal direction of cell. (e) The molecular model used in the present study. Black arrows indicate the crystal orientations of cellulose.
200
journal of the mechanical behavior of biomedical materials 42 (2015) 198 –206
behavior of wood cell wall (Adler and Buehler, 2013; Altaner and Jarvis, 2008; Burgert, 2006; Keckes et al., 2003). However, the studies of the mechanical behaviors of the constituting molecules are out of current capability of experiments. Some earlier simulation work (Adler and Buehler, 2013) explored the mechanics of wood using a coarse-grained model, but lacked atomistic scale resolution and chemical details of the constituents. This motivates us to apply molecular dynamics (MD) in the present study, as it enables us to build a bottom-up model starting from the chemical structure of the basic constituents of wood. Because studies have revealed the micro-structural design of wood cell wall material (Salmén, 2004) and the force field parameters for polysaccharides (Guvench et al., 2008; Guvench et al., 2009; Hatcher et al., 2009; Raman et al., 2010) and lignin (Petridis and Smith, 2009) have been developed, we are able to build a molecular model of wood cell wall material and perform MD simulations to achieve the goal of this study. A layered molecular model of the wood cell wall is built (Fig. 1e) and shear loading, as suggested earlier, is applied to the model. By analyzing the stress-strain response of the model and simulation trajectories during deformation, this study aims to recognize the deformation mechanisms at the molecular level.
Model A general description As shown in Fig. 1e, the molecular model in the present work includes the three major components of wood cell wall at atomistic resolution: cellulose, hemicellulose and lignin, among which cellulose has the highest mass ratio (Lewin and Goldstein, 1991; Terashima et al., 2009). Although these three components constitute the basic structure of wood cell wall material, water molecules, another factor that it not included here, act as a regulator of the material's mechanical properties. Studies have shown that the bound water within the cell wall can soften the material (Skaar, 1988). Since the bound water in the cell wall can exist at the cellulose surface and within the matrix (Skaar, 1988), the softening effects are expected to be lubricating the cellulose-matrix interface and softening the matrix. In order to comprehensively unravel water's role from these two different aspects, simulations are needed to study water's effects on matrix and fibril-matrix interface separately and it will also involve a study of water's distribution within the cell wall. These studies are out of the scope of our present modelling work. We limit our exploration in the present work on the deformation mechanisms and all the simulations are performed on the dry model as show in Fig. 1e. The cellulose fibril has a crystalline structure that has been revealed by experiments (Nishiyama et al., 2002; Nishiyama et al., 2003). The Iβ phase cellulose is chosen in the present model since it is rich in wood (Horii et al., 1987; Sugiyama et al., 1991). Linear cellulose chains, in which glucose rings are covalently connected, are held together through inter-chain hydrogen bonds (H-bonds) to form cellulose layers and these layers stack in [100] direction (Nishiyama et al., 2002). The crystalline structure and the covalent bonding along the chain
direction render the cellulose fibrils high modulus. Experiments and simulation works have reported the modulus along the chain direction ranging from 150 to 200 GPa (Diddens et al., 2008; Dri et al., 2013; Kulasinski et al., 2014; Tanaka and Iwata, 2006; Tashiro and Kobayashi, 1991). Different morphologies have been proposed for the cellulose fibrils, including square and hexagonal shape (Ding and Himmel, 2006; Fernandes et al., 2011; Terashima et al., 2009). This means that the cellulose fibrils are able to interact with matrix molecules on different crystalline planes. In the present model, the (110) crystalline surfaces (Charlier and Mazeau, 2012; Mazeau, 2011) are chosen to be the fibril surfaces. The size of cellulose fibrils ranges from several nanometers to tens of nanometers (Sinko et al., 2014), which makes it too expensive to include a full cellulose fibril in the modelling. Considering that the deformation of cellulose fibril is much smaller than that of matrix molecules because of its high stiffness, it contributes little to the overall deformation and it is not necessary to model a full cellulose fibril. As shown in Fig. 1e, only two layers (thickness is approximately 1.1 nm) of cellulose in (110) plane are chosen to represent each cellulose fibril. This two-layer setting meets the minimum requirements: it provides one layer to be added loading condition (restricted movement) and the other layer to interact with matrix molecules without external restriction. The thickness, which is beyond the cutoff for non-bond interaction, prevents the matrix molecules interact with vacuum space. Therefore, this simplification significantly reduces the computation consumption while it is still capable of exploring the deformation response. Among the two soft matrix components, hemicellulose directly attaches on the cellulose surface (Lewin and Goldstein, 1991). Although hemicellulose has similar composition as cellulose, it is in an amorphous state (Lewin and Goldstein, 1991; Terashima et al., 2009). Lignin molecules fill the remaining space between hemicellulose molecules. The thickness of matrix between cellulose fibrils is proposed to be approximately 6nm (Terashima et al., 2009). Together with the volume ratio between hemicellulose and lignin, which is proposed to be 2:3, and the density of each component (1.4 g/ cm3 and 1.36 g/cm3 respectively) (Terashima et al., 2009), the amount of both hemicellulose and lignin in the model can be calculated. Due to the limitation of computation consumption, the modelling work does not attempt to model a fibril surrounded by matrix but instead focuses on the behavior of matrix between adjacent fibrils. Periodic boundaries are applied in both directions of the cellulose layer. Thus, the present model represents actually two infinite cellulose slabs separated by an infinite layer of matrix molecules.
Building the molecular model The first step of building the model is to generate the (110) cellulose surfaces (Fig. 2a). The crystalline structure revealed by synchrotron X-ray and neutron fiber diffraction (Nishiyama et al., 2002) is used to build a cellulose block. Cellulose unit cells are copied to form a block of 4 layers of cellulose in the (110) plane. Each layer has 6 chains and each chain consists of 10 covalently bonded glucose rings. Periodic boundaries are applied in all directions and each chain is covalently connected across the boundary in the chain extending direction. Equilibration of the block is performed: 1 ns under NPT ensemble
journal of the mechanical behavior of biomedical materials 42 (2015) 198 –206
201
Fig. 2 – Overview of the model building procedure, consisting of creating (110) cellulose crystal surfaces (a), and then adding hemicellulose (b) and lignin (c) molecules. (d) A dynamic bond-forming procedure of generating linkages between lignin units is applied subsequently. The green dash line indicates the formed bond. (e) After the building of the model, it is loaded by moving the outer-most cellulose layers in opposite directions (indicated by flat green arrows). (f) The lignin units and linkages included in the present model. The reactive atoms that can form bonds are highlighted by circles. (g) The density profiles of the three components across the model's thickness direction.
(1 a.t.m. and 293 K) and another 1 ns under NVT ensemble (293 K). Once the cellulose block is equilibrated, it is cleaved into two parts, each has two (110) layers, and moved apart (Fig. 2a). In this way, two (110) cellulose surfaces facing each other are created and the space between them will be filled by matrix molecules in the following steps. The second step is to add hemicellulose molecules onto the cellulose surfaces (Fig. 2b). Although hemicellulose has various composing units (Lewin and Goldstein, 1991), we choose the Xylan type cellulose in the present study. Significant disparity is not expected for other types of cellulose, since they share similar structure. For each cellulose surface, in total 22 hemicellulose segments are added. These include 12 segments with 5 D-xylopyranosyl residues covalently connected and another 10 segments with additional L-Arabinofuranosyl (Araf) linked to the third the xylosyl residues through α-(1-3) linkage (Charlier and Mazeau, 2012). The Araf residues are used for the further connection with lignin molecules. The segments are added sequentially with random orientation and an equilibration is performed after every adding step (Charlier and Mazeau, 2012; Mazeau and Charlier, 2012). Additional equilibrations under high temperature and following low temperature are performed to obtain the final hemicellulose phase. The final step is to generate the lignin molecules in the model (Fig. 2c-e). For each hemicellulose layer, 3 of the 10 Araf residues are selected to be connected with ferulate moieties that bridging the hemicellulose with lignin (Charlier and Mazeau, 2012). Lignin molecules have various composing units
and linkages between the units (Lewin and Goldstein, 1991; Pu et al., 2008; Vanholme et al., 2010). In the present model, the two major units, G and S unit, and part of the major linkages (shown in Fig. 2f) are included. It is believed that lignin has very complex 3D structure resulting from its various types of units and linkages and possible branching structure. Experiments on the structure of lignin usually provide the ratios of different units and linkages (Pu et al., 2008). In contrast to applying pre-defined sequence of a lignin molecule (Charlier and Mazeau, 2012; Petridis and Smith, 2009), we apply a dynamic bond-forming procedure to build the amorphous lignin structure. This method follows the bio-synthesis pathway of lignin, in which the lignin units are oxidized to form phenolic radicals and the radicals can further be bonded through forming the linkages (Vanholme et al., 2010). In order to reduce the complexity, only the dominant β-O-4 linkages (Pu et al., 2008; Vanholme et al., 2010) are formed in the procedure. Other types of linkages are introduced by adding corresponding dimers (Fig. 2f) into the system. This simplification can also precisely control the ratios of different units and linkages by controlling the number of units and dimers added into the system and the number of formed β-O-4 linkages. Prior to the bond-forming procedure, the units and dimers are set such that the Cβ and O4 atoms (Fig. 2d, f) can form bonds. They are packed into the system sequentially with random positions and orientations (Abbott et al., 2013) and an equilibration with high temperature (800 K) is then performed to thoroughly distribute the units. The tool developed by L. J. Abbott et al. (Abbott et al., 2013) is applied to accomplish this
202
journal of the mechanical behavior of biomedical materials 42 (2015) 198 –206
bond-forming procedure. This tool mimics a polymerization procedure that generates polymer molecules from input monomers. It is noted that when the cellulose block is separated, the space between them is larger than the anticipated 6nm so that the lignin units can have more space to adjust their positions and orientations in order to relax during the bond-forming procedure. An energy minimization is performed after each bond-forming and once every five bonds are formed an equilibration under NVT ensemble with high temperature (800 K) is performed in order to sample the structure. Additional equilibration is performed after the forming of desired number of β-O-4 linkages and then the model is compressed to the desired thickness. In order to study the effects of different lignin structures, a bunch of models with different lignin structures: different units’ ratios (S: G¼ 0: 1, 1: 1, and 1.2: 1), different linkages’ ratios ((β-O-4): (α-O-4): (β-5): (55): (5-O-4) ranges from 50%: 12.5%: 12.5%: 12.5%: 12.5% to 76%: 6%: 6%: 6%: 6%), with and without branching in molecule, different molecule weight (average molecule size ranges from consisting of 6.5 monomers to consisting of 65 monomers), are built and studied. Once the model is built according to the above steps, it is equilibrated in order to reduce the artificial effects introduced during the building procedure. The shear loading is added by moving the outer-most cellulose layers in opposite directions (Fig. 2e). We assume that the two-layer cellulose fibril here represents a fibril of about 12 nm in thickness (Terashima et al., 2009). Thus, the equivalent distance between two adjacent cellulose fibrils is approximately 18 nm (counting half thickness of each cellulose fibril and the thickness of the matrix) and a shear strain can be calculated as dividing the relative displacement of the two moving cellulose fibrils by 18 nm. The stress in the matrix can also be recorded and further we can obtain the stress-strain response on the model. All the simulations are performed by the large-scale atomic/molecular massively parallel simulator (LAMMPS) (Plimpton, 1995). The timestep is chosen as 1 fs and the cutoff for the non-bonded interactions is 1 nm. The particle–particle particle–mesh (PPPM) method is used to compute the long-range columbic interactions and the SHAKE algorithm is applied to constrain highfrequency dynamics from hydrogen-related energy terms.
these two types of molecules. In general, the profiles here are in agreement with that in literature (Charlier and Mazeau, 2012). However, the density profile of lignin here is smoother than that suggested in the literature (Charlier and Mazeau, 2012). This is a result of applying different lignin building methods. The random packing and equilibration before the bond-forming procedure in the present study can provide more uniform distribution of lignin mass than that built by the predefined sequence method (Charlier and Mazeau, 2012). An average value of the lignin density can be calculated to be approximately 1.2 g/cm3, which is close to the experimental value of 1.36 g/cm3 (Terashima et al., 2009).
The elastic-plastic behavior Fig. 3 shows a representative stress-strain curve of the model response. The curve shows elastic behavior at the initial regime and yielding at the shear strain of approximately 0.03. This value can be applied into a theoretical model (Keckes et al., 2003) to predict the tensile yielding strain of wood cell wall materials. Keckes et al.'s study (Keckes et al., 2003) reported the stress-strain curve of a wood cell wall with initial MFA of about 46 degrees and it yielded at the tensile strain of 0.015. Our simulation result, accompanied with the theoretical model, predicts that it happens at the strain of 0.016, which agrees quite well with the experimental result. After yielding, the stress continues to gradually increase with occasional drops. In order to show the plastic deformation, unloading and reloading are performed. As shown in the figure, the strain cannot return back to the original position after unloading, indicating that irreversible deformation after yielding has occurred. During the reloading procedure, the initial elasticity recovers and once the stress reaches the plastic level, the yielding behavior restarts again. At a strain of approximately 0.18, the stress shows a significant drop, and at larger deformation the stress evolves in a “saw-tooth” shape with large amplitude fluctuations. Although this is an extension of the post-yielding behavior, the significant fluctuations indicate that the underlying deformation mechanism might be different from the one that is responsible for the deformation before the strain of 0.18 (it will be discussed in more depth in the next section).
Results and discussion The molecular model The final thickness of the matrix after equilibration is calculated as approximately 5.7 nm, which is in agreement with the desired value of 6 nm (Terashima et al., 2009). As shown in Fig. 2g, the density profiles are calculated across the model's thickness direction for each component. The curves of these components are roughly separated, showing the layered structure of the model: hemicellulose layers lying on cellulose surface and then lignin molecules fill the remaining space. There are two peaks for each cellulose part, which is a result of the crystalline layered structure. The density profile of hemicellulose also shows peaks but they are not intense, showing it is in amorphous state. The overlap regions between the profiles of hemicellulose and lignin indicate the interpenetration of
Fig. 3 – A representative stress-strain response (blue curve) of the model. Three regimes are distinguished: one initial elastic regime and two plastic regimes with different underlying mechanisms. The unloading (red curves) and reloading (green curves) tests show irreversible deformation after yielding.
journal of the mechanical behavior of biomedical materials 42 (2015) 198 –206
Since the stress level before and after the strain of 0.18 does not show significant change, it cannot show clear difference in the experimental stress-stain curves to identify when the mechanism-transition happens. Thus, this mechanismtransition strain cannot be verified by experimental data. Actually, this is also a proof that shows the advantages of molecular simulations over experiments to explore the molecular level deformation mechanisms. Unloading and reloading tests are also performed in this regime. Similarly, the non-zero crossing at the strain axis indicates irreversible deformation and the elasticity recovers during the reloading. Once the stress reaches the yielding stress, the “saw-tooth” behavior is observed again. According to the above description, the stress-strain response can be divided into three regimes: an initial elastic regime and two plastic regimes that are separated by the strain of approximately 0.03 and 0.18 (indicated by the dash lines in Fig. 3).
Deformation mechanisms In order to distinguish the underlying mechanisms for different deformation regimes, snapshots of the model under loading are analyzed. Fig. 4a-c show the snapshots at the initial state, the first plastic regime and the second plastic regime, respectively. The color contours correspond to the initial coordinates in the loading direction. As shown in the figures, the color contour inclines under loading, which indicates shear deformation is generated on the matrix. In contrast, the color contours within the cellulose parts keep vertical, indicating that the adjacent cellulose layers stick with each other firmly during deformation. This result verifies our assumption of the small deformation of cellulose and indicates that the simplification of using two-layer cellulose to
203
represent a thicker cellulose fibril is reasonable. Although the inclined color contours indicate shearing in both plastic regimes, there is a difference of the displacement continuity at the matrix-cellulose interface: while it is continuous during the first plastic regime (Fig. 4b), the displacement shows a jump at the matrix-cellulose interface during the second plastic regime (Fig. 4c). This difference is most likely a result of distinct deformation mechanisms. For the first plastic regime, the irreversible deformation after yielding should be attributed to the deformation within the matrix, since the matrix keeps sticking with the cellulose surface. A close observation of the matrix molecules shows some extent of irreversible reorganizing of matrix molecules. Fig. 4d shows an example in which one molecule squeezed between two other molecules after yielding and cannot return back to the original position when the model is unloaded. Therefore, the mechanism that is responsible for the first plastic regime is the yielding of matrix molecules. In contrast, the discontinuous displacement during the second plastic regime indicates that the matrix molecules slide along the cellulose surface. The movie of the deformation procedure shows that the sliding does not happen in a smoothly continuous way but with a “slip-stick” manner. The slip of matrix happens in a very short time that the relative displacement of the two moving cellulose and further the shear strain nearly keeps the same value. This results in the sudden drops of stress, which are nearly vertical in the stress curve with respect to shear strain (Fig. 3). After the slip event, the matrix will immediately stick with the cellulose surface at a new position and the stress starts to increase again. As mentioned in the model section, hemicellulose has similar units and structure with cellulose. The rich hydroxyl groups in both hemicellulose and the (110) cellulose surface help to
Fig. 4 – (a–c) The snapshots at the initial state, during the first plastic regime and during the second plastic regime, respectively. The color contours correspond to the initial coordinates in the loading direction (indicated by green arrows). (d) Snapshots show the reorganization of molecules within the matrix after yielding. (e) The dynamic H-bonding (the short red line) between hemicellulose and cellulose. (f) The number of H-bonds at the matrix-cellulose interface with respect to the shear strain on the model.
204
journal of the mechanical behavior of biomedical materials 42 (2015) 198 –206
construct a hydrogen bonding network, an important portion of the non-boned interactions at the matrix-cellulose interface. Therefore, the stress-strain curve can be connected with the dynamics of the interfacial H-bonds. In every slip event, the matrix molecules draw back a little bit and the former Hbonds network is broken. This is followed by a reforming of the H-bonds network in a new position once the matrix sticks with the cellulose surface again. The reformed H-bonds help the interface to regain the stiffness, and the stress thus restarts to increase. Fig. 4e shows an example of the dynamic H-bonds at the interface: a hydroxyl group of hemicellulose sequentially forms H-bonds with the acceptor oxygen atoms in cellulose. The number of H-bonds at the interface is recorded with respect to the shear strain (Fig. 4f), where the big amplitudes of “saw-tooth” indicate the dynamic breaking and reforming procedure. Since the matrix draws back, it is unloaded at the slip events and the following stick reloads the matrix again. As an indication, the slope of the “sawtooth” during the second plastic regime is close to the initial elasticity. Meanwhile, the actual deformation of the matrix is limited because of the drawing-back and thus the sliding of matrix along the cellulose surface is the dominant mechanism during the second plastic regime.
In sum, two different mechanisms are recognized for the plastic behavior of the model. Although the yielding of matrix results into the first plastic regime, the “slip-stick” mechanism is dominant during the second plastic regime. In fact, this kind of “self-healing” interface is a common strategy utilized by nature to achieve excellent mechanical properties such as the dynamic H-bonds between collagen and minerals (Qin et al., 2012). This strategy can also be applied in the interfacial functionalization for synthetic high-performance materials (Nair et al., 2012; Naraghi et al., 2013) and potentially many other areas of bio-inspired materials design.
Effect of lignin structure As mentioned in the Model section, models with different lignin structures are built and studied in order to understand its effect on the mechanical response of the model. Fig. 5 shows the stress-strain curves of 9 cases with different lignin structures. As is apparent from the data, all results share the above presented 3-regime behavior; one initial elastic regime and two plastic regimes with different underlying mechanisms, although slight difference exists between each case. In fact, the dominant interaction within the lignin matrix is the squeezing and separation between atoms and the effects of structural level difference are not apparent in the mechanical behavior.
Effect of loading rate
Fig. 5 – Stress-strain responses on models with different lignin structures. All responses reproduce the three-regime behavior. The shadow area indicates the range of the first sliding event for different models.
It is noted that the effect of high-strain rate ( 108 per second in the present model) is unavoidable due to the limitation of computational resources. Two cases with lower strain rate ( 0.2 108) are performed to reveal the effects of loading rate. The stress-strain response of one of them is shown in Fig. 6a, accompanied with the case of the same lignin structure but with higher strain rate. The results show that the 3-regime deformation procedure preserves for the case of lower strain rate. The main difference lies in the initial elasticity and the stress level during the first plastic regime.
Fig. 6 – (a) A comparison of stress-strain responses with different loading rates. (b) The average initial modulus of the nine cases with different lignin structures (the red error bar indicates the standard error) in comparison with the initial modulus of the two cases (indicated by different colors) with lower strain rate.
journal of the mechanical behavior of biomedical materials 42 (2015) 198 –206
Obviously, the matrix behaves stiffer under higher loading rate. Fig. 6b shows the average initial modulus of the 9 cases shown in Fig. 5 and the values of the two cases with lower strain rate. It has an approximate 38% drop from 2.89 GPa. The difference during the second plastic regime is not obvious, probably because of the big fluctuations. Although it is not reliable to extract the stress-strain response quantitatively, the modelings here provide the insights of the molecular-level deformation mechanisms qualitatively.
Conclusion and outlook Through the study with a molecular model of the wood cell wall material, at atomistic resolution, several deformation mechanisms at the molecular level are identified. The yielding of the matrix and sliding of matrix along cellulose surface, which start at the shear strain of approximately 0.03 and 0.18 respectively, sequentially contribute to the post-yielding deformation of wood cell wall material at molecular level. The present MD simulation work provides foundation for both the further study on water's effects and further coarse-graining modelling of wood cell wall material, which could enable us to study the effects of various factors, for example, MFA, fibril length, fibril mass ratio, defects in the material, on the mechanical properties at the larger mesoscale. Further investigation can also include experiments on the inspired wood cell wall structure created using 3D printing techniques.
Acknowledgment We acknowledge the support from the BASF-NORA program, and helpful discussions with Professor Lorna Gibson (MIT).
r e f e r e n c e s
Abbott, L.J., Hart, K.E., Colina, C.M., 2013. Polymatic: a generalized simulated polymerization algorithm for amorphous polymers. Theor. Chem. Acc. 132, 1–19. Adler, D.C., Buehler, M.J., 2013. Mesoscale mechanics of wood cell walls under axial strain. Soft Matter 9, 7138. Aizenberg, J., Weaver, J.C., Thanawala, M.S., Sundar, V.C., Morse, D.E., Fratzl, P., 2005. Skeleton of euplectella sp.: Structural hierarchy from the nanoscale to the macroscale. Science 309, 275–278. Ali, Z.M., Gibson, L.J., 2013. The structure and mechanics of nanofibrillar cellulose foams. Soft Matter 9, 1580. Altaner, C.M., Jarvis, M.C., 2008. Modelling polymer interactions of the molecular velcro type in wood under mechanical stress. J. Theor. Biol. 253, 434–445. Ashby, M.F., Gibson, L.J., Wegst, U., Olive, R., 1995. The mechanical properties of natural materials. I. Material property charts. P. Roy. Soc. A: Mat. 450, 123–140. Booker, R.E., Sell, J., 1998. The nanostructure of the cell wall of softwoods and its functions in a living tree. Holz. Roh. Werkst. 56, 1–8. Burgert, I., 2006. Exploring the micromechanical design of plant cell walls. Am. J. Bot. 93, 1391–1401. Charlier, L., Mazeau, K., 2012. Molecular modeling of the structural and dynamical properties of secondary plant cell walls: influence of lignin chemistry. J. Phys. Chem. B 116, 4163–4174.
205
Compton, B.G., Lewis, J.A., 2014. 3D-Printing of lightweight cellular composites. Adv. Mater. 26, 5930–5935. Diddens, I., Murphy, B., Krisch, M., Muller, M., 2008. Anisotropic elastic properties of cellulose measured using inelastic X-ray scattering. Macromolecules 41, 9755–9759. Ding, S.Y., Himmel, M.E., 2006. The maize primary cell wall microfibril: a new model derived from direct visualization. J. Agr. Food Chem. 54, 597–606. Dri, F.L., Hector, L.G., Moon, R.J., Zavattieri, P.D., 2013. Anisotropy of the elastic properties of crystalline cellulose Iβ from first principles density functional theory with Van der Waals interactions. Cellulose 20, 2703–2718. Eisoldt, L., Smith, A., Scheibel, T., 2011. Decoding the secrets of spider silk. Mater Today 14, 80–86. Fahlen, J., Salmen, L., 2002. On the lamellar structure of the tracheid cell wall. Plant Biol 4, 339–345. Fernandes, A.N., Thomas, L.H., Altaner, C.M., Callow, P., Forsyth, V.T., Apperley, D.C., Kennedy, C.J., Jarvis, M.C., 2011. Nanostructure of cellulose microfibrils in spruce wood. Proc. Natl. Acad. Sci. U.S.A. 108, E1195–E1203. Fratzl, P., Weinkamer, R., 2007. Nature’s hierarchical materials. Prog. Mater. Sci. 52, 1263–1334. Gibson, L.J., Ashby, M.F., 1999. Cellular solids: structure and properties, 2nd ed. Cambridge University Press, Cambridge. Gibson, L.J., Ashby, M.F., Harley, B.A., 2010. Cellular materials in nature and medicine, 1st ed. Cambridge University Press, Cambridge. Gibson, L.J., Ashby, M.F., Karam, G.N., Wegst, U., Shercliff, H.R., 1995. The mechanical properties of natural materials. II. Microstructures for mechanical efficiency. P. Roy. Soc. A: Mat. 450, 141–162. Guvench, O., Greene, S.N., Kamath, G., Brady, J.W., Venable, R.M., Pastor, R.W., Mackerell Jr., A.D., 2008. Additive empirical force field for hexopyranose monosaccharides. J. Comput. Chem. 29, 2543–2564. Guvench, O., Hatcher, E.R., Venable, R.M., Pastor, R.W., Mackerell, A.D., 2009. CHARMM additive all-atom force field for glycosidic linkages between hexopyranoses. J. Chem. Theory Comput. 5, 2353–2370. Hatcher, E., Guvench, O., MacKerell, A.D., 2009. CHARMM Additive all-atom force field for aldopentofuranoses, methylaldopentofuranosides, and fructofuranose. J. Phys. Chem. B 113, 12466–12476. Horii, F., Hirai, A., Kitamaru, R., 1987. Cp/Mas C-13 Nmr-spectra of the crystalline components of native celluloses. Macromolecules 20, 2117–2120. Keckes, J., Burgert, I., Fruhmann, K., Muller, M., Kolln, K., Hamilton, M., Burghammer, M., Roth, S.V., Stanzl-Tschegg, S., Fratzl, P., 2003. Cell-wall recovery after irreversible deformation of wood. Nat. Mater. 2, 810–814. Keten, S., Xu, Z., Ihle, B., Buehler, M.J., 2010. Nanoconfinement controls stiffness, strength and mechanical toughness of beta-sheet crystals in silk. Nat. Mater. 9, 359–367. Kulasinski, K., Keten, S., Churakov, S.V., Derome, D., Carmeliet, J., 2014. A comparative molecular dynamics study of crystalline, paracrystalline and amorphous states of cellulose. Cellulose 21, 1103–1116. Launey, M.E., Buehler, M.J., Ritchie, R.O., 2010. On the mechanistic origins of toughness in Bone. Annu. Rev. Mater. Res. 40, 25–53. Lewin, M., Goldstein, I.S., 1991. Wood structure and composition. Marcel Dekker, New York, Basel, Hong Kong. Mazeau, K., 2011. On the external morphology of native cellulose microfibrils. Carbohyd. Polym. 84, 524–532. Mazeau, K., Charlier, L., 2012. The molecular basis of the adsorption of xylans on cellulose surface. Cellulose 19, 337–349.
206
journal of the mechanical behavior of biomedical materials 42 (2015) 198 –206
Meyers, M.A., McKittrick, J., Chen, P.Y., 2013. Structural biological materials: critical mechanics-materials connections. Science 339, 773–779. Nair, A.K., Gautieri, A., Chang, S.W., Buehler, M.J., 2013. Molecular mechanics of mineralized collagen fibrils in bone. Nat. Commun. 4, 1724. Nair, A.K., Qin, Z., Buehler, M.J., 2012. Cooperative deformation of carboxyl groups in functionalized carbon nanotubes. Int. J. Solids Struct. 49, 2418–2423. Naraghi, M., Bratzel, G.H., Filleter, T., An, Z., Wei, X.D., Nguyen, S. T., Buehler, M.J., Espinosa, H.D., 2013. Atomistic investigation of load transfer between DWNT bundles crosslinked by PMMA oligomers. Adv. Funct. Mater. 23, 1883–1892. Navi, P., Rastogi, P.K., Gresse, V., Tolou, A., 1995. Micromechanics of wood subjected to axial tension. Wood Sci. Technol. 29, 411–429. Nishiyama, Y., Langan, P., Chanzy, H., 2002. Crystal structure and hydrogen-bonding system in cellulose 1 beta from synchrotron X-ray and neutron fiber diffraction. J. Am. Chem. Soc. 124, 9074–9082. Nishiyama, Y., Sugiyama, J., Chanzy, H., Langan, P., 2003. Crystal structure and hydrogen bonding system in cellulose 1(alpha), from synchrotron X-ray and neutron fiber diffraction. J. Am. Chem. Soc. 125, 14300–14306. Petridis, L., Smith, J.C., 2009. A molecular mechanics force field for lignin. J. Comput. Chem. 30, 457–467. Plimpton, S., 1995. Fast parallel algorithms for short-range molecular-dynamics. J. Comput. Phys. 117, 1–19. Pu, Y., Zhang, D., Singh, P.M., Ragauskas, A.J., 2008. The new forestry biofuels sector. Biofuel Bioprod. Bior. 2, 58–73. Qin, Z., Dimas, L., Adler, D., Bratzel, G., Buehler, M.J., 2014. Biological materials by design. J. Phys. Condens. Matter: Inst. Phys. J. 26, 073101. Qin, Z., Gautieri, A., Nair, A.K., Inbar, H., Buehler, M.J., 2012. Thickness of hydroxyapatite nanocrystal controls mechanical properties of the collagen-hydroxyapatite interface. Langmuir: ACS J. Surf. Colloid 28, 1982–1992. Raman, E.P., Guvench, O., MacKerell, A.D., 2010. CHARMM additive all-atom force field for glycosidic linkages in
carbohydrates involving furanoses. J. Phys. Chem. B 114, 12981–12994. Reiterer, A., Lichtenegger, H., Fratzl, P., Stanzl-Tschegg, S.E., 2001. Deformation and energy absorption of wood cell walls with different nanostructure under tensile loading. J. Mater. Sci. 36, 4681–4686. Reiterer, A., Lichtenegger, H., Tschegg, S., Fratzl, P., 1999. Experimental evidence for a mechanical function of the cellulose microfibril angle in wood cell walls. Philos. Mag. A 79, 2173–2184. Salme´n, L., 2004. Micromechanical understanding of the cell-wall structure. C. R. Biol. 327, 873–880. Sinko, R., Mishra, S., Ruiz, L., Brandis, N., Keten, S., 2014. Dimensions of biological cellulose nanocrystals maximize fracture strength. ACS Macromol Lett. 3, 64–69. Skaar, C., 1988. Wood-water relations, 1st ed. Springer-Verlag, Berlin, Heidelberg, New York. Sugiyama, J., Persson, J., Chanzy, H., 1991. Combined infrared and electron-diffraction study of the polymorphism of native celluloses. Macromolecules 24, 2461–2466. Tanaka, F., Iwata, T., 2006. Estimation of the elastic modulus of cellulose crystal by molecular mechanics simulation. Cellulose 13, 509–517. Tashiro, K., Kobayashi, M., 1991. Theoretical evaluation of 3-dimensional elastic-constants of native and regenerated celluloses - role of hydrogen-bonds. Polymer 32, 1516–1530. Terashima, N., Kitano, K., Kojima, M., Yoshida, M., Yamamoto, H., Westermark, U., 2009. Nanostructural assembly of cellulose, hemicellulose, and lignin in the middle layer of secondary wall of ginkgo tracheid. J. Wood. Sci. 55, 409–416. Vanholme, R., Demedts, B., Morreel, K., Ralph, J., Boerjan, W., 2010. Lignin biosynthesis and structure. Plant Physiol. 153, 895–905. Vural, M., Ravichandran, G., 2003. Microstructural aspects and modeling of failure in naturally occurring porous composites. Mech. Mater. 35, 523–536. Weiner, S., Wagner, H.D., 1998. The material bone: structure mechanical function relations. Annu. Rev. Mater. Sci. 28, 271–298.
