Friction of a driven chain: role of momentum conservation, Goldstone and radiation modes

We consider a 1D harmonic chain as shown in figure 1(a), where neighboring beads are coupled via a potential $V = \frac{K}{2}\sum_n (u_n-u_{n-1})^2$, with a coupling strength K and u_n the displacement of bead n from the ground state position. With M the mass of a bead, Newton's equation of motion for bead n is [41]

$$\begin{equation} M \ddot{u}_n = K\left(u_{n-1} + u_{n+1} - 2u_n\right) . \end{equation} \tag{ 1 }$$

The displacements $\{u_n\}$ can be Fourier transformed in space and time by introducing a lattice spacing a, yielding

$$\begin{align} \tilde{u}_q\left( \omega\right) = \int_{-\infty}^{\infty} \sum_{n = 0}^{N} e^{i\left(q n a - \omega t\right)} u_n\left(t\right) dt. \end{align} \tag{ 2 }$$

Equation (1) thus yields the familiar dispersion relation [41, 42] of the infinite length system

$$\begin{equation} \omega = 2 \Omega_0 \text{sin}\left(\frac{a q}{2}\right) \, , \quad \Omega_0 \equiv \sqrt{\frac K M} \,. \end{equation} \tag{ 3 }$$

We use Ω₀ later as the main unit for frequency.

Zoom In Zoom Out Reset image size

Next we introduce microscopic damping and consider two paradigmatic types: Langevin and DPD damping, as shown in figures 1(b) and (c).

2.2.1. Langevin damping.

In presence of Langevin damping (see figure 1(b)), bead n feels a damping force proportional to its velocity $\dot u_n$, quantified by the damping coefficient γ_L , yielding

$$\begin{equation} M \ddot{u}_n = K\left(u_{n-1} + u_{n+1} - 2u_n\right) - M \gamma_L \dot{u}_n . \end{equation} \tag{ 4 }$$

Note that γ_L carries units of frequency, which allows for easy comparison to inertial contributions. Using the Fourier transformed displacements in equation (2), the equation of motion can be written as

$$\begin{equation} \omega^2 M \tilde{u}_q\left(\omega\right) = K\left(e^{iq_L a} + e^{-iq_L a } - 2\right) \tilde{u}_q\left(\omega\right) - i \omega M \gamma_L \tilde{u}_q\left(\omega\right). \\ \end{equation} \tag{ 5 }$$

Using equation (5), the dispersion relation is then obtained to be,

$$\begin{eqnarray} &&4 \text{sin}^2\left(\frac{q_La}{2}\right) = \tilde{\omega}^2 + i \tilde{\omega}\tilde{\gamma}_L, \end{eqnarray} \tag{ 6 }$$

where we note that comparing ω and γ_L yields over- and underdamped cases, as will be discussed below. We introduced dimensionless frequency and damping, $\tilde \omega\equiv\frac{\omega}{\Omega_0}$ and $\tilde \gamma_L\equiv\frac{\gamma_L}{\Omega_0}$. For later reference, we explicitly give the complex wavenumber q_L by solving equation (6) for it,

$$\begin{equation} a q_L \equiv\tilde q_L = 2 \text{arcsin} \left( \frac{1}{2} \sqrt{\tilde \omega ^2+i \tilde \gamma_L\tilde \omega }\right). \end{equation} \tag{ 7 }$$

The displacement profile is $u_n(\omega) = C_1 (\omega) e^{i q_L(\omega) na} + C_2 (\omega) e^{-i q_L(\omega) na}$, with the constants C₁ and C₂ determined from the boundary conditions introduced in section 2.3 below.

The real and imaginary parts of the wave number in equation (7) correspond to the wavelength of modes via $\frac{\lambda_L}{a}\equiv \tilde \lambda_L = \frac{2\pi}{\textrm{Re} \tilde q_L}$, and the decay length $\frac{l_L}{a}\equiv\tilde l_L = \frac{1}{\textrm{Im} \tilde q_L}$.

Notably, this model is mathematically similar to the Rouse model for a polymer chain [37]. We will comment on the connections after equations (42) and (50) below.

2.2.2. DPD damping.

In presence of DPD damping (see figure 1(c)), the damping force acting on bead n is proportional to its relative velocity with respect to its neighbors. Considering only nearest neighbor damping, the DPD equation of motion is [5, 36],

$$\begin{equation} \begin{aligned} M \ddot{u}_n = \left(K+ M \gamma_D\frac{\partial}{\partial t}\right)\left(u_{n-1} + u_{n+1} - 2u_n\right). \end{aligned} \end{equation} \tag{ 8 }$$

In Fourier space, the equations of motion can be rewritten as

$$\begin{eqnarray} &&\omega^2 M \tilde{u}_q\left(\omega\right) = \left(K + i M \gamma_D\right)\left(e^{iq_L a} + e^{-iq_L a}- 2\right) \tilde{u}_q\left(\omega\right). \end{eqnarray} \tag{ 9 }$$

The damping coefficient γ_D carries, as γ_L , units of frequency, and we use $\tilde{\gamma}_D \equiv \frac{\gamma_D}{\Omega_0}$, dimensionless. The dispersion relation for the case of DPD damping extracted using equation (9) is,

$$\begin{equation} \text{sin}^2\left(\frac{\tilde q_D}{2}\right) = \frac{\tilde{\omega}^2} {4\left(1 - i \tilde{\omega} \tilde{\gamma}_D\right)} , \end{equation} \tag{ 10 }$$

