Erratum: Microscopic approach to the macrodynamics of matter with broken symmetries (2020 J. Stat. Mech. 103203)

. A uniﬁed set of hydrodynamic equations describing condensed phases of matter with broken continuous symmetries is derived using a generalization of the statistical-mechanical approach based on the local equilibrium distribution. The dissipativeless and dissipative parts of the current densities and the entropy production are systematically deduced in this approach by expanding in powers of the gradients of the macroﬁelds. Green–Kubo formulas are obtained for all the linear transport coeﬃcients. The consequences of microreversibility and spatial symmetries are investigated, leading to the prediction of cross eﬀects resulting from Onsager–Casimir reciprocal relations. Crystalline solids and liquid crystals are potential examples of application. The approach is clarifying the links between the microscopic Hamiltonian dynamics and the thermodynamic and transport properties at the macroscale.

as a consequence of equation (B.4) with v a = 0 and φ ab = 0. In spite of the fact that δû cγ = ∇ c δx γ , the terms involving this quantity should be included because of equation (2.22), expressing the long-range order of matter with broken symmetries.

Introduction
The spontaneous breaking of continuous symmetries is a ubiquitous phenomenon in nature. It manifests itself in the quantum vacuum of spacetime [1][2][3], in condensed phases of matter at equilibrium [4], or in dissipative structures existing far from equilibrium [5,6]. In nonrelativistic condensed matter at equilibrium, different kinds of continuous symmetries may be broken including internal gauge symmetries in superfluids or superconductors [7][8][9][10], and spatial symmetries of translations or rotations in crystals, liquid crystals, or magnetic materials [11,12]. The breaking of continuous symmetries generates long-range order, rigidity (i.e. the possibility to drag the whole condensed phase from its boundaries), as well as Nambu-Goldstone modes [8,9] in addition to the hydrodynamic modes resulting from the five fundamental conservation laws for mass, energy, and linear momentum. All these soft modes have frequencies vanishing with their wave number, since they represent perturbations with respect to equilibrium becoming slower and slower as their wave length increases. These modes may be propagative (e.g. the sound modes) or diffusive (e.g. the heat mode). In any case, they are damped because of energy dissipation coming from thermal fluctuations at positive temperature. All these effects can be described in terms of macroscopic equations ruling the time evolution of these modes. The macroscale formulation of hydrodynamics in matter characterized by broken symmetries has been achieved, in particular, for crystalline solids and liquid crystals [13,14]. A basic issue is to deduce the macroscopic equations from the underlying microscopic dynamics of atoms and molecules composing matter. For this purpose, statistical mechanics is required, not only for matter at equilibrium, but also away from equilibrium to obtain the time-dependent properties including the transport coefficients associated with energy dissipation [15][16][17][18][19][20][21]. In the regime of relaxation towards global equilibrium, linear response theory combined with projection-operator method has been much developed to deduce the transport properties of condensed phases with broken symmetries [11,22]. This approach has also been formulated for crystalline solids in [23,24]. However, nonlinear effects often arise because of advection induced by the velocity field, as it is the case in fluid turbulence. In order to deduce the Eulerian terms in hydrodynamics, a more systematic approach consists of using local equilibrium probability distributions instead of global equilibrium distributions. This method, which has been developed since the sixties [25][26][27][28][29][30][31][32][33][34][35], provides not only the microscopic expressions for the dissipativeless fluxes of Eulerian type, but also local thermodynamics and the statistical-mechanical expression for entropy production in terms of the dissipative fluxes. Furthermore, these latter fluxes can be obtained at leading order in the gradients of the macrofields together with the transport coefficients given by Green-Kubo formulas. In this framework, the consequences of microreversibility and spatial symmetries can also be investigated. However, this method has been defined and used only for normal fluids and its generalization to a system with broken symmetry is still lacking.
Our aim is here to use this systematic approach to derive a unified set of macroscopic equations applicable to systems manifesting spatial ordering, such as crystalline solids and liquid crystals. The full derivation relies on the identification of the hydrodynamic variables and local order parameters associated with the broken continuous symmetries at the microscopic level of description. These variables obey balance equations in the form of local conservation laws, which can be obtained from the microscopic Hamiltonian dynamics, as presented in section 2. On this basis, the local equilibrium probability distribution is introduced and its time evolution can be investigated using the microscopic Hamiltonian dynamics. In this way, the entropy functional can be defined, as well as the Massieu functional given by its Legendre transform. The microscopic expression for entropy production is deduced, allowing us to identify the dissipativeless and dissipative fluxes (also called current densities) as summarized in section 3. The calculations are carried out by expanding the macrofields in powers of their gradients around local equilibrium. The local thermodynamic relations for matter with broken symmetries are obtained at leading order in the gradients in section 4. The dissipativeless current densities are determined at next order in section 5, including the extra contributions arising from the broken symmetries. In section 6, the dissipative current densities are derived and the Green-Kubo formulas are found giving all the possible transport coefficients. Even if the microscopic expressions of the quantities of interest are essential for a concrete evaluation of the dissipative coefficients, the full derivation can be performed without their explicit knowledge. The applicability to crystalline solids is discussed in section 7 and the case of liquid crystals in section 8. Finally, our concluding remarks are presented in section 9.
Notations and conventions: Latin letters a, b, c, . . . = x, y, z correspond to spatial coordinates and Greek letters α, β, γ, . . . . label the hydrodynamic variables and the order parameters. The indices i, j, k, · · · = 1, 2, . . . , N label the particles, where N is the total number of particles. The microscopic expression of a field such as a densityĉ α or a current densityĴ a c α is denoted with a hat to make the distinction with respect to its expectation value denoted without the hat. Einstein's convention of summation over repeated indices is adopted.

