Hydrodynamic attractor and novel fixed points in superfluid Bjorken flow

Extending the quantum effective approach of Son and Nicolis and incorporating dissipation, we develop a MIS formalism for describing a superfluid out of equilibrium by including the Goldstone boson and the condensate together with the hydrodynamic modes as the effective degrees of freedom. We find that the evolution of the superfluid undergoing Bjorken flow is governed by the conventional hydrodynamic attractor with unbroken symmetry and an even number of novel non-dissipative fixed points with broken symmetry. If the initial temperature is super-critical, then the condensate becomes exponentially small very rapidly and the system is trapped by the hydrodynamic attractor for a long intermediate time before it reheats rapidly and switches to one of the symmetry-breaking fixed points eventually. Finally, we show that the fixed points are unstable against inhomogeneous perturbations that should lead to spinodal decomposition. We conclude that these features should be generic beyond the MIS formalism.


I. INTRODUCTION
Phase transitions in out-of-equilibrium systems is of interest for the physics of heavy-ion collisions [1,2], cold atoms [3] and the early Universe [4]. Given the success of hydrodynamic theories for understanding out-of-equilibrium matter with both weak and strong self-interactions, there is sufficient motivation to construct frameworks that can describe superfluids by incorporating Goldstone mode(s) and the order parameter(s) (condensate(s)) associated with spontaneous breaking of global internal symmetries in hydrodynamic theories. Recently the demonstration of the ubiquitous presence of hydrodynamic attractors in a wide variety of phenomenological frameworks (eg. kinetic theory, holography, etc) has given us insights into how a system hydrodynamizes far away from equilibrium such that first order hydrodynamics can describe the evolution of the energy-momentum tensor and conserved currents even in presence large pressure gradients [5][6][7][8]. It is therefore pertinent to ask if one can construct similar frameworks that incorporate Goldstone boson(s) and order parameter(s) along with hydrodynamic modes for superfluid states far from equilibrium.
A general framework for superfluids at finite density has its roots in the pioneering work of Son [9] (and generalized to finite temperature by Nicolis [10]), which utilized the quantum effective action for describing the low energy dynamics of the Goldstone bosons and superfluid vortices. An attractive feature of Son's approach is that the only microscopic input needed is simply the equation of state of the system at finite chemical potential and temperature (see [11] for a survey of related approaches and a broader discussion). Remarkably, the thermodynamic description emerges naturally from the energy-momentum tensor and conserved currents obtained from the quantum effective action and is not imposed by hand. The original context in which Son developed this approach was to describe baryonic matter at finite baryon density with spontaneous breaking of U (1) baryon symmetry relevant for physics of neutron stars. Recently, generalizations of the quantum effective action have been studied incorporating both the pions and the hydrodynamic modes at finite chemical potential and temperature [12]. Further extensions to the case of phases with spontaneous spacetime symmetry breaking such as supersolids, etc [13][14][15][16] have been investigated as well but this is beyond the scope of this work.
In this work, we include the order parameter and add dissipation by merging the Son-Nicolis framework with the Müller-Israel-Stewart (MIS) formalism [17,18] enabling a causal description. We use this to study symmetry breaking in out-of-equilibrium situations in the context of Bjorken flow. The study of evolution of many-body systems under Bjorken flow is particularly of interest not only for heavy-ion physics [7], but also for cold atomic systems in anisotropic laser traps [19,20]. We restrict ourselves here to the case of a U (1) symmetry breaking but our work can be easily generalized to non-Abelian cases.
At early time, the size of the expanding system undergoing Bjorken flow is sufficiently small so that fluctuations are significant. Therefore we should study the evolution of the system with generic initial conditions for the energymomentum tensor, conserved currents and the condensate. Remarkably, we find that as long as we start above the critical temperature, T c , the condensate decays exponentially while the phase evolves over a similar timescale to satisfy the Josephson condition. Then the full system gets trapped very close to a conventional hydrodynamic attractor over a very long period of time during which the system approaches a perfect fluid expansion with unbroken symmetries.
However, the long time physics is somewhat surprising. It turns out that the superfluid system has an even number of symmetry breaking non-dissipative fixed points in which the full system undergoes expansion at a constant temperature (determined by the equation of state and the potential) and with a constant value of the condensate. These fixed points are possible because the condensate lowers the energy with respect to the vacuum, allowing the expanding system to maintain a self-consistent constant temperature (with a non-trivial velocity profile). We find that with inital T > T c and otherwise generic initial conditions, the superfluid system switches rapidly to one of these symmetry breaking fixed points after spending a time in very close proximity to the conventional hydrodynamic attractor. The basin of attraction of these fixed points has complicated interlacing and possibly fractal boundaries. If T < T c initially, then the superfluid system generically evolves to one of the fixed points without getting trapped near the hydrodynamic attractor (unless the initial condition is close to the boundary between the basins of attraction of the fixed points).
Our work demonstrates that the evolution of superfluid matter out of equilibrium is governed by the conventional hydrodynamic attractor with unbroken symmetry (which is rather like a saddle curve in the extended phase space) and non-trivial symmetry breaking fixed points which are determined by the potential of the condensate and the equation of state. The hydrodynamic attractor governs physics at intermediate time scales if the initial temperatures are above T c , while the fixed points govern the physics at much longer time scales with a (typically) rapid crossover.
We observe that the symmetry breaking fixed points are independent of the relaxation mechanism, and are determined only by the equation of state and the potential of the condensate. Furthermore, the hydrodynamic attractor is a feature in any phenomenological framework incorporating relaxation as mentioned above. Therefore, the dynamical features of superfluid flow, especially with respect to the role of the hydrodynamic attractor and the fixed point should be universal.
However, if there is competition between various types of symmetry breaking as in neutron stars and some strongly correlated systems, we will need to incorporate a more general framework discussed in [21][22][23] that can couple multiple (super)fluids. This is beyond the scope of our present work. We also postpone a discussion regarding the relevance of our results for the quark-gluon plasma until the concluding section.
We also analyze the linearized perturbations around thermal equilibrium and show that the fixed points are unstable against inhomogeneous perturbations that should lead to spinodal decomposition. This can be studied using our MIS framework, which we postpone this to a future work.
The organization of the paper is as follows. In Section II, we discuss the construction of the MIS theory of superfluids. In Section III, we derive the equations for Bjorken flow and then present our results in Section IV. We further analyze the system by considering linearized perturbations in Section V. Finally we conclude with a discussion on the limitations of our results and necessary future work in Section VI.

