# Collective modes and quantum effects in two-dimensional nanofluidic channels

Baptiste Coquinot, Maximilian Becker, Roland R Netz, Lyderic Bocquet, Nikita Kavokine

► To cite this version:

Baptiste Coquinot, Maximilian Becker, Roland R Netz, Lyderic Bocquet, Nikita Kavokine. Collective modes and quantum effects in two-dimensional nanofluidic channels. Faraday Discussions, 2024, 249, pp.162-180. 10.1039/d3fd00115f. hal-04228462

HAL Id: hal-04228462

https://hal.science/hal-04228462v1

Submitted on 4 Oct 2023

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
Collective modes and quantum effects in two-dimensional nanofluidic channels

Baptiste Coquinot$^{a,b,c}$, Maximilian Becker$^d$, Roland R. Netz$^d$, Lydéric Bocquet$^a$, and Nikita Kavokine$^{b,c*}$

$^a$Laboratoire de Physique de l'École Normale Supérieure,

ENS, Université PSL, CNRS, Sorbonne Université,

Université Paris Cité, 24 rue Lhomond, 75005 Paris, France.

$^b$Department of Molecular Spectroscopy, Max Planck Institute

for Polymer Research, Ackermannweg 10, 55128 Mainz, Germany.

$^c$Center for Computational Quantum Physics, Flatiron Institute, 162 5^th^ Avenue, New York, NY 10010, USA.

$^d$Fachbereich Physik, Freie Universität Berlin, Arnimallee 14, 14195 Berlin, Germany.

(Dated: June 1, 2023)

Nanoscale fluid transport is typically pictured in terms of atomic-scale dynamics, as is natural in the real-space framework of molecular simulations. An alternative Fourier-space picture, that involves the collective charge fluctuation modes of both the liquid and the confining wall, has recently been successful at predicting new nanofluidic phenomena such as quantum friction and near-field heat transfer, that rely on the coupling of those fluctuations. Here, we study the charge fluctuation modes of a two-dimensional (planar) nanofluidic channel. Introducing *confined response functions* that generalize the notion of surface response function, we show that the channel walls exhibit coupled plasmon modes as soon as the confinement is comparable to the plasmon wavelength. Conversely, the water fluctuations remain remarkably bulk-like, with significant confinement effects arising only when the wall spacing is reduced to 7 Å. We apply the confined response formalism to predict the dependence of the solid-water quantum friction and thermal boundary conductance on channel width for model channel wall materials. Our results provide a general framework for Coulomb interactions of fluctuating matter in nanoscale confinement.

# I. INTRODUCTION

Fluids confined at the nanometer scale underlie many technologically important processes [1, 2], including filtration, seawater desalination [3, 4], blue energy harvesting [5, 6] and electrochemical energy storage [7]. Yet, they started to be fundamentally investigated not more than 20 years ago, and their initial theoretical description was largely inherited from macroscopic hydrodynamics, with generic walls imposing the same boundary conditions regardless of their material composition [8]. The first nanofluidic effects emerged from the realization that, at the nanoscale, one may not neglect the wall's surface charge [9], which results in coupled ion-fluid transport phenomena such as electro-osmosis and streaming currents [10]. There has been, however, accumulating evidence in recent years that surface charge is not a sufficient descriptor for the nanofluidic solid-liquid interface. From fluids near conducting surfaces [11, 12] to strongly interacting ions due to dielectric contrast [13–15], several studies pointed to the need of describing the solid walls at the level of their electronic properties.

It may indeed be expected that, close enough to a solid wall, the Coulomb potentials produced by charged particles in a liquid are screened by the dielectric response of the wall material: this effect has been termed *interaction confinement* [15]. Charged particles in a liquid are, first and foremost, ions: interaction confinement produces effective Coulomb interactions between ions in nanochannels that are modified compared to bulk Coulomb inter-

actions, leading to a wealth of correlation effects [13, 14]. But a polar liquid such as water, even though electrically neutral after time-averaging, has a molecular-level charge structure: water thus exhibits thermal charge fluctuations at terahertz frequencies and on a wide range of length scales [16] (termed *hydrons* [17]). The corresponding Coulomb fields are also subject to interaction confinement: they are dynamically screened by the thermal and quantum fluctuations of the electrons in the solid wall [18, 19]. This solid-liquid coupling has been shown to result in a "quantum" contribution to hydrodynamic friction, and in direct near-field energy transfer between the liquid and the solid's electrons [19, 20]. These effects bridge the gap between fluid dynamics and condensed matter physics, opening the way to engineering nanoscale flows with the confining walls' electronic properties [17, 21].

Fluctuation-induced effects in nanofluidics have so far been studied at the level of a single planar interface. The relevant many-body electrostatics were conveniently described in terms of surface response functions: surface analogues of the dielectric function that had been widely used, for instance, in the field of plasmonics [22, 23]. Here, we introduce *confined response functions*, that generalize surface response functions to a 2D nanochannel geometry (Fig. 1a), providing a general tool for the treatment of Coulomb interactions in 2D confinement. As an illustration, we study the confined response of water, and of solid walls described as either graphene sheets or jellium slabs [24, 25]: this allows us to predict the confinement dependence of solid-water quantum friction and thermal boundary conductance.

This paper is organised as follows. In Sec. II, we introduce confined response functions and link them to eigen-

<footnote>* nikita.kavokine@mpip-mainz.mpg.de</footnote>
modes of the Coulomb potential in the 2D nanochannel geometry. In Sec. III, we compute the confined responses of specific media. We take the examples of a graphene sheet and a semi-infinite jellium for the solid; for the liquid, we study water in the framework of both force field and ab initio molecular dynamics (MD) simulations. In Sec. IV, we describe the influence of confinement on fluctuation-induced effects, specifically quantum friction and near-field heat transfer. To do so, we carry out the field-theory derivation of quantum friction directly in the confined geometry, which leads to natural emergence of the confined response functions. Finally, Sec. V establishes our conclusions.

**Units and conventions.** We set the Boltzmann constant $k_B = 1$ (that is, we express the temperature in energy units), but otherwise use SI units throughout the text. In real space, we use the cylindrical coordinates $\mathbf{r} = (\rho, z)$. The interfaces are at $z = 0$ for a single-interface and at $z = \pm h/2$ for a confined channel. We use Fourier transforms for both $\rho$ and the time but never for the $z$-direction. We use the following convention for the $d$-dimensional Fourier transform:

