Understanding of dielectric properties of cellulose

The theoretical understanding of structural and optoelectronic properties is well-established for a range of inorganic materials; however, such a robust foundation is, in large part, absent in the case of cellulose. Existing literature reports a wide variance in experimentally observed properties for cellulose phases, which are often in contradiction to each other. Motivated by this, we perform an exhaustive rst-principles investigation into the structural and optoelectronic properties of cellulose I α and I β phases. Utilizing exchange-correlation functionals that accurately describe van der Waals interaction and leveraging state-of-the-art density functional theory methods, we strive to present a highly accurate periodic model for the cellulose phases. We integrate the framework of volume-average theory and the potential impact of water sorption to offer insights into the considerable discrepancies seen across different experimental outcomes. Thus, our study provides a reconciliatory perspective, bridging the gap between theoretical calculations and disparate experimental data.


Introduction
The development of electronic structure theory has underscored the potential to predict material properties accurately within mean-eld density functional theory (DFT), provided the atomic identities, composition, and structure are known (Zunger 2018).This is particularly proven for inorganic bulk compounds.As a re ection of this prediction ability, several open-access databases (e.g., Materials Project (Jain et al. 2013), AFLOW (Calderon et al. 2015), and Open Quantum Materials Database (OQMD) (Kirklin et al. 2015)) have been developed.These databases are a testament to the synergy of computational techniques and materials science, representing a leap in the DFT approach to material discovery, and have accumulated an impressive volume of computed properties for a vast array of materials, providing a highly valuable resource for researchers and engineers worldwide.By integrating state-of-the-art electronic structure calculations with robust and scalable data frameworks, these resources have signi cantly expanded our capacity to understand and predict material behavior, accelerating innovation across various elds.We note, however, that such databases are usually scarce on data for organic compounds.This is not surprising, as the structures of organic compounds are often very complex to guess.One of the primary hurdles in the structural characterization of such compounds is accurately determining H atom (usually one of the main elements in the organic compounds) locations.Unlike heavier atoms, which can be more easily "visualized" using X-ray crystallography due to their higher electron densities, low electron density makes H less susceptible to X-ray diffraction (Malyi et al. 2022).Moreover, many organic compounds are not periodic or represented within a partial occupation model (Tilley 2020).For instance, cellulose, one of the most abundant natural polymer organic compounds on Earth, is a fascinating material with a complex structure composed of β(1→4) linked Dglucose units.This natural polymer consists of a periodic layered structure, interconnected mainly via interchain hydrogen bonds, with van der Waals (vdW) forces playing a vital and dominant role in stacking interactions.The structural complexity of cellulose is further highlighted by its ability to exist in different crystalline forms, with cellulose I α (triclinic, P1) and I β (monoclinic, P2 1 ) being the most common forms (French 2014;Gardiner et al. 1985;Kolpak et al. 1976;Sarko et al. 1974;Sarko et al. 1976).Interestingly, despite numerous studies (French 1978;Nishiyama et al. 2002;Nishiyama et al. 2003;Woodcock et al. 1980), both these phases are described within the partial occupancy structures.For instance, the combined x-ray and neutron diffraction study depicts that the two distinct hydrogen bond arrangements could be found with partial occupancy in the range of 70-80% vs 30 − 20%, corresponding to major and minor fractions of total hydrogen bonding network in cellulose I β , respectively (Nishiyama et al. 2002).In the case of cellulose I α , the partial occupancy range of a minor fraction is higher (45%) (Nishiyama et al. 2003).Indeed, even the development of the rst-principles model that accurately identi es the representative non-partial occupancy structure and includes an accurate description of vdW forces is not present, and the available studies have a large discrepancy between DFT and experimental works.For instance, the static dielectric constant of cellulose I β calculated using DFT is 2.5-3.this raises the question: (i) can a comprehensive electronic structure theory for the primary cellulose phases be developed?and (ii) if so, can these models account for the discrepancies observed in previous studies on cellulose properties?Driven by these inquiries, we present an in-depth rst-principles analysis of the cellulose I α and I β phases, highlighting their lowest energy structures and optoelectronic characteristics.Through our study, we elucidate the in uence of vdW forces and suggest a potential reason for the signi cant variations in the dielectric constant observed in experimental studies.Importantly, we highlight different factors that can affect materials properties and explain, for instance, a possible origin of the large discrepancy in the experimentally measured and theoretically obtained dielectric constant of the cellulose phases.