II. THE MIS FORMULATION OF SUPERFLUID EFFECTIVE THEORY
The MIS formalism [17,18] enables a causal description of a relativistic hydrodynamic system. In order to generalize this approach to a relativistic superfluid, we utilize aspects of the effective theory of superfluids due to Son and Nicolis [9,10] at finite temperature. Furthermore, we include the order parameter as a dynamical variable to incorporate the possibility of out-of-equilibrium phase transition.
For simplicity, consider a complex field Σ ≡ ρe iψ charged under U (1) global symmetry which is broken spontaneously below a critical temperature T c due to the potential V (ρ, . . . ) where . . . denote other variables of the effective theory. Crucially, the velocity u µ and temperature T of the local thermal frame appear together as a new fundamental variable β µ . Note that in this approach (local) thermodynamics emerge from the effective action and is not imposed separately. We first succinctly show how this works.
Due to Lorentz invariance, the effective action in a background metric g µν and gauge field A µ (introduced to allow us to compute the energy-momentum tensor and conserved current efficiently) can depend on β µ only via the temperarture T and chemical potential µ defined by the relations T := (−β µ g µν β ν ) −1/2 and µ/T := A µ β µ , respectively, at leading order in derivatives. Similarly, due to U (1) gauge and Lorentz invariance, the dependence on Σ can be via the scalars X := (D µ Σ) † (D µ Σ) and Y := ((β · D)Σ) † ((β · D)Σ), where with D µ ψ := ∇ µ ψ − A µ . It will turn out (from absence of entropy production) that (u · D)ψ = 0, i.e. (u · ∇)ψ = µ (a.k.a. the Josephson condition) must be satisfied in equilibrium while ρ = ρ * , a constant, so that Y must vanish. On the other hand, we will find that X can assume an arbitrary non-trivial profile in equilibrium due to non-vanishing of derivatives of ψ in directions orthogonal to u µ . Since Y must vanish at equilibrium, we will drop it and retain X, T and µ as fundamental variables. Additionally we incorporate the gauge-invariant variable ρ to implement spontaneous symmetry breaking. Since derivatives of ρ (and ψ) appear in X, we must think of X and ρ as independent variables. 1 We therefore consider an effective action of the form: where F is the kinetic piece of the action and V is the symmetry breaking potential. Simple instances of F and V are in which σ, λ > 0 and are both independent of T and µ. In Son's approach where only the dynamics of the Goldstone mode are considered [9], the condensate is set to its equilibrium value. The dependence of the generalized pressure F on the chemical potential µ then also determines its X dependence. In our case, the condensate is an independent variable and we truncate our theory to two-derivatives. To relate to Son's approach, we can simply expand the generalized pressure F in derivatives of the condensate and normalize the latter to get a standard two-derivative kinetic term. This leads to our simple ansatz (2). In what follows for the rest of this section, we consider a general F (X, T, µ). The ideal energy-momentum tensor and conserved U (1) current can be readily obtain by varying the action (1) w.r.t. g µν and A µ respectively. These turn out to be respectively with ∆ µν = g µν + u µ u ν . Each of these can be rewritten as sum of a normal and a coherent superfluid component. Explicitly, where j µ ψ := ρ 2 D µ ψ and With the entropy density S defined as the relations (6) imply the standard thermodynamic identity Note that P is distinct from p which appears in (2) as it is the generalized pressure which governs the normal fluid component. The coherent part in (3) which is proportional to D µ ψD ν ψ can be reinterpreted as a second fluid with a velocity field proportional to D µ ψ and with an independent free energy. We refer the reader to [21,22,24] for more details. Note in our case we will also need a third fluid component with velocity field proportional to ∇ µ ρ as evident from (3). However, such reinterpretations will not be necessary in what follows. Varying the effective action (1) w.r.t. Σ, we obtain the equations of motion for ρ and ψ which are The above equations of motion, in conjunction with ∇ µ T µν ideal = 0 and ∇ µ j µ ideal = 0, imply the Euler equations where ∇ µ ⊥ ≡ ∆ µν ∇ ν . In order to further extend the effective theory, we need to add dissipative terms to the energy-momentum tensor and the conserved current, which we denote as π µν and q µ respectively, so that with E, P and N as defined in (6). We choose a generalization of the Landau frame to define u µ in the non-ideal case so that π µν u ν = 0. Furthermore, we modify (9) by adding dissipative sources θ 1 and θ 2 to the equations of motion of ρ and ψ respectively, so that The above in conjunction with ∇ µ T µν = 0 and ∇ µ j µ = 0 imply where D ⊥ µ = ∆ µν D ν . As a consequence, we obtain a candidate entropy current j µ s = Su µ − (µ/T )(q µ + j µ ψ ) whose divergence turns out to be The appropriate constitutive relations that lead to the positive definite divergence of the entropy current is then with κ, η, ζ, κ 1 and κ 2 positive definite functions of T and µ (and in principle of ρ too). Also is the shear-stress tensor. We note that the constitutive relation for q µ implies that the conserved U (1) current is We readily validate our earlier claim from (15) and (16) that in equilibrium where entropy production is absent, we should satisfy the Josephson condition (u · D)ψ = 0, i. e. (u · ∇)ψ = µ. Furthermore, (12) implies that at equilibrium ρ should be in the (thermal) vacuum so that ∂V /∂ρ = 0 in absence of spatial gradients of ψ. Note when ∇ ⊥µ ψ∇ ⊥ µ ψ = 0, e. g. a constant, we can have ρ taking values away from the minima of V in equilibrium at any T and µ.
Finally in order to develop a MIS-type formulation, we simply replace the constitutive relations for π µν and j µ in (16) by the dynamical equations with τ π and τ q being additional parameters that depend on T and µ (and in principle also on ρ). Clearly after the timescale of the respective relaxation times, π µν and q µ + j µ ψ will relax to their respective constitutive relations (16). This completes the construction of the MIS formalism for effective description of a relativistic superfluid. Finally we note that our construction is not quite reliable at very low temperatures where higher derivative terms can contribute also to the relaxation processes. In this context, one may use the formalism discussed in [25]. We can choose parameters such that for generic initial conditions with supercritical temperature (and also with a range of subcritical temperatures), the system never reaches sufficiently small temperatures while undergoing Bjorken flow.