$$ \hat{F}(\mathbf{q}, \omega) = \int d^d\mathbf{r} dt F(\mathbf{r}, t) e^{-i\mathbf{q} \cdot \mathbf{r} + i\omega t}, $$

$$ F(\mathbf{r}, t) = \int \frac{d^d\mathbf{q} d\omega}{(2\pi)^{d+1}} \hat{F}(\mathbf{q}, \omega) e^{i\mathbf{q} \cdot \mathbf{r} - i\omega t}. $$

The charge densities are expressed in units of $e$ and the electrical potentials include an additional factor $e$. We denote $V(\mathbf{r}) = e^2/(4\pi\epsilon_0 r)$ the Coulomb potential which becomes $V(\mathbf{q}, z) = e^2/(2\epsilon_0 q)e^{-q|z|}$ in Fourier space.

## II. ELECTRIC RESPONSE OF INTERFACIAL SYSTEMS

### A. Single-interface: surface response function

We first briefly recall the widely-used concept of surface response function [22, 23]. Consider a semi-infinite medium occupying the half-space ($z < 0$). Given the electrostatic potential $\phi_{\text{ext}}$ applied by an external source (an appropriate charge distribution) inside the medium, we wish to determine the potential $\phi_{\text{ind}}$ induced by the medium in the half-space ($z > 0$). The potential $\phi_{\text{ext}}$ must solve the Laplace equation for $z < 0$. The physically-meaningful (non-diverging) solutions are given by the evanescent plane waves

$$ \phi_{\text{ext}}(\mathbf{q}, z, \omega) = \phi_{\text{ext}}(\mathbf{q}, z = 0, \omega) F^0(q, z), \\ F^0(q, z) = e^{-q|z|}, \quad (1) $$

where we have introduced the surface weight function $F^0$: this seemingly cumbersome notation will be useful upon generalization to a confined geometry. Assuming that the medium has a linear charge density response function $\chi$, the induced potential is given by (see Fig. 1b)

$$ \phi_{\text{ind}}(\mathbf{q}, z, \omega) = -g^0(\mathbf{q}, \omega)\phi_{\text{ext}}(\mathbf{q}, z = 0, \omega)F^0(q, z) \quad (2) $$

where we have introduced the surface response function

$$ g^0(\mathbf{q}, \omega) = -\frac{e^2}{2\epsilon_0 q} \int dz dz' F^0(q, z) \chi(\mathbf{q}, z, z', \omega) F^0(q, z'). \quad (3) $$

It is worth noting that $g^0(\mathbf{q}, \omega)$ is a scalar. Since we are considering a linearly responding medium, the response to an evanescent plane wave at $(\mathbf{q}, \omega)$ is an evanescent plane wave at $(\mathbf{q}, \omega)$, so that the induced potential has the same $z$-dependence as the external potential, given by the weight function $F^0$. In other words, the evanescent plane waves form an eigenbasis for the surface response.

Let us illustrate the role of the surface response function with a simple example. We consider a static point charge $e$ at a distance $z_0$ from the interface. This charge produces an "external" potential

$$ \phi_{\text{ext}}(\rho, z) = \int \frac{d\mathbf{q}}{(2\pi)^2} \frac{e^2}{2q} F^0(q, z - z_0) e^{i\rho \cdot \mathbf{q}} \quad (4) $$

where we have introduced a Fourier decomposition into evanescent plane waves. Accounting for the response of the medium, the total potential is

$$ \phi(\mathbf{q}, z) = \frac{e^2}{2q} [F^0(q, z - z_0) - g^0(\mathbf{q})F^0(q, z_0)F^0(q, z)] \quad (5) $$

The surface response function gives the contribution of the medium's polarization to the total potential.

### B. Double interface: confined response functions

We now generalize this approach to a two-dimensional nanochannel geometry. We consider two interfaces at $z = \pm h/2$ that define the channel, with the outside medium ($|z| > h/2$) being distinct from the inside medium ($|z| < h/2$). Later, we will specify that we consider water inside the channel, but we remain general at this point.

Contrary to the previous case, there is no longer a symmetry between the two media, and we need therefore to distinguish the responses of the inner and outer medium. All the applied and induced potentials still satisfy the Laplace equation, but, in both media, the subspace of harmonic functions with wavevector $\mathbf{q}$ and frequency $\omega$ is now of dimension 2: the response function at a given $(\mathbf{q}, \omega)$ is then in principle given by a $2 \times 2$ matrix. We shall, however, express the potentials in a basis of even (symmetric) and odd (antisymmetric) harmonic functions, where the matrix turns out to be diagonal. We call these basis functions *confined weight functions* and denote them $F_{i/o}^{s/a}$, where s/a stands for symmetric/antisymmetric function and i/o for inner/outer medium (see Table 1). The amplitude of the basis functions is in principle arbitrary: it is chosen so that the confined response functions defined in the following reduce to the conventional surface response functions in the non-confined case.
<figure><img src="image_2.png" /></figure>

<figure><img src="image_3.png" /><figcaption>FIG. 1. *Electric response of interfacial systems.* a) Schematic of the confined geometry under consideration, emphasizing the role of the channel walls' electronic degrees of freedom. b-f): Electric response in different geometries. b) Single-interface case: the response function is the usual surface response function. c-f) In confined geometry there are four situations to distinguish. The responding medium can be either in the outer space (c&d) or in the inner space (e&f). The external potential can be either symmetric (c&e) or antisymmetric (d&f). We observe that the inter-solid interactions bring corrections by lowering the induced potential in the symmetric case (c) and increasing the induced potential in the antisymmetric case (d).</figcaption></figure>

<table><caption>TABLE I. Weight functions used in the definitions of the surface and confined response functions (Eq. (3) and (8)).</caption><thead><tr><th>Model</th><th>Geometry</th><th>Weight function</th></tr></thead><tbody><tr><td>Single interface</td><td>Half-space</td><td>F<sup>0</sup>(q, z) = e<sup>-q|z|</sup></td></tr><tr><td rowspan="2">Confined</td><td>Symmetric</td><td>F<sup>s</sup><sub>i</sub>(q, z) = &radic;2 cosh(qz)e<sup>-q<sup>h</sup>/2</sup> F<sup>a</sup><sub>i</sub>(q, z) = &radic;2 sinh(qz)e<sup>-q<sup>h</sup>/2</sup></td></tr><tr><td>Antisymmetric</td><td>F<sup>s</sup><sub>o</sub>(q, z) = <sup>1</sup>&frasl;<sub>&radic;2</sub>e<sup>-q(|z|-<sup>h</sup>/2)</sup> F<sup>a</sup><sub>o</sub>(q, z) = <sup>sign(z)</sup>&frasl;<sub>&radic;2</sub>e<sup>-q(|z|-<sup>h</sup>/2)</sup></td></tr></tbody></table>
Let us start with the response of the inner medium. We consider a generic external potential $\phi_{\text{ext}}(\mathbf{q}, z, \omega) = \phi^s(\mathbf{q}, \omega)F_i^s(q, z) + \phi^a(\mathbf{q}, \omega)F_i^a(q, z)$ applied on the inner medium of charge susceptibility $\chi_i$. The induced potential is then

