Abstract

The solubility of cellulose in water-based media is promoted by low temperature, which may appear counter-intuitive. An explanation to this phenomenon has been proposed that is based on a temperature-dependent orientation of the hydroxymethyl group. In this paper, this hypothesis is investigated using molecular dynamics computer simulations and NMR spectroscopy, and is discussed in conjunction with alternative explanations based on solvent–solute and solvent–solvent hydrogen bond formation respectively. It is shown that neither simulations nor experiments lend support to the proposed mechanism based on the hydroxymethyl orientation, whereas the two alternative explanations give rise to two distinct contributions to the hydration free energy of cellooligomers.

Keywords

Introduction

Cellulose is insoluble in pure water at ambient conditions. There are however water-based solvent systems for cellulose, most notably solutions of alkali salts (Budtova and Navard 2016); NaOH and LiOH, with or without the addition of a co-solvent such as urea (Cai and Zhang 2005). Interestingly, these solvent systems become more efficient as temperature is lowered, and moreover, they have an optimal dissolution temperature well below the freezing point of pure water. This may seem counterintuitive, since the solubility of polymers in general is promoted by high temperature, due to more favorable entropy contributions, −TΔS, both from the polymer conformational entropy, and from the ideal mixing part. These contributions are always working in favor of dissolution, decrease with T, and are independent of the solvent. The reason for the inverted T dependence in the case of cellulose must therefore be a consequence of cellulose-solvent interactions. However, at the same time, aqueous solubilities of short-chain oligosaccharides, from glucose up to, at least, cellopentaose (Taylor 1957) actually increase with temperature, both in pure water and, for cellobiose, also in salt solutions (Liu et al. 2016). Apparently there must be some critical degree of polymerization above which the inversion of the temperature dependence takes place.

The question then becomes: What is it that causes the unexpected temperature dependence? Suggestions that have been brought forward include polarization effects in both solute and solvent (Budtova and Navard 2016), and temperature dependent coordination numbers of the salts present in the solution (Liu et al. 2016), without any attempts to further clarify the underlying mechanisms. Lindman and Karlström have in a number of publications (Lindman and Karlström 2009; Lindman et al. 2010; Medronho et al. 2012) argued that the mechanism behind the inverted T dependence is due to the fact that the conformation of the hydroxymethyl group depends on T. The orientation of the hydroxymethyl group is defined by the rotation around the C5–C6 bond (Fig. 1) defined by the dihedral angle ω. It is usually found in one of three staggered conformations denoted tg, gt or gg, where the first letter refers to the orientation (trans or gauche) with respect to O5–C5–C6–O6, and the second refers to the orientation with respect to C4–C5–C6–O6. Due to the distribution of charge on the atoms, it is clear that the different orientations will differ in their contribution to the molecular dipole moment. Thus, Lindman and Karlström suggest, based on ideas originally developed for poly-(ethylene oxide) (PEO) (Karlström 1985), that as temperature decreases, the more polar conformations will be increasingly populated. This would lead to more favorable interactions with a polar solvent (i.e., water), thereby enhancing the solubility. In one paper (Lindman and Karlström 2009), they discuss their idea within the framework of Flory–Huggins solution theory, and compare it to two previously proposed mechanisms for reversed solubility of PEO, which were based on the disruption of, respectively, the water hydrogen bond network (Kjellander 1982), and the solute-water hydrogen bonds (Goldstein 1984). They concluded that they were not, at that time, able to determine which of these mechanisms were responsible for the observed behavior, and that the true mechanism possibly could be a combination of all three, but also that there were experimental evidence that the polarity of (ethylene oxide)-containing molecules decreased with increasing temperature, in strong support of their own hypothesis. It was obvious however that only the model based on breaking of water-water hydrogen bonds could potentially explain the similar trends in solubility found for non-polar solutes like benzene (Privalov and Gill 1989).

Molecular model of cellotetraose, without non-polar hydrogens as in the GROMOS force field

In relation to this observation, it is interesting to note that cellulose has a strong amphiphilic character, with much of its properties in aqueous environments, including solubility, thus being governed by hydrophobic effects (Bergenstråhle et al. 2010; Miyamoto et al. 2014; Lindman et al. 2010; Medronho et al. 2012). It is tempting to assume that similar mechanisms as for purely hydrophobic substances could be important for carbohydrates as well. One piece of evidence in this direction is that cellulose, just as many other hydrophobic solutes that are insoluble in water at room temperature, becomes water soluble at supercritical conditions (Tolonen et al. 2015; Wohlert et al. 2015; Adschiri et al. 1993).

