Stochastic Thermodynamics of Oscillators’ Networks

We apply the stochastic thermodynamics formalism to describe the dynamics of systems of complex Langevin and Fokker-Planck equations. We provide in particular a simple and general recipe to calculate thermodynamical currents, dissipated and propagating heat for networks of nonlinear oscillators. By using the Hodge decomposition of thermodynamical forces and fluxes, we derive a formula for entropy production that generalises the notion of non-potential forces and makes transparent the breaking of detailed balance and of time reversal symmetry for states arbitrarily far from equilibrium. Our formalism is then applied to describe the off-equilibrium thermodynamics of a few examples, notably a continuum ferromagnet, a network of classical spin-oscillators and the Frenkel-Kontorova model of nano friction.


Introduction
Dissipation and heat transfer are universal phenomena in Physics, appearing whenever a small system is coupled to the much larger environment. In this situation, it is, in practice not, possible to keep track of all the observables of the universe. Instead, some measurable macroscopic quantities, such as energy, entropy and heat flows are used to describe the evolution of the system and of its average properties.
In the presence of several environments (thermal baths or reservoirs) at different temperatures, the system reaches a non-equilibrium steady state where thermodynamical currents (such as heat, energy, spin or electrical) may flow through the system from one reservoir to the other. Close to equilibrium, those currents are proportional to the corresponding thermodynamical forces. Examples of thermodynamical forces include gradients or differences of temperature, voltage, chemical potentials and concentrations of chemical species. The language of thermodynamical forces and currents, which is nowadays the cornerstone of non equilibrium thermodynamics, was first developed by L. Onsager in the 1930s [1,2], and by R. Kubo in the 1950s [3]. The formalism can be naturally extended beyond the linear regime, as it was first observed by Schnakenberg [4] and subsequently in more recent works on stochastic thermodynamics [5][6][7].
In out-of-equilibrium setups, it is of primary importance to determine the (possibly many) currents that flow between subparts of the system, together with the corresponding forces, and the heat flow that is dissipated to the environment [8]. The first case corresponds typically to the work done on the system by the environment, while the latter is associated to the production of entropy and an increased disorder of the ensemble (the system and environment), and it is related to the efficiency of the thermodynamical process. In networks of nonlinear oscillators, the heat flow throughout the system is a coherent phenomenon that requires synchronisation of the oscillators, while the dissipated heat is incoherent [8][9][10].
The scope of this paper is twofold. On one hand, to provide a unified and concise view of many subjects that are scattered in the literature and are apparently disconnected, such as non-hermitian hamiltonians, brackets and anti-brackets formalism, canonical transformations and their connections to heat transfer in oscillators' networks. Then, and most importantly, to give a general recipe to calculate currents, dissipated heat and work in a large class of out-of-equilibrium systems. To this end, we shall adopt the formalism of stochastic thermodynamics (ST) [5] applied to the dynamics of complex-valued Langevin and Fokker-Panck equations.
Complex Langevin equations, which here represent the paradigm to treat oscillator's networks, have been long investigated and have application in a variety of physical systems, from the stochastic formulation of Quantum Mechanics [11,12] to Quantum Cromo Dynamics [13], quantum statistics [14] and lattice gauge theories [15,16]. Recently, they found application also in polymer field theory [17] and random matrices [18]. In the present paper, we shall focus on the discrete nonlinear Schrödinger equation (DNLS), a general complex Langevin oscillator equation, which has ramification in many branches of Physics and has recently attracted a certain attention.
At variance with the standard formulation of stochastic thermodynamics, which uses colloidal particles as paradigm [5], using complex Langevin equations, here we provide a systematic way to describe the non-equilibrium thermodynamics of networks of nonlinear oscillators. Following the idea of ST, we start from the stochastic trajectories of a small ensemble of oscillators coupled to a bath and we extract useful information (such as currents and entropy) from the ensemble averages of the main observables. As the system evolves, different parts of the network and the environment become statistically correlated. The currents are expressed in terms of those correlations.
The present paper is organised as follows. In Section 2 we apply the Lagrangian and Hamiltonian formalism to describe complex-valued equations of motion, following References [19][20][21][22]. In particular, we formulate the conservative and dissipative dynamics, respectively, in terms of Poisson commutators and anti-commutators for a possibly non-Hermitian Hamiltonian, a topic that has been extensively studied [19,[23][24][25][26].
In Section 3 we develop the stochastic thermodynamics formalism for complex-valued equations, and we derive a simple and general formula for entropy production, which makes transparent the breaking of detailed balance and is proportional to the heat dissipated to the environment. This section generalises previous work [8,27,28] to complex Langevin equations with multiplicative noise. Section 4 contains the formulation of the first principle of thermodynamics. This constitutes the main result of this paper and allows one to identify the heat transported and dissipated.
Section 5 provides some examples of realistic physical systems where thermodynamical currents and entropy production are calculated. We shall describe in particular the dynamics of a one dimensional continuum ferromagnet, of a network of classical magnetic spins and of the Frenkel-Kontorova model [29] for nano-friction.

Hamiltonian-Lagrange Formulation for Complex Equations of Motion
We consider here the following complex Langevin equation: where the dot indicates time derivative and ψ m = p m (t)e iφ m (t) is a complex wave function with amplitude p m = |ψ 2 m | (also referred as power) and phase φ n . The force F m is an arbitrary function of the ψ = (ψ 1 , ..., ψ N ) and their complex conjugate ψ * = (ψ * 1 , ..., ψ * M ). We assume that both the coupling between the ψs and the damping are contained in the definition of F m .
The terms ξ m , which model the stochastic baths, are complex Gaussian random processes with zero average and correlation ξ m (t)ξ * n (t ) = δ mn δ(t − t ). Here g m is an arbitrary function of the (ψ, ψ * ). Throughout the paper vectors and matrices are written in bold text, while their components are written in plain text with m and n subscripts. The quantity |g 2 m | plays the role of diffusion constant.
We assume that the latter is proportional a damping coefficient Γ m and temperature T m , according to the fluctuation-dissipation theorem. Thus, at variance with previous studies [8], we consider here the more general situation of Langevin equations with multiplicative noise. The force is given by the derivative of a complex (and possibly non-Hermitian) Hamiltonian H. Here the Wirtinger derivatives are defined as where ∂ * m is the complex conjugate and ψ m = x m + iy m . The complex conjugate equation to Equation (1) contains the forces F * m = −i∂ m H. A straightforward calculation shows that in a dissipative system with Hamiltonian H = H R + H I , where R and I are respectively the Hermitian (or reversible) and anti-Hermitian (or irreversible) components, the dynamics of an arbitrary function f of the observables (ψ, ψ * ) can be written asḟ Here the Poisson commutators (−) and anti-commutators (+) are defined respectively as We note in particular that, from Equations (2)-(4), the reversible and irreversible forces can be expressed as F R m = i{H R , ψ m } − and F I m = {H I , ψ m } + . On the other hand, the couple (ψ n , iψ * n ) are canonical conjugate variables, since one has i{ψ m , iψ * n } − = δ mn . In order to derive Equation (4), one proceeds as in the case of classical Hamiltonian mechanics, by writing the evolution equation for an arbitrary function of the observableṡ then one substitutes the equations of motions:ψ m = i∂ m (H R + H I ) and its complex conjugate equation forψ m . Note that ψ m and the conjugate variable iψ * m are considered here independent variables. We remark that the time reversal of our system corresponds to the transformation of the Hamiltonian H(t) → H * (−t). The irreversible and reversible components respectively change and do not change sign under this transformation. A typical way to make the Hamiltonian irreversible and non Hermitian is by adding a non-symmetric coupling between the oscillators. This will be discussed in the next sections, and has been treated extensively also in Reference [8].
Since the commutators and anti-commutators define respectively a symplectic and a metric structure on the space tangent to the phase space, this kind of system is called metriplectic. The formulation of dissipative dynamics in terms of anti-brackets in metriplectic structures has been extensively studied, both for classical and quantum systems [24,[30][31][32]. In those formulations, the irreversible part of the Hamiltonian is usually identified with the entropy of the system. Here we do not pursue this identification, since we will describe the irreversibility in terms of the information entropy and the associated entropy production, as it is customary in the ST formalism. Later in the paper, we shall elucidate the connection between the irreversible part of the Hamiltonian and the entropy production.
An important step here is to determine the canonical transformations, that must preserve the metriplectic structure. In practice, from the definition of force and from Equation (5) one must have thaṫ for a variable b k function of the old coordinates ψ. From the chain rule of partial derivative one has However, from the definition of commutators and anti-commutators the following equalities must also hold: The two equalities can be both satisfied only if the following holds: This means essentially that the new coordinates must be analytic functions of the old ones, since they cannot contain both a variable and its complex conjugate. Adding a complex number or performing a U(1) gauge transformation preserves the commutators [16], however note that the Bogoliubov transformations are not canonical in this case, although they are canonical transformations of the system without dissipation. Note that the system can be described using the following Lagrangian, similar to the one for the heat Equation [20]: The equations of motion for ψ * m are given by the Euler-Lagrange equations while the dynamics of ψ given by the complex conjugate equations. Equations (13) and (14) are particularly useful to determine the conserved currents of the system associated to the invariance of the Lagrangian with respect to a global U(1) transformation, according to the Noether theorem, as it will be clarified in the next sections. We remark also that one could in principle consider a stochastic Lagrangian, from which the full Langevin equation can be obtained by means of the usual Euler-Lagrange equations, as it has been noted in Reference [16].

Fokker-Planck Equation and Entropy Production
This section generalises the material presented in Reference [8] to the case of multiplicative noise. The time evolution of the probability distribution in the phase space, associated to the Langevin Equation (1), is given by the following Fokker-Planck (FP) equation: Following References [8,27], we define the reversible and irreversible probability currents as with J m = J R m + J I m and J * m the complex conjugate. By using those currents, Equation (15) assumes the usual form of a continuity equation: Thermal equilibrium corresponds to the case where the probability currents are zero, while non-equilibrium steady state corresponds to non-zero divergenceless currents, withṖ = 0.
The entropy flow Φ and entropy production Π are obtained starting from the definition of phase space entropy where · denotes the ensemble average. Computing the time derivativeṠ by means of Equation (17) we obtain:Ṡ Upon integrating by parts, and observing that the divergence of the reversible forces, We remark that the fact that the reversible forces have zero divergence is due to the expression for the force From the definition of probability currents Equation (16) one has together with the complex conjugate equation. Upon substituting the previous equation into Equation (20) givesṠ The first and second terms are respectively entropy flow and entropy production. We remark that we here consider only steady states. In this condition, assuming that the probability currents vanish at infinity [8,27], one can integrate by part the last term, which is proportional to ∂ m J + c.c. However, since the divergence of the thermodynamical currents is zero in stationary states, the last term vanishes in that case. Thus, in steady state the entropy flow Φ is minus the entropy production Π [8,27], as in the case of additive noise.
At this point we substitute integrals containing P with ensemble averages. In this way Equations (16) and (22) give which is the same expression obtained in Reference [8], with |g m | 2 playing the role of diffusion constant. As in References [8,27] we identify the quantity TΦ with the heat exchanged with the bath.
We proceed now by deriving an expression for the entropy production that makes transparent the breaking of detailed balance and the onset of irreversibility. Since in steady states one has ∑ m (∂ m J * m + ∂ m J * m ) = 0, one can apply the Hodge decomposition [33] and write the currents as where Ω is an anti-symmetric tensor and Λ a scalar. We separate the entropy flow into two components Φ 1 and Φ 2 containing respectively Ω and Λ. For the first component one has where we have used the anti-symmetry of Ω and integrated by parts discarding the boundary terms. One has that the condition of detailed balance is ∂ m F |g | 2 − ∂ F m |g m | 2 = 0, which is met when the forces are potentials and/or the temperatures are the same, |g m | 2 = |g | 2 . Note that this condition generalises the formulation of References [27,28] to the case of complex-valued forces.
However, in our system we have two coupled currents, associated respectively to the conservation of energy and of the total power p n , or "number of particles". Thus, we expect that the entropy production contains two components: one that depends on the temperature differences and one that depends on the chemical potential differences.
To see this, let us write the force as the derivative of the following Hamiltonian: where µ k is the local chemical potential. A straightforward calculation gives for Equation (25): The first term is non zero if the Hamiltonian is non-Hermitian and/or if the temperatures are different. On the other hand, the second term is non zero if the chemical potentials or the temperatures are different.
There is also another way to drive the system off equilibrium: by applying a constant chemical potential that compensates the damping [8,34]. In this case one expects that the entropy production is non zero, even if the current vanishes. To see this, we consider the second contribution to the entropy production: However, if we write the Hamiltonian as in the case of the DNLS [8,35], the term Re[∂ m ∂ * m H] is the damping of the system, Γ m . Thus one has Φ 2 ∝ iΛ(Γ m − µ m )dx. This shows that the system does not relax to equilibrium in the case where the chemical potential compensates the damping, as has been pointed out also in Reference [34].

Transported Vs Dissipated Heat
This section contains the main results of the paper. Starting from the first principle of thermodynamics, we derive the expressions for the heat dissipated and flowing through the system. To keep the notation simple, we consider the case with g m = D m ≡ αT m , with α the damping of the system. It is straightforward to generalise our discussion to the case of multiplicative noise. Following Reference [8], for a network of m = 1, ..., M oscillators, we consider the Hamiltonian where h m is the local energy, and the Hamiltonian splits into a reversible and irreversible component, Here µ m is the local chemical potential, while |ψ n | 2 plays as usual [8] the role of particle number. The first principle of thermodynamics can be expressed as a balance equation for the energy according to d dt where dx = i 2 ∑ dψ m ∧ dψ * m is the volume element of the phase space. In Reference [5] and in stochastic thermodynamics in general, the first and second terms of the previous equation are respectively identified with heat Q and work W. However, in the present case one does not have a clear distinction between heat and work. In particular, it turns out that Q is the heat dissipated to the environment, while W contains contributions both from the dissipated heat and from the heat flow that propagates between oscillators m and n in the network.
To see this, we started by calculating the value of Q. We remark that, as observed in Reference [8], only the irreversible part of the Hamiltonian enters these expressions. We use the FP equation Equations (15) and (17) and substituteṖ with the derivative of the currents J : Then, upon substituting the expression for the currents, integrating by parts and discarding boundary terms one has This corresponds to the entropy flow multiplied by the temperature, which constitutes the heat dissipated to the environment. This generalises to the complex case the results obtained in References [27,28].
From Equation (30), we calculate W as Here one can see that this expression corresponds to the work of the system, i.e., to the average forces ∂ m H along the "velocity"ψ m . Applying the substitution ∂ m H I = iF * m and its complex conjugate and substitutingψ m with the equation of motion, a straightforward calculation shows that Thus, the "work" is indeed the heat flowing through the mth oscillator, and therefore transported along the network [8]. On the other hand, Q is the heat dissipated to the bath, given by Equation (31).

Hamilton-Lagrange Description of a One Dimensional Continuum Ferromagnet
We consider here the dissipative dynamics of a continuum ferromagnet at zero temperature. In particular, we show how the symmetry of the system allows one to obtain the conserved quantities and the corresponding equations of motion.
The (linearised) magnetisation dynamics close to equilibrium of such a system is described by the following Schrödinger equation with complex potential [36] The stereographic variable ψ(x, t) = (M x + iM y )/2M s describes the precession of the magnetisation vector M = (M x , M y , M z ) in the x-y plane, around the z axis, see Figure 1a. Here M s is the saturation magnetisation. The precession frequency reads ω = γh ext , where γ is the gyromagnetic ratio and h ext is the applied field along the z direction. The quantity Γ = αω is the damping rate, proportional to the phenomenological damping parameter α. The chemical potential µ accounts for spin transfer torque, which can compensate the damping and leads to a steady state precession of the magnetisation [34]. The spin stiffness A is the strength of the exchange interaction. Note that in realistic cases one should consider a nonlinear damping [34,35] Γ(p) ≈ Γ 0 (1 + 2p), with p ≡ |ψ| 2 , that allows the system to have limit cycle oscillations when µ > Γ 0 . Apart from chemical potential term, Equation (35) can have more terms, accounting for temperature and additional time-dependent magnetic fields that drive the system out of equilibrium. In the present case however we shall first consider the linearised dynamics with zero temperature, since this is sufficient for our purpose to derive and illustrate the expressions for the spin currents. The more general case of a network at finite temperature will be discussed in the next section.
The Lagrangian density for our Schrödinger equation reads The equations of motion are given by the Euler-Lagrange equations with the dynamics for ψ * being given by δL δψ = 0. From the Lagrangian one can derive the following moments and finally a Legendre transform gives the following complex Hamiltonian density: From the Hamiltonian one obtains the equation of motion Equation (35) asψ = δH δiψ * , so that ψ and iψ * are conjugate variables. One can check that this equivalent to the derivation of the equations of motion using commutators and anti-commutators as described in Equation (4).
The conservation equation for the local spin wave power p ≡ |ψ| 2 is obtained from the invariance of the Lagrangian Equation (36) with respect to the global phase transformation ψ → e −iα ψ, with the corresponding infinitesimal transformation δψ ≈ −iαψ. The invariance of the Lagrangian with respect to such infinitesimal transformation yields δL δψ δψ + c.c. = 0, c.c. indicating the complex conjugate. A straightforward calculation gives then By using Equations (35), (36) and (40) one obtains the following conservation equation for the spin wave powerṗ where the spin current reads j p = 2AIm[ψ * ∂ x ψ], while Γ and µ act respectively as sink and source of excitations. We remark that this is precisely the same expression as the probability currents that appears in the Schrödinger equation of quantum mechanics. In the present case, it describes the transport of the z component of the magnetisation along the system. Indeed, one can check that j p is the same as the spin-wave current j = AM × ∇M written in terms of the stereographic variable ψ [37].

Entropy Production for a Network of Classical Spins
The finite-temperature dynamics of an ensemble of magnetic spins {M m }, n = 1, ..., M, inside a ferromagnet is described by the Landau-Lifshitz-Gilbert (LLG) equation of motion [38] with stochastic thermal baths. The LLG equation is a vector equation, and obtaining the associated FP equation in practice very cumbersome. A great simplification is obtained by re-writing the LLG equation in terms of the complex variable ψ m = m xm +im ym 1+m zm , where m = M/M s is the magnetisation vector normalised over the saturation magnetisation. In this way one obtains [34,36,39] where the force reads The first term corresponds to the applied field H z along the precession axis z of the magnetisation. The second term corresponds to the demagnetising field, while the last term models the coupling with the other spins. Note that the formulation of the coupling is completely general. In particular, such coupling can have different origins (exchange or dipolar interaction) depending on the coupling matrix A, which can be a function of the ψs.
The reversible and irreversible components of the forces read respectively From here one can check immediately that the divergence of the conservative forces, ∑ m (∂ m F R m + ∂ * m F R * m ) vanishes as it should, since ∂ m F R m = −∂ * m F R * m = i 1+α 2 (γH z + A mm ). This holds in the case where the coupling A is Hermitian, so that its diagonal is real.
The term g k n in Equation (42) is the strength of the noise, and models thermal fluctuations on site n. There are three components of the noise on each site, one per each direction of the magnetisation: Here D m = αk B µ 0 V m M s is the diffusion constant, with kB the Boltzmann constant, µ 0 the vacuum magnetic permeability and V m the elementary volume containing the magnetisation vector at site m, of the order of few nm 3 . T m is the temperature at site m. The ξ are Gaussian random variables with zero average and correlation ξ k m (t)ξ k m (t ) = δ kk mm δ(t − t ) The entropy production splits into the sum of two components, Φ = Φ 1 + Φ 2 , with and In the more general case where easy axis anisotropy and spin transfer torque are present, the LLG equation contains additional terms, but is still very similar to Equation (42) [36].

Entropy Production in the Frenkel-Kontorova Model
Let us consider the Frenkel-Kontorova (FK) model, which describes the motion of an oscillator chain sliding over a periodic potential in the presence of random fluctuations: where for simplicity we consider unit mass oscillators. Here η m the friction parameter, g the coupling strength between the oscillators, h the strength of the on-site potential and ξ m a real Gaussian random variable with zero average and variance ξ m (t)ξ m (t ) = δ mm δ(t − t ). The diffusion constant reads D m = 2η m k B . The last term b is the constant force applied to the chain to make it slide on the periodic potential. For our purposes, it is useful to introduce the "frequency" ω = 2g and rewrite the FK equation aṡ Then, one can use the complex coordinates ψ m = x m + i ω mẋ m and get From Equation (49) one has that the kinetic term becomes and finally one obtainsψ where c.c. indicates the complex conjugate. The complex FK equation can be obtained aṡ where the FK complex Hamiltonian reads To calculate the entropy production, one needs the irreversible (or dissipative) components of the force, given by F I = i∂ * m H I . It is straightforward to identify the dissipative component of the Hamiltonian as H I = −iη m |ψ m | 2 − 1 2 ψ 2 m − 1 2 ψ * 2 m . Thus the irreversible force is F I = η m (ψ m − ψ * m ). Then, applying Equations (22) and (23) gives We remark that, at variance with the DNLS, here the coupling is conservative and does not enter in the definition of entropy production [8,35]. Finally, we apply the transformations given in Equation (50) and go back to the real-valued variables: The last formula, which contains the particle kinetic energy, is consistent with what has been obtained in References [27,28] and is the dissipated power. Next, we compute the heat flow, defined as the correlation function between reversible and irreversible forces [8]: where we identify the heat flow to the correlator between neighbours oscillators. By substituting the expressions for the forces and changing coordinates to the real displacements gives j Q m+1 − j Q m = (x m+1 + x m−1 )ẋ m , which is the standard formulation of the heat flow for a chain of oscillators [40,41].

Conclusions
In summary, we have presented a general method, based on stochastic thermodynamics, to calculate entropy production and heat flows in complex-valued Langevin equations with multiplicative noise. The method is particularly useful to describe the off-equilibrium dynamics of oscillator networks for a variety of physical systems, as described by our examples. Possible research direction involves formulating the dynamics in terms of a master equation, following the discretisation of the Fokker-Planck equation proposed in References [27,28]. This should allow to formulate the irreversibility in terms of fluctuation theorems, relating the synchronisation of the oscillators to the propagating currents and the breaking of detailed balance.