Renormalization group crossover in the critical dynamics of field theories with mode coupling terms

Motivated by the collective behaviour of biological swarms, we study the critical dynamics of field theories with coupling between order parameter and conjugate momentum in the presence of dissipation. By performing a dynamical renormalization group calculation at one loop, we show that the violation of momentum conservation generates a crossover between a conservative yet IR-unstable fixed point, characterized by a dynamic critical exponent $z=d/2$, and a dissipative IR-stable fixed point with $z=2$. Interestingly, the two fixed points have different upper critical dimensions. The interplay between these two fixed points gives rise to a crossover in the critical dynamics of the system, characterized by a crossover exponent $\kappa=4/d$. Such crossover is regulated by a conservation length scale, $\mathcal R_0$, which is larger the smaller the dissipation: beyond $\mathcal R_0$ the dissipative fixed point dominates, while at shorter distances dynamics is ruled by the conservative fixed point and critical exponent, a behaviour which is all the more relevant in finite-size systems with weak dissipation. We run numerical simulations in three dimensions and find a crossover between the exponents $z=3/2$ and $z=2$ in the critical slowing down of the system, confirming the renormalization group results. From the biophysical point of view, our calculation indicates that in finite-size biological groups mode-coupling terms in the equation of motion can significantly change the dynamical critical exponents even in the presence of dissipation, a step towards reconciling theory with experiments in natural swarms. Moreover, our result provides the scale within which fully conservative Bose-Einstein condensation is a good approximation in systems with weak symmetry-breaking terms violating number conservation, as quantum magnets or photon gases.


I. INTRODUCTION
The success of the theory of critical phenomena is based upon a simple observation: systems with very different microscopic details behave in strikingly similar ways when correlations are sufficiently strong. This experimental fact eventually crossed over into theory with the formulation of the phenomenological scaling laws [1][2][3][4], whose key idea is that the only relevant scale ruling the spatio-temporal behaviour of a system near its critical point is the correlation length. Eventually, the great conceptual edifice of the Renormalization Group (RG) tied everything together, explaining why microscopically different systems shared so much at the macroscopic level, giving a demonstration of universality through the concept of attractive fixed points, and providing a method to calculate experimentally accessible quantities, most conspicuously the critical exponents [5][6][7][8].
Employing the same set of conceptual tools in collective biological systems could prove very helpful, given the recent massive flow of hugely diverse empirical data theory has to make sense of. In support of this strategy there is first an empirical observation regarding collective biological systems, namely systems in which a large numbers of units (cells, bacteria, insects, birds, mammals) interact locally in space and time giving rise to macroscopic patterns [9,10]: these systems often exhibit unusually strong correlations, whose spatial range is significantly larger than the microscopic scales [11][12][13][14][15]. Besides, recent experiments on natural swarms found evidence of dynamical scaling, a core mechanisms of statistical physics linking spatial correlation to temporal relaxation [3,16], whose validity in a biological context can hardly be considered a coincidence. Hence, despite the temptation, in front of the arresting complexity of biology, to confine ourselves to describing the specifics, we believe that exploring the path correlation-scaling-RG is a reasonable course of action. The hydrodynamic theory of flocking of Toner and Tu has led the way: it applied field-theoretical methods and the RG to bird flocks, namely collective biological systems in their strongly ordered phase [17][18][19]. Here, we use the RG approach to study the other side of collective behaviour, namely the near-critical disordered phase of natural swarms.
In the biophysics of collective behaviour, a prominent role is played by a class of ferromagnetic theories with continuous symmetries, both in their symmetry-broken phase (flocks), and in the near-critical disordered phase (swarms) [17,20]. When dynamics is taken into consideration, though, this universality class breaks down into smaller sub-classes, as there are different ways to implement the dynamics given the same static probability distribution of the system [21,22]. Dynamical diversity is regulated essentially by two distinct -though related -factors, namely conservation laws and symmetries. On the one hand, we have dynamical theories lacking symmetries and conservation laws (as in the classic Heisenberg model, or Model A of [21]), or in which conservation is imposed despite the absence of an explicit symmetry (as in phase separation, or Model B of [21]). On the other hand, we have theories ruled by symmetry and conservation laws, whose dynamics is characterized by the coupling between two fields, namely the order parameter and the conserved generator of the symmetry, i.e. the conjugate momentum. This second type of theories therefore have non-dissipative mode-coupling terms in the equations of motions, and were originally introduced to describe systems displaying Bose-Einstein condensation (BEC), as superfluid helium, superconductivity, and quantum magnets (Models E, F, and G of [21]). Bizarre as it may seem, recent experiments suggest that some collective biological systems, as bird flocks [23] and insect swarms [20], also have non-dissipative mode-coupling terms in their dynamical equations, and are thus akin to this second class of theories. The connection between BEC systems and flying animals reflects the great generality of the mathematical structure of collective dynamics governed by symmetry and conservation laws, whether the order parameter is the quantum phase of a condensate, or the direction of motion of a flock.
Here we will focus on this second class of theories, with the aim to study the critical dynamics of swarms. To make this introductory discussion more concrete, let us anticipate the actual dynamical field equations we are going to derive and analyze in detail in this work: with effective Hamiltonian, In the biological context the vector order parameter ψ(x, t) represents the velocity field, but it has different interpretations in BEC systems (for example, in liquid helium ψ is the expectation value of the Bose field). In all cases, though, the order parameter is coupled to its conjugate momentum, we call it spin, s(x, t), which is the generator for rotations of ψ, given the rotational symmetry of H. 1 The distinctive trait of this class of models are the mode-coupling cross terms, ∂ t ψ ∼ δ s H and ∂ t s ∼ δ ψ H, which generate a non-dissipative dynamics with the classic coordinate-momentum Hamiltonian structure; were it only for these terms, dynamics would be completely deterministic. On the other hand, the diagonal terms, ∂ t ψ ∼ δ ψ H and ∂ t s ∼ δ s H, give rise to the diffusion and transport phenomenology typical of stochastic statistical systems, and are thus complemented by the noises, θ and ζ, whose variance is proportional to the kinetic coefficients, 2Γ 0 and 2(−λ 0 ∇ 2 + η 0 ), respectively.
The crucial feature of this theory is that, in absence of dissipation, namely when the effective friction η 0 is zero, the total integral of the spin is conserved: the cross term in (2) gives rise to a continuity equation for the symmetry generator, s(x, t), prescribed by Noether's theorem, while the stochastic transport term in (2), λ 0 ∇ 2 s, is still the divergence of a current, leaving the continuity equation intact. This structure -symmetry and conservation -is a very profound feature of this class of models, as it leads to the existence of propagating hydrodynamic modes in the ordered phase, called spin waves; this mechanism give rise to 'second sound' in liquid helium [21], it is responsible for linear information propagation in bird flocks [25], and finally it explains spin-wave remnants in the near-critical phase of insect swarms [20].
Why, then, introducing in equation (2) a dissipative term, η 0 , which destroys spin conservation? In the context of biological systems the answer is quite simple: the symmetry generator, or spin, is conjugated to the velocity field; hence, by rotating the velocity, the spin is what actually makes an animal to turn. Indeed, kinematically one can prove that the spin is related to the radius of curvature of the individual trajectories [26]. Hence, at the individual level it is clear that there must be some dissipation relaxing the spin, thus making a trajectory straight in absence of external perturbations or interaction with the neighbours. On the other hand, in systems like superfluids or superconductors, the conservation law generated by the continuous symmetry of the quantum phase corresponds to number conservation and it cannot be violated. In other BEC systems, though, like quantum magnets [27], exciton condensates [28], and photon gases [29], the Hamiltonian can contain terms that weakly violate the symmetry, hence dissipating the density in the continuity equation of the momentum. The effect of weak dissipation in the ordered phase is simply to generate a damping length scale on propagating spin-waves. However, in the near-critical phase the situation is more complicated: dissipative and non-dissipative models are known to have completely different critical exponents, hence what is the effect of dissipation in this case is unclear. This question is particularly relevant for biological swarms, as experiments found a dynamical critical exponent that cannot be reconciled with the prediction of purely dissipative theories [20].
Here, by using a dynamical renormalization group approach, we study the effect of dissipation on the critical dynamics of a systems with mode-coupling terms. Our calculation shows that the dissipative term η 0 gives rise to an interesting crossover characterized by nontrivial critical exponents. The competition between conservative transport, λ 0 ∇ 2 s, and dissipative friction, −η 0 s, generates a novel conservation length scale, R 0 ; beyond R 0 the dynamics is ruled by a purely dissipative RG fixed point, so that the whole conservative (and propagating) nature of the theory is lost, whereas for distances smaller than R 0 , the conservative RG fixed point governs the dynamics, giving rise to the classic spin-wave phenomenology. We calculate the value of the dynamical critical exponents in these two regimes and of the crossover exponent, and we confirm our results through numerical simulations.
As we shall see, the conservation scale R 0 is larger the smaller the dissipation. The implications of this fact are very important in the biophysical context. The presence of dissipation in the dynamical equations of biological groups may suggest that these systems are in the same universality classes as fully dissipative models, as dissipation always wins over conservative terms in the infinitetime and infinite-distance hydrodynamic limit. However, real biological groups are of course finite-size systems (and quite moderately sized, in the case of flocks and swarms), in which dissipation has been demonstrated by experiments to be quite low [20]. Therefore, the size of these systems may actually be smaller than the crossover scale R 0 , so that, even if dissipative terms are present in the equations of motion, critical dynamics is still ruled by the symmetric and conservative structure of the equations, and therefore have critical exponents drastically different from the dissipative ones. As we shall see, for natural swarms this theoretical mechanism produces a critical dynamics whose phenomenology is remarkable similar to that found in experiments.
Although our motivation is biological, it is worthwhile to remark that our results apply to any BEC system with weak dissipation, a relevant example of which are quantum magnets [27], exciton condensates [28] and photon gases [29]. In quantum magnets Bose-Einstein condensation of magnons occurs at low temperature, due to the spontaneous breaking of the U (1) symmetry; real quantum magnets, though, contain weakly symmetrybreaking terms in their Hamiltonian, thus violating the conservation of the conjugate momentum. Another BEC system our results could be applied to is that of excitons, bosonic hole-particle excitations created by laser pumps, whose number, though, is conserved only within their lifetime (which is finite, and dependent on many factors) [28]. A similar situation arises within the context of photon gases, when a polariton condensate emerges [29]; in this case, too, depending on the polariton lifetime, one can have a violation of the number conservation symmetry, which is equivalent to an effective dissipation. In all these cases, our calculation could provide the crossover scale within which an exact BEC assumption is justified and it may describe the critical behaviour of the crossover.
Here is the plan of the paper. In Section II we will give a derivation of the microscopic dynamical equations in their biological context, whereas in Section III we will coarse-grain the microscopic equations and work out the dynamical field theory described by equations (1) and (2). In Section IV we will perform a renormalization group calculation of critical dynamics in the momentum shell context; this Section will culminate with the formulation of the RG recursive equations, while the analysis of the crossover between the two different fixed points on the critical manifold, and the corresponding crossover of the critical dynamics, will be studied in Section V. In Section VI we will give an alternative derivation of our results using the more field-theoretical Callan-Symanzik approach. In Section VII we will perform numerical simulations to validate the RG results, and finally we will present our conclusions and discuss the outlook in Section VIII. Parts of the most technical material are contained in the Appendixes. A shorter account of our results can be found in [30].

II. BIOPHYSICAL ORIGIN OF THE MICROSCOPIC MODEL
In this Section we derive the microscopic model of collective behaviour that we will use to describe the dynamics of natural swarms. Because this model was first introduced in the context of flocks, rather than swarms, we will have to take a short detour in that direction. At the end of the Section we will discuss under what approximations we will be able to perform a field-theoretical RG study of the model.

