text:
# Baryogenesis via Density Fluctuations with a Second Order Electroweak Phase Transition

## I Introduction

Baryogenesis at the electroweak phase transition in the universe is one of the most attractive models of solving the matter-antimatter problem of the universe . As compared to other models such as those based on Grand Unified Theories (GUT), primordial black holes, etc., much of the physics of electroweak theory is well understood and is accessible to laboratory experiments. Unfortunately, it is by now generally accepted that it does not seem possible to get sufficient baryon asymmetry entirely within the standard model. The requirements of a strong first order phase transition, large CP violation, and suppression of baryon violating interactions in the broken phase (so that produced asymmetry survives), seem only possible if one considers extensions of the standard model. It is also known that in the case of a second order phase transition, when the temperature decreases due to the expansion of the universe, the resulting baryon asymmetry is 8 to 9 orders of magnitude below the required value . This happens essentially because the cooling rate of the universe, and hence the rate of change of the vacuum expectation value of the Higgs field (say, in two Higgs extension of the Standard Model), is many orders of magnitude smaller than the electroweak scale . There have been some proposals such as those based on topological defects  (e.g., utilizing the symmetric core of the defect), and reheating of local regions due to evaporation of primordial black holes , where it has been argued that sufficient baryon asymmetry may arise even with a second order phase transition.

An important aspect of phase transitions in the universe, which has not received much attention, is the fact that the universe is necessarily inhomogeneous. There are density fluctuations present in the universe, which eventually lead to the formation of structure we see today. Density fluctuations, and hence temperature fluctuations, even if they are of small magnitude, can affect the phase transition dynamics in crucial ways. There has been some discussion of the effects of inhomogeneities on the dynamics of a first order phase transition, in the context of quark-hadron transition in the universe . For example, Christiansen and Madsen have discussed  heterogeneous nucleation of hadronic bubbles due to presence of impurities. It is mentioned in ref. that possible sources of such impurities could be primordial black holes, cosmic strings, magnetic monopoles, or relic fluctuations from the electroweak scale. Hadronic bubbles are expected to nucleate at these impurities with enhanced rates. Recently, Ignatius and Schwarz have proposed  that the presence of density fluctuations (those arising from inflation) at quark-hadron transition will lead to splitting of the region in hot and cold regions with cold regions converting to hadronic phase first. Baryons will then get trapped in the (initially) hotter regions. Estimates of sizes and separations of such density fluctuations were made in ref. using COBE measurements  of the temperature fluctuations in the cosmic microwave background radiation (CMBR). In an earlier work , we have considered the effect of cosmic string induced density fluctuations on quark-hadron transition and have shown that it can lead to formation of planar regions of baryon inhomogeneity which may have important effects on nucleosynthesis.

There are many possible sources of density fluctuations in the early universe such as inflation and cosmic defects (strings, monopoles, textures, domain walls). There has been extensive study of density fluctuations generated by cosmic strings from the point of view of structure formation . Though recent measurements of temperature anisotropies in the microwave background by BOOMERANG, and MAXIMA experiments  at angular scales of $`\mathrm{}`$ 200 disfavor models of structure formation based exclusively on cosmic strings , still, due to many uncertainties in the scaling models of cosmic string network evolution one can not rule them out as candidates of sources of required density fluctuations. Further, even with present models, it is not ruled out that cosmic strings may contribute to some part in the structure formation in the universe. Above all, cosmic strings generically arise in many Grand Unified Theory (GUT) models. If the GUT scale is somewhat lower than $`10^{16}`$ GeV then the resulting cosmic strings will not be relevant for structure formation (for a discussion of these issues, see ). However, they may still affect various stages of the evolution of the universe in important ways.

In this paper we will study the implications of density fluctuations produced by cosmic strings on electroweak baryogenesis for the case of a second order electroweak phase transition. We will also briefly comment on the first order transition case, a detailed discussion of that case is postponed for a future work where we will discuss the effects of inflation generated density fluctuations on electroweak baryogenesis.

It is well known that cosmic strings can produce wake like over-densities with magnitude which can be large . (As we will discuss below, this is also true for strings moving through a relativistic ideal fluid). We consider the universe to be at a temperature $`T_b`$ which is below the electroweak transition temperature $`T_c`$, and consider density fluctuations via string wakes. Typical wavelengths associated with such wakes will be much smaller than the horizon (for ultrarelativistic strings). The evolution of such density fluctuations is simple for small amplitudes ; they propagate as plane waves with speed of sound $`c_s=1/\sqrt{3}`$. Even for larger amplitudes, the propagating density disturbance generated by string wakes will typically retain short wavelength . At a given region, through which the fluctuation propagates, the density, and hence the temperature, undergoes an oscillation, with the oscillation time scale being governed by the typical wavelength of the fluctuation. As a propagating density disturbance will be like a pulse of short wavelength, the temperature in a given spatial region in its path will typically undergo one oscillation cycle with short time period.

Depending on the amplitude of temperature oscillation, and the background temperature $`T_b`$ of the universe, it may happen that the temperature temporarily exceeds the electroweak transition temperature $`T_c`$, thereby restoring the electroweak symmetry. As the temperature oscillates back to values lower than $`T_c`$, the symmetry will be broken again spontaneously. During this re-occurrence of the electroweak phase transition, baryon asymmetry will be generated. For a second order phase transition no baryons would have been generated when the whole universe went through the electroweak phase transition, due to extremely slow cooling of the universe . However, the cooling time scale, and hence the time scale of the phase transition, due to the density (temperature) fluctuation pulse will be governed by the wavelength of the fluctuation, and can be much shorter than the Hubble time scale. The baryon asymmetry generated in this case, therefore, can be enhanced by many orders of magnitude. We will show that for a $`10^{16}`$ GeV GUT scale cosmic string, the resulting density fluctuations can give the required value of the baryon to entropy ratio (in the context of a two Higgs doublet model as in , with the CP violation parameter of order 10) if the fluctuation propagates without getting significantly damped. If fluctuations are dissipated rapidly, then a volume suppression factor reduces the produced baryon asymmetry, though it is still 3-4 orders of magnitude larger than the conventional case where baryons are produced during the cooling of the universe via expansion.

The paper is organized in the following manner. In section II we briefly recall the relevant results for the electroweak baryogenesis for the case when the electroweak phase transition is of second order. Here, the baryon to entropy ratio is expressed in terms of $`\dot{\theta }`$, the time derivative of the (relative) Higgs phase, and consequently $`\dot{T}`$, where $`T`$ is the temperature. Since largest baryon asymmetry will arise from the largest value of $`\dot{\theta }`$, which translates to the requirement of largest value of $`\dot{T}`$, it is clear that the most relevant density fluctuations in our model are those with the shortest wavelength (for a given value of the amplitude of density fluctuation). In section III, we discuss such density fluctuations as expected from cosmic strings moving through a relativistic fluid. Section IV presents our results for the baryon to entropy ratio with the value of $`\dot{T}`$ governed by short wavelength density fluctuations. Section V briefly discusses the expected results when the transition is first order. Conclusions are presented in section VI.

## II Baryon to Entropy Ratio with a Second Order Electroweak Transition

We first briefly outline the results for baryon to entropy ratio (B/s) when the electroweak phase transition is of second order and the phase transition proceeds by the cooling of the universe due to its expansion. For this, we will mostly follow the discussion in ref. . As discussed in ref. , in two-Higgs doublet extensions of the standard model, there are terms in the effective action for the gauge-Higgs sector at high temperature, which are of the form $`\chi F\stackrel{~}{F}`$. Here $`\chi `$ is some combination of the Higgs fields, $`F`$ is the SU(2) field strength tensor, and $`\stackrel{~}{F}`$ is its dual. Integration by parts shows (with the use of the chiral anomaly equation) that this term is equivalent to a chemical potential term $`\dot{\chi }B`$ for the baryon number, where $`B`$ is the baryon number density  (see, also, ref. ). The value of this chemical potential $`\mu _B`$ can be written as ,

$$\mu _B\dot{\chi }7\zeta (3)(\frac{m_t}{\pi T})^2\frac{v_2^2}{v_1^2+v_2^2}\dot{\theta },$$
(1)
where $`\zeta (3)1.2`$, and $`\theta `$ is the relative phase between the two Higgs fields. $`v_1`$ and $`v_2`$ are the magnitudes of the vacuum expectation values of the two Higgs fields, and $`m_t`$ is the top quark mass. In the presence of chemical potential given in Eq.(1), one can show  that the rate at which the baryon number relaxes to equilibrium, is given by

$$\dot{B}=\mathrm{\Gamma }^{}(Bc_n\mu _BT^2),$$
(2)
where $`\mathrm{\Gamma }^{}=\frac{N_F^2\mathrm{\Gamma }}{c_nT^3}`$, and $`c_n`$ is a slowly varying parameter with value remaining close to about 0.4. $`N_F`$ is the number of quark generations, and $`\mathrm{\Gamma }`$ is the sphaleron rate per unit volume in the broken phase . Using the Green’s function method, the solution to this equation is obtained as

$$B(t)=_{t_{in}}^t\mathrm{\Gamma }^{}(t^{})c_n\mu _B(t^{})T^2(t^{})[e^{_t^{}^t\mathrm{\Gamma }^{}(t^{\prime \prime })𝑑t^{\prime \prime }}]𝑑t^{}.$$
(3)
If we make the assumption that $`_{t_{in}}^t\mathrm{\Gamma }^{}𝑑t1`$, then the upper bound of this equation is given by the maximum of the function $`c_n\mu _B(t)T^2(t)`$ which we call as $`B_0`$. The solution is then simply