which shows that, comparing $\tilde{\omega}\tilde{\gamma}_D$ to unity distinguishes over- and underdamped limits. Solving explicitly for $q_D(\omega)$,

$$\begin{equation} aq_D\equiv\tilde q_D = 2 \text{arcsin} \left( \frac{1}{2} \sqrt{\frac{\tilde{\omega}^2 }{\left(1 - i \tilde{\omega} \tilde{\gamma}_D\right)}}\right). \end{equation} \tag{ 11 }$$

Naturally, also here, the real and imaginary parts of the wave number in equation (11) correspond to the wavelength and decay length, $\frac{\lambda_D}{a}\equiv \tilde \lambda_D = \frac{2\pi}{\textrm{Re} \tilde q_D}$ and $\frac{l_D}{a}\equiv\tilde l_D = \frac{1}{\textrm{Im} \tilde q_D}$.

As mentioned in the introduction, the continuum version of equation (8) is the Kelvin Voigt Model [36]. In the absence of elastic forces, this is equivalent to the (Navier-)Stokes equation for a fluid [38]. We will comment on the connections to shearing a fluid after equations (25) and (35) below.

We introduce so called strain-controlled driving, where the displacement of the bead n = 0 is prescribed to be

$$\begin{align} u_0\left(t\right) = U_0 \text{cos}\left(\omega_0 t\right), \end{align} \tag{ 12 }$$

with ω₀ the driving frequency and U₀ the strain amplitude. The other end of the chain is held fixed, i.e. anchored so that

$$\begin{eqnarray} u_N\left(t\right) & = & 0. \end{eqnarray} \tag{ 13 }$$

For the monochromatic driving of equation (12), the displacement profile is monochromatic as well,

$$\begin{align} u_n\left(t\right)& = \textrm{Re}\left\{u_n\left(\omega_0\right) e^{i \omega_0 t}\right\}, \end{align} \tag{ 14 }$$

which introduces the Fourier amplitude $u_n(\omega_0)$ of bead n with $0\unicode{x2A7D} n\unicode{x2A7D} N$, $u_n(\tilde \omega_0)/U_0\equiv\tilde u_n(\tilde \omega_0)$. For the given boundary conditions, it is given by

$$\begin{align} \tilde u_n\left(\tilde\omega_0\right) = \frac{\text{sin}\left(\left(N-n\right)\tilde q\right)}{\text{sin}\left(N\tilde q\right)}. \end{align} \tag{ 15 }$$

For the Langevin case, this is explicitly,

$$\begin{align} \tilde u_n \left(\tilde\omega_0\right)& = \frac {\text{sin} \left(2 \left(N-n \right) \text{arcsin} \left(\frac{1}{2} \sqrt{\tilde{\omega}_0 \left(\tilde{\omega}_0 +i \tilde{\gamma}_L \right)}\right)\right)}{\text{sin} \left(2 N \text{arcsin} \left(\frac{1}{2} \sqrt{\tilde{\omega}_0 \left(\tilde{\omega}_0 +i \tilde{\gamma}_L \right)}\right)\right)}, \end{align} \tag{ 16 }$$

and for DPD,

$$\begin{align} \tilde u_n\left(\tilde\omega_0\right) = & \frac{\text{sin} \left(2 \left(N-n\right) \text{arcsin} \left(\frac{\tilde{\omega}_0 }{2 \sqrt{1-i \tilde{\gamma}_D \tilde{\omega}_0 }}\right)\right)}{\text{sin} \left(2 N \text{arcsin} \left(\frac{\tilde{\omega}_0 }{2 \sqrt{1-i \tilde{\gamma}_D \tilde{\omega}_0 }}\right)\right)}. \end{align} \tag{ 17 }$$

The dissipated power $P_{\mathrm{dis}}$ is defined as the rate of work done when driving the top bead n = 0. It takes the general form as an average over the driving period $\frac{2\pi}{\omega_0}$,

$$\begin{equation} P_{\mathrm{dis}} = \frac{\omega_0}{2\pi}\int_0^{\frac{2\pi}{\omega_0}} dt F\left(t\right) \, \dot u_0\left(t\right) . \end{equation} \tag{ 18 }$$

$\dot u_0(t) = -\omega_0U_0\text{sin}(\omega_0 t)$ is the prescribed velocity of the end bead n = 0 according to equation (12), and $F_0(t)$ is the force acting on that bead, i.e. the force required to exercise that motion. For the harmonic chain, the force is also monochromatic at ω₀, i.e. $F (t) = \textrm{Re}\{F(\omega_0) e^{i \omega_0 t}\}$ with a complex amplitude $F(\omega_0)$. $P_{\mathrm{dis}}$ may then be expressed in terms of it via

$$\begin{eqnarray} &&P_{\mathrm{dis}} = \frac{1}{2} U_0 \omega_0 \textrm{Im} \left\{F \left(\omega_0\right)\right\}.\end{eqnarray} \tag{ 19 }$$

The quantity of interest is the friction coefficient $\Gamma(\omega_0)$ of the bead n = 0. We will also refer to it as the macroscopic friction coefficient to distinguish from microscopic damping γ. Using the above definitions of Fourier amplitudes, it can be extracted as the ratio