Hamiltonian dynamics
We consider a system of N particles of mass m, all supposed of the same species. The particle i has the position r i , the velocityṙ i = dr i /dt, and the momentum p i = mṙ i .

J. Stat. Mech. (2020) 103203
Microscopic approach to the macrodynamics of matter with broken symmetries Their microscopic dynamics is ruled by the following Hamiltonian function where V is the interaction potential energy and r ij = r i − r j the distance between the particles i and j (with i, j = 1, 2, . . . , N ). The positions and the momenta are threedimensional vectors with the components r i = (r a i ) and p i = (p a i ) (for a = x, y, z). The time evolution of this system is represented in the phase space Γ = (r i , of dimension 6N , where the trajectories Γ t = Γ t (Γ 0 ) are uniquely determined by their initial condition Γ 0 . Any observable function A(Γ) defined in phase space is thus evolving in time according to dA/dt = {A, H} is the Poisson bracket between the phase-space functions A and B. For timeindependent observable functions such that ∂ t A = 0, time integration gives A t = A(Γ t ), since dΓ t /dt = {Γ t , H}. This time evolution can be expressed in terms of the Liouvillian operator LA ≡ {A, H} according to dΓ t /dt = LΓ.
The observable quantities may be global, i.e. extensive, such as the total energy given by the Hamiltonian function, or local, i.e. intensive, such as densities. In particular, the densities of energyê, massρ, and momentumĝ = (ĝ a ) are the variables associated with five fundamental conservation laws in the system. The time variations of these densities can be expressed in terms of the divergence of a current density or flux, leading to slow modes called hydrodynamic modes.
Beyond that, there may exist other densities, such as the densities of kinetic energy and potential energy, which are not directly associated with conservation laws. Their time variations are ruled by a rate, instead of a flux divergence, so that they most often generate fast modes, called kinetic modes. Nevertheless, upon continuous symmetry breaking, some of these fast local quantities may turn into slow modes, called Nambu-Goldstone modes [8,9], which thus behave as hydrodynamic modes. According to the Goldstone theorem [9], there are as many such slow modes as continuous symmetries that are broken. Therefore, the total number of hydrodynamic variables is the sum of the five conservation laws and the number of broken continuous symmetries (one, two or three depending on the type of system in consideration, i.e. a liquid crystal or a crystalline solid). For the rest of this section, we introduce the microscopic definitions of the hydrodynamic variables and their current densities.

Conserved quantities
The Hamiltonian system ruled by equation (2.1) has several conserved quantities of fundamental origin. Since the Hamiltonian function is time independent ∂ t H = 0, this function is itself a conserved quantity because {H, H} = 0, which represents the total energy E = H. Besides, the total mass M = mN is conserved. Moreover, the Hamiltonian function is invariant under spatial translations r i → r i + a with a ∈ R 3 , so that the total momentum P = i p i is conserved. Because of the symmetry of the Hamiltonian function (2.1) under rotations O ∈ SO(3), the angular momentum L = i r i × p i is also conserved.
In order to obtain the local conservation laws, we introduce the densities associated with the total mass, energy, and linear momentum according to wheren denotes the particle density. The extensive quantities are given by integrating these densities over space: M = ρ dr, E = ê dr, and P = ĝ dr. We note that a density of angular momentum could be introduced similarly. However, the density of angular momentum is not strictly local since it is defined with respect to a coordinate origin, so that its status is different from the five aforedefined densities [11]. It is known [19] that the densities (2.2)-(2.4) obey the following local conservation equations, given in terms of the following current densities or fluxes, J a ρ (r ; Γ) ≡ĝ a (r ; Γ), is the force exerted on the particle i by the particle j, and is a uniform linear density distributed on the straight line joining the positions of the particles i and j [35]. As expected, there are five hydrodynamic variables coming from the conservation laws.

Breaking of continuous symmetries
A low enough temperature, phase transitions happen from normal fluids to liquid crystals or crystalline solids. In these new phases, continuous symmetries are broken in the structure of matter at equilibrium. In nematic liquid crystals, continuous rotational symmetry is broken by the emergence of a special orientation of the molecules, while the continuous translational symmetry is broken in crystals, where only discrete translational symmetry remains. In these phases of matter, the continuous symmetries of uniform and isotropic normal fluids are thus broken. Such phenomena are not described by Gibbsian statistical distributions based on the Hamiltonian function (2.1) since this latter is invariant under continuous translations and rotations. Therefore, an external potential energy should be added to the Hamiltonian in order to break explicitly the symmetries, The symmetric Hamiltonian function (2.1) is recovered in the limit → 0. At the inverse temperature β = (k B T ) −1 and chemical potential μ, the equilibrium phase of matter can be described by the following Gibbsian grand canonical probability distribution p eq (Γ) = 1 ΞΔΓ e −β(H −μM) , (2.14) expressed in terms of the total mass M = mN and depending on the random variables Γ = (Γ N , N ) since the particle number N is also a random variable. In equation (2.14), Ξ is the partition function given by the normalization condition and ΔΓ = h 3N is the elementary phase-space volume, where h is Planck's constant. In the following, Boltzmann's constant is set equal to unity, k B = 1, except if it is explicitly written. The statistical average with respect to the probability distribution is denoted Because of the external potential V (ext) , continuous symmetries are explicitly broken. For instance in crystals, the particle density can now have an equilibrium mean value n eq (r) = n(r ; Γ) eq that is a periodic function in space, which is invariant under the discrete crystalline group of symmetry, but no longer under the continuous group of three-dimensional translations and rotations.
We note that the system can undergo a phenomenon of spontaneous symmetry breaking without the help of the external potential, i.e. for the Hamiltonian function (2.13) in the limit → 0. For = 0, all the continuous symmetries are restored for the probability distribution (2.14). However, the system manifests long-range order. For instance, in crystals, the particle density is periodic with respect to the center of mass of the crystal. Accordingly, the symmetric Hamiltonian function (2.1) should be split into the part describing the motion of the center of mass and the other part ruling the dynamics in the frame moving with the center of mass. This latter part is no longer symmetric under continuous translations and can itself define a grand canonical probability distribution with the broken symmetry. This reasoning shows that including the center of mass in the Hamiltonian function restores the continuous symmetry in Gibbsian equilibrium probability distributions. Similarly, in nematic liquid crystals, the part of the Hamiltonian ruling the rotation of the system around the orientation selected by spontaneous symmetry breaking should be separated from the rest of the Hamiltonian function in order to define equilibrium probability distributions describing the properties of the phase with broken symmetry.
Phases with broken symmetries can be characterized by local order parameters denotedx α , where the index α runs over the subset of variables originating from symmetry breaking. The decay of these variables is given by a rateĴ x α such that We may also introduce the variablesû bα ≡ ∇ bxα that obey the following equations similar to the local conservation equations (2.5)-(2.7), where the corresponding current density is defined aŝ The microscopic expression forx α depends on the phase in consideration, as discussed in sections 7 and 8 for crystalline solids and liquid crystals. However, it is possible to proceed with a general derivation of the macroscopic equations without using any explicit microscopic expression forx α . Because of long-range order, the equilibrium correlation functions of the variableŝ x α decay in space as 19) or for their Fourier transforms, In this regard, these order parameters have a singular behavior, contrary to regular conserved quantities that have correlations of medium or short range. Nevertheless, the gradientsû aα ≡ ∇ axα of the order parameters have medium-or short-ranged correlations because δû aα (q)δû bβ (−q) eq ∼ q 0 , (2.22) so that they are regular as for the conserved variables.