$$B(t)B_0_{t_{in}}^t\mathrm{\Gamma }^{}(t^{})[e^{_t^{}^t\mathrm{\Gamma }^{}(t^{\prime \prime })𝑑t^{\prime \prime }}]𝑑t^{}=B_0|c_n\mu _B(t)T^2(t)|_{max}.$$
(4)
We will use this upper limit $`B_0`$ as an estimate of the value of baryon asymmetry created. We write $`B_0`$ as follows,

$$B_0=|c_n\dot{\chi }T^2|_{max}=|c_n\frac{d\chi }{dT}\dot{T}T^2|_{max}.$$
(5)
The value of $`(\frac{d\chi }{dT})_{max}`$ will crucially depend on the details of the phase transition dynamics. Uncertainties in these estimates have been discussed in the literature . For simple estimates, one can use,

$$T\frac{d\chi }{dT}\frac{d\varphi }{dT}ϵϵ,$$
(6)
where $`ϵ`$ characterizes CP violation in the model, $`\varphi `$ being the Higgs field. Assuming order one CP violation, one can take $`T\frac{d\chi }{dT}1`$. It has been argued  that depending on the details of the specific models, this value may be larger by about two orders of magnitude. One can thus use $`ϵ`$ as a parameter, with its value ranging from 1 to about 200 (ref. ).

If the departure from equilibrium is due to the expansion of the universe, then the Friedmann equation for the radiation era can be used to relate the time derivative to the temperature derivative. In the radiation era,

$$H(T)=(\frac{8\pi ^3Gg_{}}{90})^{1/2}T^2=\frac{1}{2t}.$$
(7)
Here $`G`$ is the Newton’s constant, $`H^1`$ is the Hubble time, and we have used the plasma energy density $`\rho =\frac{\pi ^2g_{}T^4}{30}`$ with $`g_{}`$ being the number of degrees of freedom. We get,

$`\dot{T}=H(T)T,H(T_c)=H(T){\displaystyle \frac{T_c^2}{T^2}}.`$ (8)
Using these equations, we can write $`B_0`$ as,

$$B_0=c_nϵH(T)T^2.$$
(9)
The equation for the entropy density is given by

$$s=\frac{2\pi ^2}{45}g_{}T_c^3.$$
(10)
Resulting baryon to entropy ratio can now be written as,

$$\frac{B}{s}\frac{45c_nϵ}{2\pi ^2g_{}}\frac{H(T_c)}{T_c}.$$
(11)
For the electroweak case the value of $`g_{}`$ is about 100. With $`T_c=110GeV`$ and $`\frac{H}{T_c}10^{16}`$, we get,

$$\frac{B}{s}10^{18}ϵ.$$
(12)
Thus, as discussed in the literature , even with optimistic estimates of $`ϵ200`$, required baryon asymmetry cannot be generated for the case of second order transition when the transition proceeds by the cooling of the universe due to expansion.

In the next sections we will discuss the density fluctuations present in the universe. We will see that in the presence of these fluctuations, the time scale of phase transition can be smaller by several orders of magnitude, such that required baryon asymmetry may be generated with a second order transition.

## III Density Fluctuations in the Universe

Results for B/s obtained in the previous section correspond to the situation when the universe is perfectly homogeneous, and undergoes second order phase transition everywhere as the temperature falls below $`T_c`$ due to the expansion of the universe. As mentioned in the Introduction, presence of density fluctuations in the universe is unavoidable, as then only structures observed today can arise. In this section we will discuss the nature of density fluctuations produced via cosmic strings moving through the cosmic fluid. For the second order transition case, density fluctuations produced via inflation will not be helpful in generating a baryon asymmetry in our model, as we will explain at the end of Sec.IV. For simplicity, we discuss the case of straight cosmic strings which give rise to wake like density fluctuations. It is important to note that for cosmic strings, another important contribution to density fluctuations comes from small loops. However, the exact nature of density fluctuations in this case is more complicated due to the fact that loops rapidly oscillate, and the issue of time dependence is of crucial importance for our results. In a more detailed investigation their contribution must be taken into account. (Segments of very large string loops, which have sizes not much smaller than the horizon, will effectively behave as straight strings).

We first discuss the evolution of such small wavelength density fluctuations . It is well known that the density fluctuations in a relativistic ideal fluid which have wavelength of order, or larger, than the horizon, grow with the scale factor. However, when wavelength of the fluctuation is much smaller than the horizon, which is the relevant case for us, then the density fluctuation simply propagates as a plane wave with the speed of sound $`c_s=1/\sqrt{3}`$, at least when the amplitude of density fluctuation is small . Thus, such a density fluctuation will evolve as,

$$\frac{\delta \rho }{\rho }Ae^{i(𝐤.𝐫\omega t)},$$
(13)
where A is the amplitude of the fluctuation, $`𝐤`$ is the wave vector and $`\omega (=2\pi c_s/\lambda )`$ is the angular frequency of the fluctuation with wavelength $`\lambda `$. Our interest will be in the evolution of density disturbances for times shorter than the Hubble time, so the effects of universe expansion can be neglected for these fluctuations for an order of magnitude estimate. Also, the wake like density fluctuation arising from cosmic string moving through a relativistic ideal fluid (with zero chemical potential) may not lead to sustained oscillations. What one expects is that the resulting density fluctuation (hence temperature fluctuation) will propagate as a single pulse.

It is important to note here that the main factor responsible for baryon production enhancement in our model will be (as we will see below) the short time scale associated with these small scale density fluctuations compared to the Hubble time. This remains true irrespective of the fact whether these fluctuations propagate for significant distances, or they decay rapidly. If the fluctuation decays without propagating significantly, still the time scale for the temperature fluctuation will be the same as for the propagating case, as again it will be typically governed by the speed of sound $`c_s`$ and the length scale of the fluctuation. As we will see below, for the first case, the propagating disturbances can lead to required baryon production as they sweep the entire horizon volume. If fluctuations decay rapidly, then there will be a volume suppression factor which will lead to much smaller baryon production, though it will still be several orders of magnitude larger than the conventional case when baryons are produced only due to the overall cooling of the universe.

Once the nature of density fluctuations is known (depending on the source of density fluctuations), we can use Eq.(13) to determine the time variation of temperature resulting from such propagating density fluctuations. An important assumption we make in this regard is that temperature at any time and at any point of space, is determined by the local density $`\rho `$, even in the presence of short wavelength density fluctuations. This assumption of local thermal equilibrium is justified here because the wavelengths we will be considering will be about $`10^6`$ cm, with the relevant time scale of density variation being of order $`10^{16}`$ sec. This time scale, though much shorter than the Hubble time, is still many orders of magnitude larger than the time scale of relevant interactions, which keep the particles in equilibrium.

With the use of local thermal equilibrium assumption, the density $`\rho `$ is related to the temperature $`T`$ as, $`\rho =\frac{\pi ^2g_{}}{30}T^4`$, where $`g_{}100`$ is the number of degrees of freedom at the electroweak scale. From this we can relate the temperature fluctuation to the density fluctuation,

$$\frac{\delta \rho }{\rho }\frac{\rho \rho _b}{\rho _b}=\frac{T^4T_{b}^{}{}_{}{}^{4}}{T_{b}^{}{}_{}{}^{4}},$$
(14)
where $`T_b`$ is the background temperature of the plasma. Then using Eq.(13) we obtain the time - temperature relation for propagating density fluctuations.

$$\frac{T^4T_{b}^{}{}_{}{}^{4}}{T_{b}^{}{}_{}{}^{4}}=Ae^{i(𝐤.𝐫\omega t)}.$$
(15)
Given this, we take the time dependence of the temperature at a given point in space to be given by,

$$T(t)=T_b[1+A\mathrm{sin}\omega t]^{1/4}.$$
(16)
Thus, at $`t=0`$, $`T=T_b`$, the background temperature.

As discussed above, although we are using the above expression for $`T(t)`$ with a periodic time dependence, the propagating density fluctuation actually constitutes only a single pulse. Thus the above expression for $`T(t)`$ is to be used only for a single cycle. Of course that is enough for our purpose as the baryon asymmetry is to be generated during the cooling part of the cycle, irrespective of the number of cycles.

Propagation of density fluctuation given in Eq.(13), and resulting oscillation of temperature (Eq.(16)) is valid for the case when the amplitude $`A`$ is small, so the disturbance propagates as acoustic wave. For our case, density fluctuations produced in string wakes can have large amplitudes, with $`A1`$. Evolution of such density fluctuations is not as simple as given above. In fact, for a general fluctuation (i.e. with large amplitude), the propagating disturbance can lead to transfer of energy to shorter wavelengths, leading to steepening of the disturbance . In our case the initial profile of the density fluctuation will be given by the wake like shock which forms behind a cosmic string moving with supersonic velocity. Also, the baryon asymmetry produced in our model will be larger for shorter wavelengths due to shorter time scale of temperature variation. Thus, a propagating disturbance even with large amplitude will lead to baryon production at least as large as given by the initial wavelength of the disturbance. This is assuming that the disturbance does not decay away very rapidly, e.g. by subsequent shock developments etc. Though, as we have discussed above, even if fluctuation decays rapidly, the typical time scale of temperature evolution will still be given by the wavelength of the fluctuations. With this in mind, we will use the temperature variation given by Eq.(16) (again, for one cycle only as appropriate for a pulse of over-density) even for the case when the amplitude $`A`$ of density fluctuation is not too small, as in the case of a relativistic cosmic string wake.