$$\begin{eqnarray} &&\Gamma\left(\omega_0\right) = \frac{\textrm{Im} \left\{F\left(\omega_0\right)\right\}}{\omega_0 U_0} , \end{eqnarray} \tag{ 20 }$$

i.e. it is proportional to the component of F in phase with the velocity of the driven bead. Comparing equations (20) and (19), dissipated power and the macroscopic friction parameter Γ have the obvious relation,

$$\begin{eqnarray} &&P_{\mathrm{dis}}\left(\omega_0\right) = \frac{1}{2}\Gamma\left(\omega_0\right) \omega_0^2 U_0^2. \end{eqnarray} \tag{ 21 }$$

We will in the following exclusively discuss the behavior of Γ, from which $P_{\mathrm{dis}}$ may be extracted via equation (21).

We start by investigating the case of DPD damping. The expression for the force acting on the driven bead in presence of DPD damping is given by

$$\begin{align} \frac{F}{M\Omega_0^2 U_0}\equiv \tilde F_0\left(\tilde\omega_0\right) = \left(1 - i \tilde{\gamma}_D \tilde{\omega}_0 \right)\left(\tilde u_1\left(\tilde{\omega}_0\right)-\tilde u_0\left(\tilde{\omega}_0\right)\right). \end{align} \tag{ 22 }$$

It is expressed in terms of the difference of displacements of bead 0 and bead 1, and contains a contribution from both the potential (spring) as well as the damping (dashpot) connecting these particles in figure 1. As the friction force, via equation (20), is given by the imaginary part of $\tilde F$, the first term in equation (22) will pick up the imaginary part of $\tilde u_1(\tilde{\omega}_0)$, i.e. the motion out of phase with the driven bead. The second term will pick up the real part of $\tilde u_1(\tilde{\omega}_0)-\tilde u_0(\tilde{\omega}_0)$, i.e. the part in phase with the driven bead.

We will in the following consider several regimes of chain length compared to wavelength and decay length.

For a short chain, $N\ll \tilde{\lambda}_D, \tilde{l}_D$, we expand the sinus in equation (17) for small arguments, to obtain the displacement in leading orders in driving frequency,

$$\begin{align} \tilde u_n\left(\tilde{\omega}_0\right) = \left(\frac{\left(N-n\right)}{N} +\tilde{\omega}_0^2 \frac{ n\left(N-n\right) \left(2N - n \right)}{6 N} +\dots\right). \end{align} \tag{ 23 }$$

Notably, to leading orders, u is purely real: The imaginary part is of order ${\omega}_0^3$, i.e. an imaginary term of order ω₀ is missing: It cancels between numerator and denominator in equation (16). The force in equation (22) is, up to order ω₀, thus exclusively given by the first term in equation (23),

$$\begin{eqnarray} &&\tilde F_0\left(\tilde \omega_0\right) = \frac{\left(1 - i\tilde{\gamma}_D \tilde{\omega}_0\right)}{N}. \end{eqnarray} \tag{ 24 }$$

Equation (24) shows two terms. The first term is real, and yields the elastic force component. The second is imaginary, and it yields the friction coefficient Γ. In this limit, the displacement profile corresponds to quasistatic shear [38], for which the relative displacement of beads is proportional to $1/N$. Because of this, the elastic force amplitude, i.e. the elastic force acting on the top bead vanishes as $1/N$ in equation (24). This quasistatic chain profile can be interpreted as a Goldstone mode [41], which costs less and less energy the larger its wavelength. Notably, also the frictional component of force in equation (24) vanishes as $1/N$, showing that the Goldstone mode yields vanishing friction in DPD. Explicitly, $\tilde{\Gamma} \equiv\Gamma/M \Omega_0$ follows from equation (20),

$$\begin{eqnarray} &&\lim_{N\ll \tilde \lambda_D,\tilde l_D}\tilde \Gamma = \frac{\tilde\gamma_D }{N}. \end{eqnarray} \tag{ 25 }$$

Figure 2(a) shows numerical results for Γ as a function of N. The figure shows the underdamped case, i.e. $1\gg\tilde\gamma_D\tilde\omega_0$. The used parameters also fulfill $\tilde\gamma_D\tilde\omega_0\gg\tilde\omega_0^2$, for which $|\tilde q_D|\ll 1$; Wave- and decay length are thus large compared to lattice spacing a, corresponding to the continuum limit. With $\tilde q_D\approx\tilde\omega_0$ we expect equation (25) to be obeyed for $N\lesssim \tilde\omega_0^{-1}$. This is indeed confirmed by the curves in the graph, which deviate from equation (25) for $N\gtrsim \tilde{\omega}_0^{-1}$.

Zoom In Zoom Out Reset image size

We remark that equation (25) is also valid for overdamped cases, and can be mapped on the shear force felt by a simple fluid when between parallel plates of finite distance [38], see also the discussion section 6.

For a long chain, we use the exponential representation of $\text{sin}(x) = \frac{1}{2i}(e^{ix}-e^{-ix})$. If the chain is long compared to the decay length l_D , $N\gg\tilde l_D$, we can neglect the term $e^{iN a q_D}$. The dependence on N cancels in equation (17), and the result becomes independent of the boundary condition at n = N. We obtain in this limit,