Computational Methods
All the electronic structure calculations were performed using the DFT approach and projector augmented wave (PAW) method (Blöchl 1994) as implemented in Vienna Ab-initio Simulation Package (VASP) (Kresse et al. 1996a, b;Kresse et al. 1993).An energy cutoff of 550 eV was employed for the plane-wave basis set used in wavefunction expansion.The position and lattice optimization were performed until the forces were less than 0.01 eV/Å.Brillouin zones were sampled using Г-centered k-points mesh corresponding to 5000 k-points/atom and 10000 k-points/atom for structural relaxation and electronic structure calculations, respectively.To account vdW interaction accurately, the structural relaxation was performed using different levels of dispersion correction methods: optB86b-vdW (Klimeš et  2013) XC functionals.These optimized structures were further employed to accurately describe electronic properties using a revised hybrid HSE06 functional (Krukau et al. 2006).Here, due to computational cost, Г-centered k-points mesh corresponding to 1000 k-points/atom were used.The frequency-dependent dielectric properties were calculated using independent particle approximation (Gajdoš et al. 2006).The Born effective charge tensors, phonon frequencies, and dielectric permittivity tensors have been calculated using a variational approach within the framework of density-functional perturbation theory, yielding linear response functions (Baroni et al. 2001;Gonze 1997;Gonze et al. 1992;Gonze et al. 1997).

Results and Discussion
Evaluating H arrangements in cellulose I α and I β -identifying the lowest energy con gurations Using available crystallographic data and partial occupancy models (French 2014) for cellulose I α and I β phases (Figs.1a and 1b), we rst analyze different H atom arrangements to identify the lowest energy structures.Speci cally, using the primitive cells for the partial occupancy models of cellulose I α and I β , we screen different arrangements of H atoms corresponding to the C 6 H 10 O 5 chemical formula.All these structures are fully relaxed using optB86b-vdW XC functional (Klimeš et al. 2011).To account for the possibility of various H positioning at the primary alcohol O6 and secondary alcohol O2 sites labeled as in (Nishiyama et al., 2002;Nishiyama et al., 2003), we have analyzed multiple con gurations, including atomic arrangements (referred to as "single con gurated") where all O6 and O2 form a single O-H bond, as well as atomic arrangements where two H atoms are coordinated to the same O6 or O2 site, referred to as "double coordinated" sites.In total, 125 unique input con gurations for cellulose I α and 174 con gurations for cellulose I β were considered.It should be noted that the stability of the lowest-energy con gurations with respect to random atomic displacements has also been con rmed.The results, summarized in Figs.1c and 1d, show that while different atomic con gurations have different energies, the lowest-energy con gurations (given in the supplementary materials) for both phases correspond to the atomic arrangement with single con gurated O sites.This implies that the cellulose structure where the O atoms have double coordination with H atoms is energetically unfavorable.Indeed, the relaxation of such input structures often leads to the desorption of H atoms in the form of H 2 molecule.The lowest energy con guration identi ed for cellulose I α and I β agrees with that found in previous works using different methods (Chen et al. 2014;Nishiyama et al. 2008), showing the reliability of the developed rstprinciple model.