### A Density Fluctuations arising due to Straight Cosmic Strings

We now give a brief review of the structure of density fluctuations produced by a cosmic string moving through a relativistic fluid. The space-time around a straight cosmic string (along the z axis) is given by the following metric ,

$$ds^2=dt^2dz^2dr^2(14G\mu )^2r^2d\psi ^2,$$
(17)
where $`\mu `$ is the string tension. This metric describes a conical space, with a deficit angle of $`8\pi G\mu `$. This metric can be put in the form of the Minkowski metric by defining angle $`\psi ^{}=(14G\mu )\psi `$. However, now $`\psi ^{}`$ varies between 0 and $`(14G\mu )2\pi `$, that is, a wedge of opening angle $`8\pi G\mu `$ is removed from the Minkowski space, with the two boundaries of the wedge being identified. It is well known that in this space-time, two geodesics going along the opposite sides of the string, bend towards each other . This results in binary images of distant objects, and can lead to planar density fluctuations. These wakes arise as the string moves through the background medium, giving rise to velocity impulse for the particles in the direction of the surface swept by the moving string. For collisionless cold dark matter particles the resulting velocity impulse is , $`v_{impls}4\pi G\mu v_{st}\gamma _{st}`$ (where $`v_{st}`$ is the transverse velocity of the string). This leads to a wedge like region of overdensity, with the wedge angle being of order of the deficit angle, i.e. $`8\pi G\mu `$ ($`10^5`$ for 10<sup>16</sup> GUT strings). The density fluctuation in the wake is of order one. Subsequent growth of this wake by gravitational instability in the matter dominated era has been analyzed in great detail in literature .

The structure of this wake is easy to see for collisionless particles (whether non-relativistic, or relativistic). Each particle trajectory passing by the string bends by an angle of order $`4\pi G\mu `$ towards the string. In the string rest frame, take the string to be at the origin, aligned along the z axis, such that the particles are moving along the $`x`$ axis. Then it is easy to see that particles coming from positive $`x`$ axis in the upper/lower half plane will all be above/below the line making an angle $`4\pi G\mu `$ from the negative $`x`$ axis. This implies that the particles will overlap in the wedge of angle $`8\pi G\mu `$ behind the string leading to a wake with density twice of the background density. One thus expects a wake with half angle $`\theta _w`$ and an overdensity $`\delta \rho /\rho `$ where ,

$$\theta _w4\pi G\mu ,\frac{\delta \rho }{\rho }1.$$
(18)
However, the case of relevance for us is cosmic strings moving through a relativistic plasma of elementary particles at temperatures of order 100 GeV. At that stage, it is not proper to take the matter as consisting of collisionless particles. A suitable description of matter at that stage is in terms of a relativistic fluid which we will take to be an ideal fluid consisting of elementary particles. Generation of density fluctuations due to a cosmic string moving through a relativistic fluid has been analyzed in the literature . The study in ref. focused on the properties of shock formed due to supersonic motion of the string through the fluid. In the weak shock approximation, one finds a wake of overdensity behind the string. In this treatment one can not get very strong shocks with large overdensities. In refs. , a general relativistic treatment of the shock was given which is also applicable for ultra-relativistic string velocities. The treatment in ref. is more complete in the sense that the equations of motion of a relativistic fluid are solved in the string space-time (Eq.(17)), and both subsonic and supersonic flows are analyzed. One finds that for supersonic flow, a shock develops behind the string, just as in the study of ref. . In the treatment of ref. one recovers the usual wake structure of overdensity (with the wake angle being of order $`G\mu `$) as the string approaches ultra-relativistic velocities. Also the overdensity becomes of order one in this regime.

We mention here that it is not expected that the string will move with ultra-relativistic velocities in the early universe. Various simulations have shown  that rms velocity of string segments is about 0.6 for which the shock will be weak. For this case, the opening angle of the wedge $`\theta _w`$ will also be large. This will then imply a large wavelength for the density fluctuation, and consecutively, small value of $`\dot{T}`$ and small baryon to entropy ratio. It is still possible that various segments of string may move with ultra-relativistic velocities. (Note that the friction dominated regime for GUT strings ends long before these temperatures are reached .) We will take this to be the case, and will consider the situation when string produces strong shock, with a small wedge angle.

For ultra-relativistic case, we use expressions from ref.. Resulting density fluctuation in the wake of the moving string is expressed in terms of fluid and sound four velocities,

$$\frac{\delta \rho }{\rho }\frac{16\pi G\mu u_f^2(1+u_s^2)}{3u_s\sqrt{u_f^2u_s^2}},sin\theta _w\frac{u_s}{u_f},$$
(19)
where $`u_f=v_f/\sqrt{1v_f^2}`$ and $`u_s=c_s/\sqrt{1c_s^2}`$. In this case, when string velocity $`v_f`$ is ultra-relativistic, then one can get strong overdensities (of order 1) and the angle of the wake approaches the deficit angle $`8\pi G\mu `$. This is the structure of wake which is also expected for a collision-less matter as discussed above. We will assume that this is the case, and hence, take the properties of the shock as given by Eq.(18). The density fluctuations resulting from (straight) cosmic strings, as characterized by Eq.(18), have the amplitude $`A`$ of order 1, and average wavelength $`\lambda `$ of order (1/2)($`8\pi G\mu d_H`$) where $`d_H`$ is the horizon size. For a GUT scale cosmic string with $`8\pi G\mu 10^5`$, $`\lambda `$ is of order $`10^6`$ cm. Note that this value of $`\lambda `$ is much larger than the diffusion length scales $`l_{diff}`$ of various particles at these temperatures which is of order 10<sup>-14</sup> cm . Hence these fluctuations will not be dissipated by diffusion processes.

## IV Baryon to Entropy ratio due to Density Fluctuations

In the previous section, we have determined the time variation of the temperature resulting from density fluctuations. From Eqs.(5),(6), and (16), we can write,

$$B_0=(c_nϵT\dot{T})_{max}=\frac{c_nϵA\omega T_b^4}{4T^2}.$$
(20)
This estimate of maximum baryon asymmetry has been obtained with the optimistic assumption that sphaleron processes turn off when the source of baryon number has optimum value . As the sphaleron rate changes by many orders of magnitude when $`T`$ decreases below $`T_c`$, relevant value of $`T`$ in the above equation will be close to $`T_c`$. For rough estimates, we will take $`TT_b`$ (note that $`T_b`$ is of same order as $`T_c`$ in our model as we want density fluctuations to be able to restore the electroweak symmetry locally). With this, Eq.(20) is obtained with maximum value of $`cos\omega t`$ = 1. The resulting baryon to entropy ratio can now be expressed in terms of the amplitude $`A`$, and the wavelength $`\lambda `$ of the density fluctuation as,

$$\frac{B}{s}=\frac{45c_nc_s}{\pi g_{}}\frac{A}{4\lambda T_b}ϵ\frac{0.01ϵA}{\lambda T_b}.$$
(21)
Here we have used $`g_{}=100`$.

For density fluctuations generated due to cosmic string wakes, the amplitude $`A1`$, and the shortest wavelength is given by the average width of the generated wake $`\frac{1}{2}8\pi G\mu d_H`$ (since the wake has the shape of a wedge) where $`d_H`$ is the horizon size at that time. Since $`d_H`$ is of the order of 0.1 cm at the time of the electroweak scale, and $`8\pi G\mu 10^5`$ for a GUT string, we get $`\lambda _{csmc}10^6cm`$. The resulting value of baryon to entropy ratio from Eq.(21) is $`B/s10^{11}ϵ`$. We see that for the case of cosmic string, with $`ϵ`$ of order 10, one is able to get the required baryon asymmetry. (We mention here that with the amplitude of density fluctuations of order one, $`T_b`$ can be much lower than $`T_c`$ in the cosmic string case. Thus, there may be possibility of having a weaker bound on the Higgs mass which is needed to avoid sphaleron wash out of the created baryon asymmetry, see ref..)

It is important to realize that, if fluctuations propagate without getting significantly damped, then essentially all regions will participate in this generation of baryon asymmetry, so there is no volume suppression factor here. For example, the structure of density fluctuations here is in form of planar sheets (wakes). These planar density fluctuations will propagate in the direction normal to the plane, and sweep the entire horizon volume. At any given point within the horizon volume, the temperature will undergo an oscillation as given by Eq.(16). There may be some horizon volumes which will not have any strings present to generate significant density fluctuations at the relevant stage. However, the probability of occurrence of such regions will be small, and hence the order of magnitude estimate of the baryon to entropy ratio as given above will not be affected. There will be some change in the wavelengths of these fluctuations during the time they propagate across the entire horizon, but again, one will expect it not to affect the final numbers within an order of magnitude. Some level of dissipation of density fluctuations, while they propagate, may be compensated by the fact that there will be multiple sources of density fluctuations in a given horizon volume. For example, string simulations have shown  that the number of long strings in a given horizon volume is expected to be about 15. Each of these strings will lead to wakes of density fluctuations, which will sweep the horizon volume. This may also compensate for the effects of shorter wakes of a given strings (possibly getting truncated behind the string due to other strings/fluctuations). It is also possible that the baryon asymmetry generated due to one string may first get wiped out when temperature there rises above $`T_c`$ due to another string wake. In that region, the surviving baryon asymmetry will correspond to the density fluctuations which are last to pass by. Eventually, the temperature in that region will be sufficiently low so that any subsequent cosmic string wakes are not able to re-heat the region enough to destroy the already created baryon asymmetry. This will be the stage when baryon asymmetry will get frozen.

