Scaling behaviour for the water transport in nanoconfined geometries

The transport of water in nanoconfined geometries is different from bulk phase and has tremendous implications in nanotechnology and biotechnology. Here molecular dynamics is used to compute the self-diffusion coefficient D of water within nanopores, around nanoparticles, carbon nanotubes and proteins. For almost 60 different cases, D is found to scale linearly with the sole parameter θ as D(θ)=DB[1+(DC/DB−1)θ], with DB and DC the bulk and totally confined diffusion of water, respectively. The parameter θ is primarily influenced by geometry and represents the ratio between the confined and total water volumes. The D(θ) relationship is interpreted within the thermodynamics of supercooled water. As an example, such relationship is shown to accurately predict the relaxometric response of contrast agents for magnetic resonance imaging. The D(θ) relationship can help in interpreting the transport of water molecules under nanoconfined conditions and tailoring nanostructures with precise modulation of water mobility.

is obtained by the MSD calculated between 600 and 1000 ps.     Notice that in Case 10 and Case 12 NPs are initially placed randomly within the silica nanopore, whereas in all other cases NPs are initially placed on the surface of silica nanopore, where they tend to adsorb during the remaining time of computation.
Finally, the influence of temperature on has been explored, considering Case 7 as a reference and changing the temperature of the system to 350 K and 280 K, the measured are m 2 s -1 and m 2 s -1 , which implies a linear trend between and temperature in accordance with the Einstein relation. However, the anomalous dynamics of nanoconfined water at different temperatures 5,6 and the supercooled regime at low temperatures 7,8 are beyond the aims of this work and will require further investigations. where the of water is evaluated. Note that in cases where water molecules are highly confined (e.g. Case 1 in Supplementary Figure 15), the sub-diffusion regime can last tens of picoseconds, because the proximity of most of the water to a surface significantly affects the cage-breaking needed by the diffusing molecules for escaping their own hydration shell. 12,13 Finally, the convergence of in 1 ns runs is verified by evaluating it every 200 ps and by verifying that tends to an asymptotic value, as shown in Supplementary Figures 1921. For the sake of completeness, a few setups are extended up to 2 ns, in order to further assess the stability of in longer runs. In particular, two silica nanopores (Supplementary Figures 22,   23), two silica/magnetite nanoparticles (Supplementary Figures 24, 25) and two silica nanopores filled by magnetite nanoparticles (Supplementary Figures 26, 27) are considered.

Supplementary Note 2. Convergence of Molecular
Note that error bars are obtained by fitting MSD in different time intervals, whereas the dotted line indicates the value of water self-diffusion coefficient evaluated in the 600-1000 ps interval.
The extended simulations confirm that converges to a steady state value after approximately 0.6 ns, then it fluctuates (less than ±10%) around the equilibrium value.
After a few hundreds of picoseconds, the stability of the simulated protein structures is assured For the magnetite NP ( 2.0 nm), the first layer of water shows a clear peak of 1.5 , and a second layer of water is also visible with a peak of 1.1 (Supplementary Figure 1 b, e). The higher peak in the density profile for magnetite is determined by the selection of a stronger force-field.
The density distribution of water within a system composed by a silica nanopore ( 8. Note that the hydrophilic liquid-solid interaction within the first layers of adsorbed water induces a strong distortion of the hydrogen bond network, which in turn significantly affects the dynamics of confined water. 17,18 Of course, the level of hydration determines the effective pressure, namely the sum of the bulk pressure and the spreading pressure, 19 inside the low mobility region due to nanoconfinement and consequently the thermodynamic state of water in the low mobility region. Eventually the effective pressure could determine a switch from a low-density liquid to a high-density liquid. 20 Our observations here refer to the mobility of nanoconfined water, which are not greatly affected by the actual local density.

Alternative definitions of dimensionless scaling parameter
In this subsection, two alternative definitions of the scaling parameter are discussed in order to further support the one discussed in the main text. The first one is based on the notion of interfacial work; while the second makes an attempt to incorporate also the water density / hydration level.

θ accounting for interfacial work
Before considering the discussed procedure, a preliminary effort has been made to search for a proper scaling quantity. Inspired by other succesful works, where a scaling behaviour could be found by resorting to the notion of interfacial work, 21,22 our first attemp in searching for the scaling parameter was focused on the following argument.
Let us consider the schematics in Figure 2b. For an arbitrary particle , inspired by the computation of , we also define a characteristic energy : where represents the well depth of the potential energy shown in Figure 2c). Since the interfacial work between a particle and the solvent is proportional to the above , a possible guess for the scaling quantity would be: where and are the Avogadro Number and the volume of influence of particle , respectively. However, if is assumed as a unique independent variable for scaling the values, a poor correlation appears, as evident in Supplementary Figure 31. Hence, was judged not suitable for our purpose.