Structural properties of cellulose I α and I β -the role of vdW interaction
The complexity of cellulose structure arises not only from the partial occupancy discussed above but also from the signi cant effect of vdW interactions on materials properties (Malyi et al. 2018).What makes this even more complex is that while various XC functionals have been developed, there remains no universally accepted XC functional capable of accurately reproducing vdW interactions.This limitation becomes evident especially when compared to the results obtained from diffusion quantum Monte-Carlo calculations (Shulenburger et al. 2015), showing, for instance, the inability of modern XC functionals to reproduce the reference state of charge density distribution in the bilayer black phosphorene systems.Because of this, in practice, one often needs to check the sensitivity of material properties to speci c ways of treatment of vdW interactions and verify that main conclusions are not affected by the choice of the computational setup.Motivated by this, using the lowest energy structures identi ed above, we calculated structural and electronic properties for cellulose I α and I β phases using different ways to treat vdW forces.The results are summarized in Fig. 2 and, as expected, show variations in lattice constants for different XC functionals.Importantly, the lowest energy cellulose I β and I α structures have βdextrorotatory-glucose units covalently bonded along the c-direction, with weak vdW interactions along the a and a & b lattice vectors, respectively.This makes the lattice constants along the a and b directions signi cantly more sensitive to XC functional for both phases.For instance, for cellulose I β , depending on the XC functional (Fig. 2d), the lattice constants varied from 7.01-7.85,8.11-8.36,and 10.40-10.58Å for the a, b, and c lattice constants, respectively.Similar variation of lattice constants for cellulose I α structure is also observed (Fig. 2b).It is important to note that our calculations also revealed minor discrepancies in lattice angles compared to experimental ndings.This divergence arises primarily because our computational approach allows unrestricted relaxation of the internal structure and lattice vectors.Consequently, assuming that the experimentally detected structures could represent an ensemble average of different low-energy domains is plausible.Alternatively, the variations may be attributed to thermal uctuations affecting different crystal regions, leading to the occupation of diverse H-bonding arrangements.This would result in an experimentally observed global average structure.In general, however, the computed lattice vectors are comparable to the experimental one.For instance, for cellulose I β structure, the experimentally obtained lattice vectors are a = 7.64 Å, b = 8.18 Å, c = 10.37 Å, and γ = 96.5 at 15 K and a = 7.76 Å, b = 8.20 Å, c = 10.37 Å, and γ = 96.5 at 295 K (Nishiyama et al. 2008).As noted above, the accurate description of the vdW system is important not only for describing vdW interaction but also for accurate calculations of the electronic properties (Malyi et al. 2018).Indeed, as we can see, the calculated band gap energies for cellulose I β and I α structures are 5.35-5.54eV and 5.21-5.47eV, respectively.We note, however, that vdW-DF functionals are still based on soft XC functionals, which are known to underestimate the band gap energy.Moreover, there is not yet a universally accepted XC functional that can be used to describe both electronic and structural properties at the same time.In practice, such descriptions are often done by applying hybrid functional to calculate band gap energy for the frozen vdW-DF functional-relaxed structures (Yadav et al. 2023).Indeed, in many cases, this method allows for a good description of the electronic band gap energy with high accuracy.Application of such method of calculations in this case results in a range of band gap energy (HSE06) values of 7.23-7.46eV and 7.15-7.31eV for cellulose I β and I α structures, respectively.Here, it is worth mentioning that attempts have been made to quantify the experimental band gap energy.Speci cally, for cellulose I β , reported experimental band gap energy values fall within the 4.5-5.7 eV range according to multiple studies (Plermjai et al. 2019;Simao et al. 2015;Sriphan et al. 2023).However, it is important to acknowledge that, in general, the prediction of band gap energy for wide band gap insulators using standard absorption techniques is often challenging due to weak absorption, the role of point defects, etc.It should also be noted that our results are comparable to previously reported work using DFT-PBE0, which yielded values of 8.2 eV and 8.1 eV for cellulose I α and I β , respectively (Srivastava et al. 2020).