$$ \phi_{\text{ind}}(\mathbf{q}, z) = \frac{e^2}{2\epsilon_0 q} \int_{-h/2}^{h/2} dz' dz'' e^{-q|z-z'|} \dots \\ \dots \chi_i(\mathbf{q}, z', z'') \phi_{\text{ext}}(\mathbf{q}, z'') . \quad (6) $$

In the outer space $|z| > h/2$, taking advantage of the definition of the confined weight functions (see Table 1), this reduces to

$$ \phi_{\text{ind}}(\mathbf{q}, z) = -g_i^s \phi^s F_o^s(q, z) - g_i^a \phi^a F_o^a(q, z) \quad (7) $$

where we have defined a generic confined response function:

$$ g_m^c(\mathbf{q}, \omega) = -\frac{e^2}{2\epsilon_0 q} \int dz dz' F_m^c(q, z) \dots \\ \dots \chi_m(\mathbf{q}, z, z', \omega) F_m^c(q, z'). \quad (8) $$

where the interaction takes place over the domain of definition of the weight function $F_m^c$; $m = i, o$ and $c = s, a$. To summarise, the inner medium responds to a perturbation of the form $F_i^{s/a}$ by a potential of the form $F_o^{s/a}$ in the outer medium with an amplitude $-g_i^{s/a}$ (see Fig. 1c-d).

We proceed similarly for the response of the outer medium (with charge susceptibility $\chi_o$). It responds to an external potential $\phi_{\text{ext}}(\mathbf{q}, z, \omega) = \phi^s(\mathbf{q}, \omega)F_o^s(q, z) + \phi^a(\mathbf{q}, \omega)F_o^a(q, z)$ with an induced potential in the inner space $|z| < h/2$:

$$ \phi_{\text{ind}} = -g_o^s \phi^s F_i^s(z) - g_o^a \phi^a F_i^a(z). \quad (9) $$

The two components of the outer medium response are shown in Fig. 1c-f.

To summarize, we have generalized surface response functions to a 2D nanochannel geometry by identifying the new eigenmodes of the Coulomb potential, given by the confined weight functions. The response to a symmetric (antisymmetric) weight function is a symmetric (antisymmetric) weight function with amplitude given by the confined symmetric (antisymmetric) response function, in the same way that the response to an evanescent wave is an evanescent wave in the single interface case.

### C. Fluctuation-dissipation theorem and physical interpretation

To obtain a physical interpretation of the confined response functions, it is useful to relate them to equilibrium charge density fluctuations using the fluctuation-dissipation theorem (FDT). For the charge susceptibility $\chi$, the FDT reads [19]

$$ S_m(\mathbf{r}, \mathbf{r}', \omega) = f(\omega) \text{Im}[\chi_m(\mathbf{r}, \mathbf{r}', \omega)] \quad (10) $$

where $f(\omega) = 2T/\omega$ for classical dynamics and $f(\omega) = \hbar \cotanh(\hbar\omega/(2T))$ for quantum dynamics. The structure factor $S$ is defined as

$$ S_m(\mathbf{r}, \mathbf{r}', t) = \langle n_m(\mathbf{r}, t)n_m(\mathbf{r}', 0) \rangle, \quad (11) $$

where $n_m$ is the charge density of the medium. Adapted to the confined response function defined in Eq. (8), the FDT becomes

$$ S_m^c(\mathbf{q}, \omega) = f(\omega) \text{Im}[g_m^c(\mathbf{q}, \omega)] \quad (12) $$

where $S^c$ is a confined structure factor defined as