Nambu-Goldstone modes
As a consequence of the emergence of long-range order, there exist Nambu-Goldstone modes behaving as the conserved modes with vanishing dispersion relations for long enough wave length [6,8,9]. The equations ruling all the conserved and kinetic modes ψ β = ψ β could be written in the following form Let us suppose that there exists an equilibrium solution Now, we may consider small perturbations of different kinds with respect to this equilibrium solution. On the one hand, an additive perturbation may be considered giving solutions of the following form, Substituting into equation (2.23) and linearizing, we obtain the following evolution equations, Such a mode is decaying exponentially in time with a rate that is not vanishing for long enough wave length, because ∂F β /∂ψ γ eq defines a matrix with non-vanishing elements in general. On the other hand, another perturbation may be considered where the order parameters are locally varying in space and time as For this solution, we have that Now, substituting into equation (2.23) and linearizing, we find where the term with the non-vanishing coefficients ∂F β /∂ψ γ eq no longer appears. If we multiply equation (2.30) by ∂ψ β eq /∂x γ , sum over β, integrate over the volume V of the system to average out the local variations of the symmetry-breaking equilibrium solution, and relabel the quantities, we get In equation (2.31), the dots denote terms with higher spatial derivatives for x β (r, t).
Since the solution (2.28) is breaking continuous symmetries, we have that (∂ψ β eq /∂x α ) = 0 and thus the matrix N = (N αβ ) is non vanishing and can be inverted to define the quantity V= N −1 · M. If we suppose that the local order parameters behave as x (r, t) ∼ exp(ıq · r − ıωt), the frequency ω is related to the wave vector q according to showing that the dispersion relations of these solutions are vanishing with the wave number as lim q→0 ω(q) = 0, as in the case of locally conserved quantities. The eigenvalues of the matrix V · q are giving the propagation speed of the modes. The mode is diffusive if its propagation speed is equal to zero. The existence of the Nambu-Goldstone modes is thus a consequence of continuous symmetry breaking. There are as many such modes as components of the vector x = (x α ), i.e. as continuous symmetries that are broken, which is the statement of the Goldstone theorem [6,9]. In order to investigate the effects of the Nambu-Goldstone modes, we should thus consider equation (2.17) on the same footing as the equations (2.5)-(2.7) for the locally conserved quantities.

Nonequilibrium statistical mechanics
In this section, the formalism of the systematic approach previously developed for normal fluids in [25][26][27][28][29][30][31][32][33][34][35] is introduced and extended to phases with broken continuous symmetries. The identities necessary to derive the set of hydrodynamic equations are explicitly derived. We also clarify the interpretation of the approach in terms of the production of entropy.

Time evolution
First of all, let us adopt systematic notations. On the one hand, equations (2.5)-(2.7) for the locally conserved quantities, as well equation (2.17) for the gradients of the order parameters can all be written as are respectively the densities and the corresponding current densities or fluxes. At time t, the densities are given in terms of the Liouvillian operator L or the trajectories Γ t bŷ with similar expressions for the current densities.
On the other hand, any phase-space probability density p t (Γ) at time t is given by in terms of the initial probability density p 0 (Γ) and the reversed trajectory Γ −t going from the current phase-space point Γ back to the initial conditions Γ 0 of the trajectory. Consequently, the macroscopic densities can be obtained by taking the mean value of the time-independent densities over the time-evolved probability distribution p t (Γ) or, equivalently, the mean values of the time-dependent densities over the initial probability distribution p 0 (Γ), because of Liouville's theorem dΓ 0 = dΓ t , the expectation value with respect to p t (Γ) being denoted as · t . Similar results hold for the current densities and other fields.