Optical properties of cellulose I α and I β phases
The dielectric properties of a material play a pivotal role in shaping its interactions with electromagnetic waves, such as light.These properties dictate how a material polarizes in response to an external electric eld, affecting light absorption, re ection, refraction, and transmission.At the rst-principles level, these properties are often characterized using the complex dielectric function ε = ε 1 + iε 2, which, for instance, can be calculated using the independent particle approximation (Gajdoš et al. 2006).It is worth noting, however, that deriving these properties through calculations is far from straightforward.The complexity of these calculations arises from the need to use large simulation cells, the sum of a wide range of transitions from occupied to unoccupied states within a dense k-point grid, and employ su ciently accurate methods to predict the optical band gap energy.While, in principle, in selective systems, all this can be done directly, in this work, to minimize computational cost, we perform optical calculations using the scissor-correction technique, as illustrated in Fig. 3a.Speci cally, we calculate the imaginary part of the dielectric function using a given XC functional within independent particle approximation.Then, we align the electronic band gap values with those computed by the HSE06 method for the corresponding structure by shifting the imaginary part of the dielectric function (Fig. 3b).It should be noted that due to the difference in k-points density used for HSE06 and other functionals, the scissor correction is calculated using the difference in corresponding band gaps at Г point.Subsequently, we apply the Kramers-Kronig relations to calculate the real part of the dielectric function.This method has indeed been adopted in previous studies (Blaha et al. 1990;Crovetto et al. 2016;Malyi et al. 2016) and allowed accurate description of the electronic and optical properties of the materials.The results of the calculations for both cellulose phases determined using optB86b-vdW XC functional are shown in Fig. 3cf, revealing a marked degree of anisotropy in dielectric properties.Speci cally, cellulose I α exhibits a narrower dielectric anisotropy range than cellulose I β .This distinction becomes particularly evident when assessing the electronic contributions to the dielectric constant (ε ∞ ) along the xx, yy, and zz directions.
Using the optB86b-vdW XC functional with scissor correction, the calculated values are 2.59, 2.60, and 2.64 for cellulose I α , and 2.54, 2.64, and 2.73 for cellulose I β .
It is important to recognize that in actual experimental samples, the orientation and arrangement of cellulose bers can vary.As a result, it makes sense to talk about average dielectric function values along various directions.For instance, one can use an average value calculated as, , as shown in Fig. 3 for different XC functionals.Broadly speaking, the results across various methods are similar.To highlight this, we look at the directional average of the electronic contribution to the real part of the dielectric function, ε 1 .We nd that the rVV10 functional yields the highest ε ∞ values-2.75and 2.78 for cellulose I α and I β , respectively.In contrast, the vdW-DF functional gives values of 2.47 and 2.49 for cellulose I α and I β , respectively.As for the imaginary part of the dielectric function, ε 2 , its behavior under different vdW-DF methods is evident from the positions of the absorption edge and peaks, which show minor differences in energy values.Here, the vdW-DF-cx relaxed structure in both phases shows a red-shifted edge compared to other XC functionals.
When an electric eld is applied to the solid, the dielectric properties of materials are not only de ned by the electronic response, but at low frequencies, ions have enough time to move in response to the changing electric eld.Such contribution is usually discussed in a low-frequency regime (i.e., ε 1 (ω→0)), which is represented as ε 0 = ε ∞ + ε i , where ε ∞ and ε i are an electronic and ionic contribution to the dielectric constant, respectively.The latter can be obtained knowing Born effective charges and the phonon frequencies (Cockayne et al. 2000;Gajdoš et al. 2006;Vali et al. 2004).Using this methodology and the electronic contribution calculated above, we calculated the static dielectric constant for both phases, with results for different XC functionals summarized in Fig. 4c,f.The results show that I β phase exhibits substantially higher ionic contribution to the dielectric constant.Moreover, the extent to which ions contribute to the dielectric constant in both phases signi cantly depends on the speci c correction techniques used in the computations.Importantly, when using the rVV10 and vdW-DF functionals, the highest and lowest average values for the diagonal ε 0 components are found in cellulose I α (cellulose I β ), with corresponding values of 3.62 (3.72) and 3.17 (3.22), respectively.We note that the computed results for dielectric constant are in the lower range of experimentally reported data 3.11-15.00(Boutros et  To gain a deeper understanding of the variability observed in experimental results, we scrutinize the factors that in uence material properties during experimental measurements.This scrutiny is especially critical for cellulose structures, given their propensity to absorb substantial amounts of water.Such absorption can signi cantly alter the material's strength, exibility, and overall behavior.Additionally, experimental samples may consist of mixed cellulose phases in varying proportions, further complicating the interpretation of results.Consequently, these complexities could result in the experimentally measured properties of cellulose deviating from what might be expected based solely on its chemical formula. Instead, these properties might better re ect a cellulose structure incorporating water molecules or varying proportions of different cellulose phases.While the difference in dielectric constant for cellulose phases cannot explain the nearly 5-fold variation observed in the dielectric constant, the formed effect of water molecules is critical.Indeed, water molecules are polar, on the application of an electric eld, the molecule polarizes, leading to the high dielectric constant in the region with su cient density of waterwater is known to have a large static dielectric constant of 76.5-88.2(Elbaum et al. 1991;Fiedler et al. 2020;Malmberg et al. 1956).Moreover, the presence of water molecules can signi cantly affect structural property, which can also result in tuning the vdW interaction between the layers (Agarwal et al. 2017;Salmén et al. 2021).In the simpli ed picture, the effect of water sorption on the dielectric constant can be understood within the so-called volume average theory (Braun et al. 2006;Malyi et al. 2016;Navid et al. 2008), where the effective dielectric constant (ε eff ) of the composite material containing water with fraction φ and cellulose can be expressed as ε eff = (1 -φ)ε c + φε w , where ε c and ε w are dielectric constant of cellulose and water, respectively.In this way, we can see that the increase in water fraction will increase the dielectric constant of the composite material.These results can be directly correlated to experimental measurements with some assumption.Speci cally, according to Boutros and Hanna (Boutros et al. 1978), the dielectric properties of cellulose exhibit notable variations with changing moisture levels, assuming that it is directly correlated to the amount of water that can be sorbed by cellulose structure.For instance, they showed that, at frequencies of 10 Mc/sec and a temperature of 25°C, the dielectric constant of cellulose exhibits a progressive increase from 3.11 to 4.59, 5.25, 5.92, and 5.96 as relative humidity levels rise from 0-35%, 52%, 76%, and 92%, respectively.By seeing these variations, we can infer that the observed differences in the dielectric properties of the experimentally measured cellulose structures can be primarily attributable to the sorption of water molecules.Furthermore, it is reasonable to hypothesize that the lower range of experimentally observed values aligns closely with a dry cellulose structure and, hence, is in close agreement with our DFT results.

Conclusions
Utilizing rst-principles calculations, we have comprehensively investigated cellulose I α and I β , from predicting the most energetically favorable arrangements of H atoms in the partial occupancy models of both phases to analyzing their optical and electronic properties.Our results con rm that using various XC functionals, coupled with scissors correction, can reliably capture the physics and chemistry of cellulose I α and I β with high accuracy.We also reveal the strong phase-dependent anisotropy in dielectric properties.However, due to the complex structure of cellulose samples and thermal disorder, such anisotropy may not be fully detectable experimentally unless single-crystal samples are analyzed.In addition to electronic contributions (computed herein within the independent particle approximation), our work discusses the often-overlooked ionic contributions to the dielectric constants, demonstrating, for instance, the difference in static dielectric constants for the two phases.Furthermore, we bridge the gap between theoretical calculations and a wide range of experimentally obtained dielectric constants by explaining the impact of potential water sorption on the dielectric properties of cellulose phases and its connection to experimental data.Using the framework of volume-average theory, we thus explain the discrepancies observed in different experimental measurements.This reveals that the experimentally measured dielectric properties often do not correspond to an ideal, dry cellulose structure but are signi cantly in uenced by water sorption.This work thus offers a valuable benchmark for future experimental/theoretical studies and highlights the need for careful control of environmental factors (e.g., humidity), especially in cases requiring speci c dielectric properties.
-Centre of Excellence for nanophotonics, advanced materials and novel crystal growth-based technologies" project (GA No. MAB/2020/14) carried out within the International Research Agendas programme of the Foundation for Polish Science conanced by the European Union under the European Regional Development Fund and the European Union's Horizon 2020 research and innovation programme Teaming for Excellence (GA.No. 857543) for support of this work.We gratefully acknowledge Poland's high-performance computing infrastructure PLGrid (HPC Centers: ACK Cyfronet AGH) for providing computer facilities and support within computational grant no.PLG/2023/016228.Ethical Approval Not applicable Funding This work has been supported by International Research Agendas programme of the Foundation for Polish Science co-nanced by the European Union under the European Regional Development Fund and the European Union's Horizon 2020 research and innovation programme and Teaming for Excellence with the grant numbers (GA No. MAB/2020/14) and (GA.No. 857543).Con ict of interest The authors have no relevant nancial or non-nancial interests to disclose.