$$ S_m^c(\mathbf{q}, t) = \frac{1}{\mathcal{A}} \int d\mathbf{r} d\mathbf{r}' \langle n_m(\mathbf{r}, t) n_m(\mathbf{r}', 0) \rangle \dots \\ \dots e^{-i\mathbf{q} \cdot (\mathbf{r} - \mathbf{r}')} F_m^c(q, z) F_m^c(q, z') \quad (13) $$

where $\mathcal{A}$ is the area of the interface.

In Eq. (13), if the weight function $F_m^c$ is symmetric, the charge density $n_m(z)$ can be replaced with $(n_m(z) + n_m(-z))/2$: the structure factor only counts the symmetric charge fluctuations. Similarly, if the weight function is antisymmetric, the charge density $n_m(z)$ can be replaced with $(n_m(z) - n_m(-z))/2$: the structure factor only counts the antisymmetric charge fluctuations. Thus, for the channel walls, the symmetric (antisymmetric) response function accounts for in-phase (out-of-phase) coupled modes. For the inner medium, the symmetric (antisymmetric) response function accounts for monopolar (dipolar) charge fluctuations.

### D. Confined response function of the outer medium

We now provide an expression of the outer medium confined response function in terms of the usual surface response functions. We distinguish the contributions of the two solid walls to the potential induced in the inside medium: we denote them $\phi_{\text{ind}}^T$ and $\phi_{\text{ind}}^B$ for the "top" and "bottom" solids, respectively. We need to account for the solid-solid interactions: each solid responds to the external potential and to the potential induced by the other solid. Let us consider an external potential $\phi_{\text{ext}}(\mathbf{q}, z, \omega) = \phi_o(\mathbf{q}, \omega)F_o^c(q, z)$ where the weight function $F_o^c$ is either $F_o^s$ or $F_o^a$. The top solid induces a potential

$$ \phi_{\text{ind}}^{\text{T}}(z) = -g_o^0 F^0(z - h/2) \dots \\ \dots [\phi_o F_o^c(z = h/2) + \phi_{\text{ind}}^B(z = h/2)] \quad (14) $$

for $z < h/2$, while the bottom solid induces a potential

$$ \phi_{\text{ind}}^{\text{B}}(z) = -g_o^0 F^0(z + h/2) \dots \\ \dots [\phi_o F_o^c(z = -h/2) + \phi_{\text{ind}}^{\text{T}}(z = -h/2)] \quad (15) $$

for $z > -h/2$. Combining these two equations, we obtain the total induced potential as

$$ \phi_{\text{ind}}(z) = \phi_{\text{ind}}^{\text{T}}(z) + \phi_{\text{ind}}^{\text{B}}(z) = -g_o^c F_i^c(z) \phi^0 \quad (16) $$
in the inner space $|z| < h/2$. Comparing to Eq. (9), we deduce

$$g_o^s = \frac{g_o^0}{1 + g_o^0 e^{-qh}}, \quad g_o^a = \frac{g_o^0}{1 - g_o^0 e^{-qh}}. \quad (17)$$

We find that the confined response functions reduce to surface response functions at wavelengths $1/q \ll h$: the deviation from the surface response function originates from the Coulomb interaction between the two walls. As expected on physical grounds, the inter-wall interaction reduces the in-phase (symmetric) response and enhances the out-of-phase (antisymmetric) response. However, the induced potential does not exceed the applied potential in confinement if it does not for a single interface: there is no confinement-induced overscreening (see ESI.1).

### E. Interaction confinement

A first consequence of the wall electric response in a 2D channel is a modification of the effective Coulomb interactions between charged particles inside the channel: this effect has been termed *interaction confinement* [15]. Two confined charges (for example, ions in water) interact not only directly, but also indirectly, through the polarization charges induced in the channel walls. Coulomb interactions near polarizable walls have been computed in various geometries and within different models for the wall polarizability.

In ref. [15], interaction confinement in a 2D channel geometry was addressed in the framework of surface response functions, which allowed for evaluation of effective Coulomb interactions between two ions in the channel mid-plane, for various models of the channel wall material. For example, with walls described by a Thomas-Fermi model, the ion-ion interaction is reinforced with respect to bulk water at small distances, and exponentially screened at large distances.

Our confined response formalism generalizes the results of ref. [15], allowing for the evaluation of the effective Coulomb interactions between arbitrary charge distributions. Any charge distribution can indeed be decomposed into an even and an odd part. The wall response to the even (odd) part is then given by the symmetric (antisymmetric) response function. In ESI.2, we provide expressions in terms of the confined response functions for the effective Coulomb potential produced by prototypical even (point charge) and odd (dipole) charge distributions.

### III. CONFINED RESPONSE OF MODEL MEDIA

We now specialise to a 2D nanofluidic channel of height $h$ filled with water, and we investigate the effect of confinement on the electric response of water and of model channel wall materials.

### A. Outer medium: coupled plasmon modes

We consider two models for the channel wall material: a graphene monolayer and a semi-infinite jellium. The graphene surface response function is computed as detailed in ref. [19], starting from the charge susceptibility in the Dirac cone and zero-temperature approximations, at a chemical potential $\mu = 100$ meV. The semi-infinite jellium is treated in the specular reflection approximation, as detailed also in ref. [19]. In the jellium model, electrons are free in a uniform positive background, and are completely characterised by their chemical potential $\mu$ and effective mass $m$. We use $\mu = 180$ meV and $m = 0.1$ m$_e$ as a model for a doped semi-conductor: this corresponds to an electron density parameter $r_s = 5$.

For both systems, the surface response functions feature a sharply-defined surface plasmon mode and a broad particle-hole continuum (Fig. 2). The effect of confinement is most clearly visible on the plasmon mode: its energy is increases in the confined symmetric response and decreased in the confined asymmetric response. This is consistent with the physical interpretation outlined above: as the two solid walls face each other, in-phase (out-of-phase) charge density oscillations have an increased (decreased) energy cost due to the Coulomb interactions between the two solids. These interactions are significant only for charge fluctuations whose wavelength is longer than the confinement width $h$, so that the confined response functions differ from the surface response function only at small enough momenta $q$ (see Eq. (17)).

We anticipate that the formation of coupled plasmon modes between the walls of a 2D nanofluidic channel will affect transport inside the channel, and particularly fluctuation-induced effects (see Sec. IV).

### B. Inner medium: confined water spectra from simulations

We now turn to the confined response functions of water. Confinement may impact water charge fluctuations in two ways. First, the interaction between the two interfaces of the water slab is expected to result in a difference between the symmetric and antisymmetric responses. Second, the confinement-induced modifications of the water structure may intrinsically affect its charge fluctuations.

We determine the water response functions in the framework of molecular dynamics (MD) simulations. We carry out both force-field (FF) and ab intio density functional theory-based (DFT) simulations of water confined between two frozen graphene sheets, for various separations $h$ between the graphene sheets. From the simulation trajectories, we compute the charge structure factor of water, integrated along $z$ after multiplication by the weight functions summarized in Table 1, and determine the corresponding response functions through the fluctuation-dissipation theorem (Sec. II C). DFT simulations
<figure><img src="image_4.png" /><figcaption>FIG. 2. Surface and confined response functions for different solids. a-f) Interfacial response functions Im[g(q,ω)] of the solid as a function of momentum q and frequency ω. The first row (a-c) is for graphene with Fermi level µ = 100 meV. The second row (d-f) is for a jellium model with an effective mass m = 0.1 m_e and a Fermi level µ = 180 meV (electron density parameter r_s = 5), corresponding to a doped semi-conductor. The first column (a&d) corresponds to the surface response functions. The second column (b&e) corresponds to the symmetric confined response functions with a confinement of 7 Å. The third column (c&f) corresponds to the anti-symmetric confined response functions with a confinement of 7 Å.</figcaption></figure>

tions are required to capture intra-molecular modes: we use them to obtain the spectra at frequencies above 150 meV. At lower frequencies, we use the FF simulations to capture the contribution of inter-molecular modes, inac- cessible with the short simulation times of DFT. Details of the numerical parameters and procedures are given in ESI.3.

The results are presented in Fig. 3. Fig. 3b shows the water surface response function (obtained from the sim- ulation at weakest confinement) at a fixed wavevector q₀ = 0.67 Å⁻¹ and in the frequency range 0 – 200 meV. In this range, the water charge response can be decom- posed into four modes: in order of increasing frequency, the Debye mode, the hydrogen-bond stretching mode, the libration mode and the OH-bond bending mode [16, 26]. We note that the OH-stretch mode (at around 450 meV) falls outside the studied frequency range. The confine- ment dependence of these modes could in principle be analysed in terms of their molecular origin; this is, how- ever, beyond the scope of this article, and we restrict

Overall, the water response functions are remarkably