θ accounting for water density
The suggested scaling parameter in Equation (3) does not include the density of water within the analyzed configurations. In fact, the considered MD setups are characterized by a range of hydration levels (from 660 to 1080 kg m -3 ), where no heterogeneous wetting and anomalous behavior due to low water filling regimes is expected. 6,23 Towards an effort of incorporating also the hydration level, the following variable was also considered as a scaling parameter: , where is the bulk density of water (given the pressure and the temperature), is the actual water density in the setup, while is the suggested scaling variable previously defined by Equation (3).
In Supplementary Figure 32 and Figure 3 we report the results obtained when the scaling variables and are used, respectively.

Estimate of the scaling parameter for the literature data
In Figure 3 we also included results from the literature. By referring to the work of Milischuk et al., 24  , with being the volume of water within the pore. Finally, both the axial ( ) and radial ( ) diffusion coefficients are provided, hence here we report the value of: .
In reference 25 water confined within five single walled carbon nanotubes (CNTs) is simulated.
By referring to Figure 5 of reference 25 , again both the axial ( ) and radial ( ) diffusion coefficients at T=298 K are provided, hence for each configuration we are here interested in the value of: . Let and be the nominal tube diameter and the minimal distance between a carbon atom and the water molecule, respectively. Owing to the high curvature of those cases, the volume of influence has been more accurately computed as: , where is the tube length while . On the other hand, the total volume of water reads: ( ) , therefore the scaling parameter (for small CNTs) takes the simple form: In Figure 3 we reported cases with and , corresponding to chirality (8,8), (10,10), (12,12), (14,14) and (16,16), respectively, with (consistently with the value reported in reference 25 ).
Finally, we considered the results in Figure 2 of reference 26 , where the radial profile of the water diffusion coefficient around myoglobin has been reported as a function: (see also reference 27 ). Clearly, at a fixed distance from the protein ̅ , a value for the scaling parameter can be computed as: , while the corresponding (averaged) value of the diffusion coefficient reads: with denoting the location of the water molecules nearest to the protein. In Figure 3, we report five values from the curve in Figure 2 of reference 26