A. Collective behaviour and the Vicsek model
Collective behaviour in biological systems, and more specifically collective motion, is essentially a game of mutual imitation, in which each individual tries to make its own state of motion as close as possible to that of its neighbours [9]. From a physical point of view, such mechanism is clearly suggestive of a ferromagnetic-like interaction: if we focus our attention on the direction of motion of each individual, that is on the orientation of the velocity vector, such imitation game amounts to a local interaction due to which each (normalized) velocity vector tends to align to those of its neighbours, much as classical Heisenberg spins tend to align to each other. At variance with standard ferromagnets, though, in collective motion the positions of the particles change in time, as they are carried around by their own velocities, thus creating a non-equilibrium feedback between the alignment degrees of freedom and the interaction network [18,31]. The simplest yet most illuminating model describing this core mechanism of collective motion was introduced by Vicsek and co-workers [32]; it describes a set of self-propelled particles that interact with each other in a ferromagnetic way, where r i is the position of particle i, v i its velocity, and n ij (t) is the (short-ranged) adjacency matrix (who is neighbour of whom) at time t. The interaction between individuals is given by the ferromagnetic term in (4), whereĴ gives the strength of the tendency to align to each other. 2 Such alignment interaction is often called social force in the collective behaviour literature [9]. In the Vicsek model the speed is kept fixed, |v i | = 1, which is the purpose of the cross-products at the r.h.s. of (4). The term ζ i is a Gaussian white noise with variance, whereη is a dissipation coefficient and T a generalized temperature measuring the strength of the noise. The power of the Vicsek model is that it describes collective motion in its two different phases. When noise is low (or density is high, in the metric case [9]), the alignment interaction produces long-range order across the system, forming a polarized moving flock. The interesting thing is that such ordering also occurs in two dimensions, which would be forbidden by the Mermin-Wagner theorem [33] in an equilibrium ferromagnet with continuous symmetry; however, the Vicsek model has an off-equilibrium feedback between alignment and selfpropulsion promoting long-range order [18]. On the other hand, when noise is large enough (or density is low, in the metric case [9]), the system is in a disordered (paramagnetic) phase, which reproduces quite well the statistical properties of real swarms. More precisely, it has been observed that natural swarms are disordered, but highly correlated systems [34]; the velocity static correlations are reproduced (at least qualitatively) by the Vicsek model close to its ordering transition. Hence, the Vicsek model captures rather well the static correlation functions of collective motion for both flocks and swarms. Dynamics is more problematic, though, at both the qualitative and the quantitative level.

B. The Inertial Spin Model
The first hint that the Vicsek equation of collective motion required some new ingredients came from experiments on flocks, in which it was observed that disturbances in the direction of motion of the birds (that is, turns) propagate linearly, with very low dissipation [25]. Although the hydrodynamic field-theoretical description of the Vicsek model introduced by Toner and Tu [17] contains linearly propagating 'sound' modes, caused by the feedback between local density and phase fluctuations [35], experiments indicate that flocks follow a different mechanism: during the propagating event the density displays very weak fluctuations, if any; moreover, the speed of propagation of the wave has been found to be higher the higher the polarization of the group, a feature absent in the hydrodynamic theory of the Vicsek model [17] (see also the discussion in [24]). It was therefore suggested in [25] and [23] that Vicsek dynamics had to be complemented with some non-dissipative inertial couplings between order parameter (the velocity) and a conjugate momentum, in order to reproduce the structure of a conservative Hamiltonian dynamics. The resulting microscopic dynamical equations give rise to the Inertial Spin Model (ISM) of collective motion [23], where the new variable s i represents a generalized momentum conjugated to the velocities v i and it is the generator of the rotational symmetry of the interaction; it is therefore called spin, in an analogy with quantum mechanics. Associated to the momentum s i we have a generalized inertia,χ, which embodies the resistance of a particle to change its instantaneous radius of curvature [36]. One can show that, in the low noise, strongly polarized phase, the non-dissipative coupling between spin and velocity of the ISM generates linear propagating modes of the velocity fluctuations, which match quite accurately the experimental results, including the key relation between speed of propagation and polarization [37]. In absence of a dissipative term, the Hamiltonian structure of the ISM would conserve the total spin, as it happens for any generator of a symmetry. However, one can show that the spin is essentially the instantaneous curvature of the particle's trajectory [25], hence a single particle (or bird, in a flock) would maintain its radius of curvature forever, were the spin strictly conserved. This is quite unrealistic. Rather, it seems reasonable to expect curvature (and therefore spin) to be dissipated in the long run in absence of interaction or external perturbations. For this reason the ISM has also the dissipative term, −ηs i , and stochastic noise, ζ i , granting relaxation of the spin for large times. If dissipation is small, though, and the biological group has finite size, linear waves will still propagate across the system, before dissipation kicks in [23]. In other words, although in the hydrodynamic limit (infinitely large times and distances) the conservative Hamiltonian structure always becomes irrelevant, on the finite-time and finite-size scales typical of biological phenomena the interplay between velocity and spin has crucial consequences on signal propagation. Note, finally, that once dissipation is included in the equations, one can recover the Vicsek model as the over-damped limit of the ISM [23], which is quite reassuring.
The second hint that a model with non-dissipative dynamics was required came from swarms. Swarms of insects are systems apparently completely different from flocks: they show no group-scale coordination, so that their net motion is zero: swarms 'dance' above some landmark in seemingly random fashion [13]. In fact, experiments on natural swarms [13] showed that these systems have strong velocity correlations, indicating that, despite the lack of long-range order, the individuals within these groups are interacting with each other rather intensely, hence driving the system close to an ordering transition; indeed, such static correlations were qualitatively similar to those developed by the Vicsek model at its critical point [34]. More recent experiments [20] showed that swarms in their natural environment exhibit another important property of classical statistical physics, namely dynamic scaling [3]: according to this law, the dynamic correlation function of a system close to the critical point obeys the following relations, where t is time, k momentum, C 0 is the static correlation function, τ k is the relaxation time of mode k, F and f are well-behaved scaling functions, and z is the dynamic critical exponent, ruling how space and time scale with each other. The key idea of dynamic scaling is that the only relevant scale in ruling both spatial and temporal behaviour of a system close to the critical point, is the correlation length, ξ. For k = 0 we obtain, τ ∼ ξ z , a property known as critical slowing down: a system strongly correlated in space must also be strongly correlated in time [21]. Experiments showed that swarms satisfy relations (8) with a dynamic critical exponent z ≈ 1, whereas numerical simulations of the Vicsek model in d = 3 give z ≈ 2 [20]. It must be noted that z = 2 is the exact value of the dynamical critical exponent for a purely dissipative free theory (Gaussian model), and that even in the interacting case the exponent receives only very small (two loops) corrections to the value 2, if the dynamics has only dissipative terms (or even values larger than 2, as in the case of Model B [21]). On the other hand, dynamical models with non-dissipative inertial terms tend to have values of the exponent z significantly smaller than 2, as a result of the interplay between order parameter and conjugate momentum [21]. Hence, the low value of z in natural swarms was a further indication of the need for non-dissipative dynamics also in these non-polarized systems. Finally, dynamic relaxation in the Vicsek model has a classic exponential form [20], while natural swarms display a completely different shape, showing clear evidence of non-dissipative inertial behaviour for short times. More precisely, if we define the relaxation form factor [20], h =Ċ(t/τ )/C(t/τ ), in the limit t/τ → 0 we have that h → 1 for the Vicsek model (as for any exponential correlation function), while experiments showed h → 0 for natural swarms, as it would happen in a weakly damped harmonic system, where inertia dominates over dissipation [20]. Hence, the whole dynamical behaviour of swarms seems to require the existence of non-dissipative inertial terms in the equations of motion, which is exactly the extra ingredient the ISM has compared to the Vicsek model. Our plan is therefore to study the critical dynamics of the ISM in its disordered yet near-critical phase, to describe the collective behaviour of natural swarms of insects, in order to try and reproduce the experimental results of [20]. Because the ISM was originally introduced to describe the dynamics of flocks, it has been studied extensively in its deeply ordered (i.e. polarized) phase, both numerically [23], and theoretically [37], while no previous studies of the ISM in the near-critical regime have been performed.

C. Fixed network approximation
Before we proceed with the coarse-graining of the model, though, we need to decide whether to attack directly the full-fledged off-equilibrium problem, which includes the self-propelled nature of collective motion, or whether we first take on the simpler (and yet nontrivial) equilibrium problem, in which particles sit on a fixed network and thus have a time-independent interaction matrix. For a number of reasons, we will follow the second strategy. The model we want to study differs from previous known cases in two main respects: i) it contains non-dissipative terms and effective friction, the interplay of which has never been studied before, not even at equilibrium; ii) the model is a self-propelled one, hence intrinsically off-equilibrium, which may seem particularly important in the swarm phase, in which each particle changes the local neighbours quite rapidly. Our central experimental concern is to reproduce the correct dynamical critical exponent z and the correct relaxation form factor in natural swarms. The fact that the self-propelled, off-equilibrium Vicsek model in its swarm phase gives exactly the same exponent and form factor as equilibrium fully dissipative models (as the classical Heisenberg model), suggests that self-propulsion is not the primary source of the anomalous critical dynamics of swarms. Moreover, we believe that having under control the equilibrium problem puts us in a more solid position to tackle the off-equilibrium one in the future, much as knowing the physics of the equilibrium XY and Heinseberg models has been fundamental to fully understand and appreciate the Vicsek model. Hence, we will study a fixed-network version of the ISM, in which the particles belong to a lattice and the connectivity matrix does not depend on time. In this context, the order parameter no longer has the role of a physical velocity, hence we will call it ψ, the generic symbol for the order parameter, and write the microscopic model in the following way, The modulus of the order parameter is still fixed to |ψ i | 2 = 1 and the adjacency matrix n ij now corresponds to a fixed interaction network. Thanks to this approximation, dynamics can now be rewritten in Hamiltonian terms, with microscopic Hamiltonian, This is the dynamical model we now proceed to coarsegrain in order to obtain a dynamical field theory.