$$\begin{align} \tilde u_n\left(\tilde \omega_0\right) = e^{i n \tilde q_D} = e^{-\frac{n}{\tilde l_D}}\left[i\text{sin}\frac{2\pi n}{\tilde\lambda_D}+\text{cos}\frac{2\pi n}{\tilde\lambda_D}\right]. \end{align} \tag{ 26 }$$

The long chain thus shows the expected behavior of an exponentially damped profile, with the meanings of $\tilde l_D$ and $\tilde\lambda_D$ becoming explicit. This is inserted in equation (22) for the force, which then leads the following expression for friction via equation (20)

$$\begin{align} \lim_{N\gg \tilde l_D} \tilde \Gamma = \frac{e^{-\frac{1}{\tilde l_D}}}{\tilde{\omega}_0} \text{sin}\frac{2\pi}{\tilde\lambda_D}+\tilde{\gamma}_D\left(1-e^{-\frac{1}{\tilde l_D}}\text{cos}\frac{2\pi}{\tilde\lambda_D}\right). \end{align} \tag{ 27 }$$

Equation (27) gives $\tilde \Gamma$ in terms of four parameters $\tilde l_D$, $\tilde \lambda_D$, $\tilde\omega_0$ and $\tilde\gamma_D$. These are not independent of each other, and one may express via equation (10),

$$\begin{align} \frac{4}{\tilde{\omega}_0^2}& = \textrm{Re}\frac{1}{\text{sin}^2\left(\frac{\pi}{\tilde\lambda_D}+\frac{i}{2\tilde l_D}\right)} \end{align} \tag{ 28 }$$

$$\begin{align} 4\tilde{\gamma}_D& = -\textrm{Im}\frac{1}{\text{sin}^2\left(\frac{\pi}{\tilde\lambda_D}+\frac{i}{2\tilde l_D}\right)}. \end{align} \tag{ 29 }$$

This shows that $\tilde\Gamma$ is determined from $\{\tilde l_D, \tilde \lambda_D\}$, or from $\{\tilde\omega_0,\tilde\gamma_D\}$. The mixed form of equation (27) is advantageous because of its compactness.

Before discussing the limiting cases, we note that the second term of equation (27) is bound by $2\tilde\gamma_D$. It is also easy to show that the first term in equation (27) is bound, namely as

$$\begin{align} \frac{\textrm{Im}\left[e^{i \tilde q_D\left(\tilde\omega_0,\tilde\gamma_D\right)}\right]}{\tilde\omega_0}\unicode{x2A7D} 1. \end{align} \tag{ 30 }$$

The friction of the long DPD chain is thus bound by

$$\begin{align} \lim_{N\gg \tilde l_D} \tilde \Gamma\unicode{x2A7D} 1+2\tilde\gamma_D.\end{align} \tag{ 31 }$$

Regarding in more detail we start with the underdamped limit, $\tilde{\gamma}_D\tilde{\omega}_0\ll1$, for which generally the first term in equation (27) dominates. Taking formally $\tilde l_D\gg1$, the first term yields

$$\begin{eqnarray} &&\lim_{\tilde \gamma_D\tilde \omega_0\ll 1} \lim_{N\gg \tilde{l}_D} \tilde{\Gamma} = \frac{\text{sin}{\frac{2\pi}{\tilde{\lambda}_D}}}{2\text{sin} {\frac{\pi}{\tilde{\lambda}_D}}}, \end{eqnarray} \tag{ 32 }$$

which, for allowed values of $\tilde\lambda_D\unicode{x2A7E} 2$ within the Brillouin zone varies between 0 and 1. This result is notably independent of the damping coefficient $\tilde{\gamma}_D$, and corresponds to radiation of waves, i.e. in this case, energy is transported away, and dissipated far away from the source. Taking the continuum limit of $\tilde\lambda_D\gg1$,

$$\begin{eqnarray} &&\displaystyle{\lim_{\tilde\lambda_D \gg1}\lim_{\tilde\gamma_D\tilde\omega_0\ll 1} \lim_{N\gg \tilde l_D} }\tilde \Gamma = 1. \end{eqnarray} \tag{ 33 }$$

Notably, in this limit, Γ is independent of $\tilde\omega_0$, i.e. it corresponds to the limit where a simple macroscopic result for friction is found. The connection to radiation is easily seen, as $P_{\mathrm{dis}} = \frac{1}{2}M \Omega_0 \omega_0^2 U_0^2$, i.e. it is kinetic energy of a bead, $\frac{1}{2} M\omega_0^2 U_0^2$, which is excited at rate Ω₀, i.e. via a wave traveling with the speed of sound. Such radiation was also found for a small probe in a 3D system in [5], and its universality was discussed in [36].

At the end of the Brillouin zone, i.e. as $\tilde{\lambda}_D\to 2$, this term vanishes. It corresponds to $\tilde\omega \to 2$ from below. Of course, $\tilde\lambda_D$ cannot be smaller than 2 [41] (as the real part of arcsin cannot exceed $\pi/2$). For $\tilde\omega_0\gt2$, thus $\tilde\lambda_D = 2$, and the first term in equation (27) is identically zero, as no wave is emitted. $\tilde{\Gamma}$ is then given by the second term. As this term is small for underdamped cases, $\tilde\Gamma$ shows a sharp drop at $\tilde\omega_0 = 2$. For $\tilde\omega_0\gt2$, $\tilde{l}_D$ is finite and grows, even for underdamped cases, so that,