D-θ model: further details about thermodynamic insights
In this subsection, we interpret the scaled results reported in Figure 3 by a model based on experimental data of nanoconfined water.
Unfortunately, despite its fundamental importance in science and technology, the physical properties of water are far from being completely understood. 20 The liquid state of water is definitively unusual because it behaves as if there exists a singular temperature toward which its thermodynamic properties diverge. 20 Hence, the efforts of scientists from many disciplines to seek a coherent explanation for this unusual behavior make water one of the most important open questions in science today. 20 Very briefly, one could say that water in the presence of heterogeneities (e.g. in contact with foreign substances or surfaces) freezes at atmospheric pressure at the melting temperature 273.15 K. However, pure bulk water can be supercooled below its melting temperature down to the homogenous nucleation temperature 235 K, below which it inevitably crystallizes. In the range between the melting temperature and the homogenous nucleation temperature, liquid water exists in a supercooled metastable state, which is characterized by anomalies in its thermodynamic properties. 20 There is a growing interest in the prediction of properties of supercooled water. In particular, in applied atmospheric science, it is commonly accepted that the uncertainties in numerical weather prediction and climate models are mainly caused by poor understanding of properties of water in tropospheric and stratospheric clouds, where liquid water can exist in a deeply supercooled state. 28 In particular, thermodynamic properties of supercooled water diverge toward a singular temperature 228 K. 20 28 It is worth to point out that some insights into phases of liquid water can be obtained by alternative interpretations, for example based on the order-disorder transition scenario. 30 Actually, there are at least two methods to explore the no-man's-land for liquid water, namely the liquid states at temperature below . One possibility is to use molecular solutions rich in water. 30 Another one, which is receiving a growing attention because of the practical implications in nanotechnology, is to study water in nanoscopic confinement. Water near surfaces often does not crystallize upon cooling, but only recently have the properties of such water been measured. 30 For example, when confined in nanopores made of porous hydrophilic silica glass (Vycor), water does not crystallize and can be supercooled well below . 20 Nanoconfinement seems to confirm the assumption of a second critical point scenario, because it clearly reveals the existence of a peak in the specific heat capacity (Widom line) and a dynamic crossover in the self-diffusion properties of water, at 225 K, between the two dynamical behaviors known as fragile and strong, which is a consequence of a change in the hydrogen bond structure of liquid water. 20,31 Clearly this situation resembles the pseudo-critical temperature in supercritical fluids. 32 Nanoconfinement seems to damp out the peaks in heat capacity, rather than to exaggerate, the behavior in real water 30 (see, for example, Figure 1 in reference 31 ).
Some authors attribute this smoothing to the confinement length scale effects. 33 Alternatively, one could interpret the results by imaging that nanoconfinement acts as a pressure. 30 This would be consistent with the concept of spreading pressure, which is popular in the thermodynamics of adsorption and, more generally, in the non-equilibrium thermodynamics of heterogeneous systems. 19 Coming back to our scaled results reported in Figure 3, it is reasonable from the previous discussion to describe water in the region with lower molecular mobility (C) as supercooled liquid water with strong nanoconfinement (characteristic length is smaller than 0.6 nm). The thermodynamic state depends on the confinement length scale 33 but only few experimental data are available. For example, let us consider the experimental data about the specific heat capacity of ordinary water confined in silica nanopores as a function of the system temperature ( Figure 4a). 31 This thermodynamic property depends on the diameter of the nanopore.
However, for very small nanopores (<1.7 nm), the heat capacity shows no meaningful peak, which suggests that no ice is present (effectively too far from the virtual LLCP because of the surface spreading pressure). In case of strong nanoconfinement, the heat capacity seems to move gradually from a temperature-independent plateau, which starts at room temperature, to a linear function of the temperature, which is typical for ice, at very low temperature (<100 K).
In the following, we base our analysis on the experimental curve corresponding to the system under the strongest nanoconfinement (nanopore with diameter 1.7 nm), and assume that such a results is the most representative of a generic configuration with .
By numerical integration, the function ∫ can be evaluated with , and 300 K being the room temperature (Figure 4b). 31 For any reduction of the system energy , the previous function allows one to estimate Clearly, indicates the energy of adsorption. 19 The potential well is the average potential in the confined region defined by , which is comparable with the water molecule size. This means that only 1-2 water layers fit in the potential well and hence , where is the (total) potential minimum, is a good approximation. In this way, the self- which well matches both our numerical simulations and data from the literature (see Figure 3) with an excellent coefficient of determination ( 0.93), and most importantly without the need of introducing any empirical factors.

Outer sphere theory
Transversal relaxivity of contrast agents (CAs) such as SPIO nanoparticles is given by (9) being the diamagnetic contribution ( 2.8 s), and the iron concentration in solution, in mM. is the transverse relaxation time of NP, which can be modeled by the quantum-mechanical outer-sphere theory 34-36 (10) when the CA fulfills the condition (11) which indicates a motional averaging regime (MAR). In a MAR regime water protons, because of their diffusion, experience changing magnetic fields, which are effectively time-averaged.
Supplementary Equation (11) is generally fulfilled in case of ultrasmall NPs (i.e. diameters less than 10 nm). 37 is the volume fraction of NPs, which can be expressed as (12) where is the effective radius of NP (magnetic core with radius; coating layer with thickness), is the Avogadro number and is the molarity of SPIOs, in mM. is the time needed for water to diffuse across the NP surface; whereas is the angular phase shift induced by the NP at the equator line on its surface. Note that it has been proven that outer sphere theory accurately predicts of SPIOs and other superparamagnetic contrast agents. 36,37 Different definitions of and have been suggested in literature, 36 where is self-diffusion coefficient of water molecules around the CA, is the proton gyromagnetic ratio and is the Curie moment of SPIOs. In Supplementary Equation (13) the choice of is justified by arguing that the coating is not completely impermeable, and water molecules can diffuse until the crystal core, especially in the case of shorter and less dense coatings. At high magnetic field, it is possible to approximate , where is the saturation magnetization of NP. 39 Therefore, Supplementary Equation (10) may also be recast as: which can be finally simplified as (16) by introducing the volume fraction of the magnetic crystal core .