A. Equations of motion
Since we are interested in describing the large scale behavior of the system, it is convenient to pass from a microscopic description in terms of site dependent variables to a field description, where we consider smoothly varying velocity and spin fields ψ(x, t), and s(x, t), obtained by coarse graining the original variables over a small spatial volume. Upon coarse-graining, the original Hamiltonian (13) gives rise to an effective field Hamiltonian H[ψ, s] that -as in standard ferromagnetic systems -reads [38,39] where r 0 is the bare mass (negative in the ordered phase), u 0 is the bare static coupling constant and χ 0 is the effective inertia. Here, we remind, the gradient term comes from the mutual alignment interaction, which favours smoother configurations; the quadratic and quartic contributions for ψ represent a confining potential and derive from the original constraint on the ψ i and the coarsegraining entropy; while the field s remains Gaussian as its microscopic counterpart. When writing down the dynamical equation of motion for the fields, we need to take into account the presence of both the reversible and dissipative contributions present in the microscopic dynamics (11) (12), and add additional dissipative terms, which might arise upon coarse graining. Under very general assumptions [21,40], we can therefore write where the noise correlations are chosen to have a Boltzmann-like static probability distribution, i.e.
Here Γ 0 , η 0 and λ 0 are the bare kinetic coefficient of the field ψ, the bare friction coefficient and transport coefficient of the field s, respectively, while g 0 is a modecoupling constant that regulates the reversible dynamical terms and describes the symmetry properties relating the two fields: the fact that s is the infinitesimal generator of rotations of ψ is indeed specified by the Poisson commutation rules, where αβγ is the Levi-Civita antisymmetric symbol. 3 The static properties of the model only depend on the Hamiltonian H. For the field ψ they are therefore the same as in the Heisenberg model [21], with an ordering phase transition occurring for r 0 = r c . On the other hand, at the static level s is a trivial, purely massive, Gaussian field. Since there is no static coupling term between this field and the order parameter, the inertia χ 0 will not acquire any perturbative contributions; hence, in order to simplify our notation, we choose the units of s such that χ 0 = 1. The dynamic properties are ruled by the transport coefficient λ 0 , by the effective friction η 0 , by the kinetic coefficient Γ 0 and by the dynamic coupling constant g 0 ; these quantities will take perturbative contributions arising from the dynamic interaction between s and ψ, which is ruled by g 0 .
Equations (14)-(15) have two additional dissipative terms compared to the microscopic theory of Eqs. (11)- (12), namely −Γ 0 δH/δψ and λ 0 ∇ 2 δH/δs. The first term actually contains two contributions: first, a derivative of the confining potential, ψ 2 + ψ 4 , which is the coarsegrained analogue of the microscopic sharp constraint, |ψ i | 2 = 1; second, a diffusive piece, Γ 0 ∇ 2 ψ(x, t), which derives from a loss of reversibility due to the coarsegraining, and which describes the role of the fluctuations of the order parameter in the relaxation process; even though such fluctuations are negligible in the low temperature phase (where we recover the microscopic theory with Γ 0 = 0) they are crucial when considering the system close to the critical point. For this reason, even though we neglected the Laplacian of the order parameter in our previous analysis in the deeply ordered phase [37], we need to take it into account in the present study of the critical regime.
On the other hand, the origin of the spin transport term, λ 0 ∇ 2 s(x, t), is perhaps less intuitive. In the microscopic model the spin is dissipated by the friction through the term −ηs i (t), hence one might have expected just a term −η 0 s(x, t) in the coarse-grained theory. Why then are we introducing the term in λ 0 ∇ 2 s(x, t)? We will show in the following Sections that, in the context of perturbation theory and the renormalization group, such term arises naturally from the non linear interaction between the two fields once a momentum shell integration is performed. It is then necessary to include the transport coefficient λ 0 directly from the starting field equations.
We notice that for η 0 = 0, equations (14) and (15) coincide with those of Model G (antiferromagnet), or, in the planar case, of Model E (liquid helium), which have a fully conserved spin dynamics and whose critical dynamical properties have been studied long ago in a series of seminal papers [41,42]. The renormalization group analysis described in the following sections will show that in an appropriate regime (when η 0 is small), the ISM displays the same critical behaviour as these fully conservative models.

B. Free theory in Fourier space
The starting point to build the perturbative expansion of the equations of motion is the free (or non-interacting) theory, which is obtained by setting to zero the non-linear coupling constants, namely g 0 = 0 and u 0 = 0. In Fourier variables, hence using momentum k and frequency ω, the free equations of motion become, The free theory is linear and it is therefore possible to solve it exactly by merely inverting equations (18) and (19), where the free propagators (or Green functions) are the inverse of the dynamical operators in Fourier space, The propagators describe the response of the fields to noise and to external perturbations [43]. We can also define the free dynamic correlation functions, By using (20) and (21), and the noise correlators (16), we get the relations, These four quantities, propagators and correlation functions, are the building blocks of the perturbative expansion.
Calculations in the RG context are carried out in Fourier space, hence all relevant integrals are performed over the momentum, k; in the infinite size limit the lowest extreme of integration is k = 0 (otherwise, in finite-size systems, it is of order 1/L), whereas the upper extreme of integration is a momentum scale -the so-called cutoff -indicated by Λ, corresponding to the inverse of the length scale over which the coarse-graining has been performed; practically speaking, if the continuous field has been obtained by averaging the microscopic variables on a volume of linear size L, we have Λ = 1/L. In principle the coarse-graining is supposed to be performed over a scale much larger than the lattice spacing, a; in practice, though, L is still a microscopic length scale of the system, so that, broadly speaking, one often assumes that Λ is of order 1/a.
The cutoff is an arbitrary scale, which therefore appears as an extra unknown parameter of the theory.
In fact, all bare parameters in the theory, r 0 , u 0 , Γ 0 , λ 0 , η 0 , g 0 , depend on Λ, and therefore they are all equally unknown. As we shall see, the central idea of the renormalization group is to exploit constructively the arbitrariness of Λ, by studying how the bare parameters change when Λ is changed; from this flow the critical properties of the theory will emerge.

IV. RENORMALIZATION GROUP IN MOMENTUM SHELL
Broadly speaking, the renormalization group is a set of symmetry transformations that are useful to determine the scale invariance properties of a system at its critical point [5]. An RG transformation unfolds through two stages: i) integration of the short wavelength details and, ii) rescaling of length and time. The first operation amounts to integrating the fields over large values of the momentum, Λ/b < k < Λ, where b is a rescaling factor larger than -but close to -1; this integration interval is the so-called momentum shell. The effect of integration is to shift the cutoff from Λ to Λ/b; the RG idea is that the long-distance physics of a system close to the critical point, where the correlation length is large, cannot change due to an arbitrary change of the cutoff. Hence, the second stage consists in rescaling space (and consequently time) in such a way to formally restore the original cutoff Λ and to compare the newly obtained equations to the original ones. The compound effect of these two stages is to effectively change the parameters that appears in the equations of motion, hence determining a flow in the space of parameters. At the critical point, where the correlation length is infinite, the RG transformation must have left the system exactly the same, and therefore the fixed points of the RG flow provide all the important information on the large scale physical properties of the system.
The RG technique is nowadays standard and discussed in the literature both for static and dynamical problems [22,40,44,45]. In this section we adopt a momentum shell renormalization scheme [5], as this is the approach that was used in the original papers on critical dynamics [21,41,46]. This will allow us to immediately spot differences with respect to the fully conservative case [41]. In Section VI we will illustrate how the same results can be obtained using a Callan-Symanzik approach, more common in recent applications of the dynamical renormalization group [40].

A. Integration of the short-wavelength details
In the first stage of the RG we integrate out short wavelength fluctuations, namely modes with Λ/b < k < Λ. This operation (described in Appendix D) leads to a new effective theory that only depends on fields fluctuating over larger wavelengths, k < Λ/b. In the free theory, modes at different wave vectors are independent, so this operation has no practical effects. On the other hand, when non-linear interactions are present, the coupling between long and short wavelength modes makes it impossible to carry our exactly this operation, which therefore requires a perturbative expansion that we describe in detail in Appendix A. The bottom line result of shell integration is to produce additional terms in the equations of motion that effectively modify the coefficients of both the linear and the non-linear terms. We start with the linear dynamical coefficients, namely Γ 0 , λ 0 and η 0 . These parameters are contained in the free propagators, eqs. (22) and (23). Due to the shell integration, the propagators acquire some new contributions, the socalled on-shell self-energies, Σ b and Π b , so we can write, As we said, to compute the self-energies one uses perturbation theory. To carry out the expansion we follow the generating functional approach of Martin-Siggia-Rose [22,42,47], where averages of physical observables over the stochastic dynamics are rewritten as thermal averages over a functional measure. The complication is that new auxiliary fields must be introduced and the effective field theory therefore involves four fields, rather than two. The advantage is that the standard Feynman technique can be used to perform a diagrammatic expansion, and perturbation theory can be carried out in the same way as in equilibrium statistical field theory. The details can be found in Appendix A. To one loop order, the self-energies read, The self-energies modify the poles of the propagators in the frequency plane, therefore affecting the way response and correlation functions decay in time. In particular, the k → 0 and ω → 0 expansion of the self-energies and of their derivatives modify the kinetic and transport coefficients, so that we can define their renormalized values, First of all, we notice an important point: from (31) we immediately see that Π b (k = 0) = 0, and therefore we conclude that the effective friction η 0 has no perturbative corrections. As discussed in Appendix A, from the diagrammatic point of view this is a consequence of the structure of the vertex, which makes all perturbative contributions to Π b (k = 0) equal to zero; therefore this result is valid to all orders in perturbation theory. Physically, this fact is a consequence of the rotational symmetry of the theory: even though η 0 breaks the conservation law of the spin, the symmetry is still at work, implying that it is impossible for the conservative theory to produce a non-conservative friction through coarse-graining. The accuracy of the perturbation expansion is increased by substituting the bare mass r 0 with its renormalized value r [45], which represents the inverse static susceptibility and goes to zero when the systems approaches the critical temperature. Since we are interested in the critical behavior, from now on we will evaluate all integrals at r = 0, namely at the critical point; given that integrals are on the shell there are no infrared singularities and the self-energies are finite. We thus have, These equations show that there is a great difference in the role of the two parameters λ 0 and η 0 . If we start from a frictionless model that has η 0 = 0, the coarse-graining of the RG will not generate a friction, η R = 0, through shell integration. On the contrary, even if we start from a model without spin transport coefficient, λ 0 = 0, integration over short wavelengths inevitably generates a transport term λ R = 0. In other words, the interaction between the spin s and the primary field ψ generates a non-zero transport coefficient λ R even if λ 0 = 0 in the original microscopic theory. For this reason, we included from the outset the parameter λ 0 in the coarse-grained field equations 4 .
To make further progress we must address a rather crucial algebraic detail. While the integral in (36) is straightforward, the one in (35) requires some care. The cutoff Λ is large, while the rescaling factor b is close to 1; hence, the shell integration is performed over large values of the internal momentum p. If the effective friction η 0 is finite (or zero), then the term of order p 2 dominates over η 0 at the denominator, so that the overall integrand will have a 1/p 4 behaviour for large momentum. As we shall see later on, the hypothesis that η 0 < ∞ is by no means harmless, and we shall need to return over this point. Yet, for now we will work under this hypothesis, and recognise that it is therefore convenient to rewrite the integral as, We can now change variables, defining p = Λx, and obtain, It is now convenient to introduce a set of effective parameters, through which we can express these integrals in a simpler way (so to speak), where K d is the unit sphere volume in dimension d. We can thus finally write the perturbative expression of all three kinetic parameters after shell integration, where we have introduced the dimensionless crossover parameter X 0 , and where to compute the integrals we have exploited the fact that in the limit b → 1, that is for an infinitesimal RG transformation, the shell becomes infinitesimal, so we have written the integrals as the shell thickness, 1−1/b ∼ log b, times the integrand evaluated at x = 1.
The dimensionless parameter w 0 is rather harmless, and it will play only a moderate role in what follows; quite conveniently, it will remain finite in all fixed points we will find. On the other hand, f 0 is crucial, as it acquires the role of the effective dynamical coupling constant: the perturbative expansion, which naively one would think as a series in powers of g 0 , is in fact a series in powers of f 0 . From the dimensional form of f 0 , and in particular from the fact that it contains a term Λ d−4 , RG connoisseurs can already deduce that the dynamical upper critical dimension of the theory will be d c = 4, the same as the static one. This will be made explicit once we will have solved the recursive RG equation for f 0 further on. We recall that this is a consequence of the fact that both integrands go like 1/p 4 for large momenta, hence giving a logarithmic behaviour at d = 4, and that, in turns, this is a consequence of having assumed that η 0 is finite in the integral of equation (35). We will return on this hypothesis later on.
The second important effective parameter emerging from the equations is the length scale, R 0 , given by the ratio between the transport coefficient and the effective friction of the spin. Because of its definition, we can intuitively expect that if R 0 is large (η 0 small) the dynamics of the spin is ruled by a conservative diffusion mechanism. On the contrary, if R 0 is small (η 0 large), we expect the dynamics of the spin to be ruled by a dissipation mechanism. It is worth noticing that for R 0 = ∞, namely η 0 = 0, the crossover parameter is equal to 1, and one correctly gets the same equations as the fully conservative Model G. On the other hand, when the conservation length scale is very small, R 0 ∼ 0, which happens for η 0 λ 0 , one gets X 0 ∼ 0, so that Γ 0 receives very weak perturbative corrections at one loop. We will return on this crucial point later on and we will see that this interplay between non-conservative dissipation and conservative transport coefficient of the spin plays a key role, giving rise to a non-trivial crossover between two different RG fixed points with different dynamic critical exponents.
Up to now we focused on the coefficients of the linear terms in the equations of motion, Γ 0 , λ 0 , η 0 . However, in general shell integration produces corrections to all terms, including the non-linear ones. Therefore, the dynamical coupling constant, g 0 , could in principle get a perturbative correction from shell integration, in particular from the renormalized vertex. However, it can be shown that -due to the structure of the interaction ver-tices -there are no perturbative corrections at all orders (see Appendix B about vertex corrections), As in the case of the lack of corrections to η 0 , this result is a consequence of the symmetry properties of the system. Indeed, we remind that the field s is the generator of the rotational symmetry of the ψ field. Even though the global spin is not conserved in our case, the symmetry still generates some Ward identities that protect g 0 at all orders (see Appendix C).
Finally, let us note that the static coupling constant, u 0 , is renormalized as usual in the standard equilibrium theory [45]. However, the lowest order corrections to the dynamical coefficients due to the vertex u 0 are at two loops [21], whereas we perform here a one loop calculation. Therefore we do not need to address static renormalization any further in what follows. This is actually a nice feature of all theories with mode-coupling terms: the dynamical vertex is triple (see Appendix A), hence one obtains very sizeable corrections to the critical exponents already at one loop order, without using the two-loop corrections of the static vertex.