From above discussion we see that the most important factor responsible for the enhancement of baryon to entropy ratio in our model is the short time scale of these density fluctuations. In addition, propagation of these over-density pulses leads to baryon production in almost entire region (say in a given horizon volume). That is, there is no volume suppression factor here. We now discuss the possibility if these fluctuations do not propagate for large distances, and decay very rapidly. In such a situation, to make a conservative estimate, we can take the over-dense wake of each string as a fixed region in space where the overdensity will dissipate typically with speed of sound $`c_s`$. This means that the time scale of temperature variation in this wake will again be the same as discussed above. From this we conclude that the resulting baryon to entropy ratio in this case will be suppressed by a volume factor $`f_s`$ compared to the value given in Eq.(21). $`f_s`$ will be given by the ratio of the total volume inside the wakes of all strings in a given horizon volume to the horizon volume. Taking account of about 10 long strings in a horizon volume, with typical wake thickness of order 10$`{}_{}{}^{5}d_{H}^{}`$, we get,

$$f_s=\frac{10\times 10^5d_H\times d_H^2}{d_H^3}10^4$$
(22)
With various factors taken as earlier, we get net $`B/s10^{15}ϵ`$ for this case. With optimistic value of $`ϵ100`$, resulting value of $`B/s`$ is $`10^{13}`$, about three orders of magnitude smaller than the required one. There may be a possibility of further enhancement in this value by changing various parameters of cosmic string network e.g. string scale and string velocity (which determine the parameters of string wakes), as well as by properly taking into account the effect of density fluctuations produced by all strings at that stage, i.e. including contributions of all string loops as well, integrated over a suitable duration of time. The main point we would like to make is that even in this case where density fluctuations do not propagate, our model leads to the enhancement of baryon to entropy ratio by at least 3 - 4 orders of magnitude compared to the case when baryogenesis happens only due to the overall cooling of the universe.

An important assumption in deriving equations such as Eq.(5),(6) (using which Eq.(20) has been derived) is that at any stage during the variation of the temperature, the value of the Higgs field is essentially given by the vacuum expectation value ($`vev`$) at that temperature. This assumption requires that the time scale of temperature variation should be much larger than the typical time scale of the Higgs field evolution. More importantly, as discussed above, as the temperature rises above $`T_c`$ (Eq.(16)), the Higgs field should settle down to value zero so that the electroweak symmetry is restored. This is very important as during the rise of temperature, one will get production of antibaryons via above equations (as sign of $`\dot{\chi }`$ will be reversed during this stage). However, if the electroweak symmetry is restored after that, then all these net antibaryon asymmetry will be erased due to unsuppressed sphaleron processes in the symmetric phase. Thus, only the net baryon asymmetry will remain which will be generated during the stage when the temperature is decreasing.

For checking this requirement on the time scale of temperature variation, we take a simple effective potential to model the second order phase transition (as in ref.),

$$V(\varphi ,T)=D(T^2T_c^2)\varphi ^2+\frac{\lambda }{4}\varphi ^4.$$
(23)
We take the values of the parameters as suitable for the Standard model, $`T_c=110`$ GeV, $`D0.18`$ and $`\lambda 0.1`$. The evolution equation of $`\varphi `$ (neglecting spatial dependence) is,

$$\ddot{\varphi }+\frac{3}{2}\frac{\dot{\varphi }}{t}+V^{}(\varphi ,T(t))=0,$$
(24)
where the time dependence of $`T`$ in $`V`$ is given by Eq.(16). Using this equation, we have checked the evolution of $`\varphi `$ for one oscillation period, starting from the initial time $`t=0`$ in Eq.(16), which corresponds to $`T=T_b`$ (the background temperature) initially. Top figure in Fig.1 shows the time variation of the temperature resulting from cosmic string wakes. Here we have taken $`T_b=100`$ GeV, $`A`$ = 0.8 in Eq.(16), and $`\omega `$ corresponds to $`\lambda =4\pi G\mu d_H10^6`$ cm. The values of $`T_b`$ and $`A`$ are taken as sample values such that the maximum value of $`T`$ (Eq.(16)) reaches slightly above $`T_c=110GeV`$. Bottom figure in Fig.1 shows plot of the evolution of the Higgs field $`\varphi `$ (solid line), and the $`vev`$ of $`\varphi `$ (shown by solid dots at few values of $`t`$). Note that the plot of $`\frac{\mathrm{\Delta }T}{T}`$ deviates from purely sinusoidal behavior due to large value of $`A`$.

From Fig.(1) we see that $`\varphi `$ traces the evolution of the $`vev`$ faithfully. For time duration when $`T>T_c`$, $`\varphi `$ settles down to the value zero, thereby restoring the electroweak symmetry. This restoration of symmetry happens for the values of $`t`$ between $`0.9\times 10^7`$ to $`3.5\times 10^7`$ GeV<sup>-1</sup>. As we mentioned above, this temporary symmetry restoration is necessary to wipe out the antibaryons which would have been generated during the rising phase of the temperature. It is also important to note from these figures that the time duration for which the symmetric phase lasts (i.e. when $`T>T_c`$) is of the order of $`10^7GeV^1`$, which is much larger than the typical time scale of a sphaleron transition in the symmetric phase which is roughly of order $`(\alpha _w^2T)^110GeV^1`$. (Spatial extent of regions with density fluctuation is also many orders of magnitude larger than the sphaleron size $`(\alpha _wT)^1`$.) Thus, equilibrium calculations used in this paper make sense, sphaleron transitions during the time when $`T>T_c`$, should be able to erase any antibaryon asymmetry generated earlier during the rise of temperature.

We now comment on the density fluctuations generated via inflation. Generic prediction of most inflationary models  is that the resulting density fluctuations are scale invariant, with a magnitude of about $`10^4`$ (when the fluctuation re-enters the horizon), as constrained by COBE observations of temperature fluctuations $`\frac{\delta T}{T}10^5`$. In the context of our model, inflationary density fluctuations may contribute to baryon production only if such small temperature increase can lead to significant change in the baryon violation rates via sphaleron transitions. This requires that the background temperature of the universe $`T_b`$ be very close to the electroweak transition temperature $`T_c`$. However, in such a situations the produced baryon asymmetry will be erased as the sphaleron transitions at $`T=T_b`$ will not be much suppressed. Basically, for the scales involved, a fractional temperature difference of order $`10^5`$ is unlikely to lead to large change in sphaleron rates. However, the situation is entirely different or the case of a first order transition, where even very small temperature differences can lead to significant changes in the bubble nucleation probability. We briefly discuss it in the next section.

## V First Order Transition Case

Our conclusions about the second order transition should be applicable in the case of a very weak first order transition (or a cross-over). For example, after a rapid production of bubbles of the broken phase, the bubbles will collide and the transition will be completed. Presence of propagating density fluctuations will again reheat these regions temporarily, restoring the electroweak symmetry. Subsequent cooling will then produce a net baryon asymmetry as described above, with the rate of change of temperature governed by the oscillation frequency of the density fluctuation.

If the first order transition is not very weak, then the effect of density fluctuations should manifest in entirely different manner. Now, the presence of density (and hence) temperature inhomogeneities will split the region into cold regions and hot regions (at a given time). Geometry of such regions will depend on the source of density fluctuations. For example, for inflationary density fluctuations, such regions can be taken as localized overdense, or underdense regions, as discussed in ref. . Phase transition (say, via bubble nucleation) will first happen in the colder regions, forming large interfaces which then extend into the hotter regions. Thus, even for the case of inflationary density fluctuations, the process of electroweak baryogenesis may get affected significantly. For example, it is possible that even for the case of a weak first order transition, different temperature zones may lead to propagating interfaces emanating from the colder regions and extending into the hotter regions. This in effect may lead to a picture similar to that of a strong first order transition which is one of the requirements for a successful baryogenesis model.

Also, the propagation of the interface will reheat the surrounding region. Due to temperature variations, density of baryons produced in this case may then vary from cold regions to the hot regions. Similarly, for the case of cosmic strings, wedge shaped regions will form where the temperature will be higher than the surrounding (as discussed in ref.). Bubble nucleation and coalescence may happen in the colder regions first, with the interface then extending into the planar hotter regions. This may lead to variations in the baryon density of planar nature. It will be very interesting to find whether any such baryon inhomogeneities can survive until nucleosynthesis to affect elemental abundances.

## VI Conclusion