III. BJORKEN FLOW
In order to study Bjorken flow we need to consider the Milne background and set u µ = (1, 0). Furthermore, we should consider all other variables T , µ, ρ, ψ, π µν and q µ to be dependent only on the proper time τ . For a concrete example, let us consider F and V to be of the form given in (2). We further set µ = 0 and p = T 4 as in a conformal equation of state. Similarly, we assume that where C η , etc are constants and we set ζ = 0. Alternatively, we can parameterize κ 1 = C κ1 T c and κ 2 = C κ2 T 3 c . None of the qualitative features that will feature in our discussion will depend on the details of such parametrizations. Denoting τ −derivatives via a prime, the equations of motion (9) now take the form Furthermore, we consider a diagonal form of π µν : It is convenient to use the dimensionless pressure anisotropy χ := 3Π/4T 4 . The conservation of the energy-momentum tensor T µν given in (11) provides the equation for the evolution of the temperature Assuming p to be conformal, it is natural to modify the MIS equation (18) by replacing u · ∇ on the LHS with a Weyl-covariant derivative as in BRSSS formalism [26] so that the evolution of Π is given by As such, we have a five-dimensional phase space given by T , χ, ψ , ρ and ρ when we set µ = 0 (for simplicity).

IV. RESULTS
Superfluid Bjorken flow has formidable complexity, but the results are best understood by considering possible attractor or saddle surfaces/curves/points first. It is easy to note from our equations that we can always consistently set Σ = 0. The evolution of T and χ are then exactly same as the case of ordinary Bjorken flow discussed in [5]. It is well known that the system flows to a unique hydrodynamic attractor curve given by χ att (τ T ). However, we find this curve is actually a saddle curve away from the Σ = 0 surface. However, the system has two fixed points! To see this we set T = t * T c , χ = 0 and note that the equations of motion have solutions with ψ = 0 and since ∂V /∂ρ = 0 at ρ = ±ρ * . Note setting ψ constant implies that we are moving on a line in the complex Σ plane, so ρ can assume any sign. The equation for T given by (24) determines t * via in which we take 0 ≤ t * < 1. In Fig. 1, we see that irrespective of our choice of constants (σ, λ, T c ), we are always able to find a real solution within this range. Choosing T c = σ = λ = 1, we find that t * ≈ 0.551847 and ρ * ≈ 0.669442. The reader may readily note that (27) orginates from a specific choice of the equation of state i.e. p(T ). Also the fixed point is determined only by the quantum effective action (the ideal component) and is independent of the relaxation terms. The existence of such fixed points should therefore be generic in the Bjorken flow equations of the superfluid. The reason that such fixed points with constant temperature can exist in this expanding superfluid droplet is simply because the energy gained by non-reduction of temperature is compensated by the lower energy of the expanding vacuum with broken symmetry. Generally, of course the system will have even number of such fixed points since the symmetry breaking vacuua are doubly degenerate for fixed ψ and for any acceptable solution of t * . At these symmetry breaking fixed points there is no entropy production because ρ = ψ = π µν = q µ = 0. Since these are fixed points of an expanding system (with non-trivial velocity profile), these are not thermal.
We are now ready to report our results. The system of equations for the simple case (2) which give the evolution of (T , χ, ρ, ψ) under Bjorken flow with vanishing µ are (21), (22), (24) and (25). We recall that the phase space is five-dimensional when µ = 0 with coordinates T , χ, ρ, ρ and ψ . The characteristics of the evolution of the system, however, depend only on whether we initialize with T > T c or T < T c .
For the purpose of generating plots, we choose C η = 1/(4π) and C τ π = (2 − ln 2)/(2π) as obtained by matching MIS with the holographic description of N = 4 SYM theory at infinite 't Hooft coupling and for infinite number of colors [5]. Also we choose C κ1 = C κ2 = 1, and σ = λ = 1 in units T c = 1. Likewise, we will count our proper time, τ , in units of T c = 1.
a. Case 1 -initially T > T c : In this case, the typical evolution of the system is as shown in the left panel of Fig. 2. Early on, both ρ and ψ vanish very rapidly, particularly ρ ≈ τ −1/2 e − √ σTcτ . The system approaches the Σ = 0 surface (with ρ = ρ = ψ = 0) where the evolution is purely hydrodynamic. As a result, the system gets trapped near the hydrodynamic attractor curve χ att (τ T ) where there is no symmetry breaking (ρ remains exponentially small). The entrapment occurs for a long time during which the system approaches a perfect fluid type expansion (T ≈ τ −1/3 ) due to its evolution close to the attractor curve. However, at late time the system suddenly reheats quickly and switches to one of the symmetry breaking fixed points with T = t * T c and ρ = ±ρ * (with ψ = 0 and χ = 0). Note that during the entire evolution of the system, the temperature remains above t * T c , including during the switching between the hydrodynamic attractor behavior and the final approach to its asymptotic value.
In Fig. 3, we further illustrate the above by fixing the initial temperature to 2T c and initializing with (ρ , ψ , χ) = (0, 0, 1), but varying the initial value of ρ. As mentioned above, the qualitative features do not change if we choose initial values of T , ρ , ψ and χ differently as long as T > T c initially. For any initial value of ρ, the system gets trapped close to the vertical ρ = 0 line, representing the hydrodynamic attractor curve where we have perfect-fluid like expansion T ≈ τ −1/3 for a considerable time (see the inset) until it switches to one of the symmetry breaking fixed points marked with red and cyan colors. Varying initial values of ρ, we find alternating intervals in which the system 954 (corresponding to the solid and dashed lines, respectively), at early times, the temperature (red lines) quickly goes to perfect-fluid type expansion T ≈ τ −1/3 (shown as a gray dashed line) before switching to one of the symmetry breaking points at late time. This switching time corresponds to the rapid growth of the condensate, ρ, shown in black. At late times, both T and ρ asymptote to their respective final values, t and ±ρ , which we denote in both plots as blue and orange lines, respectively. In the right panel corresponding to initial values of ρ as 2 and 1 (solid and dashed lines, respectively), with initial T < Tc, the condensate rapidly approaches its asymptotic value. The evolution of the temperature is non-monotonic even at early times. In these plots t * ≈ 0.551847 and ρ * ≈ 0.669442 and we start with ρ = ψ = 0 and χ = 1.
chooses one of these fixed points at late times. The switching time depends on the initial condition and increases with the inital value of T very rapidly.
To study the basin of attraction of the fixed points, we fix the initial values (ρ , ψ , χ) = (0, 0, 1) and vary initial values of T and ρ. We readily observe from Fig. 4 that for initial values of T and |ρ| greater than about 5T c and 7T c , respectively, the basin of attraction of these two fixed points (marked with yellow and blue for ρ * and −ρ * , respectively) are interlaced in a complex manner, making long term prediction difficult. It is likely that the basin boundaries are fractal in these regimes.
We emphasize that even if the system is initialized with non-vanishing ρ and ψ , it rapidly approaches the hydrodynamic attractor curve with vanishing Σ and lingers there for a long duration of time before rapidly switching to one of the fixed points. The basin of attraction of the fixed points has similar features.
If we start from the boundary of the basin of attraction of the two fixed points initially, then the system should get trapped by the hydrodynamic attractor curve on the vanishing Σ surface (with ρ = ρ = ψ = 0) and would asymptotically approach it without switching to any of the fixed points.
b. Case 2 -starting with T ≤ T c : The evolution is qualitatively different to the previous case, although the system settles to one of the fixed points for a generic initial condition. As shown in Fig. 2, the temperature behaves non-monotonically at initial time even if we start with T ≈ t * T c . The condensate ρ approaches ±ρ * while ψ and χ vanishes over the same timescale. Thus, the system approaches one of the fixed points without getting trapped near Σ = 0, unless the initial condition is close to the boundary of the basin of attraction as demonstrated in Fig. 3. In the latter case, the system evolves close to the hydrodynamic attractor for a certain period of time. This is expected because if we start from the boundary of the basins of attraction, then the system should approach the hydrodynamic attractor curve in the long run instead of evolving to one of the fixed points. We readily observe from Fig. 4 that the projections of the basin of attraction on T − ρ plane separates into two simply connected regions.
In the following section, we analyze the MIS theory about thermal equilibrium and argue that the fixed points are actually unstable against inhomogeneous perturbations. The fate of this instability can be studied using our MIS framework.