B. Rescaling of space and time
After integration over the shell, we are left with a theory which has new renormalized parameters, and also a new, smaller cutoff, Λ/b. In order to compare the new theory with the old one, and therefore to be able to write a set of recursive equations for the parameters, we rescale space, and therefore momentum k, by a factor b, in such a way to formally restore the original cutoff, Λ/b → Λ. It must be noted that frequency does not have a similar cutoff, hence in principle we would have no formal need to rescale ω; however, this is deceiving: in order to reabsorb all powers of b in the novel equations of motions one can see that it is necessary also to rescale frequency [21]. Physically, this means that to the rescaling of space and time cannot proceed independently: space and time are tied together by the -yet unknowndynamic critical exponent z, through the following rescaling relations, Due to the non-linear form of the equations of motion, we do not a priori know how the spatial integration affects the dynamics and we therefore allow for a generic dynamic critical exponent, z. As we shall see, its value determines how the order parameter relaxes close to criticality.
Rescaling momentum and frequency actually means changing physical units, hence each parameter will also rescale according to its naive physical dimensions, which can be expressed in powers of b in the following way, Because the perturbative contribution of the shell integrals are given in terms of the effective parameters, f 0 , w 0 , R 0 , it is also necessary to write their corresponding rescaling laws, For d > 4 the naive dimension of the effective coupling constant f 0 becomes negative; as we shall see this implies that at its non-trivial fixed point the effective coupling constant will be of order , with, as it happens to the static coupling constant, u 0 [8], confirming the fact that the dynamical upper critical dimension is d c = 4. We notice that, in general, the rescaling of k and ω implies also a rescaling of the fields, which also carry some physical dimensions. However, once again we remark that the current calculation is at one-loop level, whereas the fields renormalize at two-loop level; hence, in the following, for all practical purposes all will proceed as if the fields do not renormalize. 5

C. Renormalization group recursive equations
The two stages described above, shell integration and rescaling, must now be put together to define one step of the RG transformation; in this step, a generic parameter P is brought from its initial bare value, P 0 , to a new value, P 1 through an RG equation with the structure, where the power of b comes from the rescaling step, so that D P is the physical dimension of P, whereas the term in the bracket comes from the shell integration. For an infinitesimal RG transformation b ∼ 1, hence we can write 1 + δ P log b = b δ P , so that δ P is an effective correction to the naive scaling dimension D P of the parameter. Of course, δ P , which is the result of the shell integration, will depend on all the other parameters of the theory. One can then iterate this step l times, giving rise to a recursive RG equation for P, where we emphasise that all integrals that appear at the r.h.s. through the factor δ P l , must be evaluated at the running value of the parameters, namely at their value at the RG step l, whereas the naive physical dimension D P is fixed once and for all. By using this procedure for the dynamical parameters of our theory, we obtain the following RG recursive relations, and we recall that we are working at T = T c , namely on the critical manifold. From these equations we can finally write a closed set of recursive relations for the effective coupling constant f l , the dimensionless parameter w l , and the conservation length scale R l , where X l depends on R l and w l through equation (45). We note that the full scaling dimension of the conservation length scale R is determined by its naive dimension, b −1 , plus a perturbative contribution, 1+ 1 4 f l log b = b 1 4 f l , hence developing an anomalous scaling dimension that will be crucial in ruling the crossover.
The derivatives of f , w and R with respect to (− log b) are called beta-functions, and measure how the parameters change when performing an infinitesimal RG transformation, 6 The zeros of these functions define the fixed points of the RG flow and thus have a crucial role in the theory. The beta-functions also will provide a link between the momentum shell RG approach followed so far and the Callan-Symanzik approach described in Section VI.

D. Fixed points and dynamic critical exponent
The fixed point values of the RG equations (that we are going to indicate with an asterisk) rule the critical behaviour of the system. The exponent z can be found by requiring that the fixed point value of the kinetic coefficient of ψ, namely Γ * , is finite [21,22]. This condition, as we shall see, is what we need to investigate the relaxation behavior of the field ψ close to criticality. From equation (53) we obtain z as: The dynamic critical exponent is therefore given by the fixed point values of the parameters f , w and X.
From the corresponding recursion equation (57) it is evident that R can have two fixed points, namely: Since the fixed point of f is expected to be of order (see Eqs. (58)), the scaling dimension of R is negative. Therefore, the R * = 0 fixed point is IR-stable while the R * = ∞ fixed point is IR-unstable: any large, but finite, initial value of R 0 , decreases under the RG equation (57), driving the systems to the R * = 0 fixed point. Inserting back the possible values of R * in the other equations, we therefore find two fixed points for the global set of parameters.
1. The IR-unstable conservative fixed point The first fixed point with R * = ∞ and X * = 1, which we call IR-unstable (or conservative), is: This fixed point describes a dynamics with z = d/2, typical of conservative models such as Model G and Model E [41]. Indeed, dissipation becomes irrelevant (η * = 0), and the conservation law expressed by the symmetries of the Hamiltonian, rules the dynamics at all scales. If the system has R 0 = ∞ (i.e. η 0 = 0) the RG flow will converge to this fixed point, the only stable one for zero dissipation. However, as mentioned above, any other value of R 0 will cause the flow to eventually converge to the other fixed point.

The IR-stable dissipative fixed point
The second fixed point is characterized by R * = 0, or equivalently X * = 0, and we call it IR-stable (or dissipative): In this case, dissipation takes over (η * = ∞) and the dynamic critical exponent that we obtain is z = 2, which is common for models with a completely dissipative dynamics [21]. What we have depicted here is a scenario that includes the presence of two fixed points with different dynamical behaviors and different dynamic critical exponents, namely z = d/2 (conservative dynamics) and z = 2 (dissipative dynamics). Even though one of such fixed points is unstable along one direction, it is stable along the others and -as it will be discussed in the next section -it can rule the RG flow at intermediate iterations. In other terms, there is a crossover in the RG flow in parameter space and, as a consequence, also in the behavior of physical observables.

A. RG flow on the critical manifold
To investigate the dynamic crossover we studied the RG flow from a numerical point of view. In the limit of infinitesimal RG transformation (b → 1), the recursion relations (57) become a system of coupled differential equations. We introduce the continuous variable x = l log b; Eqs. (57) can then be rewritten in the continuum limit (replacing, for instance, f (x) = f l ): where the prime stands for a derivative with respect to x.
The set of Equations (61) can be studied numerically, for any given initial condition. In Fig.1 (upper panel), we show the resulting flow in the (X, f ) plane, each line corresponding to a different set of initial values of η 0 , λ 0 , Γ 0 (and therefore of f 0 and X 0 ). The flow always proceeds from the conservative to the dissipative fixed point, as expected. However, how fast it does so, depends on the initial condition X 0 . When this value is close to 1, which means that the friction η 0 is small, the flow of parameters approaches the z = d/2 fixed point and remains close to it for many RG iterations. Then, it eventually moves towards the stable fixed point with z = 2. In the lower panel of FIG.1 we show the same dynamic crossover in terms of z, of the coupling constant f and of X. From this figure we clearly see that there is a well defined intermediate regime where the flow is regulated by the unstable fixed point, the parameter X driving the dynamic critical exponent from one value (z = d/2) to the other (z = 2).

B. Crossover in the critical dynamics
What we have discussed so far is the crossover taking place along the RG flow in parameters space. This crossover has important observable consequences in the relaxational behavior of the system. In the previous section we showed that there is a parameter with the dimensions of a lenght-scale, R, which plays a crucial role in the RG flow. As we shall see, it is precisely the interplay between R and the relevant physical lenght-scales in the system, to determine the way it relaxes.
1. Crossover in k at ξ = ∞ Let us start by considering the system at the critical point (r = 0). In this case the correlation length is infinite and the only physically relevant length scale is the one at which we are observing the system, namely 1/k. To study relaxation at this scale, we can measure the characteristic frequency of the order parameter, which is directly connected to the dynamic critical exponent through the dynamic scaling hypothesis (8), which, for ξ = ∞ and for small k, reads, where we have introduce the characteristic frequency, ω c , as the inverse of the relaxation time. Since there are two possible values of z one can wonder at this point which is the one to consider in this relationship. It turns out that this depends on k vs. R. To see this, we remind that the characteristic frequency is the pole of the propagator of the field ψ and we can therefore study its infrared behavior by looking at G −1 ψ at small k. This is very convenient because -by construction -propagators along the RG flow are related to each other, and we can link what happens in parameters space to the physical behavior of the system.
At every step l of the RG, the physical propagator verifies the relation, where with P l we indicate the set of the parameters after l steps of RG (and P = P 0 in the l.h.s.). The scaling factor on the right side of (63) is just the scaling dimension of the propagator. What we are doing is to consider an initial point in parameters space corresponding to our physical system (l.h.s.), and then follow the RG flow in the critical manifold starting at that point. As l increases, the propagator on the r.h.s. is evaluated at farther points along the RG line. Since we know that there is a crossover along the RG flow we are writing this expression with a dynamic critical exponent, which explicitly depends on the recursion step l. If we choose l such that b l = Λ/k, i.e. the maximum possible value, the inverse of the propagator satisfies: Here we have evaluated the function on the right side at the fixed point values of the parameters P. This is justified if l is large enough to approach the vicinity of a fixed point (i.e. small k). Which one of the two fixed points is reached -and therefore the value of z * above -depends on the starting point (i.e. the set P) and on the precise number of iterations. More precisely, the condition that discriminates between the two possible fixed points is: because it determines the value of the variable X l in expression (59). Let us consider the situation, which interests us more, where the starting point of the flow is close to the IR-unstable fixed point. The initial value of R is therefore large, corresponding to a system with low η 0 . If R l Λ −1 holds for all the iterations, the flow will explore only the neighborhoods of the unstable fixed point and the values of P * in (64) are the ones of the conservative dynamics. Therefore, in this case: However, it may happen that, even starting at the same point in parameters space, the number of iterations is so large that eventually the condition R l Λ −1 becomes satisfied, and the flow approaches the stable fixed point corresponding to z = 2. In this case: Since the number of iterations is fixed by the value of the wave-number k (b l = Λ/k), the condition R l Λ −1 can be translated into a condition on k. The recursion relation for the conservation length scale gives R l = R 0 b l(−1+f * /4) ; since we are considering a flow starting close to the conservative fixed point, we can set f * = , which gives the anomalous scaling dimension of R at the conservative fixed point, from which we see that the length scale R has a scaling dimension equal to its naive dimension at the upper critical dimension, d = 4, as expected. From the relation R l = R 0 b l(−d/4) we can finally identify a threshold value k c marking the limit between the two different scenarios described above, namely: To summarize, we therefore find that at criticality the relaxation behavior of the order parameter -as captured by the critical exponent z -depends on the relation between the scale at which we observe the system and the value of the length-scale R 0 , i.e.
We therefore have found the third non-trivial critical exponent of the theory, namely the crossover exponent [22], which, as we have seen, is intimately related to the anomalous dimension of the conservation length scale.
In Fig.2 we show the regions corresponding to the two dynamical behaviors in the (k −1 , R 0 ) plane.
2. Crossover in ξ at k = 0 In many cases, and in particular when looking at experimental data, real systems are not exactly at the critical point. For all practical purposes we need to extract information on the critical behavior of the system also when its correlation length ξ is finite, even if large. Predictions can be obtained following a reasoning much similar to the one above, but taking explicitly into account the dependence of the propagator on temperature, i.e. on the correlation length. Besides, since there is a characteristic What we are doing is, again, to consider a point in parameters space corresponding to our physical system, and to relate the physical propagator with propagators of models along an RG line starting at that point. The difference with the previous case is that now the RG flow takes place off the critical manifold, therefore not only the parameters change upon iteration, but also the correlation length, i.e.
with ξ 0 = ξ (i.e. the correlation length of the physical system). We can choose the number l of iterations such that b l = ξΛ. If ξ is large enough that the system comes close to a fixed point, then the inverse propagator satisfies the relation: Since the pole of the propagator for k = 0 is the global characteristic frequency of the system, we immediately get the relaxation behavior as ω c (ξ) ∼ ξ −z . As before, the value of z depends on which one of the two fixed points is approached at the end of the RG flow after l iterations. The discriminating condition is always R l Λ −1 . For R l Λ −1 the fixed point is characterized by z = d/2, then the characteristic frequency diverges as: For R l Λ −1 , the other fixed point is reached (z = 2) and we have: Since the number of iterations is fixed by ξ (i.e. b l = ξΛ), the discriminating condition R l Λ −1 now identifies a threshold value ξ c for the correlation length that can be obtained using the recursion relations of both R and ξ: To conclude, we therefore find that critical slowing down is ruled by two different critical exponents depending on how large the correlation length is (i.e. how close the system is to the critical point) with respect to the conservation length-scale R 0 , i.e.
thus giving the same crossover exponent as in the k description. A graphical representation of the different critical regimes can be found in Fig.2. To summarize, it is therefore the interplay between the correlation length ξ and the conservation length-scale R 0 that defines what kind of critical dynamical behavior is observed. We also note that -due to the non-trivial recursion relation for R (see Eqs. (57)) -these two lenghtscales rescale differently upon RG transformations. As a consequence, the region corresponding to the conservative critical dynamics is larger than in the case of naive scaling.