$$\begin{eqnarray} &&\displaystyle{\lim_{\tilde\omega_0\gg 2}\lim_{\tilde{\gamma}_D\tilde{\omega}_0\ll1} \lim_{N\gg \tilde{l}_D} }\tilde{\Gamma} = \tilde{\gamma}_D. \end{eqnarray} \tag{ 34 }$$

This is the friction coefficient of a chain where only the driven bead moves, and the rest of the chain is at rest. This case, where no wave is emitted, is reminiscent of evanescent waves [15] found in 3D.

In the overdamped limit, $\tilde{\gamma}_D\tilde{\omega}_0\gg1$, the second term in equation (27) dominates. In the sound wave limit, $|q_D|\ll 1$, we obtain, using $\tilde{\lambda}_D/(2\pi) = \tilde{l}_D = \sqrt{2\tilde \gamma_D/\tilde\omega_0}$,

$$\begin{align} \displaystyle{\lim_{|q_D| \ll 1}\lim_{\tilde{\gamma}_D\tilde{\omega}_0\gg1} \lim_{N\gg \tilde{l}_D} }\tilde{\Gamma} = \sqrt{\frac{\tilde{\omega}_0\tilde{\gamma}_D}{2}}. \end{align} \tag{ 35 }$$

In this limit, the friction coefficient behaves as $\sqrt{\tilde{\omega}_0}$, i.e. in the overdamped limit, there is no obvious regime of simple frequency independent macroscopic friction coefficient.

Interestingly, $\tilde\Gamma$ in equation (35) equals the underdamped result multiplied by $\sqrt{\frac{\tilde{\omega}_0\tilde{\gamma}_D}{2}}$. As this term, by construction, is large in the overdamped limit, the overdamped limit generally has a large friction coefficient compared to the underdamped case.

As a side note, we remark again that the overdamped limit can be mapped on a viscous fluid. Indeed, the result of (35) corresponds to the (shear) friction felt by a wall moving laterally in an incompressible fluid with shear viscosity $\tilde{\gamma}_D$ [38] (so called vorticity diffusion).

If, in the overdamped limit, we reach the end of the Brillouin zone at $\tilde\omega\approx\tilde\gamma_D$, we find

$$\begin{align} \displaystyle{\lim_{\tilde\omega_0 \gg \tilde\gamma_D} \lim_{\tilde{\gamma}_D\tilde{\omega}_0\gg 1} \lim_{N\gg \tilde{l}_D} }\tilde{\Gamma} = \tilde{\gamma}_D, \end{align} \tag{ 36 }$$

which is again the friction of the bead when the second bead is not moving. The crossover between equations (35) and (36) occurs when $\tilde\Gamma_D$ approaches $\tilde\gamma_D$ from below.

Figure 2(b) shows $\tilde\Gamma$ as a function of $\tilde\omega_0\tilde\gamma_D$ for various chain lengths N. With the chosen value $\tilde\gamma_D = 80$, we have $\tilde\omega_0\ll1$ and $\tilde\omega_0\ll\tilde\gamma_D$ for the regimes shown in the graph, i.e. the curves correspond to the continuum limit. For the largest N, we both observe the radiation behavior of $\tilde\Gamma = 1$ of equation (33) as well as the scaling of equation (35). The oscillations seen in the curves result from interference effects with reflected waves. They are not contained in the analytical formulas for a long chain, and will be discussed in the next subsection.

The different friction scaling regimes are further elucidated in figure 2(c) where we have plotted $\tilde{\Gamma}$ as a function of $\tilde{\omega}_0 \tilde{\gamma}_D$ for a fixed chain length $N = 1 0^4$ and different values of $\tilde{\gamma}_D$. In the underdamped limit, the oscillations are observed followed by the phonon radiation mode with $\tilde{\Gamma} = 1$ independent of the damping $\tilde{\gamma}_D$. As discussed in equation (34), $\tilde{\Gamma}$ shows a sharp drop and settles to the value $\tilde{\Gamma} = \tilde{\gamma}_D$. In the overdamped regime, $\tilde{\Gamma}$ follows $\sqrt{\tilde{\omega}_0 \tilde{\gamma}_D}/2$ in the continuum limit following equation (35). Finally, the friction coefficient $\tilde{\Gamma}$ saturates to $\tilde{\gamma}_D$ (see equation (36)) for large values of $\tilde{\omega}_0 \tilde{\gamma}_D$ when the decay length is $\tilde l_D \ll 1$.

Equations (37) and (38) below hold both for DPD and Langevin systems so that we omit indices D or L for the first part of this subsection.

The amplitude $\tilde{u}_n(\tilde{\omega}_0)$ in equation (15) encodes resonances [43], occurring for $N\textrm{Re}[\tilde{q}] = p \pi$ with integer p, corresponding to $\lambda = 2\frac{N}{p}$, i.e. to positive interference of original and reflected waves. More precisely, resonances also require underdamped conditions with $N\textrm{Im}[\tilde{q}]\ll 1$. Being an underdamped phenomenon, this effect is strong for small damping, and it yields a dominant contribution from the term $\tilde F_0(\tilde\omega_0) = (\tilde u_1(\tilde\omega_0)-\tilde u_0(\tilde\omega_0))$ in both equation (22) as well as equation (39). It is the imaginary part of $\tilde u_1$ which contributes to friction in this case, which with the above resonance condition, and for $N\textrm{Im}[\tilde{q}]\ll1$ (i.e. the wave reflects many times), is