V. LINEARIZED FLUCTUATIONS
Here, we compute the linearized fluctuations around thermal equilibrium at zero chemical potential. Note that in this section we are not considering an expanding fluid and we use the standard Minkowski coordinates. We show that the broken phases with T ≤ t * T c < T c (with t * T c corresponding to the value of the temperature at the fixed point of The color coding is yellow for +ρ and blue for −ρ . Note that below a certain value of temperature, there is a clear separation between a choice of initial conditions determining which ρ the system will evolve to. This is not the case for higher temperatures and higher |ρ|, where predictability from the initial conditions becomes difficult. The boundary of the basins of attraction is possibly a fractal at higher tempertures. the Bjorken flow) are unstable in our MIS theory. While the non-hydrodynamic relaxation poles remain on the lower half-plane of complexified frequency, the speed of sound, the diffusion constant and the sound attenuation coefficient change sign at T = t * T c and thus are negative for T < t * T c . This hydrodynamic instability should lead to spinodal decomposition [27,28].
Without loss of generality, we take the momentum of the fluctuation in the z-direction. The equations of motion, (21) and (22), energy momentum tensor conservation (24) and MIS equation (25) are linear equations of V ≡ (δT, δρ, δψ, δu, δπ) in Fourier space. We can recast the equations in matrix form Q(ω, k) · V = 0, where Q is a 5 × 5 matrix. We find the dispersion relations by setting the determinant of Q to zero. The perturbations of the hydrodynamic variables are given by shear: π zx ∼ δπ zx e −iωt+ikz and u x ∼ δu x e −iωt+ikz , sound: T ∼ T 0 + δT e −iωt+ikz , π zz ∼ δπ zz e −iωt+ikz and u z ∼ δu z e −iωt+ikz .
The fluctuations in the scalar sector are We distinguish between the broken and unbroken phase by the equilibrium value of the condensate. When ρ 0 = 0, we are in the unbroken phase, T > T c . In this case, the determinant of Q factors neatly into the phase, condensate and hydrodynamic dispersions. Since ρ 0 vanishes, the dispersion of the phase is trivial. We find the expected sound and shear modes of a conformal system (with η = C η T 3 0 ) as well as the non-hydrodynamic modes: one from the MIS equation and two from the scalar field respectively. Note that the superfluid mode is overdamped when κ 1 > 2σ(T 0 − T c ).
Next, we turn our attention to the broken phase, T < T c , and ρ 0 = ± σ(Tc−T0) λ . The dispersion relation of the phase is now no longer trivial, indicating the presence of a Goldstone mode Clearly we get two modes which at low k take the form Thus one behaves as a diffusion mode and the other as a relaxation mode in both sound and shear channels. The shear sector in the broken phase with perturbations (28) has the following determinant where we introduced the mass of the condensate, m 2 σ := 3λρ 2 0 + σ(T c − T 0 ) = 4σ(T c − T 0 ). The determinant factorizes into the condensate and hydrodynamic modes. For small k, we see that the shear modes are as well as the MIS mode (32). We note that the diffusion constant has a pole in T 0 precisely when T 0 = t * T c and its sign changes as we further lower T 0 . A negative diffusion constant leads to what is known as uphill diffusion (against the concentration gradient) enhancing inhomogeneities and thus spinodal decomposition. As a side remark, the two scalar (condensate) modes become one diffusive mode in the limit m 2 σ → 0, (see the Goldstone sector discussed above) which in the context of the chiral phase transition can be thought of as the transition from pion propagating modes to the diffusive quark mode [29]. Of course, in the present context, we are working with a complex scalar field and not attempting to capture the O(4) dynamics relevant for the chiral phase transition.
Finally, in the sound channel with perturbations (29), the dispersion takes the following form We see that in the small k limit, we have the following modes Above, the first mode is clearly the sound mode but the square of the speed of sound c 2 s becomes negative for T 0 < t * T c triggering spinodal decomposition. In Fig. 5, we plot c 2 s and the sound attenuation coefficient Γ which behaves similarly as the diffusion coefficient in the shear sector, diverging at T 0 = t * T c and becoming negative as we further lower T 0 .
We note that the instabilities occur only in the hydrodynamic sound and shear modes only for T 0 ≤ t * T c . This indicates that the fixed points of the Bjorken flow could be unstable against inhomogeneous perturbations. However, when we start with supercritical temperatures (and a range of subcritical temperatures), the temperature remains above t * T c for its entire evolution for a large class of initial conditions. We would like to investigate the fate of the instabilities of the fixed points using an appropriately modified MIS framework in the future.