C. A new upper critical dimension
So far we have been studying the RG flow in the vicinity of the conservative, z = d/2, fixed point. Our original motivation was indeed to describe experimental findings on swarms of insects, where a low-dissipation critical dynamics has been observed. As we have seen, in the neighbourhood of this fixed point we have an upper critical dimension d c = 4 and the effective dynamic coupling constant is the parameter f 0 . However, we also showed that the conservative fixed point is unstable, hence the RG flow inevitably brings the system to the dissipative fixed point, z = 2. The problem is that, in the vicinity of this fixed point, the on-shell self-energy Σ b has to be treated quite differently form the previous case, and f 0 does not play the role of the effective dynamic coupling constant anymore. Let us see this in more detail.
In proximity of the IR-stable fixed point, the running effective friction η l becomes very large, eventually diverging. In this regime, our previous assumption to have a mild, finite value of the friction η l in equation (35) must be revised, and the integral must be rearranged differently, We see that, as the running friction η l goes to infinity, approaching the stable fixed point, the large p behaviour of the integrand turns from 1/p 4 to 1/p 2 , thus giving, From this last equation it is clear that, in the proximity of the IR-stable fixed point, the actual effective coupling constant in the perturbative expansion of Γ l , is no longer f l , but q l = g 2 l Λ 2−d /Γ l η l , whose naive scaling dimension is d − 2, not d − 4; accordingly, the integral now has a logarithmic UV divergence at d = 2. We conclude that the upper critical dimension for this fixed point is no longer 4, butd c = 2, and that the actual parameter of expansion is,˜ In d = 3, which is the case of interest for us, the dimensions of q is negative, which is equivalent to say that the only stable fixed point is q * = 0 (this can also be seen explicitly by writing the RG recursive equations for q l ). Therefore, the self-energy contribution in (80) vanishes and the kinetic coefficient has no perturbative contributions (at one loop), thus giving, so that the only way to keep finite the kinetic coefficient at its fixed point is to have, in agreement with the previous result. In this regime ψ behaves dynamically as an independent field, i.e. its relaxation has no contributions from the mode-couplings term in the equations of motion. This very non-trivial crossover between two different upper critical dimensions will be made more explicit in the Callan-Symanzik approach, which we describe the following Section.