We have recently published a study in which the hydroxymethyl orientation in methyl β-cellobioside and cellopentaose in water was studied using solution state NMR and Molecular Dynamics (MD) simulations (Angles d’Ortoli et al. 2015) in which it was found that the distribution of dihedral angles around the C5–C6 bond exhibited only small variations with temperature. However, a possible relation between ω and the aqueous solubility of the carbohydrate was not fully addressed. In the present study, molecular dynamics (MD) simulations and NMR spectroscopy are used to test Lindman and Karlström’s hypothesis. Specifically, we seek answers to the following questions:

1.

Can simulations reproduce the temperature dependence of cellulose solubility?

2.

Does the conformation of the hydroxymethyl group depend on temperature, and, specifically, does the population in the tg conformation increase with T?

3.

Does the dipole moment of the hydroxymethyl side chain decrease with increasing temperature?

4.

Can other explanations for the anomalous temperature dependence be ruled out? Specifically, what is the role of hydrophobic solvation forces?

To that end, simulations of cellotetraose in water were used both to calculate the hydration free energy in pure water as a function of temperature, and to study the conformation of the hydroxymethyl group. Two different carbohydrate force fields were employed in order to as far as possible avoid model dependent effects. In addition, novel NMR experiments were performed to validate the simulations. To answer the first question (i) above, we show, in the Theory section, that although their solubility actually increases with T, results from short chain oligosaccharides can be extrapolated to cellulose.

Theory

Any dissolution process can be modeled as a two-step process; a first step in which the solute is brought from its solid (crystalline) state into its vapor phase, and subsequently a step where it becomes solvated. Thus, one can write the total free energy change upon dissolution as \(\Delta G = \Delta_{\text{sub}} G + \Delta_{\text{hyd}} G,\) where the first term represents the free energy of sublimation, and the second one the free energy of hydration. Following Wohlert (2014), the sublimation term consists of one enthalpic part, and one entropic part, which, in turn, has one ideal and one conformational contribution. The ideal contribution can be obtained using the Sackur-Tetrode equation for the translational degrees of freedom, and a rigid-rotor approximation for the rotational degrees of freedom. This contribution will depend on the chain length n as log(n). The conformational contribution should, on the other hand, be more or less linear in n since it is a function of the number of internal degrees of freedom, although cross-correlations between vibrational modes do play a significant role for the internal entropy (Baron et al. 2009). The enthalpic part of \(\Delta_{\text{sub}} G\) is, to a good approximation, linear in n but independent of temperature (Wohlert 2014). The free energy of sublimation thus contains contributions that scale differently with chain length, and also some that depend on temperature, and some that do not.