VI. DISCUSSION
In this work, we setup the MIS formalism for superfluids, extending the quantum effective action approach and study the Bjorken flow. We find that the dynamics of the system is governed by the existence of the conventional hydrodynamic attractor with unbroken symmetries and new fixed points determined by the potential for the condensate and the equation of state of the system. Particularly, if the initial temperature is above T c , we find that the Josephson condition is satisfied by the phase very rapidly while the condensate becomes exponentially small over a similar timescale and undergoes slow-roll during which the system gets trapped in the vicinity of the conventional hydrodynamic attractor. This persists for a long time until the system reheats and switches over to one of the symmetry-breaking fixed points. If the initial temperature is below T c , the system does not generically get trapped by the conventional hydrodynamic attractor, but rather evolves to one of the fixed points. In both cases, the hydrodynamic attractor traps the system forever if we start with initial conditions at the border of the basin of attraction of the fixed points.
As discussed before, we expect that our key result that the superfluid Bjorken flow is governed by a combination of hydrodynamic attractor and an even number of non-dissipative fixed points should be a generic phenomenological feature. The fixed points are determined by the ideal component of the dynamics obtained from the quantum effective action with the relaxation mechanism playing no role while the hydrodynamic attractor appears in any strong coupling or weak coupling description of thermalization.
The analysis of fluctuations about thermal equilibrium indicate that the system has unstable hydrodynamic modes when the temperature corresponds to values lower than that at the fixed points. This indicates that the fixed points are unstable against spinodal decomposition. In the future, we would like to investigate this phenomenon by introducing inhomogeneities in the initial conditions and using our MIS theory with necessary improvements.
It would be important to study the Bjorken expansion of the superfluid in a UV-complete setup at strong coupling, particularly in a holographic model. Of particular interest would be to understand if the approach to the hydrodynamic attractor and subsequent switching to the symmetry breaking fixed points occurs rapidly as in the MIS setup.
Furthermore, to be of relevance to the quark-gluon plasma we will need to consider interacting multi-component systems as in [24]. Our model is essentially that of the quark sector. The gluonic degrees of freedom should be introduced as a separate fluid. Furthermore, it may be necessary to include both weakly coupled and strongly coupled degrees of freedom. Such hybrid fluid models with a higher dimensional hydrodynamic attractor has been discussed in [30,31]. In case of a multicomponent system, the reheating transition of the superfluid to the fixed point may not be possible or may not occur irreversibly because of the transfer of energy to the other components of the full system.
Also, it would be interesting to explore more general forms of expanding flows, particularly those with spherical symmetry (see [32] for a relevant model). As the superfluid condensate is naturally driven to slow roll on a hydrodynamic attractor with a natural mechanism of exit from the slow-roll, while reheating of the system upon switching to one of the fixed points, this raises the possibility that one can use such a superfluid as an inflaton in cosmology with a natural preheating mechanism in place.