VI. CALLAN-SYMANZIK APPROACH
In this section we derive the RG results within a different renormalization approach, in which the large scale properties of the system are deduced from a differential equation (called Callan-Symanzik equation or renormalization group equation). This equation in turn follows, as we explain below, from a reparametrization invariance of the renormalized dynamic theory, which is introduced to deal with the strong cutoff (Λ) dependence of the original theory (which leads to divergences in the Λ → ∞ limit). This approach is complementary to the momentum-shell renormalization developed in Sec. IV. Its principles are described in several texts, e.g. [39,[48][49][50]. Our treatment of the ISM under the Callan-Symanzik approach follows the lines of the dynamic renormalization study of model E by De Dominicis and Peliti [42] (see also [50]). The CS approach involves the following steps: 1. Write a renormalized theory, i.e. reparametrize the original dynamic functional in a way that all Λdependence (equivalently, divergencies that appear for Λ → ∞) of physical observables is absorbed into a finite set of constants. This is done at an arbitrary momentum scale µ.
2. Using the fact that the renormalization can be done at arbitrary values of µ, write a differential equation describing how relevant renormalized observables (in our case the response and correlation functions) change as µ = |µ| is varied. This is the RG equation, sometimes called CS equation. Combining this with dimensional analysis, one finally obtains a differential equation that describes the change of the renormalized observable as the observation scale (external momentum) is varied. The coefficients of this equation are the β-functions, which are computed from the renormalization constants at a given order in perturbation theory. The equation is solved by the standard method of characteristics.
3. The solution by the method of characteristics shows that the behavior of the response and correlations functions at large scales can be obtained by studying the response/correlation at a reference observation scale but with scale-dependent coupling constants. How the coupling constants change when increasing the observation scale is ruled by the β functions, and the trajectories in parameter space induced by a change of scale are called RG flow. Thus one finally studies the RG flow, with particular attention to fixed points, which will lead to scaling behavior of the response.
We describe these steps in the following subsections. Many aspects of the calculation are identical to Models E and G and for these we refer to the article by De Dominicis and Peliti [42]. We only describe in detail the aspects that are novel in the ISM.

A. Renormalized theory and renormalization factors
The diagrammatic expansion of the dynamical action involves integrals in momentum space that are divergent (in the space dimension of interest) for large integration momenta (ultraviolet divergences) unless some regularization procedure is adopted (like the cutoff Λ for large momentum we used in the momentum-shell calculation, see Appendix A). To construct a renormalized theory means to reparameterize the functional in terms of a different set of coupling constants and fields in such a way that the divergences (or equivalently the details of the regularization procedure) are confined in a finite set of constants.
Instead of using a cutoff, here we renormalize according to the dimensional regularization plus minimal subtraction prescription [50]: diverging integrals are evaluated in a dimension low enough that they are convergent, then analytically extended to non-integer dimension. The original divergences then show up as poles in the dimension variable. The minimal subtraction procedure consists in introducing the renormalized parameters so that they absorb only those poles.
Renormalization thus starts with the identification of all the ultraviolet divergences of the theory and with the definition of the renormalized constants to absorb them. Looking at the perturbative expansion, we see that the introduction of η 0 leaves the free propagator of the ψ field (22) unchanged with respect to the model G case, while in the free propagator of the s field G 0,s (k, ω) a kindependent term is added, so that the k → ∞ behavior of the free propagators is unchanged. Then, since the structure of the diagrams is identical to that of models G and E (because the interacting part is the same), the divergences in ISM arise in the same diagrams. Then from ref. [42] we know that the theory is renormalizable in d = 4 (which is the upper critical dimension of the theory). The divergent diagrams relevant to the dynamic renormalization arise in the expansion of G ψ (k, ω) and G s (k, ω), in particular in the derivatives Both divergences are logarithmic in d = 4. There are two additional divergences in G ψ (k, ω) that we do not need to consider. One is the quadratic divergence in G −1 ψ (k = 0, ω = 0) that is absorbed into a renormalized mass (susceptibility) in the static theory. Since we work here at the critical point defined by r = 0, in practice this means setting r 0 = 0 in all the diagrams we consider. There is also a logarithmic divergence in ∂G −1 ψ (k, ω)/∂ω that however does not arise at the 1-loop level (and which is related to field renormalization).
The divergences are taken care in the following way: we consider the relevant divergent quantities (e.g. the derivatives in Eq. (84)) and evaluate them at ω = 0 and at a given value of the momentum k = µ (so as to eliminate infrared divergences). We then replace the original kinetic/transport coefficients and coupling constants by renormalized counterparts that absorb the divergences, in a such a way that -once expressed in terms of the new parameters -the quantities of interest are finite. The renormalized parameters are defined through multiplication by Z-factors; when considering the derivatives in (84) this amounts to introduce renormalized kinetic coefficients The two remaining dynamic couplings, g 0 and η 0 do not pick up perturbative renormalization. In the case of g 0 this is a consequence of a Ward identity deriving from the fact that s generates the rotations of ψ (Appendix B and C) [42]. In the case of η 0 the reason is that it is not involved in absorbing divergences due to the fact that G −1 s (k = 0, ω = 0) is finite (see next section). We introduce however η and g as adimensional counterparts of η 0 , and g 0 , where K d = 2π d/2 (2π) −d /Γ(d/2) is introduced for convenience and µ is the arbitrary momentum scale used to evaluate the propagators during renormalization. (Note that in this section we choose the frequency units so that Γ 0 and λ 0 are adimensional, i.e. [ω] = [k 2 ]). The Z-factors now have to be determined at a given order in perturbation theory so that all the renormalized propagators (and in consequence correlation and response functions) are finite, i.e. the Z-factors are divergent in a such way that all observable quantities (expressed as averages with the renormalized theory) are finite. One can then in principle determine all the renormalized parameters of the model in terms of a finite number of observations (at fixed wavenumber and frequency). Finally, since at one-loop, as mentioned above, the fields are not renormalized, the relation between original and renormalized propagators and correlations is G ψ (k, ω, u 0 , g 0 , Γ 0 , λ 0 , η 0 ) = G R ψ (k, ω, u, g, Γ, λ, η; µ) C ψ (k, ω, u 0 , g 0 , Γ 0 , λ 0 , η 0 ) = C R ψ (k, ω, u, g, Γ, λ, η; µ) We discuss the determination of the Z-factors in below. These play a leading role in the CS procedure since they give the nontrivial contributions to the RG equation coefficients (sec. VI B).

B. RG equation and dynamic critical exponent
In this approach the dynamic critical exponent is identified after finding the dynamic scaling form of the correlation functions. The procedure is the standard one, which we briefly recall. From the perturbation expansion at ω = 0, one can see that the correlation function C ψ can be written in terms of the static coupling u 0 , the effective dynamic couplings (f 0 , w 0 , R 0 ) and Γ 0 , where Γ 0 and ω always appear in the combination ω/Γ 0 . So we can rewrite Eq.(88) as (89) Since the lhs is independent of the arbitrary scale µ, deriving with respect to log µ one obtains the RG equation: where l = u, f, w, R and the β and ν functions are The only dimensional arguments are k, ω and µ, and the dimension of C R ψ is 2. Then dimensional analysis leads to an Euler equation which can be used to eliminate the µ derivative: Restricting ourselves to changes in the scale of k, i.e. k = µ/b, we have that k·∇ k = −b∂/∂b. Then combining (90) and (92) we get (we have omitted terms that only appear beyond one loop). We solve (93) using as initial condition b = 1, i.e. the value of the correlation at a reference k = µ, at some frequency ω, and at the physical values of the couplings Γ, u, f , w, R. The solution, found by the method of characteristics, is and where ν Γ depends on b through the couplings f, w, R. The dependence of these on b is given by the functionŝ u(b) etc. (the running coupling constants), which are the solution of the system where the β functions must be computed perturbatively from the relation between the bare and renormalized couplings (85). At one loop, the flow of the static coupling constant u is completely uncoupled from the dynamic couplings, so we do not take it into account in the following.
The meaning of (94) is that the correlation function at the physical values of the couplings u ≡ (u, f, w, R) and at a rescaled wave vector µ/b is equal to the response function at the original scale µ but evaluated for different couplingsˆ u(b). Fixed points thus are sets of coupling values u * such that all β functions vanish simultaneously: it is clear that if the flow starts at such a point, or approaches it for some large value of b, it will stay there for all larger b. If in addition the function C R ψ (k, ω/Γ, u) is continuous at u = u * , then all the k-dependence at large b (small k) is contained in ωb 2 /Γ(b): equation (94) is then the scaling law we seek, and we can read off the scaling behavior from its second argument even if we don't know the form of C R ψ . For example, if the flow is near a fixed point for i.e. the value of η Γ at the fixed point gives the correction to the naive dynamic critical exponent. So we proceed next (secs. VI C and VI D) to deterine the Z-factors that furnish the β-functions, and then (sec. VI E) to find the fixed points of the flow (96) and their infrared (i.e. b → ∞) stability. The infrared stable fixed points will rule the scaling behavior at large lengthscales. Unstable fixed points may, depending on initial conditions, lead to transient scaling laws observable in certain regimes.

C. Determination of Z-factors
We must determine the two dynamic Z-factors Z λ and Z Γ (in a two or higher loops calculation a third factor, related to field renormalization, would arise, but we do not need it here). First Z λ is fixed by requiring that ∂G −1 s /∂k 2 | ω=0,k=µ be finite. We have where Π is the same self-energy 7 as in equation (31). From Eq.(87) we then have From this equation two conclusions follow: the first is that the λ 0 k 2 term cannot be left out from a renormalizable theory. This is a consequence of the fact that, in an expansion of Π in the external wavevector, it is the k 2 coefficient that is divergent, not that of k 0 (in fact Π(k = 0) = 0). Thus if λ 0 is absent, it is useless to define η = Z η η 0 and try to absorb the pole of Π into Z η , because η drops out from (99). This is equivalent to the finding, in the momentum-shell scheme, that the renormalization transformation generates a λ coefficient even if it is absent in the original theory. The second conclusion is that Z λ is determined solely by the behavior of Π, and, since at one loop this self energy is unchanged with respect to model G, we can without further discussion write it from the model G result [42]: At one loop, the differences between model G and ISM are only found in Z Γ , which we proceed to compute now. The propagator of the ψ field is where Σ is the self-energy (30), which we recall here for convenience: We must now consider the renormalized derivative (103) 7 Notice, however, that from now on all the integrals in k in the self-energies will no longer be performed on-shell, but rather between 0 and ∞, and for this reason we drop the subscript b from the self-energy symbols. and choose Z Γ so that it is finite. In the dimensional regularization procedure this means that the Z-factor cancel the poles that appear for d = 4 (i.e. terms proportional to 1/ , = 4 − d), so that (103) is free of poles. Computing the derivative from (102) for general dimension d = 4 − at the critical point and ignoring the contribution from a convergent integral one has where R 2 = λ/η = µ 2 Z λ λ 0 /η 0 ,μ = µ/µ and, (105) The renormalized counterpart can therefore be written as, The one-loop term has a pole in d = 4 (from the second Γ function): this is how the original divergence of the integral in d = 4 manifests itself in dimensional regularization. The minimal renormalization prescription stipulates that this pole be identified so that an equivalent pole but with opposite residue can be added to 1/Z Γ , thus making the renormalized vertex finite. So we expand the second term: setting I(R, w) As long as R = 0, all singular behavior is contained in the pole at d = 4, i.e. the term proportional to 1/ in the last line of (107). Then defining Z Γ = 1 + (2/ )f /(1 + w) renders G R ψ finite. Thus naively one finds that Z Γ is independent of R, and, since Z λ is also independent of R at one loop, this leads to β-functions for the parameters f and w that are independent of R, and thus to flow equations identical to model G for f and w, uncoupled to the flow of R. However, this is wrong: we have already seen, in the momentum-shell scheme, that the presence of η 0 profoundly affects the flow of f and w, with the notable macroscopic consequence of a change in the dynamic critical exponent.
Even without the insight we have from the momentumshell calculation, one could guess that the naive expectation cannot be right: since R 0 has the dimensions of a length, one expects that its stable fixed point is 0, and indeed below we shall find from the β-function (115) that R ∼ b −1+f * /4 for b → ∞, with f * of order . In the above equations, one sees that the limit R → 0 requires special treatment: in (107) a logarithmic divergence appears in the expansion of I(R, w), and even before expanding one sees that (105) is problematic because I ψ (R → 0, w) vanishes for d < 4, while lim R→0 lim d→4 I ψ (R, w) = 1.
The difficulty here is that the length scale R occurs in the Gaussian part of the dynamical functional, and that it appears in the loop integrals in such a way that their convergence properties change at one of the fixed points of R. To deal with this we use a generalized minimal substraction as discussed by Frey and others [51,52]. This method involves incorporating the singular R dependence into the renormalization Z-factors (and consequently into the β-functions): we thus "enrich" the pole with a crossover factor X(R, w) extracted from I(R, w) which absorbs the singular R → 0 behavior. Let us rewrite (106) as where X(R, w) is such that X −1 (R, w)I ψ (R, w) (and in consequence all coefficients of its -expansion) is wellbehaved for all values of R (including R = 0, R = ∞). We discuss in the next subsection how to fix this factor, but before let us write the renormalization factors including the as-yet unknown X(R, w): We conclude this subsection writing ν Γ , ν λ , and the β functions for the couplings f , w and for R. These functions determine the RG flow and the asymptotic scaling properties of the observables (sec. VI B). Recalling the definition (91) of β functions and ν exponents, and developing the to first order in f we get so that the β-functions (to be compared to (58)) are the following, D. Determination of the crossover factor X To determine X(R, w), the idea is to to fix it in such a way that the renormalization factor Z Γ contains all the singularities near both fixed points R = 0 and R = ∞. When R is nonzero, the only singularity in (106) is the pole at d = 4 originating in the Gamma function Γ(2 − d/2). Thus the first condition we impose is that lim R→∞ X(R, w) = 1.
To find another condition, we must study G −1 ψ for vansishing R. To do this let's define a new parameter and rewrite (104) as where we have introduced a constant m 2 to avoid an infrared divergence in d = 2. We can now set R = 0 in the integral to find We now find a pole at d = 2, corresponding to the fact that for R = 0 the integral in (118) diverges for d ≥ 2: the critical dimension for q is d c = 2, not 4. We then expand around d c = 2, and find where,˜ We note that β q above is correct up to second order in˜ , i.e. around the new critical dimension d c = 2. Equation (121) can give us the condition to impose on X(R, w) near the other fixed point: from (113) and (115) we can write so that to recover (121) we impose A simple choice for X(R, w) is then . ( with a(d = 2) = 1. Referring to (108), we find so that X(R, w) defined as in (125) with a = /2 fulfills the three conditions i) X(R → ∞) = 1, ii) X(R → 0, d → 2) = (1 + w)R 2 , and iii) makes X −1 (y, w)I(y, w) finite and non-vanishing for all values of R, so that all terms not included in Z Γ (i.e terms of order 0 or higher) are regular and don't cause further trouble. In particular lim R→0 lim d→4 X −1 I ψ = lim d→4 lim R→0 X −1 I ψ = 1. A subtle but important point however remains to be made. At face value, our choice of X(R, w) recovers the flow (121) only at d exactly equal to 2, while (121) has actually been obtained in general dimension. Thus it seems that one would want a = 1, which however accounts for the R → 0 behavior of I(R, w) only at d = 2. The way out of this seeming inconsistence is to remember that (121) is valid in general dimension but only up to second order in˜ = 2 − d = − 2. This means that one should actually write this exponent as a = 1 +˜ /2, and (125) as It is then clear that (since the nonzero fixed point of q will be of order˜ ) (123) indeed recovers (121) in the limit R → 0 for general dimension up to order˜ 2 . But then it also becomes clear that (121) cannot fix X(R, w) at order˜ and beyond. So a = 1 +˜ /2 is fine near d = 2 and also near d = 4 (where (106) must be expanded).
From these considerations, in what follows, we will simply set a = 1 when writing the β functions also in d = 3, which is the dimension of interest here. We do so because this is the simplest choice that gives a correct description of the flow at the one-loop level: when R → 0 the O(˜ ) contributions to X(R, w) cannot be fixed without going to two loops so we may as well omit them, and on the other hand near the R → ∞ fixed point X(R, w) → 1 independently of the exponent. So the final form of β functions we find within the CS scheme is identical to equations (58) we previously obtained with the momentum-shell technique.

E. Fixed points
Equations (115) and (96d) show that R * = 0 and R * = ∞ are fixed points of R. The flow can be solved formally asR Since we expect (and confirm below) that the fixed point of f will be of order , we see that the integrand within the exponential is negative (for large b) at least, and that R * = 0 is IR-stable while R * = ∞ is IR-unstable.
The R * = ∞ fixed point corresponds to X(R, w) = 1. Since it is unstable, it is only relevant for b → ∞ when the system starts at 1/R = 0. This corresponds to a very important special case, namely η 0 = 0, i.e. model G (equivalent to model E for what concerns the scaling properties). It is also relevant at moderately large scales for 1/R very small, when the flow stays near X = 1 long enough that the other couplings approach the model E fixed point before R becomes so small that X is significantly different from 1 (see secs. V.D and VI F).
Equations (96b) and (96c) for X = 1 were studied by DeDominics and Peliti [42], who considered model E at two loops, and we refer to them for the analysis of the fixed points and their stability. In summary the relevant (IR-stable in the (f, w) subspace) fixed point is implying, It is interesting to remark that the result for z, although obtained here at one loop, has to be valid at all orders in perturbation theory as long as w * and f * are different from 0: if w is non-null, (114) implies ν Γ = ν λ , and (113) then gives 2ν Γ = − .

R * = 0 -dissipative fixed point
When R = 0, X(R, w) = 0 which gives immediately regardless of w and f (as long as they are finite). Setting X(R, w) = 0 in (96b) and (96c) one finds two solutions: and The stability of the two fixed points can be studied linearizing the flow (96) around the fixed point u * : where W is the Jacobian matrix evaluated at the fixed point. The fixed point is stable if W is positive definite, i.e. its eigenvalues are all positive. We find for the two cases above which has eigenvalues 0, 1, and − , while with eigenvalues , /2 and 1 − 2 . So the only IR stable fixed point (at 1 loop, near d = 4) is f * = 2 , w * = 0, R * = 0, which implies as we have seen that the critical exponent is z = 2.

F. Crossover
We have just concluded that the only IR-stable fixed point gives z = 2, so that for large observation scales (k → 0) the critical dynamics is like that of a purely dissipative model (like model A) when the starting (physical) value of η is nonzero. However, for non zero but small η, such that the starting R is very large, the initial value of X(R, w) is very close to 1 and will stay so until R(b) is of order one (e.g. for w = 3, X is larger than 0.99 for R > 5). So one can expect that f and w will at first move as if X = 1, i.e. towards the conservative (model G) fixed point, staying in its neighborhood until R decreases significantly, and in effect the numerical study of the flux (Fig. 1) confirms this expectation. Then experimentally one will observe model G critical behavior (z = d/2) for moderate (i.e. not too small) values of k, possibly lasting a rather wide interval, until at some point for k → 0 the asymptotic z = 2 exponent will be seen. We show here how to obtain the scaling of the wavevector k c (marking the end of the model G behavior) with the physical value of R (cf. sec. V E).
Assume then that the physical value of R isR(b = 1) ≡ R 1 1 so that X(R 1 , w(1)) ≈ 1. Assume also that one is observing at a scale k = µ/b such that the flow has already reach the neighborhood of the z = d/2 fixed point (which in particular implies ν Γ ≈ − /2). We ask how small we must make k so that the system moves away from this fixed point and the scaling law changes. For this to happen, X must be significantly less than 1, so let's impose that X < c, with c = 0.99 say. This requires R(k) < R c = c/[4(1 − c)]. Now since we are near the conservative fixed point we can use (96d) with Now the crossover wavevector will be such that R(k c ) = R c , so that Hence, we find the same crossover exponent as in momentum shell, namely, Let us notice that the crossover exponent κ is nontrivial: from naive dimensional analysis one would have guessed k c ∼ R −1 1 . However, the renormalized R is dimensionless, and the RG result is actually taking into account the nontrivial effects of the hidden microscopic lengthscale.
Finally, let us note that the crossover exponent derives its value from ν λ at the model G fixed point which, as we have mentioned before, takes the value − /2 at all orders in perturbation theory [42], and thus so must the crossover exponent.

VII. NUMERICAL SIMULATIONS
To test our results, we performed numerical simulations of the microscopic ISM model on a fixed lattice in d = 3. We implemented the dynamical equations (9) (10) using a generalized Verlet integrator [53][54][55] for second order equations. Details of the algorithm can be found in [56]. The lattice spacing is Λ −1 = 1. We fixed the parametersĴ = 1,χ = 1, and performed simulations at several values of the temperature T and of the friction coefficientη. Since the temperature sets the correlation length and the friction regulates the conservation length scale R 0 we can in this way explore the (ξ, R 0 ) plane of Fig.2. For all values of T andη considered, we computed the correlation length ξ and the relaxation time τ , and inferred the exponent z from the scaling relation Eq. (8) between them. In this way, we could investigate the dynamical critical behavior and compare results with the predictions of the RG computation. Before illustrating the results, let us briefly explain the procedure followed to compute the main quantities required for our analysis, namely ξ and τ .

A. Static behaviour and determination of ξ
Since we are interested in the critical behavior of the system, we need first of all to locate the transition temperature and characterize the static critical properties of the system. To do so, we perform numerical simulations of Eqs.(9) (10) in the stationary regime, and use decorrelated dynamical configurations to compute equal time equilibrium averages (from now on indicated with · · · ). From the static point view the ISM on a lattice is completely equivalent to a standard ferromagnetic model, we therefore expect static properties to reproduce the well known results of the Heisenberg model. For a system of N velocities/vectors, it is possible to define the polarization as:  (N = 512, 1000, 2197, 4096, 8000). An ordering transition occurs at approximately Tc 1.5. b) Susceptibility as a function of temperature, same sizes as in panel a); the maximum of each curve is located at a temperature that decreases with increasing the size of the system, and approaches the critical temperature Tc in the thermodynamic limit. c) Finite size scaling of the susceptibility. Curves at N = 2197, 4096, 8000 satisfy finite size scaling with exponents ν = 0.707 and γ/ν = 1.973, as predicted by the theory of the Heisenberg model [39].
measuring the degree of global alignment, and its modulus, the scalar polarization φ. The average value of this quantity is plotted in Fig.3a as a function of temperature, and clearly shows the occurrence of an ordering transition. The critical temperature can be conveniently located by looking at the fluctuations of φ , namely the susceptibility where β is the inverse of the temperature. We analyzed these quantities for a wide range of temperatures (0.1 ≤ T ≤ 5.0) and sizes (N = 512, 2197, 4096, 8000). From Fig.3a,b we can see that the critical temperature is located approximately at T c 1.5. The critical point moves towards lower temperatures as the linear size of the system increases, in according to finite size scaling.
To measure the correlation length ξ we first computed the static connected correlation function C(r): where r ij is the distance between two sites i and j, and δψ i = ψ i − ψ i . Since we are mostly interested in the paramagnetic phase of the model, the relevant one to describe experimental data of insect swarms, we focused on temperatures T > T c , approaching the critical point from above. The behaviour of the correlation function (not displayed) is as expected for a Heisenberg model, we therefore computed the correlation length from the expression rC(r) = exp(−r/ξ), exploiting the fact that the anomalous dimension is small [39]. We combined this information with data on the susceptibility to obtain an estimate of the ratio between critical exponents γ/ν for sizes N = 2197, 4096, 8000. Simulations at N = 8000 give γ/ν = 1.905, in agreement with the literature [39]. We therefore used this size of the system for all following analysis. Finally, to further test the equivalence of the static properties of ISM with the Heisenberg model we performed a finite size scaling analysis on the susceptibility, as displayed in Fig.3c.