We have studied the implications of small wavelength density fluctuations produced by cosmic strings on electroweak baryogenesis. We have shown that due to rapid changes in the temperature arising from short wavelength density fluctuations, electroweak symmetry may get restored locally, depending on the amplitude of temperature oscillation and the background temperature $`T_b`$ of the universe. As the temperature oscillates back to a value below the transition temperature, electroweak symmetry will be broken again spontaneously. During this re-occurrence of the electroweak phase transition, baryon asymmetry will be generated. The cooling time scale, and hence the time scale of the phase transition, due to such density (temperature) fluctuation waves will be governed by the wavelength of the fluctuation, and can be much shorter than the Hubble time scale. The baryon asymmetry generated in this case is therefore enhanced by many orders of magnitude. If these short wavelength density fluctuations propagate for significant distances without decaying, then baryogenesis happens in every region of space, as different regions are swept by the propagating over-density (hence temperature) pulse. We have shown that for a $`10^{16}`$ GeV GUT scale cosmic string, the resulting density fluctuations can give the required value of the baryon to entropy ratio (in the context of a two Higgs doublet model as in , with the CP violation parameter of about 10). On the other hand, if density fluctuations decay away rapidly, then there is a volume suppression factor, and the resulting baryon asymmetry is smaller by a factor of order 10<sup>-4</sup> (though it is still larger by at least 3 - 4 orders of magnitude compared to the conventional case where baryons are produced when the universe cools via expansion).

Many aspects of our model need to be worked out in more detail. Exact nature of propagation and decay of these density perturbations needs to be worked out. We have used order of magnitude estimates for the chemical potential for baryons (Eq.(6)), more careful estimates need to be worked out. Further, the baryon asymmetry should be calculated by taking account of several propagating density fluctuation waves. Net asymmetry should be calculated by detailed accounting for the generation of antibaryon asymmetry during the period when $`T`$ rises, its diffusion outside the region of density fluctuation, and the decay of this antibaryon asymmetry when $`TT_c`$ so that the sphaleron transitions are unsuppressed. Finally, the first order transition case may show many interesting possibilities, especially for the case of inflationary density fluctuations.

Similar calculations, as reported in this paper, can be done for other sources of density fluctuations, such as those generated by textures etc. We expect similar results for the baryon asymmetry in those cases. For example, for the case of textures, collapsing textures (at any stage) will give rise to small scale density fluctuations, which will then propagate as density waves. The amplitude of fluctuations may be high for the case of textures also since the relevant fluctuations will be generated during the last stages of texture collapse, though the volume factor will lead to an overall suppression of baryon production.

ACKNOWLEDGEMENTS

We are very thankful to Sanatan Digal, Sayan Kar, Rajarshi Ray, and Supratim Sengupta for useful comments. We especially thank Raghavan Rangarajan for detailed comments and many useful suggestions.


format: text

text:
# Thermodynamic approach to dense granular matter: a numerical realization of a decisive experiment.

Years ago Edwards proposed a thermodynamic description of dense granular matter, in which the grains (the ‘atoms’ of the system) interact with inelastic forces. The approach is intriguing but is not justified from first principles, and hence, in the absence of conclusive tests of its validity, it has not been widely accepted. We perform a numerical experiment with a realistic granular matter model specially conceived to be reproducible in the laboratory. The results strongly support the thermodynamic picture.

The similarities between a driven granular system and a system of molecules in a fluid have prompted the use of statistical ideas to describe granular matter , such as for example kinetic theories of gases to study rapid flows of low-density systems of inelastic particles. In this context, it is customary to associate the mean kinetic energy of the particles with a ‘granular temperature’, having no real thermodynamic meaning.

For the opposite limit of dense, slowly flowing granular matter, a more radical train of ideas was initiated more than a decade ago by Edwards and coworkers: the proposal of a statistical ensemble , and through it, thermodynamic notions such as entropy and temperature (the latter unrelated to the ‘granular temperature’ above). Although the idea of a thermodynamic description of granular matter was recognized as attractive, it was not universally accepted because there is no known first principle justification of Edwards’ statistical ensemble, such as there is for ordinary statistical mechanics of liquids or gases (Liouville’s theorem).

Recently, analytical developments and numerical work on schematic models originally devised for glasses have given a new perspective for the dynamics of dense granular matter . The theoretical support provided by this framework for a statistical approach has spurred a renewed interest in Edwards’ thermodynamics.

On the experimental side, years ago Nowak et al.  studied in detail the density fluctuations in a vibrated granular material, and proposed that these fluctuations should reflect an underlying thermodynamics, much as they do in an ordinary thermal system. However, as far as the actual verification of Edwards’ approach, the evidence they found was if anything rather negative: we shall discuss below some possible reasons for this. On the other hand, Edwards’ approach has never been tested with simulations of realistic models of granular materials, characterized by the fact that energy is supplied by external driving (via shear or tapping) and dissipated by inelastic collisions and slippage between the grains.

In this work we describe a numerical experiment specially conceived to be reproducible in the laboratory, using a realistic model of sheared granular matter. Considering particles of different sizes inside the medium, we compute a dynamical temperature from an Einstein relation relating random diffusion and mobility (transverse to the shear direction) at long time scales. If there is an underlying thermodynamics, it will impose that the dynamical temperature be independent of the tracer’s shape: a very strong condition which we show to hold. We then perform an explicit computation to show that the temperature arising from Edwards’ thermodynamic description and the dynamic one measured from Einstein’s relation coincide. This last step cannot be performed in the laboratory, so the numerical simulation provides the missing link between thermodynamic ideas and diffusion-mobility checks.

Consider a ‘tracer’ body of arbitrary shape immersed in a liquid in equilibrium at temperature $`T`$. As a consequence of the irregular bombardment by the particles of the surrounding liquid, the tracer performs a diffusive, fluctuating ‘Brownian’ motion. The motion is unbiased, and for large times the average square of the displacement goes as $`|x(t)x(0)|^2=2Dt`$, where $`D`$ is the diffusivity. On the other hand, if we pull gently from the tracer with a constant force $`f`$, the liquid responds with a viscous, dissipative force. The averaged displacement after a large time is $`[x(t)x(0)]=f\chi t`$, where $`\chi `$ is the mobility. Clearly, the same liquid molecules are responsible both for the fluctuations and for the dissipation. Although both $`D`$ and $`\chi `$ strongly depend on the shape and size of the tracer, they turn out to be always related by the Einstein relation $`D/\chi =T`$, where $`T`$ is the temperature of the liquid.

Conversely, if in a fluid of unknown properties we find that several tracers having different diffusivities and mobilities yield the same ratio $`D/\chi =T`$, we may take this as a strong evidence for thermalisation at temperature $`T`$. Indeed, recent analytic schemes for out of equilibrium glassy dynamics have suggested the existence of such a temperature $`T`$ governing the slow components of fluctuations and responses of all observables .

In order to test the existence of this temperature for dense granular matter, we perform, with a realistic numerical model, a diffusion-mobility experiment in conditions that can be reproduced in the laboratory. In the model, deformable spherical grains interact with one another via non-linear elastic Hertz normal forces, and non-linear elastic and path-dependent Mindlin transverse forces . The normal force, $`F_n`$, has the typical 3/2 power law dependence on the overlap between two spheres in contact, while the transverse force, $`F_t`$, depends on both the shear and normal displacements between the spheres. For two spherical grains with radii $`R_1`$ and $`R_2`$: $`F_n=\frac{2}{3}k_nR^{1/2}w^{3/2},`$ $`\mathrm{\Delta }F_t=k_t(Rw)^{1/2}\mathrm{\Delta }s.`$ Here $`R=2R_1R_2/(R_1+R_2)`$, the normal overlap is $`w=(1/2)[(R_1+R_2)|\stackrel{}{x}_1\stackrel{}{x}_2|]>0`$, and $`\stackrel{}{x}_1`$, $`\stackrel{}{x}_2`$ are the positions of the grain centers. The normal force acts only in compression, $`F_n=0`$ when $`w<0`$. The variable $`s`$ is defined such that the relative shear displacement between the two grain centers is $`2s`$. The prefactors $`k_n=4G/(1\nu )`$ and $`k_t=8G/(2\nu )`$ are defined in terms of the shear modulus $`G`$ and the Poisson’s ratio $`\nu `$ of the material from which the grains are made (typically $`G=29`$ GPa and $`\nu =0.2`$, for spherical glass beads). We assume a distribution of grain radii in which $`R_1=0.105`$ mm for half the grains and $`R_2=0.075`$ mm for the other half. The observables are measured in reduced units: length in units of $`R`$, force in units of $`GR^2`$, time in units of $`\sqrt{\rho R^2/G}`$, where $`\rho `$ is the density of the particles. Internal dissipation is included in two ways: (1) via a viscous damping term proportional to the relative normal and tangential velocities, and (2) via sliding friction: when $`F_t`$ exceeds the Coulomb threshold, $`\mu F_n`$, the grains slide and $`F_t=\mu F_n`$, where $`\mu `$ is the friction coefficient between the spheres (typically $`\mu =0.3`$).