robust to confinement. Down to *h* = 1.4 nm, the lower- frequency intermolecular modes remain unaffected. We observe, however, a slight red shift, and an increased os- cillator strength for the bending mode. A significant effect on the intermolecular modes is visible only at 7 Å confinement. In the symmetric response, the Debye and hydrogen-bond-stretch modes are amplified, while the libration mode is suppressed; in the antisymmetric response, the libration peak is strongly amplified while the other modes are unaffected. As illustrated in Figs. 3e-f, in the symmetric case, the perpendicular component of the applied electric field changes direction across the channel, while it maintains a constant sign in the anti- symmetric case. This likely indicates that the libration mode is mostly excited by the electric field perpendicu- lar to the interface, while the lower frequency modes are excited by the parallel component.
<figure><img src="image_5.png" /><figcaption>FIG. 3. Water response functions from simulation. a) Snapshot of a FF molecular dynamics simulation with 7 Å confinement b-d) Response functions Im[g(q, ω)] of water as a function of the frequency ω obtained from the simulations, at fixed wavevector q = 0.67 Å⁻¹ for FF and q = 1 Å⁻¹ for DFT. Different confinements are used: 7 Å, 14 Å (only FF), 18 Å (only DFT), 34 Å and 60 Å (only FF). b) Surface response function computed from simulations at weak confinement (FF at h =60 Å and DFT at h =34 Å). c) Symmetric response function for different confinements. The dashed black line is the non-confined surface response function. d) Antisymmetric response function for different confinements. The dashed black line is the non-confined surface response function. e) Schematic of a channel containing a single water monolayer on which a symmetric potential is applied. f) Same as (e) with application of an antisymmetric potential.</figcaption></figure>

## IV. CONFINED QUANTUM FRICTION AND HEAT TRANSFER

In this section, we discuss the effect of 2D confinement on fluctuation-induced interfacial effects – solid-liquid quantum friction and near-field radiative heat transfer – making use of the confined response function formalism developed above.

### A. Single-interface quantum friction

We start by briefly summarising the physics of quantum hydrodynamic friction. The classical friction between a liquid and a solid is usually determined by the solid's surface roughness [27, 28]. However, it was recently shown that the classical contribution is supplemented by a fluctuation-induced or "quantum" contribution, due to the coupling of water charged fluctuations (termed *hydrons*) to electronic excitations within the solid [19]. When undergoing quantum friction, a liq-
<figure><img src="image_6.png" /><figcaption>FIG. 4. Fluctuation-induced effects. (a) Schematic of the interactions involved in single-interface quantum friction. Friction results from the Coulomb coupling of liquid fluctuations (hydrons) and electronic fluctuations. (b) Schematic of the interactions involved in confined quantum friction. Compared to the non-confined case, there is an additional Coulomb interaction between the charge fluctuations in the two solid walls. (c) Dyson equation for the renormalization of the bottom wall's charge density response function (polarization "bubble") by the intra-wall Coulomb interactions. The bare bubble is empty and the renormalized bubble is filled with gray. (d) Dyson equation for the electron-water susceptibility $\chi_{ew}$. (e) Definition of the exchange term $\Pi$ (see text). (f) Dyson equation for the renormalization of the bottom wall's charge density response function by the inter-wall Coulomb interactions. The stripes indicate a fully-renormalized bubble. (g) Dyson equation for the inter-wall charge density response function, which vanishes in the absence of inter-wall Coulomb interactions ($\Pi = 0$).</figcaption></figure>

uid transfers momentum directly to the solid's electrons. Similarly, a liquid may transfer energy directly to the solid's electrons: this is near-field radiative heat transfer [20, 29].

In the case of a single solid-liquid interface, if the liquid is flowing at velocity $\mathbf{v}$, it is subject to a quantum friction force $\mathbf{F} = -\lambda A \mathbf{v}$ with the quantum friction coefficient $\lambda$ given by [19]

$$ \lambda = \frac{\hbar^2}{8\pi^2 T} \int_0^\infty d\omega dq \frac{q^3}{\sinh^2\left(\frac{\hbar\omega}{2T}\right)} \Delta\gamma[F^0]. \quad (18) $$

Anticipating the generalization to the confined case, we

have introduced the notation

$$ \Delta\gamma[F] = \frac{\mathrm{Im}[g_e(q, \omega)] \mathrm{Im}[g_w(q, \omega)]}{|1 - g_e(q, \omega)g_w(q, \omega)|^2}, \quad (19) $$

where the $g$'s ($e$ corresponds to the solid, $w$ to the liquid) are generalized response functions computed with the weight function $F$ (see Table 1). Similarly, if there is a temperature difference $\Delta T$ between the solid and the liquid, the solid-liquid heat flux is given by $Q = \kappa A \Delta T$, where the thermal boundary conductance $\kappa$ reads [20]

$$ \kappa = \frac{\hbar^2}{4\pi^2 T^2} \int_0^\infty d\omega dq \frac{q\omega^2}{\sinh^2\left(\frac{\hbar\omega}{2T}\right)} \Delta\gamma[F^0]. \quad (20) $$

## B. Confined quantum friction

We now generalize the above results to the 2D confined geometry presented in Fig. 4b. We wish to compute the total quantum friction force applied by the solid walls on the flowing liquid, and the total heat transfer rate between the liquid and the two solid walls. The liquid and the solid's electrons are described by their fluctuating charge densities $n_w$ and $n_e$, which have a Coulomb interaction of the form

$$ \mathcal{H}_{\text{h-e}}(t) = \int d\mathbf{r}_w d\mathbf{r}_e n_w(\mathbf{r}_w, t) V(\mathbf{r}_w - \mathbf{r}_e) n_e(\mathbf{r}_e, t). \quad (21) $$

The dynamics of the systems are governed by the interaction Hamiltonian that also comprises the electron-electron Coulomb interactions: $\mathcal{H}_{\text{int}} = \mathcal{H}_{\text{e-e}} + \mathcal{H}_{\text{h-e}}$.

The friction force and heat transfer rate are given by

$$ \langle \mathbf{F} \rangle = - \int d\mathbf{r}_w d\mathbf{r}_e \nabla V(\mathbf{r}_w - \mathbf{r}_e) \langle n_w(\mathbf{r}_w - \mathbf{v}t, t) n_e(\mathbf{r}_e, t) \rangle, \quad (22) $$

$$ \langle Q \rangle = \int d\mathbf{r}_w d\mathbf{r}_e V(\mathbf{r}_w - \mathbf{r}_e) \partial_t \langle n_w(\mathbf{r}_w, t) n_e(\mathbf{r}_e, t) \rangle, \quad (23) $$