Molecular Dynamics simulations of nanopores and nanoparticles
A first set of geometries analyzed by MD simulations are cylindrical nanopores and spherical nanoparticles (NPs). Nanopores are made out silica; NPs are made out magnetite or silica. Subsequently, a nanopore in the above brick is obtained by removing all atoms whose distance from the brick center is smaller than a fixed length (i.e. nanopore radius). For the sake of simplicity, since our main concern is to investigate on the trend of the self-diffusion coefficient of water molecules at different SPIOs concentrations, here amorphization of silica is neglected.
Upon the creation of the pore (see the rightmost picture in Supplementary Figure 5), all silicon atoms located along the "cut surface" with only one bonded oxygen atom are removed. In addition, one hydrogen atom is attached to all oxygen atoms that are missing one bond with silicon (surface oxygen). This is achieved by imposing that the angle formed by silicon, oxygen and hydrogen is 128.8 degrees (elevation), 42  130° (SPASIBA, empirical force field) for Fe 3+ . 44 In our simulations, we found that the numerical results are not very sensitive to this angle and hence we used the same angle of the silanol group (129°) for the sake of simplicity in the geometry preparation, with a random azimuth angle (see Supplementary Figure 7b). Such an approach of modeling the surface aims at mimicking real SPIO particles without coating, 45 and it represents a simple technique already validated in other works for iron nanoparticles in aqueous environment. 46 On the other hand, here the effect due to complex coatings will be indirectly taken into account by a sensitivity analysis of the nonbonded interactions (Lennard-Jones and partial charges parameterizations).
Towards the end of minimizing the effects due to the residual electrical dipole induced by the charges located at the NPs surface, in all the studied setups NPs are inserted in pairs within the nanopore, where the first particle (which undergoes a random rigid rotation on each Cartesian axis) is initially opposed to its mirror image with respect to the midpoint of the line segment connecting the centers of the two particles (see Supplementary Figure 8). In the initial configuration of our simulations, all NPs pairs have the latter segment parallel to the pore axis.

Molecular Dynamics simulations of proteins
The Note that crystal water molecules are removed from original PDB files, in order to fully solvate the protein by means of the GROMACS' tool genbox.

Molecular Dynamics simulations of carbon nanotubes (CNTs)
The (12351 water molecules) nm 3 SPC/E water boxes, respectively. Notice that water density is 1000 kg m -3 in all simulated CNT setups, which corresponds to fully hydrated CNT surfaces.