$$\begin{align} \textrm{Im}\left[\tilde{u}_1\right] = \frac{\text{sin}\left(\frac{N -1}{N}p\pi\right)}{N \textrm{Im}\left[\tilde{q}\right]}. \end{align} \tag{ 37 }$$

This yields for the friction coefficient, for both Langevin and DPD damping, in this limit,

$$\begin{align} \tilde\Gamma = \frac{\text{sin}\left(\frac{N -1}{N}p\pi\right)}{\tilde\omega_0 N \textrm{Im}\left[\tilde{q}\right]}. \end{align} \tag{ 38 }$$

This relation shows that the friction coefficient $\tilde\Gamma$, for a given finite N, can grow to arbitrary large values. Notably, it diverges with $\textrm{Im}[\tilde q]\to0$, providing the counter intuitive observation that smaller damping can lead larger friction. This originates from the fact that as the radiated waves are reflected a large number of times, each damping unit (bond) encounters multiple radiation waves and consequently leads to larger damping.

On the contrary, for $(N-1)\textrm{Re}[\tilde{q}] = p \pi$ with integer p, the numerator of $\tilde{u}_1(\tilde{\omega}_0)$ gets small, so that, for underdamped chains $\tilde\Gamma$ shows minima there. This situation has been termed anti-resonance [39, 40].

Specifically, for DPD, regarding equation (10), the overdamped regime is found for $\tilde \omega\ll \tilde\gamma_D^{-1}$. Hence, for any value of damping $\tilde\gamma_D$, a sufficiently small value of $\tilde\omega$ renders the chain underdamped. The resonance occurs then for a sufficiently large N.

Figure 2(d) shows the variation of $\tilde{\Gamma}$ with $\tilde{\gamma}_D\ll 1$. For the graph, we fix $\tilde{\omega}_0$ such that $N \text{arcsin}[\frac{\tilde{\omega}_0}{N}] = \pi$ for the different values of N. This corresponds to $N\textrm{Re} [\tilde{q}_D] = \pi+\mathcal{O}(\tilde\gamma_D)$. As expected, $\tilde{\Gamma}$ diverges as $\tilde{\gamma}_D \to 0$, following equation (38).

In presence of Langevin damping, the Fourier amplitude of the force F₀ acting on the driven particle is given by

$$\begin{eqnarray} &&\tilde F_0\left(\tilde\omega_0\right) = \left(\tilde u_1\left(\tilde\omega_0\right)-\tilde u_0\left(\tilde\omega_0\right)\right) +i \tilde\gamma_L\tilde\omega_0 \tilde u_0\left(\tilde\omega_0\right). \end{eqnarray} \tag{ 39 }$$

This force has two components, it results from the interaction with the second particle (position u₁) as well as the trivial part, i.e. the Langevin damping of the driven bead (via coefficient γ_L ), the second term in equation (39). The contribution of the first term to friction is found by inserting the expressions for $u_0(\omega_0)$ and $u_1(\omega_0)$ from equation (16) in equation (39). We again consider several regimes of chain length compared to wavenumber in equation (7).

As for the DPD case, we start with the limit of small $\tilde{\omega}_0$, as this corresponds to large wave- and decay lengths. In this limit

$$\begin{align} \tilde u_n\left(\tilde\omega_0\right) = \frac{N- n}{N}+i\tilde\gamma_L\tilde\omega_0\frac{n \left(n^2-3 nN+2N^2\right)}{6 N }+\dots. \end{align} \tag{ 40 }$$

The real part is to leading order, as expected, identical to the DPD case. In contrast to that, here a term linear in $\tilde{\omega}_0$ exists, and it will yield the bead's friction coefficient.

For the force amplitude, we obtain,

$$\begin{align} \tilde F_0\left(\tilde{\omega}_0\right) = \frac{1}{N} + i \tilde{\gamma}_L \tilde{\omega}_0\left(\frac{1}{3} N +\frac{1 }{6 N} + \frac{1}{2} \right) +\dots. \end{align} \tag{ 41 }$$

Again, the real part is the anticipated Goldstone result, vanishing as $1/N$ for large N. In contrast, the Goldstone mode here leads to an imaginary part which grows with N. Following equation (20) the macroscopic friction coefficient Γ is given by

$$\begin{equation} \lim_{N\ll \tilde{\lambda}_L, \tilde{l}_L} \tilde{\Gamma} = \tilde{\gamma}_L \left(\frac{1}{3} N +\frac{1 }{6 N} + \frac{1}{2} \right). \end{equation} \tag{ 42 }$$

As already mentioned, this result if fundamentally different from the DPD result, equation (25), as it grows with N. This is because for Langevin damping, the excited Goldstone mode yields the dissipation of a 'massive' mode as damping breaks translational symmetry: Every bead dissipates with respect to a background medium.

Interestingly, assuming a simple (quasistatic) profile and motion of the chain, one may naively expect a prefactor of $\frac{1}{2}$ for large N in equation (42); the simple profile suggests that the average velocity of the beads is half the velocity of the driven bead. However, already to order ω₀ (which contributes to friction) u_n deviates from that profile. Notably, the factor of $\frac{1}{3}$ can be understood by integrating the square of the simple profile, via $\int dx (1-x)^2 = \frac{1}{3}$. Thus, the dissipated energy can be naively understood in this limit as being additive, while the friction force is not.