where the integration over $\mathbf{r}_w$ (resp. $\mathbf{r}_e$) runs over the space occupied by the liquid (resp. solid). The correlation functions appearing in Eqs. (22) and (23) may be computed in perturbation theory with respect to $\mathcal{H}_{\text{int}}$, as has been detailed in ref. [19]. Since the system is subject either to liquid flow or to a temperature gradient, the perturbative expansion needs to be carried out in the non-equilibrium Keldysh formalism [30]. Ultimately, both the friction force and the heat transfer rate can be obtained in terms of the solid-liquid charge density correlation function $\chi_{ew}$, which, upon resummation of the perturbation series, is found to satisfy the following Dyson equation:

$$ \chi_{ew} = \chi_e * \chi_w + \chi_e * \chi_w * \chi_{ew}. \quad (24) $$

Here * stands for convolution in space and time, multiplication by the Coulomb potential, and contraction of the
Keldysh indices. A Feynman diagram representation of Eq. (24) is given in Fig. 4d. The Keldysh indices carried by the $\chi$'s are not important for the geometrical discussion that follows. We will therefore not write them out explicitly and we refer the reader to ref. [19] for further details: once the space-time convolutions have been dealt with, the computation is completely analogous to ref. [19].

In the single interface case, upon Fourier transformation in time and in space parallel to the interface, Eq. (24) immediately becomes a scalar equation for surface response functions. In the channel geometry, we need to introduce confined response functions, first as $2 \times 2$ matrices in the indices $\nu, \xi = \text{T,B}$:

$$ (g_{ab})^{\nu\xi}(\mathbf{q}, \omega) = -\frac{e^2}{2\epsilon_0 q} \int_{\mathcal{Z}_a^{\nu}} \int_{\mathcal{Z}_b^{\xi}} dz_a dz_b \chi_{ab}(z_a, z_b, \mathbf{q}, \omega) \dots \\ \dots e^{q[\epsilon^a(\sigma^\nu z_a - h/2) + \epsilon^b(\sigma^\xi z_b - h/2)]} \quad (25) $$

where a,b=e,w, $\epsilon^w = +$, $\epsilon^e = -$, $\sigma^T = +$ and $\sigma^B = -$. The space has been divided into three regions: the central region $\mathcal{Z}_w$, and the bottom (B) and top (T) wall regions: $\mathcal{Z}_e = \mathcal{Z}_e^B \cup \mathcal{Z}_e^T$. In the convolution over $z$ in Eq. (24), summing over the top and bottom solids corresponds to summing over the indices B,T. Thus, in terms of the confined response functions, Eq. (24) becomes

$$ g_{ew} = -g_e \cdot g_w + g_e \cdot g_w \cdot g_{ew}, \quad (26) $$

where the dot represents the matrix product, and $g_a \equiv g_{aa}$. Neglecting the off-diagonal terms in Eq. (26), we would recover two copies of the single interface Dyson equation, one for the top and one for the bottom interface. The off-diagonal terms represent cross-talk between the walls (for example, the top solid wall responding to a water fluctuation near the bottom wall), which is expected to vanish for weak confinement.

We assume in the following that the top and bottom wall materials are the same, so that the system is symmetric under the mirror transformation $z \to -z$. As a consequence, $g^{TT} = g^{BB}$ and $g^{TB} = g^{BT}$. Therefore, all the confined response matrices can be diagonalised in the form

$$ P g P^{-1} = \begin{pmatrix} g^{TT} + g^{TB} & 0 \\ 0 & g^{TT} - g^{TB} \end{pmatrix} = \begin{pmatrix} g^s & 0 \\ 0 & g^a \end{pmatrix} \quad (27) $$

where we have used

$$ P = \frac{1}{\sqrt{2}} \begin{pmatrix} 1 & 1 \\ 1 & -1 \end{pmatrix}, \quad (28) $$

which satisfies $P^{-1} = P$. We thus recover the definition of the confined response function in its eigenbasis, as introduced in Sec. 2. Multiplying Eq. (26) by $P$ on the left and on the right, we obtain two scalar equations for the symmetric and antisymmetric response functions:

$$ \left\{ \begin{aligned} g_{ew}^{s}(\mathbf{q}, \omega) &= -g_{e}^{s} g_{w}^{s} + g_{e}^{s} g_{w}^{s} g_{ew}^{s} \\ g_{ew}^{a}(\mathbf{q}, \omega) &= -g_{a}^{s} g_{w}^{a} + g_{e}^{a} g_{w}^{a} g_{ew}^{a} \end{aligned} \right. \quad (29) $$

Once the Dyson equation is reduced to a scalar equation, we may follow the steps of refs. [19, 20] to obtain the friction coefficient and thermal boundary conductance in the confined geometry:

$$ \lambda = \frac{\hbar^2}{8\pi^2 T} \int_0^\infty d\omega dq \frac{q^3}{\sinh^2\left(\frac{\hbar\omega}{2T}\right)} (\Delta\gamma[F^s] + \Delta\gamma[F^a]), \quad (30) $$

$$ \kappa = \frac{\hbar^2}{4\pi^2 T^2} \int_0^\infty d\omega dq \frac{q\omega^2}{\sinh^2\left(\frac{\hbar\omega}{2T}\right)} (\Delta\gamma[F^s] + \Delta\gamma[F^a]). \quad (31) $$

The confined response functions thus emerge naturally in the theory of fluctuation-induced effects in a 2D channel.

### C. Diagrammatic approach to confined response function

Using the diagrammatic approach developed for the fluctuation-induced effects, we may interpret the confined response functions of the solid walls as surface response functions that have been renormalised by the inter-wall interactions within the random phase approximation (RPA).

We start from the intra-wall response function $\chi_e(z, z')$, that has been renormalized by the intra-wall Coulomb interactions at the RPA level, according to the Dyson equation shown diagrammatically in Fig. 4c. These are identical in the top and bottom walls ($\chi_e^{BB} = \chi_e^{TT}$), and, at this stage, there is no inter-wall response $\chi_e^{BT}$, since we do not allow for electron tunneling between the walls.

In the presence of Coulomb interactions between the walls, we may introduce the exchange term $\Pi = \chi_e^{BB} * \chi_e^{TT}$ (Fig. 4e). Still at the RPA level, the intra-wall response function is then renormalised according to (Fig. 4f)

$$ \tilde{\chi}_e^{BB} = \chi_e^{BB} + \Pi * \tilde{\chi}_e^{BB}. \quad (32) $$

The inter-wall response is no longer vanishing, and satisfies (Fig. 4g)

$$ \tilde{\chi}_e^{BT} = \Pi + \Pi * \tilde{\chi}_e^{BT}. \quad (33) $$

Fourier-transforming these Dyson equations as detailed above, we obtain relations between confined and surface response functions. Using in particular that

$$ \frac{e^2}{2\epsilon_0 q} \int_{\mathcal{Z}_B^B} dz \int_{\mathcal{Z}_T^T} dz' \Pi(z, z') e^{-q(z-z'-h)} = (g_e^0)^2 e^{-qh}, \quad (34) $$

we obtain

$$ g_e^{BB} = \frac{g_e^0}{1 - (g_e^0)^2 e^{-2qh}}, \quad g_e^{BT} = - \frac{(g_e^0)^2 e^{-qh}}{1 - (g_e^0)^2 e^{-2qh}}. \quad (35) $$
<figure><img src="image_7.png" /><figcaption>FIG. 5. Effect of confinement on the fluctuation-induced phenomena. a) Quantum friction coefficient λ normalised by the single interface quantum friction coefficient λ^0 as a function of confinement between graphene walls. b) Thermal boundary conductance κ normalised by the single interface thermal boundary conductance κ^0 as a function of confinement between graphene walls. c) Quantum friction coefficient λ normalised by the single interface quantum friction coefficient λ^0 as a function of confinement between jellium walls with an effective mass m = 0.1 m_e and a Fermi level µ = 180 meV (electron density parameter r_s = 5) d) Thermal boundary conductance κ normalised by the single interface thermal boundary conductance κ^0 as a function of confinement between jellium walls with an effective mass m = 0.1 m_e and a Fermi level µ = 180 meV (electron density parameter r_s = 5). In all panels, the continuous lines are obtained using the single interface surface response function for water, while the crosses are obtained with the confined response function of water at the relevant confinement, as obtained from molecular dynamics simulations.</figcaption></figure>