B. Dynamic behaviour and determination of τ
To investigate the dynamical behavior of the system one has to look at time dependent quantities. In particular, the characteristic time scale τ is by definition the scale over which fluctuations of the order parameter become decorrelated. To compute it, we introduce the spatio-temporal correlation function, that is: with T max the length of the simulation. The number of operations needed to calculate this quantity is in general ∼ T max N 2 ; however what we actually need for the scaling analysis is the correlation function at k = 0 (see previous section), which is numerically less demanding: From this quantity, we computed the characteristic time scale τ from the condition, This condition corresponds to requiring that half of the total integrated area of the dynamic correlation function in the frequency domain comes from the interval −ω c < ω < ω c , with ω c = 1/τ . This definition of τ has the advantage of capturing the relevant time-scale both when relaxation is dissipative, and when propagating modes are present, and it is the standard definition adopted in the literature on dynamic critical phenomena [16].

C. Dynamic crossover
Our primary objective is to observe the crossover in the dynamic critical behavior predicted by the RG computation. The simplest protocol to do that would seem to fix the value of the dissipation coefficient (and therefore R 0 ), and extensively vary the correlation length by tuning T . In the (ξ, R 0 ) plane of Fig.2 it corresponds to a straight horizontal line crossing from the red conservative region on the left to the dissipative green one on the right. In numerical simulations, when plotting τ vs. ξ, we should then observe two different power laws, one with exponent z = d/2 for small ξ, and another one with z = 2 at large ξ. The problem with this protocol is that to see a power-law crossover one should span several orders of magnitudes in ξ; three decades is the very minimum, but L = 10 3 , gives N = 10 9 in d = 3, which is quite awful, considering that the largest relaxation time would be of order, τ ∼ ξ 2 ∼ L 2 ∼ 10 6 . This is not possible, and the maximum size we used is well below (L ≤ 20). In other terms, as illustrated in Fig.4, we can in practice explore only a finite horizontal interval in ξ included between the lattice spacing and the maximum achievable ξ, which is not long enough to appropriately measure the exponents in the two dynamical regimes (segment b in Fig.4). Interestingly, the same figure shows that for many values of R 0 the accessible interval does not even intercept the crossover line, but entirely lies within the same dynamic regime (segments a and c in Fig.4). In this case, in numerical simulations at fixedη we should expect that only one power-law is observed in the τ vs. ξ plot. For this reason, the best way to capture the presence of the crossover is to run simulations at different values of the effective friction. If the picture above is correct, when η is small, corresponding to large R 0 , we should measure z = 3/2 (segment a), while for large enoughη we should measure z = 2 (segment c).
The numerical findings fully confirm this scenario. In Fig.5 we show results for three different sets of simulations, respectively forη = 1, 2, 4. We cannot use larger values for the effective dissipation, because the maximum relaxation time becomes too long to equilibrate the system. For the smaller values (η = 1, 2), the data are in good agreement with a dynamic critical exponent z = 3/2, while forη = 4 the characteristic time scales with the correlation length with an exponent z = 2. We therefore conclude that the ISM exhibits a dynamic crossover in critical behavior, as predicted by the RG approach.
To further support the existence of two distinct regimes in dynamical behavior, we tested the full dynamic scaling hypothesis (8) on the dynamic correlation functions. In Fig. 6, upper panels, we display the normalized C(k = 0, t) for all the temperatures that we analyzed, and for two different values ofη. In the lower panels, we report the same curves but plotted as a function of the rescaled variable t/ξ z , where we used the values of z obtained from the previous analysis. The figure shows that dynamic scaling is nicely verified, but with different exponents (z = 3/2 and z = 2, respectively) at small and large values of the friction coefficient.

D. Natural swarms and inertial dynamics
Both theoretical computations and numerical simulations describe a dynamic crossover between two different critical regimes, which is ruled by the interplay of the correlation length ξ vs. the conservation length scale R 0 . Our analysis has important consequences when considering systems of finite size. What we have shown is that, even in presence of dissipation, the critical behavior of the system can be ruled by a conservative critical dynamics with exponent z = 3/2 (in d = 3) in an extensive region of parameters (case R 0 > ξ 3/4 ; red region in Fig.4). This result is particularly relevant if we think back at the biological motivation of our study: explaining experimental data in natural swarms of insects. As discussed in Section II, swarms exhibit dynamic scaling, but with an exponent smaller than the one predicted by models of collective motion with a purely dissipative dynamics. This is why we considered the ISM in the first place: to put back inertial terms in the dynamical equations and to understand whether they can produce a z < 2 on the collective scale of living groups. The answer to this question is therefore yes. The exponent that we get in the conservative region, z = 1.5 is not yet the value observed in the data (z = 1), but it is a big step forward as compared to the prediction of the Vicsek model (z = 2). This strongly indicates that the ISM captures an important ingredient -inertia -absent in previous models. Numerical simulations of the ISM also reproduce another feature measured in natural swarms, which is not reproduced by previous models. Experimental correlation functions in natural swarms display a concave shape at short times, incompatible with the exponential relaxation predicted by the Vicsek model [20]. The ISM on the other hand displays the same kind of behavior as in the swarms data. To show this, in Fig.7a we compare the dynamical relaxation of natural swarms with simulations of the ISM and of the Vicsek model in the paramagnetic phase. We can see that ISM reproduces the curvature of the experimental correlation for t → 0, contrary to the Vicsek model. The consistency between ISM and natural swarms becomes even more striking when we compute the relaxation form factor [20], The limit of this function for t → 0 is equal to 1 if the  The relaxation form factor h(t/τ ) ≡Ċ(t/τ )/C(t/τ ), goes to 1 for overdamped exponential relaxation, while it goes to 0 for inertial relaxation [20]. The fixed-network ISM reproduces the correlation form of real swarms in a rather compelling way.
dynamics is purely exponential, as in the case of the Vicsek model. On the other hand, an inertial dynamics approaches zero for small times: this is the case of ISM, and of natural swarms. We therefore conclude that the ISM in the paramagnetic phase qualitatively well describes the inertial dynamics of natural swarms.