The numerically exact result for Γ, as a function of N, is plotted in figure 3(a) with $\tilde\gamma_L = 1/100$. As $\tilde\gamma_L\gg\tilde\omega_0$ for the shown curves, the graph shows the overdamped continuum limit. In this case, $q_L\approx \sqrt{\tilde\gamma_L \tilde\omega_0} (1+i)/\sqrt{2}$. We expect equation (42) to hold for $N\ll 1/|\tilde q_L|$, and indeed, the curves in the graph deviate at around $N \approx \sqrt{1/\tilde\omega_0\tilde\gamma_L}$.

Zoom In Zoom Out Reset image size

Finally, we remark that the friction of equation (42) is related to the center of mass diffusion of a Rouse polymer [37]. The additional factor of $1/3$ is due to the anchoring, and is absent for a chain with a free end.

If the chain is long, $N\gg\tilde l_L$, we repeat the steps of section 4.3, and expand the sinus for large imaginary argument,

$$\begin{align} \tilde u_n\left(\tilde \omega_0\right) = e^{i n \tilde{q}_L} = e^{-\frac{n}{\tilde{l}_L}}\left[i\text{sin}\frac{2\pi n}{\tilde{\lambda}_L}+\text{cos}\frac{2\pi n}{\tilde{\lambda}_L}\right]. \end{align} \tag{ 43 }$$

This is inserted in equation (39) for the force. Here, only the imaginary part of u₁ contributes, so that, via equation (20)

$$\begin{align} \lim_{N\gg \tilde{l}_L} \tilde{\Gamma} = \tilde{\gamma}_L+ \frac{e^{-\frac{1}{\tilde{l}_L}}}{\tilde{\omega}_0} \text{sin}\frac{2\pi}{\tilde{\lambda}_L}. \end{align} \tag{ 44 }$$

Let us be reminded that the first term in equation (44) is the trivial contribution of the damping of the driven bead. We may express the amplitude of the second term via equation (6) as

$$\begin{align} \frac{4}{\tilde{\omega}_0^2}& = \frac{1}{\textrm{Re}\left[\text{sin}^2\left(\frac{\pi}{\tilde{\lambda}_L}+\frac{i}{2\tilde{l}_L}\right)\right]}, \end{align} \tag{ 45 }$$

showing that $\tilde\Gamma$ is fully determined from either $\tilde l_L$ and $\tilde \lambda_L$, or from $\tilde\omega_0$ and $\tilde\gamma_L$.

Before discussing specific limits, let us ask whether $\tilde\Gamma$ is bound as for the DPD case in equation (31). With equation (45), we have $\textrm{Re}[\text{sin}\tilde{q}_L/2] \unicode{x2A7E} \textrm{Im}[\text{sin}\tilde{q}_L/2] \unicode{x2A7E} 0$ and $\tilde\omega_0^2 = 4(\textrm{Re}[\text{sin}(\tilde q_L/2)]^2- \textrm{Im}[\text{sin}(\tilde q_L/2)]^2)$, and we obtain,

$$\begin{align} \Gamma = \gamma_L +\frac{\textrm{Im}\left[e^{i\tilde q_L}\right]}{2\sqrt{\textrm{Re}\left[\text{sin}\left(\tilde q_L/2\right)\right]^2 - \textrm{Im}\left[\text{sin}\left(\tilde q_L/2\right)\right]^2}}. \end{align} \tag{ 46 }$$

We see that the denominator can get arbitrarily close to zero for $\textrm{Re}[\text{sin}(\tilde q_L/2)]^2 \to \textrm{Im}[\text{sin}(\tilde q_L/2)]^2$. Because the numerator stays finite in this limit, it corresponds to a pole, i.e. $\tilde\Gamma$ can get arbitrarily large. This is a pronounced difference to the DPD case, where $\tilde\Gamma$ for the long chain is bound according to equation (31). There is however a (trivial) bound also here, but it is a bound from below

$$\begin{align} \lim_{N\gg \tilde l_L}\tilde\Gamma\unicode{x2A7E} \gamma_L.\end{align} \tag{ 47 }$$

As for DPD, we start with the underdamped case, which is here given for $\tilde{\omega}_0\gg \tilde{\gamma}_L$. Using formally $\tilde l_L \gg 1$, we have

$$\begin{eqnarray} &&\displaystyle{\lim_{\tilde{\omega}_0\gg\tilde{\gamma}_L} \lim_{N\gg \tilde{l}_L}} \tilde{\Gamma} = \tilde{\gamma}_L+ \frac{\text{sin}{\frac{2\pi}{\tilde{\lambda}_L}}}{2\text{sin}{\frac{\pi}{\tilde{\lambda}_L}}}, \end{eqnarray} \tag{ 48 }$$

which also agrees with the DPD case, equation (34), up to the first term. As for DPD, in the sound wave limit, i.e. $|\tilde \lambda_L|\gg 1$, and

$$\begin{align} \lim_{|\tilde q_L|\ll 1}\lim_{\tilde{\omega}_0\gg\tilde{\gamma}_L}\lim_{N\gg \tilde{l}_L} \tilde{\Gamma} = \tilde\gamma_L +1. \end{align} \tag{ 49 }$$