Local equilibrium distribution
The key assumption of the formalism is the initial condition being the local equilibrium distribution where λ = (λ α ) = (λ c α ) are inhomogeneous fields conjugated to the density fieldsĉ = (ĉ α ), and the asterisk * corresponds to the integration over space The normalization condition (2.15) for the local equilibrium distribution (3.6) gives the functional The expectation value with respect to the local equilibrium distribution (3.6) is denoted by · leq,λ . In this formalism, the expectation values of the densities can be obtained by taking the functional derivative of the functional (3.8) with respect to the conjugated fields as follows, The entropy is defined as leading for the local equilibrium distribution (3.6) to the entropy functional which is the Legendre transform of the previously introduced functional (3.8). The conjugated fields are thus given by the following functional derivatives, (3.12) this relation being called the second identity in [35]. Vice versa, the Legendre transform of the entropy functional (3.11) gives back the functional (3.8).

Time evolution of the local equilibrium distribution
The basic idea of the formalism is that the time-evolved probability density (3.4) should remain close to the local equilibrium distribution (3.6) with the conjugated fields λ t considered at time t. In this regard, these latter should be determined by the conditions ĉ α (r ; Γ) t = ĉ α (r ; Γ) leq,λ t ≡ c α (r, t), (3.13) according to which the expectation values (3.5) of the densities with respect to the probability distribution p t (Γ) are equal to their expectation values with respect to the local equilibrium distribution with the conjugated fields λ t at time t. The conditions (3.13) are thus defining the macroscopic densities c α (r, t), which are given by the functional derivatives (3.9) of the functional (3.8) with respect to the conjugated fields λ = λ t at time t. Now, we consider the time evolution of the probability density starting from the initial condition given by the local equilibrium distribution (3.6) with λ = λ 0 , Since the normalization of this probability distribution should be preserved during the time evolution, we must have that (d/dt) p t (Γ)dΓ = 0. The calculation using equation (3.1) leads to the following relation, which should hold for any conjugated field λ 0 that may thus be replaced by λ to get This relation has been obtained notably in [31,35] and is called the first identity in [35].
Remarkably, we have that with the quantity as shown in [25,35]. Therefore, the expectation value of any observable A(Γ) with respect to the time-evolved probability distribution p t (Γ) is thus given in terms of the expectation value with respect to the local equilibrium distribution according to which is called the third identity in [35]. In particular, the conditions (3.13) are equivalent to the following relations,

Entropy production and dissipative current densities
The identity (3.18) is a universal relation, which is reminiscent of the nonequilibrium work and integral fluctuation theorems [35][36][37][38][39][40][41][42][43]. Choosing A(Γ) = e −Σ t (Γ) in this third identity gives [35] e −Σ t (Γ) t = 1, (3.20) which implies by Jensen's inequality [44] that In open systems, the entropy S changes in time due to the exchange d e S with the environment and its production d i S inside the system: dS = d e S + d i S. Since the system is here isolated, there is no exchange with the environment d e S = 0, so that the change in time of the entropy is equal to the entropy production dS = d i S. In this regard, the result (3.21) may be interpreted as the non-negativity of the entropy production, in agreement with the second law of thermodynamics. Using equations (3.8) and (3.9), we have that According to the definition (3.11) of the entropy functional, the relation (3.13), equation (3.1), and integrations by parts, the entropy production is thus given by (3.24) Now, using the identity (3.18) with A taken asĴ a c α , the expectation values of the current densities with respect to the phase-space probability distribution (3.4) can be decomposed as [25] and The entropy production is thus given by where the first term vanishes because of the identity (3.15). This term is thus expressing the conservation of entropy in adiabatic (isoentropic) processes induced by the dissipativeless current densities defined by equation (3.26). The second term in equation (3.28) is in general non vanishing and related to the production of entropy, leading to the definition of the dissipative current densities by equation (3.27). Accordingly, the entropy production can be expressed as in terms of the dissipative current densities (3.27) and the gradients of the conjugated fields, which play the role of thermodynamic forces, also called the affinities, which is in accordance with macroscopic nonequilibrium thermodynamics [45][46][47][48][49]. As we will show explicitly below, the three identities Let us now proceed with the explicit computation for generic systems with spatial ordering. We will then show in sections 7 and 8 the applicability of the general approach to specific examples such as crystalline solids and liquid crystals.

Local thermodynamics
The systematic approach establishes a clear connection between microscopic expressions and phenomenological quantities. In this section, the local thermodynamics for matter with broken symmetries is obtained.

The local Euler, Gibbs, and Gibbs-Duhem relations
The identity (3.15) shows that the dissipativeless current densities (3.26) are leaving the entropy constant in time. All the processes involved by the dissipativeless current densities may thus be considered as reversible (i.e. adiabatic or isoentropic). The approach of nonequilibrium statistical mechanics based on the local equilibrium distribution (3.6) is therefore providing the local thermodynamic relations in every element of matter, as shown here below.
At leading order in the expansion in the gradients, the functionals giving the local relations up to terms of second order in the gradients. In equation (4.3), the first relation can be identified as the local Euler relation, the second as the local Gibbs relation for the entropy density, and the third as the associated Gibbs-Duhem relation. At this stage, a comparison becomes possible with previous works on the thermodynamics of matter with broken symmetry [13,14], in which the relevant thermodynamic relations are given in the laboratory frame by Eulerrelation: e = T s + μρ + v a g a + φ aα u aα − p, (4.4) Gibbs relation: de = T ds + μdρ + v a dg a + φ aα du aα , (4.5) Gibbs-Duhem relation: dp = sdT + ρdμ + g a dv a + u aα dφ aα , (4.6) where e ≡ E/V is the mean energy density, T the temperature, s ≡ S/V the entropy density, μ the chemical potential, ρ ≡ M/V the mean mass density, v a the velocity, g a = ρv a the mean momentum density, u aα the gradients of the mean order parameters x α , φ aα the fields thermodynamically conjugated to u aα , and p the hydrostatic pressure.