VIII. CONCLUSIONS
We have performed a one-loop RG calculation of the critical dynamics of a statistical system with inertial nondissipative couplings, in presence of a dissipative term which violates the conservation law of the symmetry generator. Our calculation was motivated by recent experiments on the collective dynamics of natural swarms of insects [20], although the dynamical field equations we studied are relevant also for BEC systems with terms weakly violating the symmetry in the Hamiltonian [27][28][29]. We find that the RG flow has two fixed points, a conservative yet unstable one, and a dissipative stable fixed point, associated to the dynamical critical exponents z = d/2 and z = 2, respectively. The crossover between the two fixed points is regulated by a conservation length scale, R 0 : for scales much larger than R 0 , the dynamics is ruled by the dissipative fixed point, while for scales smaller than R 0 critical slowing down is governed by the conservative fixed point. Numerical simulations on the microscopic model confirm our results.
The crossover length scale, R 0 , is determined by the ratio between the transport coefficient, λ 0 , and the effective friction, η 0 , of the spin field. If the coarse-grained parameter η 0 is certainly connected to its microscopic counterpart,η, in the original model, the same cannot be said for the transport coefficient, as in the original microscopic model there is no transport term. The interesting fact, then, is that the conservative transport term, λ 0 ∇ 2 s, is explicitly generated by the renormalization group, through the spin self-energy Π at one loop. Therefore, we are in one of those rare cases in which a crucial length scale of the system, i.e. R 0 , cannot be guessed purely on the basis of dimensional analysis of the microscopic equations of motion (possibly with some renormalized anomalous dimensions). Of course, one could have guessed (admittedly rather smartly) that the presence of a symmetry and conservation law, albeit violated by η 0 , should require a conservative transport term. But in case our intuition were not so good, the RG would require by itself the existence of such term, and therefore the emergence of a crossover length scale, thus confirming its power in dictating what is relevant and what is not in strongly correlated systems.
The fact that the crossover length scale R 0 is larger the smaller the dissipation has important consequences for biological systems. Real biological groups always have finite size, hence in order to study their behaviour we cannot just take for granted the hydrodynamic limit (infinitely large times and distances), but we have to cope with the actual size of the system. In both flocks and swarms, experiments have shown that dissipative terms are rather weak [20,25], hence suggesting that the conservation length scale R 0 is quite large. Under these circumstances one may have a conservation length scale that is larger than the system's size, R 0 > L. In this case, one expects to find a dynamical critical exponent equal to that of the fully conservative RG fixed point, namely z = 3/2 in d = 3, and a dynamic correlation function with strong signature of non-exponential inertial relaxation. Thanks to this finite-size critical crossover, the fully conservative phenomenology should hold at all practically attainable values of the correlation length, which is always limited by the system's size.
From the point of view of the comparison between theory and experiments in natural swarms, our calculation therefore puts us in a semi-satisfactory situation. Certainly we can say that the form of the dynamical correlation functions of natural swarms, and in particular the non-exponential inertial nature of the short times dynamics, is rendered by the ISM in a much more compelling way than the Vicsek model (Fig.7); actually, our simulations show almost no quantitative difference between theory and experiments in this respect. Concerning the dynamical critical exponent, z, the situation is still open, although we would say that the result of the present calculation -namely z = 3/2 in finite-size weakly-damped 3d swarms -definitely goes in the right direction. Experiments give z ≈ 1, even though values up to z = 1.2 would probably be acceptable, given the noise in the data [20]. On the other hand, the Vicsek model, and in fact any model dominated at short times by purely dissipative dynamics, gives z ≈ 2. This is quite understandable, as all these models belong (at equilibrium) to the same dynamical universality class as classic Heisenberg (Model A of [21]), which has z = 2 at one loop level, with very small corrections at two loops [21]. Moreover, when off-equilibrium (self-propelled) effects are taken into account, numerical simulations performed over time and space scales comparable to real swarms still give z ≈ 2 [20], completely incompatible with the data. The present calculation, on the other hand, shows that, once nondissipative terms are introduced in the dynamics, and provided that dissipation is not too strong, the dynamical critical exponent changes already at one loop, giving z = 3/2 in three dimensions. This is a value significantly closer to the experimental exponent than that of purely dissipative models. Hence, it seems to us that nondissipative terms are important to reproduce the correct critical dynamics of real swarms.
Of course, one must now ask how to bridge the gap between the one-loop RG exponent, z = 3/2, and the experimental value, z ≈ 1. There are several possibilities. First, one should try to have more statistics in the experiments, possibly with larger swarms, to check whether or not the data are really inconsistent with z = 3/2; work in this direction requires considerable technical effort on the experimental side (in particular, higher definition and faster acquisition systems). Secondly, one may hope that a two-loop calculation improves things. We are not very optimistic in this respect, though. Normally, two-loop corrections to the exponents are quite small, so it seems hard to bridge the gap between 1.5 and 1 in this way; furthermore, in the non-dissipative case the value z = 3/2 is actually valid at all order of the perturbative series, courtesy of the Ward identities generated by the symme-try [41]. Although in our case there is dissipation, we suspect that, as long as the system is in the proximity of the conservative RG fixed point, z = 3/2 will resist any attempt to be perturbatively changed.
Finally, there is the third and most promising source of corrections to z, namely off-equilibrium effects due to the self-propulsion of the individuals. Even though these are not sufficient to change the critical exponent in the Vicsek model, it could be that the compound effect of having non-dissipative inertial couplings and a self-propelled dynamics, further shifts the exponent in the correct direction. Studying this case from the theoretical point of view (i.e. by using RG), will be quite non-trivial, as one needs to use the approach of Toner and Tu [17], including in the theory one extra field, the density, coupled to velocity and spin, much as it has been done in [37] for the low temperature phase. However, at low temperature one could exploit the spin-wave expansion to linearize the equations, while close to T c , which is the case of interest for swarms, one needs to fully take into account the non-linearities through the RG. Performing even a one loop calculation with three fields (which become six once we use the Martin-Siggia-Rose representation) really does not look like a piece of cake. Still, one should try. In the meanwhile, numerical simulations of the full-fledged selfpropelled ISM close to criticality should be performed, to see from the data if there is case for hope.
Actually, we have some reasons to be optimistic. The fact that models with non-dissipative terms have dynamical critical exponent z significantly smaller than the purely dissipative value 2 may be interpreted as a critical counterpart of the linear spin-wave behaviour at low T : in this regime, 'second sound' modes propagate linearly, with dispersion relation ω = ck [36]. Naively, this relation would suggest z = 1 for these systems, but this is not the case, because close to T c parameters renormalize, so that that the second sound speed, c, goes to zero as some function of k; this RG-induced k dependence changes the exponent from the trivial 1, to the final z = 3/2 in this kind of models [21]. Despite this correction, though, the exponent remains significantly lower than the purely dissipative 2, as a relic of the low-temperature spin-wave dynamics. We may hope that a similar mechanism will be at work when self-propulsion will be taken into account. The first obvious effect of self-propulsion on a system with non-dissipative mode-coupling terms is to produce ballistic (i.e. linear) motion of each individual, even in the disordered collective phase. The dynamic critical exponent does not measure the motion of the individuals, of course, but rather the relaxation law of the velocity fluctuations; however, similarly to what happens with the renormalization of linear spin-waves, one may hope that some relic of the ballistic regime creeps into the critical phase calculation of z, thus lowering it below the static equilibrium value, 3/2, eventually bringing it closer to the experimental value. Further experimental, numerical and theoretical effort will tell whether this educated guess is just wishful thinking or not.

Perturbation expansion at one loop
To compute average quantities with the measure (A4) one proceeds as usual to develop the exponential contribution due to the interaction action, being left with a perturbation expansion where only free propagators and free correlations appear, connected to each other through the various interaction vertices. When building the full averages in such a way, we have to take into account that both ψ ψ 0 and ψψ 0 are non zero. To graphically distinguish between them, we will represent the propagators with an arrow and the correlation functions with a line, since propagators are time ordered while correlation functions are not. We will use the same rules also for propagators and correlation functions of s, but using wavy lines. It is more convenient to write down the perturbative expansion of G using the Dyson equation [58]: For which, we use the following diagrammatic notation: G 0,ψ α,β = C 0,ψ α,β = (A17) G 0,s α,β = C 0,ψ α,β = (A18) Here the blob indicates the sum of all 1PI diagrams with an incoming ψ (or s) field and an out comingψ (orŝ) field and with amputated external legs: namely, the selfenergies Σ αβ and Π αβ . The diagrammatic expressions for the self-energies of s and ψ at one loop are: where external legs are amputated. It is possible to translate these diagrams into integrals using standard Feynman diagrams rules: G 0,ψ (p, ω )C 0,s (k − p, ω − ω ) + (k 2 − p 2 )C 0,ψ (p, ω )G 0,s (k − p, ω − ω ) Performing the frequency integration, we get the following expressions for the self energies: We thus find that the self-energies only have, as expected, a diagonal non-zero contribution for α = β. We shall therefore drop the coordinate index and simply indicate them as Σ and Π, as in Eqs. (30) (31) in the main text.
When the perturbative corrections are calculated integrating over the shell, as in the RG approach we use in Section IV, the integrals are performed between Λ/b and Λ. On the the other hand, in the Callan-Symanzik ap-proach of Section VI all p integrals are performed between 0 and ∞.
Both sides in this expression implicitly also depend on h(x, t). Let us then derive with respect to this last field and then set it to zero. We get: where we relabelled integration variables in the second integral for future convenience. Here, R h γδ (x, t; x 1 , t 1 ) = ∂ ψ γ (x, t) /∂h δ (x 1 , t 1 )| h,H=0 is the linear response of the order parameter to its conjugate field. Using response theory, the l.h.s. of (C5) can also be written as where now in the r.h.s R hH αβδ (x, t; x 1 , t 1 ; x , t ) = ∂ 2 ψ α (x, t) /(∂h δ (x 1 , t 1 )∂H β (x , t ))| h,H=0 is the nonlinear quadratic response. Equating the r.h.s. of (C5)(C6) we finally get dx R hH αβδ (x, t; x 1 , t 1 ; x , t ) = g 0 αβγ R h γδ (x, t ; x 1 , t 1 ) −η 0 t t dt e −η0(t −t ) R h γδ (x, t ; x 1 , t 1 ) (C7) with t 1 < t < t. For η 0 = 0 this relation corresponds to the Ward identity reported for model E in [41].

Appendix D: Shell integration
To perform a RGT, as described in the main text, we need to implement two different steps: integration of short wavelength fluctuations, and rescaling. To this end, once fixed the coarse-graining factor b, it is convenient to rewrite the fields as the sum of two distinct components, one fluctuating on short wavelengths Λ/b < k < Λ and and the other on larger ones 0 < k < Λ/b, i.e. ψ(k, ω) = ψ < (k, ω) + ψ > (k, ω) (D1) At this point, one integrates out explicitly from Eq.(A4) the ψ > fields, to remain with a measure and a new effective action that only depend on the ψ < fields. To perform this integration one proceeds, again, using perturbation theory. The basic ingredients of this perturbation expansion (free propagators and vertices) are the same as the ones discussed in the previous sections, the difference being that they refer to ψ > fields only, while the ψ < are kept fixed as external sources. The perturbation series therefore consists in diagrams with external ψ < legs and internal loops integrated over > propagators. It can be recasted in exponential form, as usual, by only retaining one particle irreducible diagrams. These diagrams, that have external ψ < fields attached, will therefore modify the original terms appearing in the action. For example, for the Gaussian part of the action we get where Σ b has the same expressions as in Eq.(A23), but where integrals are performed only in the shell Λ/b < k < Λ. From this expression we immediately see that the behavior of the self-energy Σ b at small k effectively modifies the coefficient of k 2 . We are then left with a free part of the action similar to the original one, but where integrals run only up to Λ/b. The second step of the RGT, namely the rescaling of k, ω, ψ andψ, has the purpose of reinstating momentum integrals over their original integration range. At one loop the renormalization of the field is trivial (i.e. related to its physical dimensions) and has therefore not been addressed explicitly in the main text. The result is a new action formally of the same kind as the original one but with a new renormalized kinetic coefficient Γ b . A similar procedure can be applied also to the free action of the field s, and to the interacting part. All the coefficients and coupling constants will get renormalized by the shell integration and rescaling. If we call P the set of all parameters entering the action, i.e. P ≡ {r 0 , u 0 , Γ 0 , η 0 , λ 0 , g 0 }, a RGT will therefore imply P −→ P b S(P) −→ S b = S(P b ) Multiple iterations of the RGT therefore define a flow in the space of parameters, i.e. in the space of the statistical models defined by the action (A5).