Molecular orientation distribution of regenerated cellulose fibers investigated with rotor synchronized solid state NMR spectroscopy

A regenerated cellulose fiber is, in contrast to cotton, a man-made fiber. In the fiber production, the cellulose polymer is subject to various processing steps, affecting the underlying molecular orientation distribution, which is a determining factor for mechanical properties of the fiber. In this work, the molecular orientation distribution was determined in a 13C natural abundance Lyocell regenerated cellulose fiber bundle using rotor synchronized magic angle spinning NMR spectroscopy (ROSMAS) to investigate the chemical shift anisotropy (CSA). The recorded signal intensities were compared with an analytical model of the experiment to find the order parameters reflecting the orientation of the fiber. The CSA tensor was calculated using density functional theory for the crystalline cellulose II structure, commonly found in regenerated cellulose, and is required as an input parameter. The expected order parameter values were only found when approximating the glycosidic bond and its CSA tensor as being parallel to the molecular frame with the order parameter P2=0.45±0.02\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$P_{2}=0.45 \pm 0.02$$\end{document} compared to P2=0.46±0.02\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$P_{2}=0.46 \pm 0.02$$\end{document} obtained with wide angle X-ray scattering on a fiber bundle. To make this method accessible to the community, we distribute the Matlab script for the simulation of spectra obtained by the ROSMAS experiment at github.com/LeoSvenningsson/ROSMAS.


Introduction
The regenerated cellulose fibers produced by the Lyocell process or the more recent ionic liquid dissolution have shown to be a possible environmentally friendly alternative to traditional cotton production (Kosan et al. 2008;Sixta et al. 2015;Wanasekara et al. 2016) and the widely used viscose process. These methods produce cellulose materials for textiles, packaging films and non-woven products from naturally rich cellulose sources, such as wood. Regardless of production method, a cellulose textile fiber is characterized by a significant portion of both crystalline and amorphous contributions with each crystallite having dimensions from ten to a few hundred nanometers (Klemm et al. 2011;Krässig 1993). The crystalline cellulose II structure is formed, among other methods, using the Lyocell or the ionic liquid process and comprises two anti-parallel 180 helical chains, that are known to have different bond angles and chemical shift anisotropy (CSA) as commonly detected by X-ray scattering (Langan et al. 2005;Krässig 1993) and NMR spectroscopy (Idström et al. 2016), respectively. In addition, the cellulose contains a significant fraction of anhydroglucose units exposed as fibrilar surfaces (Idström et al. 2016). The relatively large surface area affects properties, such as absorption and adsorption of water and other compounds, of which water content is known to alter various characteristics of a fiber.
This work focuses on the quantitative aspects of the ROSMAS NMR spectroscopy experiment, which can investigate uniaxial orientation distributions in bundles of fibers by studying the effects of the CSA tensor. Evaluation of the molecular orientation distribution using ROSMAS typically requires the CSA, and its orientation in the molecular frame, to be known a priori. Therefore, density functional theory is used to study the crystalline carbon CSA for the two dissimilar chains that constitute cellulose II, shown in Fig. 1.
The overall molecular orientation distribution in a fiber is only dependent on the angle h between the molecular chain segment and the fiber axis. Therefore, the orientation distribution function (ODF) can be represented as linear contributions of Legendre polynomials weighted by their respective order parameters P ' h i. The molecular chain segments are completely ordered when P 2 h i ¼ 1 and disordered when P 2 h i¼ 0. Rotor synchronization has previously been used to determine up to P 8 h i on poly(p-phenyleneterephthalamide) (Wilhelm et al. 1993). In this work, a solidstate NMR experiment is demonstrated to obtain the P 2 h i order parameter for a naturally abundant 13 C Lyocell style regenerated cellulose fiber.