The conjugated fields
The Euler relation (4.4) can be written to give the entropy density instead of the energy density. After substitution into the entropy functional (3.11) and taking the functional derivatives (3.12), the conjugated fields are obtained as in terms of the local inverse temperature β = (k B T ) −1 . Furthermore, the comparison between the phenomenological Euler relation (4.4) and the theoretical expression given by the third relation in equation (4.3) allows us to identify the Legendre transform of the entropy density as the local thermodynamic potential ω = βp, proportional to the hydrostatic pressure p and referred to as a Massieu function (here, per unit volume) [49]. Now, the Gibbs relation (4.5) written for the entropy density gives the relations

Dissipativeless current densities
Expression (3.29) for the entropy production shows that entropy is conserved if the dissipative parts of the current densities are vanishing. This is the case for the leading parts of the current densities giving the dissipativeless current densities (3.26). In normal fluids, these parts lead to Euler's equations. Our purpose is now to obtain their expressions in phases with broken symmetries. With this purpose, we first consider Galilean transformations from the laboratory frame to the frame moving at the velocity v of the system. In this way, the expectation values of the current densities with respect to the local equilibrium distribution (3.6) can be calculated according to the definition (3.26).
In phases with broken symmetries, extra contributions are expected, which need to be determined, in particular, using the microscopic expressions of the corresponding decay ratesĴ x α introduced in equation (2.16). These microscopic expressions will be given in sections 7 and 8 for crystals and liquid crystals. However, the most general form of their expectation value is known on the basis of time-reversal symmetry [13]. Once this general form is fixed, the contributions of broken symmetry to the dissipativeless current densities of momentum and energy can be obtained using the identity (3.15) expressing the conservation of entropy, which finds its origin in the adiabaticity of the reversible processes ruled by the dissipativeless current densities (3.26).

Galilean transformations
The particle momenta p i in the laboratory frame and those p i0 in the frame moving at the local velocity v(r) are related to each other according to p i = p i0 + mv(r i ) by the principle of Galilean relativity. Carrying out this change of variables in the microscopic expressions for the energy densityê, the momentum densityĝ a , and their current densities, we find

2)
J a e =Ĵ a e0 +ê 0 v a +Ĵ a where the quantities with the subscript 0 are those with the momenta p i0 replacing p i , and Similar expressions hold for the current densitiesĴ a u bα of the gradients of order parameters.
In equation (5.3), the contribution (5.5) would vanish if the velocity field v was uniform. Furthermore, we note thatΔ a goes as the square of gradients, in particular, the square of the gradients of the velocity field, so that this term can be dropped if the corrections of O(∇ 2 ) are neglected. Now, the terms with odd powers of the momenta p i0 are vanishing after averaging over local equilibrium in the frame moving with the element of matter. Consequently, we get ĝ a 0 leq = 0, Ĵ a e0 leq = 0, and ρ leq = ρ. Moreover, the internal energy density is given by ê 0 leq = e 0 . The velocity field is thus defined as the ratio of the momentum to the mass densities: v a ≡ ĝ a leq / ρ leq . As aforementioned, we need to consider the possibility that symmetry breaking may introduce extra contributions in the local equilibrium expectation values of the currents densities of energy, momentum, and the gradients u aα of the order parameter, so that we obtain the following forms for these expectation values, where e = e 0 + ρv 2 /2, the hydrostatic pressure p is separated from the contributions J a e | BS andJ a g b | BS of broken symmetries, and J a u bα BS = δ abJ x α . (5.11) Since the order parametersx α are usually even under time-reversal symmetry, their ratê J x α should be odd. Therefore, the most general form of the dissipativeless mean rates is given byJ where the coefficients A aα and B abα are deduced from the microscopic dynamics in sections 7 and 8, instead of being postulated in the macroscopic phenomenological theory of [13]. Here, we note that B abα = 0 in crystals and A aα = 0 in nematic liquid crystals.
In equations (5.6)-(5.10), the terms vanishing with the velocity are due to the advection of the corresponding quantity by the motion of the element of matter at the velocity v = (v a ).
The next issue is to determine the contributionsJ a e | BS andJ a g b | BS to the dissipativeless current densities of energy and momentum by using the identity (3.15).

Dissipativeless current densities of energy and momentum
The dissipativeless currents are satisfying the identity (3.15). As shown in equation (3.28), this identity is equivalent to the requirement that the dissipativeless terms conserve the entropy. This identity can be used to deduce the energy and momentum current densities from the expression (5.12), as follows.
Using the leading-order contributions of the conjugate fields λ α , derived in equations (4.7)-(4.10), and writing the dissipativeless currents asJ a c α ≡ Ĵ a c α leq,λ , the identity (3.15) is giving up to higher-order corrections. We note thatJ a ρ = ρv a because of equation (5.7). Using the Gibbs-Duhem relation (4.18) with the differential d replaced by the gradient ∇ a , we get (5.14) Besides, the scalar product of the term ∇ a (βp) with the velocity v a can be transformed according to after using equation (5.11) with the mean decay rate of the order parameter given by equation (5.12). In order to satisfy this identity, we first require that ∇ a v bJ a g b | BS + ∇ a φ aαJ x α can be expressed as a divergence, which leads tō and the conditions Now, the divergence ∇ a (φ aα A bα v b ) resulting from the last two terms in equation (5.18) is combined with the factor β in front of them to give an extra contribution to the first term in ∇ a β, so that we finally obtain The dissipativeless parts of the momentum and energy current densities are thus determined as with the reversible stress tensor defined by These results are consistent with those obtained in [13,14].
To summarize, the dissipativeless current densities are obtained from the condition that they should conserve the entropy, which is equivalent to the identity (3.15) of the formalism. Galilean invariance and the behavior of the variables under time-reversal symmetry provide the identification of the most general form of the dissipativeless current densities. This approach, based on the conservation of entropy, does not require the knowledge of the explicit microscopic expressions for the variablesx α .