Molecular Dynamics force-fields
Two types of interactions are considered in the MD simulations: i) Bonded interactions, among the atoms forming nanopores, nanoparticles, CNTs or proteins; ii) nonbonded interactions, between the water molecules and the solid surfaces, described via van der Waals and Coulomb potentials.
First, bonded and nonbonded interactions of proteins are modeled using the GROMOS96 43a2 force field, 72 which has been widely used for studying water dynamics in the proximity of protein surfaces. 73,74 Second, in the silica structure, the bonded interactions are modeled by means of two harmonic terms, adopted to describe the silicon-oxygen and oxygen-hydrogen interactions. That is, a bond stretching potential between two bonded atoms and at a distance (around the equilibrium distance ), and a bending angle potential between the two pairs of bonded atoms and (around the equilibrium angle ) are considered as follows: with parameters reported in Supplementary Table 12.
Here and hydrogen atoms, respectively.
In the magnetite structure, sufficiently high values for the force constants in Supplementary Equation (17) are assumed for all bonded interactions (rigid particle assumption), namely 400000 kJ mol -1 nm -2 and 400 kJ mol -1 rad -2 . In fact, our main concern is to investigate water self-diffusion coefficient of water in nanoconfined geometries, which is affected by nonbonded interactions. Therefore, here it is not of interest to accurately describe the fast dynamics within the magnetite NPs, thus the rigid particle assumption does not affect the measurements of self-diffusion coefficient of water, as confirmed by preliminary sensitivity analyses with respect to and (results not shown).
In the CNT structure, the carbon-carbon bonded interactions are also modeled by two harmonic terms (Supplementary Equation (17)), where 478900 kJ mol -1 nm -2 , 0.142 nm, 562.2 kJ mol -1 rad -2 and 120°. [75][76][77] Nonbonded interactions among silica atoms (consisting in both van der Waals and electrostatic interactions) are also taken into account through: (i) a Lennard-Jones term with mixed parameters consistently chosen according to the following Lorentz-Berthelot combination rules (ii) a Coulomb term (20) with being the permittivity in a vacuum, and , the partial charge of atoms and , respectively. Non-zero partial charges are assigned only to atoms at the surface of the pore and belonging to a silanol group, whereas all other atoms (bulk of silica) are considered neutral (i.e. zero partial charge =0). More specifically, partial charges in a silanol group are assigned following the criterion of the overall neutral charge for the entire system, with nonbonded parameters for silica reported in Supplementary Table 13.
Similarly to silica, partial charges are assigned only to atoms at the magnetite surface and belonging to the Fe-O-H groups, whereas zero partial charges (i.e. =0) are imposed at all other bulk atoms. Moreover, partial charges are assigned within the Fe-O-H groups in order to ensure neutrality of the whole system. Both the adopted parameterization for Lennard-Jones potentials and the partial charges are reported in Supplementary Table 14.

Algorithm for automatic evaluation of θ
By referring to the flow-chart in Supplementary Figure 30, we report the main steps involved in the computation of the scaling parameter . For the sake of completeness and without loss of generality, a few GROMACS commands are also given as example. We note in fact that, the described procedure can be properly rearranged by a well educated user of other MD software packages. 10. Compute .

Details of experimental procedures and materials
For the hydrophilic SPIOs loaded into the SiMPs, the samples as provided by the vendors (Sigma-Aldrich) presented several aggregates, and the following procedure was performed in order to purify the original solution and select the SPIOs with the higher stability in solution.
Upon sonication (~ 15 min -Bransonic Ultrasonic Cleaner), the sample were centrifuged (6 minutes, at 12,000 rpm) and the supernatant was collected. This step was repeated twice. The SiMPs were fabricated by using previously reported protocols. 78 The fabrication process consists of three major steps: Formation of porous silicon films; photolithographic patterning of particles; and Reactive Ion Etch (RIE). The porous structure was tailored by electrochemical etching while the particle sizes were precisely defined by photolithography. Since the porous structure and the particle size are controlled independently, a wide range of sizes, shapes and pore morphologies can be obtained using such an approach. In this work, particles with 1,000 nm in diameter and 400 nm in thickness were used. These particles have a mean pore size of 40 nm and a porosity of about 60%. The fabrication process is briefly described as follows: Starting with heavily doped P-type (100) wafer with resistivity of 0.005 ohm-cm as the substrate (Silicon Quest, Inc, Santa Clara, CA), the wafer was assembled on a home-made anodizing cell with the polished surface immersing in 1:3 HF(49%):ethanol solution. An etching current of 6 mA/cm 2 was applied for 125 sec to generate 400 nm porous silicon film. Then a high electrical current with current density ~76 mA/cm 2 was applied for 8 sec to form the instable release layer. An 80 nm low-temperature oxide (LTO) was deposited on the porous silicon film in a LPCVD furnace.
A standard photolithography process was used to pattern the 1,000 nm circles on the film using a As per the loading of SPIOs into SiMPs, the stock solutions of SPIOs as purchased (5 nm from Sigma) was first purified to select particles with higher aqueous stability. Then, SiMPs were lyophilized to dryness for 8 hours and 2×10 8 SiMPs were exposed to 100 μL of the purified