Rotor synchronized MAS
A complete description of the rotor synchronized MAS experiment was published by Harbison et al. (1986), of which only a brief description of the experiment is given here. The molecular orientation is determined from the CSA of crystalline cellulose II in a regenerated cellulose fiber. The molecular anisotropy is measured by applying a CP-MAS pulse sequence, as shown in Fig. 2, which is synchronized with starting position of the rotor and sampled n times over an entire rotation. The pulse sequence depicted in induces a phase dependency in the indirect dimension. When a bundle of uniformly oriented fibers is placed at an off angle from the rotor axis, it is possible to extract information about the molecular orientation. Fourier transform is performed in both the t 1 and t 2 direction. Molecular orientation is expressed in x 1 dimension separated by Mx R , where x R is the rotation speed. Chemical shift sidebands are separated by Nx R in the x 2 dimension. The signal intensity for a specific peak position is denoted I M;N for integers M and N. The M ¼ 0 contributions contain only positive signals while M 6 ¼ 0 contains both positive and negative signals. In general an ordered material has strong spinning sidebands in M 6 ¼ 0 while a disordered material only has spinning sidebands in M ¼ 0 since only the ordered material is phase dependent.
For a quantitative measurement, experimental data is evaluated with a theoretical model which describes the ROSMAS spectrum with an assumed orientation distribution of the CSA. Harbison et al. based their model on Herzfeld and Berger powder pattern analysis (Herzfeld and Berger 1980), incorporating the rotor synchronization and the ODF as a sum of Legendre polynomials as in Eq. (1).
The ODF assumes a cylindrical symmetry of the fiber and rotational symmetry around the chain axis of anhydroglucose, thus accommodating the twofold helical structure of cellulose II. The order parameters are typically calculated by a linear fitting of Eq. (2). I M;N À Á exp is the experimentally acquired integrated intensity and I ';M;N À Á calc is the calculated intensity using the method of Harbison et al. (1986).
CSA tensor from DFT calculations The theoretical intensity I ';M;N can be calculated if the CSA tensor r is known, i.e. the principal axis components of cellulose II and their orientation in relation to the molecular frame. The principal axis CSA components have previously been analyzed (Hesse and Jäger 2005;Witter et al. 2002) but, less is reported about the shift spatial orientation. To evaluate the CSA tensor, we rely on Kohn-Sham density functional theory and the gauge independent atomic orbitals (GIAO) method for the chemical shielding anisotropy tensor. All electronic structure calculations were performed using the Gaussian09 program (Frisch et al. 2016). The chosen representation of the CSA tensors assumes a right-handed coordinate system where the x-axis is along unit cell vector a and the z-axis is along c. The sequence of rotations defining the Euler angles is: a rotation by a around the z-axis, a rotation by b around the (new) xaxis, followed by a rotation by c around the (new) zaxis. There are two non-equivalent cellulose chains in cellulose II, termed the 'origin' and 'center' chain according to the position in the unit cell as reported in the crystallographic study from which the starting structure was taken (Langan et al. 2005). The unit cell contains a twofold helical axis such that the asymmetric unit contains only one anhydroglucose unit and the coordinates of the other one are obtained by this symmetry operation. This has the consequence that the shift tensors of carbon atoms 1 and 7, 2 and 8, etc are related by an 180 incrementation of a as in Fig. 3. The chemical shielding tensor, r, was converted to the shift tensor to allow direct comparison to experiment. The matrix representation of the CSA, d, is related to the shielding tensor as where 1 is the unit matrix and r ref is the isotropic shift of a standard reference. Determining which value for r ref is appropriate for use with a simulated r is not trivial. In the remarkably useful 'multi-standards' approach, r ref is computed by the same methodology as the system of interest for a chemically similar compound with an experimentally known chemical shift relative to the standard reference (Sarotti and Pellegrinet 2009). In this way, the compounding of errors is minimized.
High resolution X-ray structures of cellulose have been reported, both for mercerized and regenerated cellulose (Langan et al. 2001(Langan et al. , 2005. No significant differences between cellulose II obtained in these two ways have been identified. We used the heavy atom coordinates from Langan et al. (2005) as the starting structure for geometry optimization. Even though X-ray scattering is an accurate structure determination technique, geometry optimization tends to improve the quality of predicted NMR properties (Holmes et al. 2015). Hydrogen atoms, not present in the X-ray structure, were added according to two bonding patterns (Chen et al. 2015). In both H-bond patterns, there are continuous chains of hydrogen bonds with donation in the repeating sequence O6 (origin chain) to O6 (center chain) to O2 (center chain) to O2 (origin chain) for pattern A and reversed for pattern B. Note that the possibility of two similar but distinct hydrogen bonding patterns may explain any temperature dependence in the NMR properties as the relative populations would change with temperature if they are close enough in energy, but this is unlikely to be the case here. The positions of the hydrogen atoms was optimized by a DFT calculation using periodic boundary conditions corresponding to the lattice vectors reported by Langan et al. (2005). After an initial, arbitrary, placement of the hydrogen atoms, a preliminary geometry optimization with all heavy atoms held fixed was performed using a minimal basis set. Basis functions were progressively added and the structure re-optimized until the 6-311G(d,p) level was reached. Then, the positions of all, including heavy, atoms were optimized with this basis set with the unit cell parameters held constant at the experimental values. The orientation of set of the unit cell could not be held constant but changed by about one degree during the calculations, corresponding to a small collective re-arrangement. Results are presented for the original orientation. The HSE exchange-correlation functional was used for all geometry optimizations (Henderson et al. 2009). This hybrid functional was constructed specifically to perform well for solidstate calculations under periodic boundary conditions. The 'exact' exchange energy contains long-ranged contributions and hybrid functionals are therefore typically inefficient to the point of inapplicability under periodic boundary conditions. This problem is ameliorated in the HSE by making the functional 'hybrid' only for short distances and treating larger distances by 'pure' DFT. The accuracy of HSE is essentially equivalent to the global-hybrid exchangecorrelation functional PBE0 (a.k.a. PBEh), on which it is based.
As the gauge independent atomic orbitals (GIAO) method does not appear to be implemented for solidstate calculations under periodic boundary conditions in the Gaussian09 program, these calculations were performed on finite clusters. Given sufficient cluster size, this gives results equivalent to calculations with the full solid-state environment taken into account (Holmes et al. 2015). Cellotriose motifs were cut from the crystal lattice and hydrogen atoms were added to satisfy the valency of broken bonds and form cellotriose molecules. The positions of the terminal hydrogens were optimized as above, except without periodic boundaries, with all other atoms held fixed. The cellotriose molecules were then assembled into clusters consisting of one central molecule surrounded by six additional ones, such that all the nearestneighbor glucose residues were present in the cluster. The shielding tensors reported refer to the six carbon atoms of the middle glucose residue of the central cellotriose molecule, corresponding to the asymmetric unit. One cluster was constructed for each of the origin and center chains for each of the hydrogen bonding patterns, for a total of four clusters. An example is shown in Fig. 4.
The B97-2 exchange-correlation functional was used for the calculation of chemical shielding tensors (Wilson et al. 2001). This functional was originally constructed with shielding tensor calculations in mind, and a recent benchmark study confirms that it gives better predictions of 13 C chemical shifts than is typical for DFT (Flaig et al. 2014). The pcSseg family of basis sets was employed, (Jensen 2015) which is specifically designed for the calculation of NMR properties. It consists of five members named pcSseg-n with n ¼ 0; 1; 2; 3; or 4 for basis sets of ascending size. Basis functions were assigned to atoms according the locally dense scheme suggested in Reid et al. (2014). Namely, pcSseg-3 was used for the six carbon atoms of the central glucose residue and the hydrogen atoms bonded to them; pcSseg-2 for all heavy atoms bonded to those carbon atoms and hydrogen atoms bonded to those atoms; and pcSseg-1 for all other atoms. In this way, a large basis set could be employed for a detailed description of the central carbon atoms for which the shielding tensors are calculated, while allowing a system size sufficiently large that the interactions with the nearest-neighbor molecules could be included. It would have been possible to carry out the chemical shielding calculations in a periodic system using any of a number of programs based on plane-wave basis sets, but this would have precluded the use of the highperforming pcS-family basis sets.

NMR spectroscopy
A bundle of 17 lm Lyocell fibers was spun on a flat spool with a thin layer of polystyrene glue applied between wounds. Strips were cut out at 60 angle from Fig. 4 Cluster for shielding tensor calculation for origin chain with H-bond pattern B. Larger atom size indicates larger basis set the fiber direction and stacked on top of each other so that it fitted a 4 mm rotor, as in Fig. 5. Any additional space was filled with silicon dioxide powder for rotational stability. The small rotor size enables stronger proton decoupling, which reduces the cellulose signal linewidth, compared to a 7 mm rotor.
NMR measurements were carried out using a Bruker 500 MHz Avance III, operating at 125.8 MHz for 13 C. The ROSMAS experiment was conducted at ambient temperature and 1500 Hz (AE1) rotor spinning speed to sufficiently separate a few C1 sidebands from other signals. The rotational frequency should also comply to x R d 11 À d iso and x R d iso À d 33 to avoid signal averaging. The t 1 rotor period was divided in 10 equal steps before initiating a ramped cross polarization sequence with a constant 13C rf strength of 60 kHz while 1H rf is ramped from 45 to 90 kHz with SPINAL64 83 kHz proton decoupling averaged over 5000 scans for each phase step. The experiment is repeated to investigate potential signal obscuring effects, thereafter averaged with the first experiment. A recycle delay of 3 s, cross-polarization time of 1.5 ms and acquisition time of 20.48 ms were applied with the 13 C chemical shift externally referenced against adamantane.
Wide angle X-ray scattering A wide angle X-ray scattering measurement was performed on a Lyocell fiber bundle at Chalmers University of Technology with a 0.9 mm beam diameter Rigaku 003? high brilliance microfocus Cu-radiation source at 130 mm distance with the sample irradiated over 20 min.

ROSMAS analysis of oriented regenerated cellulose fibers
The 1500 Hz M ¼ 0 ROSMAS spectrum (Fig. 6) has wiggles, clearly seen between 140 and 200 ppm, originating from the aromatic ring in the polystyrene adhesive, centered around 125-130 ppm. The aromatic sidebands overlap with cellulose C1 peak in the M ¼ 0 spectrum, but not in the M 6 ¼ 0 spectra (Fig. 7) since the polystyrene adhesive is disordered. The x 2 shaped baseline is also effectively filtered out from the M 6 ¼ 0 spectra since it is independent of the rotor phase. The direct x 2 dimension is processed without window function to minimize effects of line broadening affecting spectral deconvolution. Idström et al. (2016) have previously determined the isotropic shift for the five unique components of the C1 signal in regenerated cellulose. The model explaining these five signals, is based on there being two crystalline signals, two signals of tentative surface origin and one amorphous signal. The two crystalline signals are what is referred to the origin chain and center chain in Fig. 3. Consequently, DFT calculations are performed for the crystalline part, and only this contribution will represent the molecular orientation. Spectral deconvolution was performed on both the M ¼ 0 C1 peak near 106 ppm and the M ¼ 1 C1 peak near 118 ppm. Deconvolution of the C1 peaks is not straightforward since it, in absence of complimentary C4 and C6 information due to signal overlapping, lacks sufficient unique features to accommodate all of the 5 signal sources. For this reason, special conditions are enforced such as: the two crystalline sources having the same orientation distribution and a pure Lorentzian lineshape while the other signals are purely Gaussian. In addition, a linear baseline is used to approximate the signal overlap from the polystyrene aromatic peak overlap and other potential sources in the M ¼ 0 band. There are interfering C2,3,4,5 sidebands located around 100 ppm and 112 ppm at M ¼ 0 and M ¼ 1 , respectively, approximated by a Gaussian lineshape. An expected surface signal, located around 104 ppm and 116 ppm at M ¼ 0 and M ¼ 1, respectively, is too weak and therefore omitted from the deconvolution. Two deconvoluted ROSMAS C1 peaks in M ¼ 0 and M ¼ 1 are shown in Fig. 8. Notably, there is only minor difference in the relative intensities between the four deconvoluted signals between M ¼ 0 and M ¼ 1, respectively. The results indicate that most of the crystalline and amorphous cellulose is ordered in the fiber direction.

The CSA tensor
The calculated chemical shift tensors, visualized in Fig. 3, are presented in Table 1. The accuracy of the calculations is expected to be within a few ppm. However, differences in shielding may be significantly more accurate than this, which is why so many digits have been retained. Only the relative shifts of the principal axis d components are important for the ROSMAS intensity calculation, therefore the much larger isotropic error can be ignored. The C1 CSA principal axis components in Table 1 and previously reported experimental data on regenerated cellulose II are within a few ppm (Hesse and Jäger 2005). However, greater variation are shown between previously reported experimental C2,3,4,5,6 CSA principal axis components and our DFT results in Table 1.

The molecular orientation of regenerated cellulose fibers
The order parameters were calculated by Eq. (2) for a C1 shielding tensor in Table 1. A stretched fiber model was also explored, where the glycosidic COC bond and the magnetic shielding d 33 was assumed to be parallel with the molecular frame. The stretched fiber model was successfully used to establish good correlation between orientation parameters experimentally received data from polarized Raman spectroscopy and wide angle X-ray scattering on the same regenerated cellulose fibers published elsewhere (Svenningsson et al. 2019). An additional assumption is that the principal axis CSA tensor from DFT remains unchanged even though the ð1 ! 4Þ bond angles are altered. Figure 9 shows the azimuthal angle of the 002 reflection from a Lyocell regenerated cellulose bundle with a wrapped Lorentzian curve fit. The conversion to order parameters is performed using Eq. (4), with f h ð Þ as the WAXS wrapped Lorentzian curve fit function.
The curve fit is performed using experimental data between AE 45 to reduce obscurities from the 021 reflection. An X-ray scattering measurement offers a more direct method of measuring crystalline molecular orientation in a fiber. Table 2 contains the calculated order parameters with DFT and stretched fiber model using deconvolution and peak height as intensity. The results are shown together with X-ray scattering and polarized Raman experiments of the same fiber, published elsewhere (Svenningsson et al. 2019).
The expected order parameter P 2 h i for a fiber bundle should lie below the order parameter for a single fiber since a fiber bundle increases the molecular disorder. The ROSMAS bundle P 2 h i order parameter introduces additional disorder from sample preparation and should be lower than X-ray scattering of a bundle, which occurs only for the stretched fiber model. The results suggest that the crystalline structure in a fiber changes after mechanical stretching so that the ð1 ! 4Þ linkage is aligned to the z axis in the molecular frame. A critical implication of this result is that the common understanding of the cellulose II crystalline structure, typically found in regenerated cellulose, might be different in an exited stretched state typically found in fibers.  By inferring the stretched fiber model for amorphous cellulose, its molecular order is similar to the molecular order of crystalline cellulose. This is further motivated by the polarized Raman study showing similar molecular order as X-ray scattering on a single fiber since the Raman vibrational mode contains both crystalline and amorphous contributions according to Wan et al. (2017).
Because of the long experimental time of 84 hours and rather large error originating from deconvolution, originating from the complexity of the fitting many signals with few unique features, it may be preferred to perform a simple peak height estimation. Such estimation is valid since the ROSMAS M ¼ 0 and M ¼ 1 signals have equal line width and the overlapping signals change uniformly with molecular order as established with the stretched polymer condition. The uniform signal change is only predicted when the CSA tensor is parallel to the molecular frame and have the same orientation distribution. Since, amorphous peak height and orientation is qualitatively similar to the crystalline origin chain signal and the amorphous tail is only weakly affecting the origin chain signal, potential obscuring of the P 2 h i order parameter by a peak height approximation should remain less than 1%. The order parameter can therefore be calculated to P 2 h i ¼ 0:45 AE 0:02 with a signal to noise of 147 (M ¼ 0) and 29 (M ¼ 1) from our respective origin chain peaks, a great improvement compared to deconvolution. The signal to noise error is similar to approximated errors from other sources, thus further reduction in signal to noise is therefore unnecessary. Since linewidth is less problematic for peak height measurements, compared to deconvolution, a 7 mm rotor would be favored. Expecting an approximated double signal-to-noise ratio from a larger sample volume, the total number of scans could be reduced by approximately 16 times if a 50% reduction of signal to noise is allowed, thus enabling the experiment to be performed around 6-8 h. The ROSMAS experiment offers an excellent opportunity of signal separation when investigating chemically modified, composites or colored celluloses when X-ray and polarized Raman spectroscopy data may be obscured.
The author have shown that the orientation distribution in a single fiber fits well with wrapped Lorentzian function published elsewhere (Svenningsson et al. 2019). In this work we show that it is also a reasonable approximation for a fiber bundle, as seen in Fig. 9. Therefore, a wrapped Lorentzian ODF can be reconstructed from only information of the P 2 h i order parameter. The same approximation can be used to exchange the Legendre polynomial function to a wrapped Lorentzian ODF and directly compute the ODF by iterative methods. The more computationally heavy method, requiring each iteration to solve a fivedimensional integral, is however rendered obsolete when ODF reconstruction from P 2 h i is available. For completeness, the order parameter is P 2 h i ¼ 0:44 AE 0:01 and P 4 h i ¼ 0:25 AE 0:01 with the wrapped Lorentzian based ROSMAS approach. Again, the error from other other sources are likely greater than the pure signal to noise based error used fr error estimation.
As an additional note, it is theoretically possible to determine the entire CSA tensor with ROSMAS, as done by Gabriëlse et al. (1995). This method however, is only realistically performed with isotope-labeled samples and is further more hindered by the cellulose II complex crystalline structure compared to Table 2 Order parameters calculated from the ROSMAS experiments on a fiber bundle by deconvolution and peak height intensity along with X-ray scattering on a bundle of Lyocell fibers. Additional results are included based on polarized Raman spectroscopy and X-ray scattering on a single fiber published elsewhere (Svenningsson et (Song et al. 1996;Antzutkin et al. 1994), enabled by the more resent DNP signal enhancement methods (Takahashi et al. 2012;Gutmann et al. 2017;Kirui et al. 2019). A 13 C T 1q relaxation filtering to remove amorphous signal was also explored, however it was unsuccessful due to extremely short cellulose T 1q relaxation times.

Conclusion
The ROSMAS NMR spectroscopy method was demonstrated on regenerated cellulose fibers. The result reveals a difference in the crystalline cellulose II CSA in a stretched cellulose fiber, compared to cellulose without a mechanical bias. Amorphous cellulose is similarly aligned as crystalline cellulose when a stretched fiber model is implied. A more simple peak height measurement in combination with the stretched fiber model of the CSA tensor is the preferred approach to analyze molecular orientation in regenerated cellulose fibers. The molecular order parameter P 2 h i ¼ 0:45 AE 0:02 is in agreement with X-ray scattering P 2 h i ¼ 0:46 AE 0:02 for a bundle of fibers. ROSMAS NMR spectroscopy in combination with X-ray and polarized Raman spectroscopy are powerful tools to quantitativevly analyze anisotropy in cellulose and other polymers. The script for calculating the analytical spectra of the ROSMAS experiment has been made public at github.com/LeoSvenningsson/ROSMAS.