The hydration part must be approximately linear in n, since both excluded volume and solvent-exposed surface area of a stiff chain scales linearly with the number of segments. Furthermore, we will assume that the free energy of hydration per segment, \(\Delta_{\text{hyd}} G'\), is a function of T, \(\Delta_{\text{hyd}} G^{\prime} = \Delta_{\text{hyd}} G^{\prime}\left( T \right)\). Thus, combining \(\Delta_{\text{sub}} G\) and \(\Delta_{\text{hyd}} G\), keeping only those terms that depend on both n and T, the dissolution free energy can be written:

This means that even if the first term dominates for small n, it will become less important when n is large. Thus, even if \(d\left( {\Delta G} \right)/dT < 0\) for short chain cellooligomers, this trend becomes weaker with larger n, and, as long as \(\left| b \right| < \left| {d\left( {\Delta_{\text{hyd}} G'} \right)/dT} \right|\) and \(d\left( {\Delta_{\text{hyd}} G'} \right)/dT > 0\), it is inevitably reversed when n becomes large enough.

For the purpose of studying cellulose solubility as a function of temperature, it is thus sufficient to show that \(\Delta_{\text{hyd}} G\) increases with temperature in the relevant temperature range, and moreover that it is sufficiently large to overcome any opposing contribution from conformational entropy of the oligosaccharide chain.

Another approach that was used by Lindman and Karlström (2009) is based on Flory–Huggins solution theory, in which the phase behavior of the cellulose/water system is determined by the free energy of mixing:

$$\Delta G_{\text{mix}} = U_{\text{mix}} - TS_{\text{mix}}$$

(3)

Here, the entropy of mixing, at constant number of molecules (note that Lindman and Karlström (2009) use constant number of monomers), will scale as log(n), meaning that we can identify it as the ideal contribution above. This term always favors mixing, and thus we need a mixing energy that, for large enough values of n, has a temperature dependence that is stronger than the entropy term. The mixing energy was formulated using an effective interaction parameter which they denoted \(w_{ 1 2}^{\text{eff}}\), with the demand that it was temperature dependent. This interaction parameter is obviously related to the hydration free energy as formulated above, and it is clear that it contains not only direct water-cellulose interactions, but also indirect contributions originating from the perturbation of water-water interactions. As such, it may include both energetic and entropic contributions.

Materials and methods

Molecular dynamics simulations

Molecular dynamics simulations of β-cellotetraose in water were performed using GROMACS 5 (Páll et al. 2015; Pronk et al. 2013). Aqueous cellotetraose was modeled using both the GROMOS 56 carbohydrate force field (Hansen and Hünenberger 2010) in combination with the SPC/E model (Berendsen et al. 1987) for water, and the CHARMM C35 parameter set (Guvench et al. 2008; Guvench et al. 2009) together with TIP3P water (Jorgensen et al. 1983). Water molecules were kept completely rigid using SETTLE (Miyamoto and Kollman 1992), and all covalent bonds in cellotetraose were kept at their equilibrium distances using P-LINCS (Hess 2008). The equations of motion were integrated using a Langevin leap-frog algorithm (Goga et al. 2012) with a basic time step of 2 fs, a friction constant of 0.5 ps−1 and several reference temperatures between 263 K and 450 K. Pressure was kept at 1 atm using a Parrinello-Rahman barostat (Parrinello and Rahman 1981), with a coupling constant of 20 ps.

The cut-off for non-bonded interactions was 1.4 nm, and long-range electrostatic interactions were included, in for case of GROMOS using Reaction-Field electrostatics with a relative permittivity of 61, and for the case of CHARMM, using PME (Darden et al. 1993; Essmann et al. 1995).

Two sets of simulations were performed. One consisted of 50 ns at each temperature, which was used to extract equilibrium properties such as hydroxymethyl populations and J coupling constants. The other set of simulations employed a thermodynamics integration scheme, in which cellotetraose-water interactions were decoupled using a coupling constant,λ. The decoupling was performed as two separate parts, one for the Coulomb interactions, and one for the Lennard-Jones interactions. Water-water interactions and carbohydrate-carbohydrate interactions were left unmodified, meaning that the completely decoupled case corresponded to a box of pure water, and one cellotetraose molecule in the gas phase. For Coulomb interactions, 13 values for λ between 0 (unmodified) and 1 (decoupled) were used: λ = [0 0.2 0.4 0.5 0.6 0.7 0.8 0.90 0.95 0.97 0.99 0.997 1.0], and for Lennard-Jones interactions 11 values were used: λ = [0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1]. For each λ, the instantaneous value of dG/dλ was calculated at each step of a 5 ns long simulation. The total free energy difference was subsequently computed using Bennet’s acceptance ratio (Bennet 1976) as implemented in GROMACS.

To enhance sampling at the lower temperatures, the simulations were set up as a replica exchange scheme (Sugita and Okamoto 1999) using 12 different temperatures between 263 and 329 K at 6 K intervals. Exchange between neighboring replicas was attempted every 100 steps. Simulations at higher temperatures were run separately. Replica exchange was thus used for both the equilibrium simulations and for the free energy calculations.

Ring conformations were predominantly found in 4C1, with only a minor fraction found in 1C4 (<1 % at 299 K, <2 % at 329 K), as determined by the pseudo-dihedral angles introduced by Strauss and Picket (Strauss and Pickett 1970; Rao et al. 1998).

Extraction of NMR coupling constants

In order to obtain experimental J values of H5–H6pro-R and H5–H6pro-S in the glucosyl residues for conformational analysis of the ω torsion angles the 1H NMR spectrum of cellotetraose were to be assigned as far as possible, aided by the NMR chemical shift prediction software CASPER (Lundborg and Widmalm 2011). Like in our previous study on methyl cellobioside and cellopentaose (Angles d’Ortoli et al. 2015) we utilized 1D 1H NMR and 1D CSSF-TOCSY spectra to obtain some of the pertinent \(^{3} J_{\text{HH}}\) values. A 2D 1H,13C-HSQC NMR spectrum revealed that up to five 13C resonances were possible to resolve in the spectral region of the hydroxymethyl groups, i.e., at 60.46–61.16 ppm. It would therefore be possible to extract additional or complementary \(^{3} J_{\text{HH}}\) information via a direct peak-maxima analysis using a 2D perfect-HSQC experiment (Castañar et al. 2015). That this is possible was validated using methyl cellobioside for which the \(^{3} J_{\text{HH}}\) had been extracted via NMR spin-simulation using the PERCH NMR software (Angles d’Ortoli et al. 2015). The direct peak-maxima analysis from 1D traces of the perfect-HSQC spectrum of methyl cellobioside resulted in excellent agreement to those determined by the NMR spin-simulation approach, namely, |ΔJ| ≤0.1 Hz. Additional \(^{3} J_{\text{HH}}\) values were thus obtained for cellotetraose from a direct peak-maxima analysis of its perfect-HSQC spectrum. Since not all \(^{3} J_{\text{HH}}\) values could be determined unambiguously for the eight spin-systems using the three different types of experiments, namely, 1H NMR, CSSF-TOCSY and perfect-HSQC, we decided to use the information gathered and to perform an iteration between tentative NMR assignments, J couplings and the 1H NMR spectrum per se using the PERCH NMR spin-simulation software, employing a methodology previously carried out for a large number of oligosaccharides(Roslund et al. 2011; Rönnols et al. 2013).

Results and discussion

Lindman and Karlström (2009), Lindman et al. (2010 and Medronho et al. (2012) have argued that a possible mechanism that makes ΔhydG of cellulose increase with T is related to the orientation of the hydroxymethyl group, which, according to a previously developed model(Karlström 1985), would attain less polar orientations at high temperatures, thus disfavoring dissolution. To test their hypothesis, the present paper seeks answers to a set of questions that were stated in the Introduction. These questions are treated in turn below.

Calculated free energy of hydration

As discussed above, a necessary condition for a decreasing solubility with increasing T is that the free energy of hydration is increasing with T, \(d\left( {\Delta _{\text{hyd}} G} \right)/dT > 0\), since the other contribution to the total free energy of dissolution, the free energy of sublimation, is always decreasing with T. To that end, the free energy of hydration, \(\Delta_{\text{hyd}} G\), for cellotetraose was calculated using MD simulations. It was found to increase with increasing temperature by approximately 30 kJ mol−1 over the temperature interval considered here, for both the GROMOS and CHARMM parameter sets (Fig. 2). Although the absolute values differ substantially between the two force fields (Table 1), the increase with T becomes almost surprisingly similar. From this, the increase in free energy per glucose unit was obtained as \(d\left( {\Delta_{\text{hyd}} G} \right)/dT =\) 0.11 kJ mol−1 K−1. The difference in absolute entropy between crystalline and gas-phase cellobiose has been calculated from computer simulations (Wohlert 2014). It was found to be completely dominated by the ideal contributions, with the configurational part being on the order 0.01 kJ mol−1 K−1 per glucose unit around room temperature. This is significantly smaller than \(d\left( {\Delta_{\text{hyd}} G} \right)/dT\), which means that the two conditions stated in the Theory section are met: the hydration free energy change is both positive and larger than the contribution from the configurational entropy change upon sublimation. Thus, simulations of short-chain oligosaccharides can reproduce the trend in solubility found experimentally for cellulose. In the following sections, possible reasons leading to \(d\left( {\Delta_{hyd} G} \right)/dT > 0\) will be examined in more detail.

Calculated hydration free energy of cellotetraose in water, using two different force fields; GROMOS (circles) and CHARMM (triangles). The total hydration free energy (blue) is obtained as the sum of a polar (black) and a non-polar (red) contribution. (Color figure online)

Table 1

Calculated free energy of hydration as a function of temperature

T [K]

ΔhydG (kJ mol−1)

GROMOSa

CHARMMa

263

−195

−269

269

−188

−267

275

−185

−263

281

−185

−261

287

−183

−258

293

−179

−255

299

−179

−252

305

−177

−250

311

−175

−247

317

−169

−244

323

−167

−243

329

−169

−240

aEstimated errors are on the order of 2 kJ mol−1

Rotation of the hydroxymethyl group

Here, the possibility that the increase in hydration free energy with temperature, \(d\left( {\Delta_{\text{hyd}} G} \right)/dT > 0\), is caused by a simultaneous change in conformational distribution of the rotamers of the ω torsion angles is investigated by analysis of results from NMR spectroscopy and MD simulations.

Results from NMR spectroscopy

The 1H NMR spectrum of cellotetraose shows severe spectral overlap of the resonances since the all the glucosyl residues are β-(1→4)-linked. In addition, the reducing end residue is present in an α/β-mixture. Thus, there are eight spin-systems of seven protons each (H1 to H6pro-S) denoted by Gα and Gβ at the reducing, followed by G′α and G′β, G″α and G″β, and finally G′″α and G′″β at the terminal end of the oligosaccharide (Fig. 1). The spectral region where the H6pro-R and H6pro-S resonances in cellotetraose reside is shown in Fig. 3a and the closely similar spectral appearance in Fig. 3b, based on NMR spin-simulation, shows that the details of the experimental spectrum have been captured in the NMR spin-simulation; thus, both 1H NMR chemical shifts and J couplings were determined for the spin-systems in cellotetraose and the pertinent ones for subsequent hydroxymethyl conformational analysis are compiled in Table 2. The populations of the ω torsion angles derived from an analysis of JH5,H6pro-R and JH5,H6pro-S coupling constants (Stenutz et al. 2002) in cellotetraose are presented in Table 3.

Experimentally determined and calculated populations (%) for the rotamers of the ω torsion angles in cellotetraose in D2O at 300 K

Pop. ω

G′′′β

G′′′α

G′′β

G′′α

G′β

G′α

Gβ

Gα

NMR

gt

52.72

52.14

43.13

41.42

43.07

45.79

44.08

36.98

gg

37.58

40.25

45.64

51.20

48.30

43.67

47.28

53.84

tg

9.70

7.61

11.23

7.38

8.64

10.54

8.64

9.18

MD

gt

54.7

–

34.0

–

33.6

–

40.0

–

gg

40.8

–

61.0

–

62.2

–

54.6

–

tg

4.4

–

5.0

–

4.2

–

5.4

–

Results from equilibrium MD simulations

A population analysis of ω at different temperatures between 263 and 329 K using the GROMOS force field is shown in Fig. 4. For the reducing and terminal residues (G and G′″) the experimental populations at room temperature are very well reproduced compared to the NMR experiments; for the G′ and G″ the deviation is larger but, most importantly, the internal order is correct (tg < gt < gg).

Rotamer populations; tg (black), gg (red), and gt (green), of the hydroxymethyl side chain in cellotetraose as a function of temperature calculated using MD with the GROMOS force field. Experimental NMR data at 300 K temperature is shown by asterisks of the corresponding color. (Color figure online)

J coupling constants are readily calculated from a simulations trajectory. They have the advantage that they can be directly compared to experimental data without having to assume a model for the population statistics. As shown in Table 3, calculated J coupling constants are comparable to experimental data, especially for the two outer residues.

Populations calculated from simulations using GROMOS show only marginal temperature dependence for the two outer residues, whereas changes in the two inner segments are larger. They are however much smaller than has previously ben determined using other parameter sets: Population statistics for ω as a function of temperature using the CHARMM force field was recently published for cellotetraose. As opposed to the GROMOS force field, it shows strong temperature dependence for both inner (G′,G″) and outer (G,G′″) residues (see Table 3 in Angles d’Ortoli et al. (2015)). Temperature dependence of similar magnitude was also noted using GLYCAM04 (Shen et al. 2009). In this work, simulations give an increase in tg of between 2 and 4 % for the present temperature interval, whereas changes in gt and gg generally are larger. The experimental temperature dependence of ω was not determined here, but it was previously determined for cellopentaose to be on the order of 1 % between 293 and 335 K (Angles d’Ortoli et al. 2015).

The dipole moment of the hydroxymethyl side chain

The assumption that less polar orientations of the hydroxymethyl group are increasingly populated as the temperature increases is central to the hypothesis of Lindman and Karlström. From the molecular structure of the glucose unit, it is easy to see that the dipole moment of the sequence O5–C5–C6–O6 should follow the order tg < gg ≈ gt (Fig. 1), although the orientation of the hydroxyl proton is likely to be important for the overall dipole moment as well. A sufficient increase of the tg conformation with increasing T would therefore lead to a molecule with smaller dipole moment at higher temperature and could thus potentially explain the T dependence of cellulose solubility.

Molecular dipole moments of the central glucose units were calculated from the trajectories from the GROMOS simulations at each temperature, both for the glucose units as a whole, and for the hydroxymethyl regions separately, and subsequently sorted into three groups pertaining to each of the three hydroxymethyl orientations. Here, the hydroxymethyl region was defined by the sequence O5, C4, C5, C6, O6, and H6, which constitute a group whose dipole moment is heavily dependent on ω, but is overall charge neutral.

The dipole moments followed the order tg < gg < gt, both for the glucose units as a whole, and for the hydroxymethyl regions (Table 4). The calculated values per hydroxymethyl orientation did not change with temperature and, more importantly, neither did the total dipole moments, averaged over all hydroxymethyl orientations (Fig. 5). Apparently, the small increase of the tg orientation was not sufficient to affect the overall polarity of the hydroxymethyl side chain.

Table 4

Calculated dipole moments in Debye for different hydroxymethyl orientations

Calculated average dipole moments as a function of temperature using the GROMOS force field for both the central glucose units as a whole, and the hydroxymethyl region. Error bars are smaller than the symbols

Thermodynamics of hydration from free energy MD simulations

Both the polar and the non-polar contributions to the total hydration free energy, which are obtained by decoupling electrostatic and van der Waal’s interactions separately in the free energy calculation, are increasing with T: the polar contribution by 10 kJ mol−1 and the non-polar by 15–20 kJ mol−1 over the present temperature interval (Fig. 2). These two respective contributions to the total free energy increase obviously have different sources.

One large contributing factor to the energy in this system comes from hydrogen bonds; from bonds between solvent molecules and the oligosaccharide and from solvent–solvent hydrogen bonds. Interestingly, both the average number of solvent–solute and solvent–solvent hydrogen bonds decreases significantly with increasing temperature (Fig. 6). The average (absolute) number of hydrogen bonds between the oligosaccharide and water molecules decreases by 2.7 over the temperature interval. By assuming a free energy of 9 kJ mol−1 per hydrogen bond (van der Spoel et al. 2006), this amounts to an increase in energy of 24 kJ mol−1 over this temperature range. Some of this will naturally be recovered by forming hydrogen bonds to other solvent molecules instead, but that amount will also decrease with temperature since the number of solvent–solvent hydrogen bonds decreases simultaneously by 0.3 bonds per water molecule. Hydrogen bonding energies are thus more than sufficient to explain the polar contribution to the increase in hydration free energy. This result can be related to one of the alternative models discussed in Lindman and Karlström (2009) based on solute solvent hydrogen bond interactions (Goldstein 1984), where the increase in free energy with temperature (expressed as an effective Flory–Huggins interaction parameter) was a result of a decreased probability of having low-energy, hydrogen-bonded, interactions at the expense of more non-hydrogen bonded interactions, which were assumed to be higher in energy.

Average number of solute–solvent (left) and solvent–solvent (right) hydrogen bonds as a function of temperature. Solvent–solvent bonds are normalized to the number of water molecules

The non-polar contribution to \(\Delta_{\text{hyd}} G\) is intriguing, since its major component represents the free energy cost of creating a cavity in the solvent, large enough to accommodate the solute. For a non-polar species, this would of course be the only contribution and therefore, when water is the solvent, it represents the hydrophobic part of the solvation.

Free energies of solvation of small hydrophobic solutes, like the noble gases and small hydrocarbons, are known to increase with temperature around room temperature. However, at higher temperatures the energy levels off, and eventually starts to decrease, leading to a parabolic shape for the hydration free energy curve, which can be said to be a hallmark of hydrophobic solvation (Garde et al. 1996; Privalov and Gill 1989; Chandler 2005). Free energy calculations using the CHARMM force field were extended to 490 K. The resulting free energy curve exhibits a parabolic shape; it levels off at around 400 K after which it is decreasing with T (Fig. 7). The accepted interpretation of such behavior in the case of non-polar solutes is that around room temperature, the creation of a cavity requires specific spatial organization of the hydrogen-bonding pattern of water (Chandler 2005), which comes at a large entropic cost, hence the increase of \(\Delta_{\text{hyd}} G\) with T. As temperature increases however, hydrogen bonding becomes less important (see Fig. 6), and is eventually overtaken by other contributions. This is obviously analogous to the model proposed by Kjellander (1982) and formulated mathematically by Lindman and Karlström (2009), based on that the effective interaction parameter between solvent and solute depends on whether the solvent is hydrogen bonded or not, and that the probability of solvent–solvent hydrogen bonds decreases with temperature. At the point where \(d\left( {\Delta_{\text{hyd}} G} \right)/dT = 0\), the entropy of hydration \(\Delta_{\text{hyd}} S = 0\). It has been shown that the temperature at which this happens is universal for the transfer of hydrophobic solutes into water, independently of the phase from which it originates, being it a liquid, gas, or even the interior of a protein (Murphy 1994). This temperature of entropy convergence is approximately 400 K (Privalov and Gill 1989; Murphy 1994; Garde et al. 1996). Incidentally, this seems to be valid for the non-polar contribution of carbohydrates as well.

Free energy of hydration of cellotetraose at high temperatures, from MD simulations using the CHARMM force field. The total contribution (blue) is obtained as a non-polar (red) and polar (black) part. (Color figure online)

Conclusions

The aim of the present study was to investigate mechanisms behind the fact that the solubility of cellulose in aqueous solvent is promoted by low temperatures, and specifically to seek answers to the following questions:

Can simulations reproduce the temperature dependence of cellulose solubility?

MD simulations give a free energy of hydration for cellotetraose that increases with temperature between 263 K and 329 K. Since ∆hydG scales linearly with the chain length n, while ∆subG scales as log(n) (Eq. 1), the solubility of cellulose above some critical degree of polymerization will increase when the temperature is decreased.

Does the conformation of the hydroxymethyl group depend on temperature, and, specifically, does the population in the tg conformation increase with T?

The orientation of the hydroxymethyl group is relatively insensitive to temperature. This means that the mechanism leading to d∆hydG/dT > 0 proposed by Lindman and Karlström (Lindman and Karlström 2009; Lindman et al. 2010; Medronho et al. 2012) based on that the equilibrium population of ω is shifted towards more polar conformations as the temperature is lowered, is likely to be of minor importance. However, the fact that the tg conformer is increasing with temperature means that there could be a small contribution, although it is counteracted by the concomitant increase of the gt conformer.

Does the dipole moment of the hydroxymethyl side chain decrease with increasing temperature?

Both the dipole moment of the hydroxymethyl side chain and of the glucose unit as a whole remain constant in this temperature interval. The slight increase of the tg conformer is too small to have any appreciable effect on the dipole moment.

Can other explanations for the anomalous temperature dependence be ruled out?

By separating ∆hydG into one polar and one non-polar contribution it was found that the polar contribution to d∆hydG/dT was about half the non-polar one around room temperature. This means that in pure water at low temperatures hydrophobic solvation forces dominate the temperature dependence of the solubility of cellooligomers. However, at elevated temperatures, the non-polar contribution levels off before eventually starting to decrease. This parabolic temperature behavior is analogous to what have been seen before for small hydrophobic solutes. The polar contribution on the other hand continues to increase even at high temperatures. This becomes a very clear manifestation of the amphiphilic character of carbohydrates and underlines the importance of accounting for hydrophobic effects when considering physical properties of cellulose.

Hydrogen bonding is important for the solution behavior of cellooligomers. Both solvent–solvent and solute–solvent hydrogen bonds decrease with increasing temperature. Whereas the latter is directly responsible for the increase of the polar part of the hydration free energy with temperature, the former can be related to the non-polar contribution, analogous to the hydration of purely hydrophobic solutes.

Notes

Acknowledgments

This work was supported by grants from the Swedish Research Council (GW), The Wallenberg Wood Science Center funded by The Knut and Alice Wallenberg Foundation (JW), and the Swedish Foundation for Strategic Research (SSF) through contract ICA 10-0086 (MB-W). Computer simulations were run at the PDC Center for High Performance Computing at KTH (PDC-2015-36).

Copyright information

Open AccessThis article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made.