Dissipativeless local conservation equations
Now, the expressions obtained for the dissipativeless current densities can be substituted back into the local conservation equation (3.30), here neglecting the dissipative part of the current densities. We introduce the total time derivative of any field f(r, t) along a streamline of matter advected by the velocity v = (v a ) according to Expanding the Eulerian local conservation equations, we obtain their following Lagrangian forms, 29) in terms of the reversible stress tensor (5.24) and up to higher-order corrections. These equations of motion thus describe adiabatic processes leaving constant the entropy.

Dissipativeless equations for the conjugated fields
Since the conjugate fields β, βμ 0 , and βφ aα are functions of the variables e 0 , ρ, and u aα , their time evolution can be obtained from the dissipativeless equations (5.26)-(5.29).
As shown in detail in appendix A using the Maxwell relations (4.14)-(4.17), they obey the following equations, which is obtained from equation (5.28) with the reversible stress tensor (5.24), using equation (4.20) and the fact that u aα = ∇ a x α = O(∇) by definition.

Dissipative current densities
The systematic approach makes it possible to identify the dissipative current densities and to provide a microscopic basis for the transport coefficients in terms of Green-Kubo expressions. In this section, we derive these quantities explicitly, using the identity (3.18) of the formalism. The full set of hydrodynamic equations is obtained in this way.

Heat current density
Having identified the dissipativeless current densities in section 5, we turn to the derivation of the dissipative current densities defined by equation (3.27) and contributing to the entropy production (3.29). This latter can be expressed in terms of the gradients of the conjugated field given by equations (4.7)-(4.10) at leading order in the gradients. Using the fact that the dissipative current density of mass is equal to zero J a ρ = 0, the relation J a u bα = δ ab J x α , and the expansions ∇ a (βf ) = f ∇ a β + β∇ a f for any field f, equation (3.29) can be written as in terms of the heat current density defined as This result shows that the dissipative current density of energy would reduce to the heat current density under the conditions where v a = 0 and φ aα = 0, i.e. in the frame moving with matter and in the absence of the contribution from φ aα to the stress.

Deduction at leading order
Let us proceed with the calculation of the dissipative current densities (3.27). At the leading order in the gradient expansion, they are given by Therefore, their derivation relies on the evaluation of the quantity (3.17), which can be transformed into in terms of the following deviations of the densities and currents densities with respect to their local equilibrium expectation values, The detailed calculations to get the integrand of equation (6.4) are carried out in appendix B. Under the conditions v a = 0 and φ aα = 0 where the dissipative current density of energy coincides with the heat current density according to equation (6.2), we find at first order in the gradients that with the following definitions, δĴ a e ≡ δĴ a e − ρ −1 (e 0 + p)δĝ a , (6.8) δĴ a g b ≡ δĴ a g b + ∂σ ab ∂e 0 ρ,u δê + ∂σ ab ∂ρ e 0 ,u δρ, (6.9) Substituting equation (6.7) back into equation (6.4), the dissipative current densities (6.3) are given at leading order by up to higher-order corrections. Assuming that the characteristic spatiotemporal scales of the conjugated fields λ are much larger than the correlation length and time of the current densities (6.8)-(6.10), we can replace ∇ c λ γ (r , τ ) by ∇ c λ γ (r, t) in the previous equation to find dτ dr δĴ a c α (r, 0)δĴ x γ (r , τ − t) leq,t . (6.12) Since the conjugated fields evolve on a longer time scale than the correlation time of the current densities, the local equilibrium distribution at time t may be considered as the equilibrium distribution at the local values of the conjugated fields in the frame where v a = 0 and φ aα = 0, i.e. for given local values of temperature and chemical potential. Since the equilibrium distribution is stationary, we have that δâ(r, 0) δb(r , τ − t) = δâ(r, t − τ ) δb(r , 0) . Replacing t − τ by τ , the integral over time becomes t 0 dτ δâ(r, 0) δb(r , τ − t) = t 0 dτ δâ(r, τ ) δb(r , 0) , where the limit t → ∞ can be taken since the time scale t of the conjugated fields is longer than the correlation time of the current densities.
Furthermore, the material properties over spatial scales larger than the microscopic structure of the phase can be defined by averaging over space. The microscopic global currents are thus introduced aŝ Microscopic approach to the macrodynamics of matter with broken symmetries δĴ a g b (t) = δĴ a g b (t) + ∂σ ab ∂e 0 ρ,u δÊ(t) + ∂σ ab ∂ρ e 0 ,u δM (t), (6.15) where δP a = δĝ a dr, δÊ = δê dr, and δM = δρ dr are the deviations in the total momentum, energy, and mass with respect to their local equilibrium value. Since the total momentum, energy, and mass are conserved, they are constants of motion (i.e. they do not fluctuate in time, but they are still random variables because the initial conditions in phase space may give different values to these constants of motion, for instance, in the grand canonical ensemble of distribution (2.14) with = 0). They may thus be added into the equilibrium time correlation functions because δĴ a c α eq = 0. Consequently, the unprime quantities can be replaced by the prime ones in the time correlation functions. Therefore, equation (6.12) becomes which holds for the dissipative current densities (J a c α ) = (J a e , J a g b , J x β ) under the conditions v a = 0 and φ aα = 0, where the dissipative energy current density is equal to the heat current density J a e = J a q .