This result, corresponding to the radiation mode, is also identical to the result from DPD damping in the underdamped limit (up to the trivial term). As the second term in equation (49) is independent of friction coefficient, it should indeed be reached by either case. Notably, only for $\tilde\gamma_L\ll1$, the radiation term dominates friction, which is a condition that is additional to the one of being underdamped.

The end of the Brillouin zone, i.e. $\tilde \lambda_L = 2$ is reached for $\tilde\omega_0\approx 2$. At this point, the second term in equation (48) vanishes, and $\tilde{\Gamma} = \tilde{\gamma}_L$ for $\tilde\omega_0\gtrsim 2$.

The overdamped limit corresponds to $\tilde{\omega}_0\ll\tilde{\gamma}_L$. In the continuum limit, $\tilde{\lambda}_L/(2\pi) = \tilde{l}_L = \sqrt{2/\tilde{\gamma}_L \tilde{\omega}_0}$, and we have,

$$\begin{align} \lim_{1\ll\tilde{\lambda}_L,\tilde{l}_L}\lim_{\tilde{\omega}_0\ll\tilde{\gamma}_L} \lim_{N\gg \tilde{l}_L} \tilde{\Gamma} = \sqrt{\frac{\tilde{\gamma}_L}{2\tilde{\omega}_0}} = \frac{1}{2}\tilde\gamma_L \tilde l_L. \end{align} \tag{ 50 }$$

$\tilde{\Gamma}$ in equation (50) shows non-trivial behavior, as it grows as an inverse square root of $\tilde{\omega}_0$. As anticipated around equation (46), $\tilde\Gamma$ can grow without bounds, and we see now how the mentioned pole is approached. As long as $\tilde l_L\ll N$, $\tilde\Gamma$ grows with $\tilde l_L$, as the number of particles contributing to dissipation grows (second equality in equation (50)). As $\tilde l_L \sim 1/\sqrt{\tilde \omega_0}$ in the considered case, the friction thus diverges as frequency goes to zero in the seen fashion.

We remark that the result of equation (50) can be translated to the short time diffusion of a segment in a Rouse polymer [37]. In time domain, $\tilde\Gamma$ of equation (50) yields $\tilde{\Gamma}(t) \sim t^{3/2}$ and the mean square displacement (MSD) $\sim t^{1/2}$. We omit here a detailed translation from friction to diffusion kernel [44].

In the overdamped limit, the end of Brillouin zone is reached for $\tilde{\omega}_0\tilde\gamma_L\approx 8$, and the second term in equation (44) vanishes, and

$$\begin{align} \lim_{\tilde{\omega}_0\ll\tilde{\gamma}_L} \lim_{N\gg \tilde{l}_L} \tilde{\Gamma} = \tilde{\gamma}_L, \end{align} \tag{ 51 }$$

which, as for the DPD chain, is dominated by the bare friction of the driven bead.

Figure 3(b) shows numerically exact results for $\tilde{\Gamma}$ as a function of $\tilde{\omega}_0/\tilde{\gamma}_L$, for different values of N. Notably, for a given N, equation (44) holds for $\tilde\omega_0\gg \sqrt{2}/N$. It deviates for smaller frequencies, as the decay length $\tilde l_L$ then exceeds N.

The different regimes are summarized in figure 3(c) which shows $\tilde{\Gamma}$ as a function of $\tilde{\omega}_0 /\tilde{\gamma}_L$ for a fixed chain length $N = 1 0^4$ and different values of damping parameter $\tilde{\gamma}_L$. When $\tilde{\omega}_0 \ll \tilde{\gamma}_L$ (the overdamped limit), $\tilde{\Gamma}$ exhibits a regime following equation (44) followed by saturation as one decreases $\tilde{\omega}_0$ beyond $l_L \approx N$. The radiation mode ($\tilde{\Gamma} = 1$) is observed at the onset of the underdamped limit for $\tilde\gamma_L\ll 1$ followed by the regime corresponding to the end of the Brillouin zone, where the friction coefficient is dominated by the bare friction of the driven bead, $\tilde{\Gamma} = \tilde{\gamma}_L$.

The discussion of resonant damping in section 4.4 is analogous here, i.e. equations (37) and (38) apply. Resonances are not observed in figure 3 due to choice of parameters. In figure 3, the regime of intermediate chain lengths coincides with the overdamped regime.

The resonant behavior requires $\tilde{\gamma}_L \ll \tilde{\omega}_0\approx\pi/N$, i.e. a sufficiently small value of $\tilde\gamma_L$. Figure 3(d) shows $\tilde{\Gamma}$ as a function of $\tilde{\gamma}_L$ in the resonant chain regime. Similar to the DPD case, $N\textrm{Re} [\tilde{q}_L] = \pi+\mathcal{O}(\tilde\gamma_L)$ is kept fixed by choosing $\tilde{\omega}_0$ such that $N \text{arcsin}[\frac{\tilde{\omega}_0}{N}] = \pi$. $\tilde{\Gamma}$ is observed to diverge as one approaches $\tilde{\gamma}_L \to 0$ in the resonant regime.

Interestingly, in contrast, for DPD, occurrence of (anti)resonances requires $\tilde{\gamma}_D\tilde\omega \ll 1$ and $\tilde{\omega}_0\approx\pi/N$, i.e. they can be observed for any value of damping parameter γ_D , for sufficiently large chain length N and small driving frequency $\tilde{\omega}_0$.