and deduce the symmetric and antisymmetric components of the confined response function

$$g_e^s = \frac{g_e^0}{1 + g_e^0 e^{-q\hbar}}, \quad g_e^a = \frac{g_e^0}{1 - g_e^0 e^{-q\hbar}}. \qquad (36)$$

We thus recover Eq. (17), that we previously obtained from purely electrodynamic considerations.

### D. Effect of confinement on the fluctuation-induced effects

We now investigate the effect of confinement on the quantum friction coefficient and thermal boundary conductance of water in a 2D channel, with walls made of either graphene or a semi-infinite jellium, with the parameters detailed in Sec. III. The results are presented in Fig. 5.

We first focus on the effect of inter-solid interactions and thus evaluate the fluctuation-induced effects using the single-interface surface response function for water
(continuous lines in Fig. 5). Interestingly, we observe opposite trends for graphene and for jellium walls. In the case of graphene, for both friction and thermal conductance, the antisymmetric contribution is enhanced and the symmetric contribution is reduced with confinement. Indeed, both effects are governed by the coupled plasmon modes of the walls, and the out-of-phase mode has lower energy than the in-phase mode, thus making a larger contribution. For our jellium model, the plasmon is well above the thermal energy (around 300 meV), and the fluctuation-induced effects are governed by single-particle excitations: we find that, in this case, the confinement enhances the symmetric contribution and reduces the antisymmetric contribution. The antisymmetric contribution dominates the behaviour of the total friction and thermal conductance, but the overall confinement effect remains lower than 10%, except for the thermal conductance with graphene walls, where it reaches 50%. In general, the confinement effect is stronger for graphene walls than for jellium walls because the electronic fluctuations that mediate quantum friction and near-field heat transfer have a longer wavelength (smaller momentum $q$) in the case of graphene.

We now turn to the effect of the confinement-induced changes in the water fluctuations. Our simulations have shown that these changes become significant only at 7 Å confinement (Fig. 3), with an amplification of the Debye peak in the symmetric response and of the libration peak in the antisymmetric response. This translates into an enhancement of the friction coefficient and thermal conductance for both solid models, by up to a factor of 2 in the case of graphene. The study of fluctuation-induced effects specifically in 7 Å confinement is thus of particular interest. For instance, the thermal conductance of the interface between graphene and nanoconfined water may be probed with optical pump - terahertz probe spectroscopy: the optically-excited graphene electrons would be expected to cool faster than in the non-confined case [20].

V. CONCLUSIONS

In this paper, we have developed a theoretical framework for studying confinement- and fluctuation-induced effects in two-dimensional nanofluidic channels: effects of (fluctuating) Coulomb interactions between the liquid and the solid. The key element of our framework is the description of the response of the solid walls and of the confined liquid to the Coulomb potentials that they apply to each other. In the case of a single solid-liquid interface, the surface response function – the reflection coefficient for evanescent plane waves – was found to be the most convenient descriptor. This convenience was due to the evanescent plane waves being eigenmodes of the Coulomb potential: the response to an evanescent wave is an evanescent wave. Generalizing this idea to the 2D channel geometry, we have introduced *confined response*

functions, that play the role of reflection coefficients for the potential eigenmodes of the confined system. Our approach is systematic, and potentially extendable to more complex geometries.

The confined response functions reveal electrodynamic cross-talk between the walls of a 2D nanochannel. Investigating model materials that exhibit a surface plasmon mode, we found that the plasmons of the two walls couple as soon as the confinement is comparable to the plasmon wavelength. While the coupling of collective modes through Coulomb interactions is in principle a well-known phenomenon [31], our framework allows for the investigation of its effect on nanoscale fluid transport. From the fluid side, we have investigated confined water through molecular dynamics simulations. We found that the water confined response remains essentially bulk-like in the thermal frequency range down to 1.4 nm confinement, but undergoes significant changes when the confinement reaches 7 Å.

As an application of our framework, we have investigated quantum friction and near-field radiative heat transfer between water and the walls of a 2D nanochannel. We have generalized the derivation of refs. [19, 20] to a confined geometry and found that the confined response functions naturally emerge. In channels wider than 7 Å, the friction coefficient and thermal boundary conductance are modified compared to their bulk values when the confinement is comparable to the typical wavelength of the relevant charge fluctuations. At 7 Å confinement, a significant enhancement occurs for both effects, due to the drastic modification of the water response. Observing a confinement-induced modification of quantum friction or heat transfer thus appears most promising for systems where charge fluctuations are on longer wavelengths, yet these are also the systems where the effects are the weakest. Nevertheless, it has been shown that graphene-water heat transfer can be measured with ultrasfast spectroscopy [20], and a small quantum friction coefficient can still result in a large quantum friction force (termed quantum feedback) if the electrons are driven by a phonon wind [17]. Our results thus provide guidelines for engineering fluctuation-induced effects in nanoscale fluid transport.