Green-Kubo formulas for the transport coefficients
As a consequence, equation (6.17) leads to the following expressions for the dissipative current densities of heat, momentum, and order parameters, with the transport coefficients given by the following Green-Kubo formulas, dt δĴ a e (t)δĴ x α (0) eq , (6.23)

J. Stat. Mech. (2020) 103203
Microscopic approach to the macrodynamics of matter with broken symmetries Here, the limit V → ∞ is taken at constant chemical potential, in order to define the bulk transport properties in arbitrarily large systems, removing in this way possible finite-size effects encountered in molecular dynamics simulation [50]. We note that, in isotropic phases of matter, the rank-three tensors χ abc and θ abα are vanishing, so that they are expected to be small in general.

Time-reversal symmetry
Since the Hamiltonian function (2.1) has the symmetry H(ΘΓ) = H(Γ) under the timereversal transformation Θ(r i , p i ) = (r i , −p i ), the equilibrium probability distribution (2.14) is also symmetric and we have the Onsager-Casimir reciprocal relations [51][52][53] ∞ 0 dt δĴ α (t)δĴ β (0) eq = α β ∞ 0 dt δĴ β (t)δĴ α (0) eq , (6.27) where α = ±1 if the current δĴ α is even or odd under time reversal (and there is no Einstein summation here). Sinceê andx α are even, their currentsĴ a e andĴ x α are odd. Besides, the currentsĴ a g b are even becauseĝ b is odd. Consequently, we have the Onsager-Casimir reciprocal relations κ ab = κ ba , η abcd = η cdab , ξ aα = ξ αa , and ζ αβ = ζ βα . This latter relation explains how the coefficient in front of ∇ b φ bα in equation (6.18) is the same as the one in front of ∇ a T /T in equation (6.20). Moreover, we have the Onsager-Casimir reciprocal relations χ abc = −χ bca and θ abα = −θ αab , which explains the changes of sign in front of ∇ c T /T in equation (6.19) and in front of ∇ a v b in equation (6.20) with respect to the corresponding terms with the same coefficients. The coefficients χ abc generate a coupling between momentum transport and temperature gradients in a similar way as at the interface between two phases [54], which is the mechanism inducing the phenomenon of thermophoresis [55].

Entropy production
Because of the antisymmetry of the coefficients χ abc and θ abα , the associated terms do not contribute to entropy production and thus dissipation. In order to confirm this result, the entropy production (6.1) is calculated from equations (6.18)-(6.20), giving where, indeed, the terms with the coefficients χ abc and θ abα do not appear. In this regard, these terms may be considered as dissipativeless contributions to the current densities. In particular, the terms with the coefficients θ abα in equation (6.20) are similar to those with the coefficients B abα in equation (5.12). The non-negativity of the entropy production results in particular from the conditions η abab 0, κ aa 0, ζ αα 0, and κ aa ζ αα (ξ aα ) 2 /T [47].

Macroscopic equations
Finally, the macroscopic equations read ∂ t ρ + ∇ a (ρv a ) = 0, (6.29) ∂ t e + ∇ a (J a e + J a e ) = 0, (6.30) with the dissipativeless and dissipative current densities and rates given by equation (5.12), (5.22), (5.23), (6.2) and (6.18)-(6.20), as obtained with the expansion in powers of the gradients. The system of macroscopic equations can be closed using the thermodynamic relations. Keeping the terms that are linear in the gradients and the velocity, we find ∂ t ρ −ρ ∇ a v a , (6.33) which reduces to equation (6.1) of [13] in the case where χ abc = 0 and θ abα = 0. The coefficients κ ab can be interpreted as the heat conductivities and η abcd as the viscosities. Thanks to the derivation of these macroscopic equations from the microscopic Hamiltonian dynamics, the coefficients now have a microscopic foundation since they are given by the Green-Kubo formulas (6.21)-(6.26). Moreover, this derivation reveals the possibility that the coefficients χ abc and θ abα can be non vanishing. Finally, note that the microscopic expressions for the order parameters and its associated current densities were not used in the derivation. However, these expressions are essential in practice and in order to evaluate the transport coefficients with the Green-Kubo formulas. To show the validity of the approach, let us consider here below two explicit examples, namely, the crystalline solids and the liquid crystals.

Crystalline solids
The method applies in particular to crystals where the continuous symmetry under the group of spatial translations is broken into the discrete symmetry of one of the 230 crystallographic space groups. Accordingly, crystals have eight hydrodynamic modes with dispersion relations vanishing with the wave number according to the Goldstone theorem. These modes are the two longitudinal sound modes, the four transverse sound modes, the mode of heat conduction, and the mode of vacancy diffusion. These modes are damped because of energy dissipation and the present results give the Green-Kubo formulas for the coefficients ruling their damping. A comprehensive analysis of the hydrodynamic modes in a crystalline solid is however beyond the scope of this work. The aim of this section is to show that the systematic approach presented here is consistent with and complements earlier results on the hydrodynamic equations in crystals [13,14,23,24].