We perform molecular dynamics (MD) simulations for a binary system of 200 large and small spheres in a periodic 3D cell. Our calculation begins with a numerical protocol designed to mimic the experimental procedure used to prepare dense packed granular materials. The simulations start with a gas of spherical particles located at random positions in a periodically repeated cubic cell of side $`L`$. The system is then compressed and extended slowly until a stationary situation with a specified value of the pressure and volume fraction— above the random close packing fraction— is achieved. We then apply a gentle shear in the $`y`$-direction at constant volume by moving the periodic images at the top and bottom of the cell with velocities $`\dot{\gamma }L/2`$, where $`\dot{\gamma }`$ is the shear rate (Lees-Edwards boundary conditions). In this simple shear flow the gradient of the velocity along $`z`$ is uniform (no shear bands). We apply a shear rate of $`\dot{\gamma }=10^4`$, which allows the system to be in the quasi-static regime. We focus our study on the slow shear rate, quasi-static limit where the system is always close to jamming, but moves just barely enough to avoid stick-slip motion . Here the external pressure is large enough ($`10`$ MPa in a typical simulation) and the average coordination number is high enough (typically $`7`$) that deformation and elasticity of the particles play the dominant role, as opposed to the collision dominated rapid flow regime described by kinetic theories. We checked that shear induced segregation is absent at the times scales of our simulations. The contacts between particles are enduring and the internal stresses in the system are transmitted via a network of ‘force chains’ . independent of the shear rate.

After a transient of the order of the inverse shear rate, we start measuring the spontaneous $`|x(t)x(0)|^2`$ and force-induced displacements $`[x(t)x(0)]/f`$ along the $`x`$-direction for the two types of particles with different sizes. The results are shown in Fig. 1. We notice that the diffusivities and the mobilities are different for the two type of particles, since they have different sizes. However, when we draw the parametric plot of $`|x(t)x(0)|^2`$ versus $`[x(t)x(0)]/f`$ (Fig. 2) we find parallel straight lines for large time scales, implying an extended Einstein relation:

$$|x(t)x(t_w)|^2=2T_{dyn}\frac{[x(t)x(t_w)]}{f},$$
(1)
valid for both particles with the same $`T_{dyn}`$ for widely separated time scales $`tt_w`$. This suggests that $`T_{dyn}`$ can be considered to be the temperature of the slow modes.

For small time scales (fast rearrangements, near the origin in Fig. 1) the fluctuation-dissipation plot shows that there is not a well defined Einstein-relation temperature. We expect this result since the fast motion of the grains depends on the microscopic interactions dominated by inelastic collisions . We also calculate the ‘granular kinetic temperature’ defined in terms of the velocity fluctuations in the $`x`$-direction. Unlike $`T_{dyn}`$, we find that this temperature is different for the two types of particles, and their values are two orders of magnitude smaller than $`T_{dyn}`$. We expect this result since the kinetic granular temperature is dominated by the fast (high frequency) modes, which are not thermalised to $`T_{dyn}`$.

We also repeat the numerical experiment for a system of Hertz spheres without transverse forces ($`\mu =0`$: experimental realizations of elastic spheres with viscous forces but without sliding friction are foams and compressed emulsions ) and find that $`T_{dyn}`$ is well defined at long time scales for this case as well (see Fig. 1). Thus, our results suggest that the validity of a structural temperature for long-scale displacements (larger than a fraction of the particle size) holds in the presence of viscous forces between grains or even of a sliding threshold (Coulomb’s law).

The existence of a single temperature is an instance of the zero-th law of thermodynamics, for which we find positive evidence here, at variance with the experimental result in Nowak et al . Three possible reasons for the apparent violation of the ‘zero-th law’ in their experiments are i) the effect of an unknown height-dependent pressure (as pointed out by authors), ii) a rather high tapping amplitude (Edwards ensemble is in principle only valid in the high compaction limit ) and iii) the fact that the density fluctuations considered are integrated over all frequencies, thus including also ‘fast’ relaxations.

Next, we treat the question whether it is possible to relate the dynamical temperature obtained above to the thermodynamic construction of slowly driven out of equilibrium systems proposed by Edwards. Whereas in the Gibbs construction for equilibrium statistical mechanics one assumes that the physical quantities are obtained as average over all possible configurations, Edwards ensemble consists of only the blocked configurations (static or jammed) at the appropriate energy and volume. The strong ‘ergodic’ hypothesis is that all blocked configurations of given volume and energy can be taken to have equal statistical weights. This formulation leads to an entropy $`S_{Edw}(E,V)`$, and the corresponding temperature $`T_{Edw}^1=S_{Edw}/E`$ and compactivity $`X_{Edw}^1=S_{Edw}/V`$.

In order to calculate $`T_{Edw}`$ and compare with the obtained $`T_{dyn}`$ we need to count the number of blocked configurations at a given energy and volume. (For this calculation we concentrate in the case without tangential forces and sliding friction, in order to avoid path dependency which would lead to an ambiguity in the definition of blocked configurations—see below). To do this in practice, we extend the ‘auxiliary model’ method  to the case of deformable grains. We define an auxiliary model composed of the true deformation energy ($`E`$, the Hertzian energy of deformation of the particles) and of a term that vanishes in (and only in) the blocked configurations: $`E_{block}_a\left|\stackrel{}{F}_a\right|^2`$. (Here $`\stackrel{}{F}_a`$ is the total contact force exerted on particle $`a`$ by its neighbours). We perform equilibrium MD with two auxiliary temperatures $`(T^{},T_{aux})`$, corresponding to the partition function $`\mathrm{exp}[E/T^{}E_{block}/T_{aux}]`$. Annealing $`T_{aux}`$ to zero selects the blocked configurations ($`E_{block}=0`$), while $`T^{}`$ fixes the energy $`E`$. We start by equilibrating the system at high temperatures ($`T_{aux}`$ and $`T^{}`$ $`\mathrm{}`$) and anneal slowly the value $`T_{aux}`$ to zero and tune $`T^{}`$ so as to reach the observed value of $`E`$ during shear. At the end of the process ($`T_{aux}0`$), $`T^{}(E)=T_{Edw}(E)`$, since in this limit we are sampling the configurations with vanishing fraction of moving particles at a given $`E`$.

Figure 3 shows various annealing protocols. We plot the value of the true deformation energy as a function of $`T_{aux}`$ during annealing using different values of $`T^{}`$ as indicated in the figure. At the end of the annealing process the system stabilizes at a given deformation energy. Only when $`T^{}`$ is equal to $`T_{dyn}`$ obtained through the Einstein relation during shear, we find that the final $`E`$ coincides, within the accuracy of the simulations, with the mean energy of the system under shear. For other values of $`T^{}`$ the annealed system equilibrates at energies which are above of below the distribution of deformation energies during shear as seen in Fig. 3. This proves the equality between Edwards’ and the dynamic temperature.

We conclude with some remarks: i) Since the blocked configurations are the same whatever the inter-grain dissipation coefficient, Edwards’ ensemble (and hence its temperature) are insensitive to viscous dissipation.

ii) On the contrary, tangential forces and sliding friction block certain configurations, depending on how they are accessed. The ensemble of blocked configurations is then ill-defined. We have not tried to construct any ensemble, but content ourselves with the observation that $`T_{dyn}`$ is also in this case independent of the particle size (Fig. 2)— we suspect that thermodynamic concepts apply, but that the relevant ensemble goes beyond Edwards’ construction as it stands.

iii) We have tested the validity of the thermodynamics in an ideal homogeneous system with periodic boundary conditions by explicitly avoiding structural features of dense granular flows such as shear bands and segregation of the species. Even though it remains to be seen whether the thermodynamic picture can account for these inhomogeneous effects, our ideal system may prove to be useful in deriving constitutive relations to be used in macroscopic theories of slow granular flows.

To summarize: we have performed numerically an experiment with dense granular systems specially conceived to be a ‘dress rehearsal’ for the real laboratory one. The independence of the Einstein-relation temperature on the tracer provides a strong test for the thermodynamic ideas. We have also showed that this temperature, obtained from measurable quantities, indeed matches Edwards’ one. If, as now seems likely, this result is confirmed in the laboratory, this will give an experimental foundation for the use of the powerful tools of statistical mechanics as the framework to study this kind of far-from-equilibrium dissipative dynamical system.


text:
# Field induced transitions in a kagomé antiferromagnet

\[

## Abstract

The thermal order by disorder effect in magnetic field is studied for a classical Heisenberg antiferromagnet on the kagomé lattice. Using analytical arguments we predict a unique $`H`$$`T`$ phase diagram for this strongly frustrated magnet: states with a coplanar and a uniaxial triatic order parameters respectively at low and high magnetic fields and an incompressible collinear spin-liquid state at a one-third of the saturation field. We also present the Monte Carlo data which confirm existence of these phases.

\]

Geometrical frustration in lattice spin models is responsible for complete suppression of conventional magnetic order and appearance of nonmagnetic spin-liquids or states with exotic order parameters . Applied magnetic field can further enhance frustration. A well known example is a weakly frustrated triangular lattice antiferromagnet, which acquires an additional continuous degeneracy in external field . Investigation of high-field effects in strongly frustrated magnets poses a new challenge for experimental and theoretical studies .

The Heisenberg antiferromagnet on the kagomé lattice (see Fig. 1) is a strongly frustrated spin model, which is approximately realized in a number of insulating layered magnets: SrCr8Ga4O19  and Ba2Sn2Ga3ZnCr7O22  (both with $`S=3/2`$), KFe3(OH)6(SO4)2 ($`S=5/2`$) . Gd3Ga5O12 ($`S=7/2`$) is another frustrated magnet on a related three-dimensional garnet lattice of corner-sharing triangles, which is often called a hyper-kagomé lattice. This magnet has a weak exchange constant $`J1`$ K and a peculiar unexplained phase diagram in magnetic field . Motivated by the above materials with large values of spin we investigate in this Letter the finite-temperature magnetization process of a classical antiferromagnet on the kagomé lattice. We predict three distinctive field regimes below the saturation field $`H_{\mathrm{sat}}`$, where exotic spin phases are stabilized.