AUTHOR CONTRIBUTIONS

N.K. and L.B. designed the project. B.C. carried out the theoretical analysis. M.B. performed the molecular dynamics simulations, with supervision from R.N. B.C. and N.K. co-wrote the paper, with input from all authors. All authors discussed the results and commented on the manuscript.

CONFLICTS OF INTEREST

There are no conflicts to declare.
ACKNOWLEDGEMENTS

The Flatiron Institute is a division of the Simons Foun-
dation. B.C. acknowledges funding from a J.-P. Aguilar
grant of the CFM Foundation, and thanks the Flat-
iron Institute and the Max Planck Institute for Poly-
mer Research for their hospitality. M.B and R.R.N ac-

knowledge funding by the Deutsche Forschungsgemein-
schaft (DFG, German Research Foundation) - Project-
ID 387284271 - SFB 1349. L.B. acknowledges funding
from the EU H2020 Framework Programme/ERC Syn-
ergy Grant agreement number 101071937 n-AQUA. N.K.
acknowledges support from a Humboldt fellowship. We
thank Lucy Reading-Ikkanda (Simons Foundation) for
help with figure preparation.

[1] S. Faucher, N. Aluru, M. Z. Bazant, D. Blankschtein, A. H. Brozena, J. Cumings, J. Pedro de Souza, M. Elimelech, R. Epsztein, J. T. Fourkas, A. G. Rajan, H. J. Kulik, A. Levy, A. Majumdar, C. Martin, M. McEldrew, R. P. Misra, A. Noy, T. A. Pham, M. Reed, E. Schwegler, Z. Siwy, Y. Wang and M. Strano, *The Journal of Physical Chemistry C*, 2019, **123**, 21309–21326.

[2] L. Bocquet, *Nature materials*, 2020, **19**, 254–256.

[3] W. A. Phillip and M. Elimelech, *Science*, 2011, **333**, 712–717.

[4] J. R. Werber, A. Deshmukh and M. Elimelech, *Environmental Science and Technology Letters*, 2016, **3**, 112–120.

[5] B. E. Logan and M. Elimelech, *Nature*, 2012, **488**, 313–319.

[6] A. Siria, M.-L. Bocquet and L. Bocquet, *Nature Reviews Chemistry*, 2017, **1**, 0091.

[7] J. Chmiola, G. Yushin, Y. Gogotsi, C. Portet, P. Simon and P. L. Taberna, *Science*, 2006, **313**, 1760–1763.

[8] N. Kavokine, R. R. Netz and L. Bocquet, *Annual Review of Fluid Mechanics*, 2021, **53**, 377–410.

[9] R. B. Schoch and P. Renaud, *Applied Physics Letters*, 2005, **86**, 1–3.

[10] S. Marbach and L. Bocquet, *Chem. Soc. Rev.*, 2019, **48**, 3102–3144.

[11] A. Marin-Laflèche, M. Haefele, L. Scalfi, A. Coretti, T. Dufils, G. Jeanmairet, S. Reed, A. Serva, R. Berthin, C. Bacon, S. Bonella, B. Rotenberg, P. Madden and M. Salanne, *Journal of Open Source Software*, 2020, **5**, 2373.

[12] A. Schlaich, D. Jin, L. Bocquet and B. Coasne, *Nature Materials*, 2022, **21**, 237–245.

[13] N. Kavokine, S. Marbach, A. Siria and L. Bocquet, *Nat. Nanotechnol.*, 2019, **14**, 573–578.

[14] P. Robin, N. Kavokine and L. Bocquet, *Science*, 2021, **373**, 687–691.

[15] N. Kavokine, P. Robin and L. Bocquet, *The Journal of Chemical Physics*, 2022, **157**, 114703.

[16] S. Carlson, F. N. Brünig, P. Loche, D. J. Bonthuis and R. R. Netz, *The Journal of Physical Chemistry A*, 2020, **124**, 5599–5605.

124, 5599-5605.

[17] B. Coquinot, L. Bocquet and N. Kavokine, *Phys. Rev. X*, 2023, **13**, 011019.

[18] A. I. Volokitin and B. N. J. Persson, *Reviews of Modern Physics*, 2007, **79**, 1291–1329.

[19] N. Kavokine, M.-L. Bocquet and L. Bocquet, *Nature*, 2022, **602**, 84–90.

[20] X. Yu, A. Principi, K.-J. Tielrooij, M. Bonn and N. Kavokine, *Nat. Nanotechnol. (accepted)*, 2023.

[21] M. Lizée, A. Marcotte, B. Coquinot, N. Kavokine, K. Sobnath, C. Barraud, A. Bhardwaj, B. Radha, A. Niguès, L. Bocquet and A. Siria, *Phys. Rev. X*, 2023, **13**, 011020.

[22] A. Liebsch, in *Electronic excitations at metal surfaces*, Springer, 1995, ch. 3, pp. 49–105.

[23] J. M. Pitarke, V. M. Silkin, E. V. Chulkov and P. M. Echenique, *Reports on Progress in Physics*, 2007, **70**, 1–87.

[24] N. D. Lang and W. Kohn, *Physical Review B*, 1970, **1**, 4555–4568.

[25] A. Griffin and J. Harris, *Canadian Journal of Physics*, 1976, **54**, 1396–1408.

[26] D. J. Bonthuis, S. Gekle and R. R. Netz, *Langmuir*, 2012, **28**, 7679–7694.

[27] J.-L. L. Barrat and L. L. Bocquet, *Faraday Discussions*, 1999, **112**, 119–127.

[28] D. M. Huang, C. Sendner, D. Horinek, R. R. Netz and L. Bocquet, *Physical Review Letters*, 2008, **101**, 1–4.

[29] S. A. Biehs, R. Messina, P. S. Venkataram, A. W. Rodriguez, J. C. Cuevas and P. Ben-Abdallah, *Reviews of Modern Physics*, 2021, **93**, 025009.

[30] J. Rammer, *Quantum Field Theory of Non-equilibrium States*, Cambridge University Press, 2007.

[31] E. H. Hwang, R. Sensarma and S. D. Sarma, *Physical Review B*, 2010, **82**, 195406.