Order parameters
For crystals, the variablesx α associated with continuous symmetry breaking are the componentsû a of the displacement field. Since the broken symmetry is a spatial symmetry in three dimensions, we can replace the Greek index α by a Latin index a. The microscopic expression for the displacement field is known [23,24,56,57]. In cubic crystals, it readŝ u a (r ; Γ) ≡ − 1 N BZ dk (2π) 3 e ık·r dr e −ık·r ∂n eq (r ) ∂r an (r ; Γ), (7.1) given in terms of the particle densityn =ρ/m, the equilibrium particle density n eq (r), an integral over the Brillouin zone (BZ), and a normalization factor N . The equilibrium particle density is a periodic function of space, which has the symmetry of the crystallographic group, so that it can be expanded over the reciprocal lattice as follows, n eq (r) = G n eq,G e ıG·r (7.2) with G being a reciprocal lattice vector. In a uniform phase where ∂n eq /∂r a = 0, the displacement field would be vanishing, as expected for order parameters. The associated rate defined by equation ( In crystals, the microscopic strain tensor is defined as the symmetric rank-two tensor the possibility of a non-vanishing piezoelectric tensor, which is also of rank three [58]. Nevertheless, the cross effects described by the coefficients χ abc and θ abc often play a negligible role. Finally, the hydrodynamic modes can be obtained from the macroscopic equations (6.33)-(6.36), which are consistent with earlier results [13,14,23,24]. The present approach has the further advantage of establishing relationships between the microscopic quantities and thermodynamics.

Liquid crystals
Liquid crystals are composed of nonspherical molecules, which interact with different types of intermolecular forces. Rotational or translational symmetries may be broken in liquid crystals, because of the emergence of a privileged orientation, e.g. in nematics, or two-dimensional columnar order, e.g. in some phases of discotic liquid crystals [12,59].
For rotational symmetry breaking, we note that the total angular momentum can be decomposed as L = L 0 + k L k in terms of the angular momentum L 0 with respect to the origin of the laboratory frame and the angular momenta L k of the molecules {k} with respect to their center of mass (or any other property). Such angular momenta L k allow us to carry out local rotations in the system. The corresponding Nambu-Goldstone modes can be defined as soft modes associated with such local rotations, as discussed in subsection 2.4.

Order parameters
For apolar nematogens (i.e. nematic molecules), an external electric field E a (r) will explicitly break the rotation symmetry. In this case, the total external potential energy in the Hamiltonian function (2.13) is given by V (ext) tot = − 1 2 E a (r)q ab (r)E b (r)dr (8.1) with the local traceless polarizability tensorq ab (r). In general, this local order parameter can be taken as the quadrupolar contribution to the density of some property associated with the nematogens, q ab (r ; Γ) = k i∈k where q i is the relevant property attached to the atom i ∈ k in the molecule k, and r k = (r a k ) is a position at the center of the molecule k [12]. The order parameter obeys the evolution equation ∂ tq ab (r ; Γ) = −Ĵ q ab (r ; Γ) (8.3)

Green-Kubo formulas
Again, the transport coefficients will be given by the Green-Kubo formulas (6.21)-(6.26). The coefficients χ abc and θ abα describing the dissipativeless cross effects may be expected to be equal to zero in most phases of liquid crystals. However, since piezoelectricity also exists in some liquid crystals [60,61], it is possible that the cross effects due to χ abc and θ abα exist here as well, although being small.

Conclusion
In this paper, we have shown that the macroscopic equations ruling the time evolution of matter with broken continuous symmetries can be derived in a unified microscopic approach based on the local equilibrium distribution, extending to crystalline solids and liquid crystals results previously obtained for normal fluids.
In the presence of broken continuous symmetries, the description should be extended to include the microscopic expressions of the order parameters beside the microscopic densities of mass, energy, and momentum, which are locally conserved. The time evolution of these variables is generated at the fundamental level of description by the underlying Hamiltonian microdynamics. The manifestation of spontaneous symmetry breaking is the emergence of as many Nambu-Goldstone modes as there are broken continuous symmetries. Those modes have frequencies ω(q) vanishing with their wave number q, as for the hydrodynamic modes associated with the local conservation laws. All these modes are damped because of their interaction with the other degrees of freedom. This damping is determined by the transport coefficients responsible for energy dissipation and entropy production.
Here, we have deduced these properties using the microscopic approach based on the local equilibrium distribution and its time evolution ruled by the Hamiltonian microdynamics. For every density, we have systematically obtained the dissipativeless and dissipative parts of the corresponding current density using expansions in powers of the gradients of the macrofields. With this approach, the dissipativeless part of each current density can be inferred including their nonlinear dependence on the velocity field, and their dissipative part can be identified in direct relation with entropy production. In this way, Green-Kubo formulas have been derived for all the possible transport coefficients, giving them a microscopic foundation. The unprecedented achievement is that the approach is systematically providing the microscopic expressions for all the possible macroscopic thermodynamic and transport properties at leading order in the gradients for condensed phases with spatial ordering.
The symmetries under time reversal and the point group of unbroken spatial rotations have been used to reduce the number of transport coefficients. Time-reversal symmetry leads to the Onsager-Casimir reciprocal relations between the coefficients coupling different transport processes. These latter may contribute to entropy production if their coupling is symmetric under time reversal, or be dissipativeless if their coupling is antisymmetric. Among the former, there are the heat conductivities, the viscosities, the coefficients associated with the order parameters, and the cross effects between heat transport and the order parameters. Among the latter, cross effects are furthermore predicted between heat and momentum transport, and between momentum transport and the order parameters. In isotropic phases of matter, these latter cross effects are absent according to Curie's principle based on the full continuous group of three-dimensional rotations. However, such cross effects become possible in the presence of anisotropies, such as planar interfaces between two isotropic phases [54]. Here, we have found that, in some classes of crystalline solids and liquid crystals, properties may be coupled together with a rank-three tensor of non-vanishing coefficients because of anisotropy, as it is the case for piezoelectricity [58]. The microscopic expressions of all these transport coefficients are here given by Green-Kubo formulas. For crystalline solids, previously obtained results are recovered [23,24], showing the validity of the methods presented in this work. Moreover, the unified approach also provides the microscopic expressions of transport properties in liquid crystals, depending on their broken continuous symmetries.
We note that the approach can be extended to quantum systems [17,18,[26][27][28]30]. In future work, we hope to use the methods developed in the present paper, in particular, to investigate the hydrodynamic and Nambu-Goldstone modes and their damping in crystalline solids and liquid crystals.