The Hamiltonian of a nearest-neighbor Heisenberg antiferromagnet on the kagomé lattice with classical unit spins can be written up to a constant term as:

$$\widehat{}=\frac{1}{2}\underset{}{\overset{N_{}}{}}\left(J𝐒_{}^2𝐇𝐒_{}\right),$$
(1)
where the sum runs over all triangles, $`𝐒_{}`$ is the total spin of a triangular plaquette and $`N_{}=\frac{2}{3}N`$ is the number of plaquettes on an $`N`$-site lattice. The zero-field classical constraint $`𝐒_{}0`$ fails to define a unique ground state. The ground state coplanar configurations are constructed by fixing spins on a first plaquette to a $`120^{}`$-structure with left or right chirality for the triad $`\widehat{𝐞}_a`$, $`\widehat{𝐞}_b`$, and $`\widehat{𝐞}_c`$ and, then, tiling this triad over the whole lattice in a way that three spins on every plaquette are different . The number of all such states for the kagomé lattice is known exactly: $`U^N`$ with $`U=1.1347\mathrm{}`$ . Nonplanar ground states are constructed from planar configurations by identifying so-called weathervane defects .

The thermal order by disorder effect can appear because of a different entropy of short-wavelength fluctuations above degenerate configurations . All coplanar states for a kagomé antiferromagnet show identical harmonic spectra. They have one flat zero-energy branch with $`N_4=\frac{1}{3}N`$ modes, which corresponds to anharmonic quartic excitations . A nonzero harmonic mode described by a classical coordinate $`y`$ has energy increase $`\mathrm{\Delta }E_2Jy^2`$ and contributes $`\frac{1}{2}T\mathrm{ln}(J/T)`$ to the thermodynamic potential, whereas a soft quartic mode with $`\mathrm{\Delta }E_4Jy^4`$ makes a reduced contribution of $`\frac{1}{4}T\mathrm{ln}(J/T)`$. Coplanar configurations have the largest number of soft modes and are, therefore, selected by thermal fluctuations . Soft modes in coplanar states correspond to alternate tilting of spins out of the ground state spin plane around elementary hexagons. There are $`\frac{1}{3}N`$ hexagons on the kagomé lattice. The counting of soft modes from such a geometrical point of view, thus, agrees with the spin-wave analysis.

Harmonic fluctuations do not select between various planar configurations, though the highest statistical weight corresponds to a so-called $`\sqrt{3}\times \sqrt{3}`$ structure . The low-temperature phase has nematic correlations of the chirality vectors defined on every plaquette . Besides, the broken symmetries include selection of the spin triad. The residual symmetry group is determined by the two elements: $`C_6\widehat{𝒯}`$ and $`U_2`$, $`C_6`$ is a spin rotation by $`\pi /3`$ about normal to the spin plane, $`\widehat{𝒯}`$ is the time-reversal, and $`U_2`$ is a rotation by $`\pi `$ about one of the three in-plane axes. The order parameter for low-temperature phases of the kagomé antiferromagnet is a three-point correlation function taken on one plaquette:

$$S^\alpha (𝐫_1)S^\beta (𝐫_2)S^\gamma (𝐫_3)=S^{\alpha \beta \gamma },$$
(2)
which has only a spin part described by a fully symmetric traceless tensor $`S^{\alpha \beta \gamma }`$ . For coplanar states at $`H=0`$ this tensor is parameterized as

$$S^{\alpha \beta \gamma }=\underset{abc}{}\widehat{e}_a^\alpha \widehat{e}_b^\beta \widehat{e}_c^\gamma .$$
(3)
Note, that the three-spin (triatic) order parameter breaks explicitly the time-reversal symmetry, which determines its nontrivial coupling to an applied field.

In a finite magnetic field the classical energy Eq. (1) reaches the minimum provided that

$$𝐒_{}=𝐇/(2J)$$
(4)
for every plaquette. Since $`|𝐒_{}|3`$, above the saturation field $`H_{\mathrm{sat}}=6J`$ all spins are aligned parallel to $`𝐇`$. We are interested in the field range $`0<H<H_{\mathrm{sat}}`$, where the classical ground state possesses an infinite degeneracy. The degeneracy is further enhanced by a field because the minimum energy constraint does not require spins on one plaquette to lie in the same plane. Shender and Holdsworth made the only attempt to understand the finite field behavior of the kagomé antiferromagnet . They noticed that spin coplanarity responsible for the increased number of soft modes is preserved if the spin plane is parallel to an external field. Hence, thermal fluctuations create magnetic anisotropy, which orients at low temperatures the spin plane.

There is an additional effect missed in the above consideration and related to the triatic nature of the spin order parameter. Once selection of the spin plane takes place, an extra macroscopic degree of freedom appears: orientation of the spin triad inside the plane relative to $`𝐇`$. At $`T=0`$ all orientations have the same classical energy, though the spin triad distorts differently, when it forms different angles $`\phi `$ with the magnetic field. This distortion modifies harmonic fluctuations inside the spin plane. In such cases with a one-parameter continuous degeneracy thermal fluctuations always produce the order by disorder effect and select a homogeneous state for one particular value of the parameter . The orientational field effect for the triatic order parameter is described by an invariant $`E_{\mathrm{an}}H^\alpha H^\beta H^\gamma S^{\alpha \beta \gamma }H^3\mathrm{cos}3\phi `$, which is present in the Landau functional because of a broken time-reversal symmetry in the triatic phase. The minimum energy is reached for $`\phi =0`$ or $`\pi `$ depending on a sign of the prefactor, i.e., one of the triad vectors has to be either parallel or antiparallel to the applied field.

To find the equilibrium orientation of the spin triad in external field, we performed spin-wave calculations. The harmonic spectrum of a coplanar state at $`H0`$ depends not only on the triad orientation but also on a way the triad is tiled over the lattice. We used, therefore, the $`\sqrt{3}\times \sqrt{3}`$ and the $`q=0`$ structures  as two example cases. For both structures $`\phi =\pi `$ is selected by fluctuations with $`E_{\mathrm{an}}=5.0T(H/H_{\mathrm{sat}})^3N\mathrm{cos}3\phi `$ and $`1.3T(H/H_{\mathrm{sat}})^3N\mathrm{cos}3\phi `$, respectively. Thus, thermal fluctuations in a weak magnetic field create anisotropy for the spin plane as well as inside the plane. The kagomé antiferromagnet at low fields is described by the coplanar triatic order parameter (3) with one of the spin triad vectors being antiparallel to the field direction, Fig. 1.

Generally, classical ground states at an arbitrary field $`0<H<H_{\mathrm{sat}}`$ are all noncollinear. Collinear configurations appear among the ground state manifold only at special rational values of the applied field. The up-up-down ($`uud`$) structure has $`|𝐒_{}|=1`$ and, hence, is stable only at $`H_c=2J=\frac{1}{3}H_{\mathrm{sat}}`$. The problem of calculating the total number of the $`uud`$ states for the kagomé lattice at $`H=H_c`$ can be mapped onto the problem of dimer coverings of a dual hexagonal lattice, which is solved exactly . The number of $`uud`$ states scales as $`V^N`$, with $`V=1.1137`$… A special role of such collinear spin configurations in magnetic field for various frustrated models has been recently emphasized in . Standard harmonic analysis for the kagomé antiferromagnet at $`H=\frac{1}{3}H_{\mathrm{sat}}`$ yields three two-fold degenerate excitation branches for an arbitrary $`uud`$ state:

$$\omega _𝐤^1=0,\omega _𝐤^{2,3}=3\pm \sqrt{3(1+2\gamma _𝐤)},$$
(5)
with $`\gamma _𝐤=\frac{1}{3}(\mathrm{cos}k_x+2\mathrm{cos}\frac{1}{2}k_x\mathrm{cos}\frac{\sqrt{3}}{2}k_y)`$. Therefore, the collinear states have $`N_4=\frac{2}{3}N`$ soft quartic modes, twice more than any coplanar state. This fact has a simple geometrical origin. Local soft modes correspond to alternate tilting of spins around elementary hexagons. For collinear states such distortions have two polarizations in two directions perpendicular to the field, while coplanar states have soft modes only in the polarization transverse to the spin plane.

Thermal fluctuations reduce the free energy of the $`uud`$ states compared to other noncollinear classical ground states at $`H=\frac{1}{3}H_{\mathrm{sat}}`$. Further lifting of degeneracy within the discrete subset of the $`uud`$ states does not occur due to the same mechanism as for coplanar states at $`H=0`$ . Every collinear state has special lines of weathervane defects, which contain alternating sequences of $`u`$-$`d`$-$`u`$-$`d`$ spins. Flipping all spins along such a line by $`180^{}`$ produces another collinear state at the free energy cost $`\mathrm{\Delta }FT`$. A macroscopic temperature independent number of such defects is generated at low $`T`$’s leaving all $`uud`$ states accessible for the magnet and creating a collinear spin liquid. Enhanced number of soft modes is also reflected in the specific heat $`C`$: out of $`2N`$ degrees of freedom for $`N`$ classical spins, every quadratic mode contributes $`\frac{1}{2}k_B`$ to the specific heat, while a quartic mode does only $`\frac{1}{4}k_B`$ . Hence, $`C=\left[1N_4/(4N)\right]k_B`$ per one spin and we expect $`C=\frac{5}{6}k_B`$ for the collinear spin liquid compared to $`C=\frac{11}{12}k_B`$ for coplanar states. Specific heat measurement in numerical simulations provides a useful way to distinguish different spin states.

The free energy for the collinear spin liquid is lower by $`NT\mathrm{ln}(J/T)`$ compared to any noncollinear state. But a small deviation of field away from $`H=\frac{1}{3}H_{\mathrm{sat}}`$ produces only a little change in the classical energy. The collinear spin liquid is, therefore, stabilized by anharmonic interactions to a finite field range in a way similar to stabilization of the incompressible quantum fluid states in the Fractional Quantum Hall Effect. The collinear spin liquid phase yields a weak magnetization plateau at $`1/3`$ of the saturated magnetization: spins in the $`uud`$ states cannot tilt towards the field and a weak variation of the magnetization for the collinear spin liquid is produced by thermal excitations alone. The above arguments also suggest the first-order transition from the collinear spin liquid state to the low- and the high-field phases at low temperatures.

In the field range $`\frac{1}{3}H_{\mathrm{sat}}<H<H_{\mathrm{sat}}`$ the ground state configurations are again all noncollinear. There appear, however, new quasicollinear states with two coinciding spins of the basis triad $`\widehat{𝐞}_b=\widehat{𝐞}_c`$, which have zero chirality on every plaquette (Fig. 1). These configurations belong to the subset of coplanar states. Hence, a quasicollinear state has all the soft excitations of coplanar states: $`\frac{1}{3}N`$ quartic modes in the polarization transverse to the spin plane. In addition, they acquire extra in-plane soft modes. For a single plaquette, such a clapping-type mode corresponds to simultaneous tilting of the two parallel spins in opposite directions inside the spin plane (Fig. 1). To preserve a weak anharmonic energy change the clapping deformation has to be extended on the lattice along a loop which contains only the $`b`$-type spins. The shortest such loop is a perimeter of a hexagon, while number of appropriate hexagons is maximized if a quasicollinear state is tiled in the $`\sqrt{3}\times \sqrt{3}`$ structure (Fig. 1). Thus, the $`\sqrt{3}\times \sqrt{3}`$ quasicollinear state has the largest number of soft modes $`N_4=\frac{4}{9}N`$ at $`\frac{1}{3}H_{\mathrm{sat}}<H<H_{\mathrm{sat}}`$.

Nevertheless, the true breaking of the translational symmetry does not take place. The disordering mechanism is again related to the presence of loops of alternating spin orientations or weathervane defects . Loop flipping $`a`$-$`b`$-$`a`$-$`bb`$-$`a`$-$`b`$-$`a`$ at zero energy cost contributes $`T\mathrm{ln}(J/T)`$ to the thermodynamic potential because of a destruction of a few soft modes. Such losses are outweighed at small concentrations of defects by the entropy gain. At finite temperatures there will be a macroscopic number of defects, which destroy the long-range $`\sqrt{3}\times \sqrt{3}`$ translational symmetry and create instead a phase with a uniaxial triatic order parameter:

$$S^{\alpha \beta \gamma }=n^\alpha n^\beta n^\gamma \frac{1}{4}\left(n^\alpha \delta _{\beta \gamma }+n^\beta \delta _{\gamma \alpha }+n^\gamma \delta _{\alpha \beta }\right),$$
(6)
where indices run over $`x,y`$ and $`𝐧=\widehat{𝐞}_a(\widehat{𝐞}_a\widehat{𝐳})`$. Higher free energy costs make the number of defects in the uniaxial triatic state to be temperature dependent: $`N_{\mathrm{def}}TN`$ with diverging magnetic correlation length $`\xi 1/T`$. Lacking a conventional magnetic order the high-field triatic phase asymptotically approaches the ordered $`\sqrt{3}\times \sqrt{3}`$ quasicollinear structure. This behavior is reflected in the temperature dependence of the specific heat which should tend to $`C=\frac{8}{9}k_B`$ at low $`T`$.

At the saturation field $`H_{\mathrm{sat}}=6J`$ the ground state becomes nondegenerate and corresponds to parallel alignment of spins. The harmonic excitations in this collinear $`uuu`$ phase are given by the same expressions (5) as for the $`uud`$ states. There are $`N_4=\frac{2}{3}N`$ soft modes at this magnetic field and the specific heat again reaches the value $`C=\frac{5}{6}k_B`$. Soft modes also produce a universal magnetic behavior of the frustrated spin system at the saturation. For $`H>H_{\mathrm{sat}}`$ all three branches (5) acquire an additional field dependent shift $`\delta \omega =(HH_{\mathrm{sat}})`$. Taking into account only the most singular contribution of the first flat branch, we derive that at $`HH_{\mathrm{sat}}`$ the classical partition function scales as

$$Z(H,T)=T^{1/4}f(u),u=(HH_{\mathrm{sat}})/\sqrt{T}.$$
(7)
The susceptibility $`\chi =dM/dH`$ determined from (7) is $`\chi =[f^{\prime \prime }(u)/f(u)][f^{}(u)/f(u)]^2`$. At $`u=0`$ ($`H=H_{\mathrm{sat}}`$) the susceptibility becomes temperature independent and all the curves $`\chi (H,T)`$ cross at one point, see Fig. 2.

So far we considered only the short wave-length degrees of freedom and completely disregarded long wave-length excitations, which determine the nature of phase transitions in two dimensions. The non-Abelian topological defects discussed previously for the triatic order parameter  are suppressed by the field. The two triatic phases break only $`U(1)`$ symmetry for rotations about the field direction. Berezinskii-Kosterlitz-Thouless type transitions separate them from the high temperature paramagnetic phase, see Fig. 1. In contrast, the first-order transition to the collinear spin-liquid state survives and corresponds to a gas-liquid type transition on the boundary with the paramagnetic phase. The collinear $`uuu`$ phase also has an enhanced entropy contribution and is stabilized below $`H_{\mathrm{sat}}`$. The corresponding line of first-order transitions ends at a critical point as sketched in Fig. 1.

We have also performed Monte Carlo (MC) simulations for the model (1). The standard Metropolis algorithm was used with up to $`10^7`$ MC steps per every point. Simulations were done for periodic $`3L^2`$-site lattices with $`L=6`$, 12, 18; the data presented in Fig. 2 are for a 972-spin cluster. The magnetization curve at $`T=0`$ is a straight line with a slope $`\chi =dM/dH=\frac{1}{6}`$. At low temperatures $`\chi (H)`$ develops a dip near $`H=\frac{1}{3}H_{\mathrm{sat}}`$ corresponding to a weak plateau at $`\frac{1}{3}M_{\mathrm{sat}}`$. Two peaks, which surround the dip, indicate first-order transitions to the plateau phase. The peaks become rounded and, then, completely disappear at higher temperatures. We estimated $`T^{}0.025J`$ as the highest temperature, where the collinear spin-liquid still exists (Fig. 1). Above this temperature the dip corresponds to a smooth crossover. The field dependence of the specific heat also follows the above predictions. At $`T=0.002J`$, the specific heat starts near $`\frac{11}{12}k_B`$, corresponding to the coplanar triatic state, at low fields. As the field approaches $`\frac{1}{3}H_{\mathrm{sat}}`$, $`C(H)`$ goes down to $`\frac{5}{6}k_B`$, corresponding to the collinear spin liquid, and, then, recovers back to $`\frac{8}{9}k_B`$, corresponding to the uniaxial triatic state. The specific heat drops again near the saturation with $`C|_{H_{\mathrm{sat}}}=\frac{5}{6}k_B`$, while the minimum value is reached at somewhat lower fields. Such a behavior signals stabilization of the $`uuu`$ phase beyond its classical boundary: bare negative modes below $`H_{\mathrm{sat}}`$ are renormalized into soft modes, which further decrease $`C(H)`$. The specific heat also exhibits a peak at $`H4.8J`$ ($`T=0.002J`$) and $`H3J`$ ($`T=0.01J`$) with a significant finite-size dependence. We attribute this peak to a phase transition from the collinear $`uuu`$ state into a noncollinear uniaxial triatic state (6). To check this we have calculated field dependence of the squared triatic order parameter (2) at $`T=0.002J`$ for three cluster sizes and extrapolated it to the thermodynamic limit. Results are shown in the inset in Fig. 2. The uniaxial triatic phase disappears at about the same field as the position of the peak in the specific heat. The present data do not resolve clearly the order of the phase transition, which must be determined in a detailed numerical study.

We have shown that the geometrical approach based on soft mode counting is a powerful tool for investigation of the field behavior of a classical kagomé lattice antiferromagnet. Our preliminary analysis shows that a Heisenberg antiferromagnet on the related garnet lattice has the same type of the phase diagram, Fig. 1. In particular, the field range of an asymptotically ordered quasicollinear state coincides roughly with a region for a field-induced ordering observed in Gd3Ga5O12 .

I am grateful to J. T. Chalker, E. I. Kats, V. I. Marchenko, O. A. Petrenko, and J. D. M. Champion for useful discussions and to A. A. Snigireva for helpful remarks.