Global flow structure and exact formal transseries of the Gubser flow in kinetic theory

The topological properties of the global flow structure of the Gubser flow for the Israel-Stewart theory and a kinetic model for the Boltzmann equation are investigated by employing Morse-Smale theory. We present a complete classification of the invariant submanifolds of the flow based on the skew product flow formalism of nonautonomous dynamical systems and determine all the possible flow lines connecting any pair of UV/IR fixed points. We also construct the formal transseries solutions to the Gubser dynamical system around the early-time (UV) and late-time (IR) fixed points. It is proven that these solutions are purely perturbative (or power-law asymptotic) series with a finite radius of convergence. Based on these analyses, we find that Gubser-like expanding kinetic systems do not hydrodynamize owing to the failure of the hydrodynamization process which heavily relies on the classification of (non)hydrodynamic modes in the IR regime. This is in contrast to longitudinal boost-invariant plasmas where the asymptotic dynamics is described by a few terms of the hydrodynamic gradient expansion. We discuss the generic conditions for the existence of a hydrodynamic attractor that is an invariant manifold determined by the long-wavelength modes of the distribution function. We finally compare our results for both Bjorken and Gubser kinetic models.


Introduction and summary
It is often argued in science and engineering textbooks that the regime of validity and applicability of hydrodynamics is linked to the existence of a hierarchy between the mean free path λ m.f.p. and the system size L. The existence of such a hierarchy implies that the gradients of the fluid hydrodynamic fields are small compared to their corresponding thermal equilibrium configurations. In a modern approach, this also gives rise to the interpretation of hydrodynamics as an effective field theory for the slowly varying long-wavelength modes of many-body systems in the form of a gradient expansion.
This viewpoint, however, has been questioned recently in light of the phenomenological success of fluid dynamical models to describe rapidly expanding fireball of nuclear matter created in Ultra-relativistic Heavy Ion Collisions (cf. Refs. [1][2][3] and references therein). As a matter of fact, different toy models in the strong and weakly coupled regimes have hinted at a possibility of hydrodynamics being able to extend further to describe dynamical aspects of extremely non-thermal equilibrium configurations [4][5][6][7][8][9][10][11][12][13][14]. Furthermore, it was shown recently in holographic systems [15,16] that the radius of convergence of the dispersion relations is orders of magnitude larger than the standard small-frequency limit. Similar conclusions were drawn in some kinetic theory models [17]. Therefore, we need a consistent fluid dynamics formulation to explain hydrodynamics far from equilibrium using first principles.
There are many unknown elements of this new formulation despite a few important lessons learned in the recent years. For instance, different numerical studies indicate that the equations of motion for a class of fluids undergoing Bjorken flow [18][19][20][21][22][23][24][25][26][27][28][29][30][31][32][33][34][35] can be dimensionally reduced to a single equation for the inverse Reynolds number Re −1 as a function of the Knudsen number Kn 1 . Being an explicitly time-dependent (also known as "nonautonomous") dynamical system with an attracting IR fixed point 2 at late times Kn 1, the model enjoys solutions that approach this fixed point which naturally merge before equilibriating. Roughly speaking, the set of points shared by the solutions to the model, aka "flow lines" in a dynamical system interpretation, is dubbed as hydrodynamic attractor or simply attractor. Since there are two UV fixed points in the system sourcing the flow of matter at early times Kn 1, the main question is whether the attractor can be continued back to one of these fixed points 3 . In other words, is there a complete flow line initiated at a UV fixed point that connects to the hydrodynamic attractor?
Finding an answer to the question might ultimately boil down to considering the fastslow decomposition of the geometry of distribution for a given system [37,38]. Although a distribution obeying Boltzmann equation includes wider class of solutions in general, it can directly correspond to hydrodynamics through a lifting/projection process once a particular set of solutions are chosen. Above all else lies the fact that the existence of a slow manifold (or hydrodynamic attractor) is a necessary and sufficient indicator for hydrodynamization. However, the existence proof is a nontrivial problem for which one has to figure out both the global and local structures of the underlying dynamical system.
A definite answer to this difficult problem can at least be given for certain kinetic models of expanding plasmas undergoing Bjorken flow. The longitudinal boost-invariant systems are often described in terms of the longitudinal proper time τ = (x 0 ) 2 − (x 3 ) 2 and they present an unresolvable singularity in the Boltzmann equation at τ = 0 [36]. It has been known for a while that certain interesting features of the hydrodynamic expansion of the Bjorken flow are better understood if we use the variable w = τ T (T being the temperature) [11]. However, in this new variable the UV structure of flow lines of the original dynamical system of ODEs is completely altered topologically 4 . This becomes a 1 For the Bjorken flow, the inverse Reynolds number is π/p in which π is the independent shear viscous component and p is the equilibrium pressure, whereas the Knudsen number is Kn ∼ (τ T ) −1 . We should also note that Kn is sometimes denoted as w. 2 We admit that the use of word fixed point in a time-dependent system is loose as it does not really exist for all times. However, at late times, we can take an open set localized at the equilibrium values for the configuration-space variables except time and consider that the flow remains oblivious to the time steps, making it a limit fixed point. 3 It is possible for this continuation to not generate a unique critical line in case both the UV fixed points are "saddle", a fact that does not obviously hold in the configuration space of the Bjorken flow. Note that a critical line is really a complete flow line that basically lies on the boundary of basin of attraction [36]. 4 The variable w does not capture any UV information from the original system and rather serves as a toy model for the Bjorken flow in τ . A simple pathology to think of is the fact that writing the equations in w takes away the singularity at τ = 0 and transforms it into an unphysical singularity in another configurationspace component, the normalized shear viscous tensor. This makes w a topologically forbidden coordinate crucial matter if we are to study the UV completion of the flow lines. That is why it is sufficient to focus on the topology of the invariant manifold which contains all the stable flow lines that connect the physical information of the UV to that of IR. Thanks to the symmetries of the dynamical system, it is an extremely important point to press on that the topology of the invariant manifold has to survive any coordinate transformation respecting those symmetries, a requirement that is not satisfied for the w variable.
In a more general class of hydrodynamic theories, the search for solutions depends on the existence of an "invariant manifold" which substantiates many properties of the original kinetics. Since a reduced description of hydrodynamics only needs parametrization of microscopic variables projected as macroscopic fields v, p, ε, a "lifting procedure" is required to acquire knowledge about a complete microscopic picture. This procedure is well-defined on a "slow" invariant manifold -or simply slow manifold 5 , where the hydrodynamic flow meets the kinetic motion of a fluid as the kinetics becomes well-behaved at a small Knudsen number and hydrodynamic field becomes bounded close to the local equilibrium [37][38][39]. If the slow manifold exists, a complete understanding of hydrodynamics can be obtained based on these macroscopic fields. Therefore, an important task in this formulation is to find the slow manifold.
Another misconception that we have to shed some light on is that in nonlinear nonautonomous systems, the dynamical system itself evolves with time so speaking of a global attractor that servers as a UV-IR complete flow line as in autonomous systems is not generally an option. Let us look at what this means with an example. Consider the Boltzmann equation with some collision kernel C[f (X(τ ), τ )] 6 which explicitly depends on time τ and some configuration-space variable X. Now, solving C[f ] = 0 puts a condition on the distribution function. We can take this condition to be f (X * (τ ), τ ) = f * where f * is a real constant. Here, X * (τ ) determines the zeros of the equation f (X * (τ ), τ ) − f * = 0 and suppose for simplicity that there is only one zero and C[f (X(τ ))] is decreasing as we tend to this point. Normally, we would have considered this zero to be an attracting fixed point of the Boltzmann equation, but now obviously it is evolving with time. However, the evolution of X(τ ) with appropriate initial values is going to get closer to X * (τ * ) at every time step 0 < τ ≤ τ * . If we repeat the same exercise of obtaining the evolution of individual variables until τ * for a new τ * , we can observe a line across all values of τ * . This line is deemed as an "attractor solution" since it naturally attracts the points falling in the basin of attraction of X * (τ ) until every time step τ * , even though as we can see, it is nothing but the trace of the zero of collision kernel across all times which really becomes the thermal equilibrium of the system at τ * → ∞. In other words, the system did not know about its equilibrium state in the past since there was none to begin with. Only it is at τ * → ∞ that the system reaches change as far as the UV limit is concerned [36]. We have listed two other reasons at the end of Sec. 5. 5 In non-relativistic kinetic theory it is usually considered that if the Navier-Stokes limit is achievable, the slow manifold is then defined by this limit where long-wavelength modes attenuate at the some decay rate. The fast-slow decomposition connects the initial configurations of macroscopic fields (initial condition for v, p, ε) to the hydrodynamization in the long-wavelength limit where the flow of hydrodynamics can be consistently defined from these macroscopic fields. 6 Normally, we should have included temperature in f , which is obviously ignored here, since it does not affect the conclusion.
a fixed point X * (τ * ) and an equilibrium distribution f * . Therefore, we need to really talk about attractors in the past or future of an evolving fluid with explicitly time-dependent dynamics; see Fig. 1. Figure 1: Phase portrait of some nonautonomous system. At the point ρ 0 the equilibrium state in the future has not been formed yet. If we wait infinitely long enough, this state will start to form at ρ. As this is happening, the flow lines at late times become further dense around the equilibrium point. This dense set is the future (forward) attractor. Notice that there is no attraction in the past of ρ due to the non-existence of an equilibrium state (fixed point), a statement that does hold in autonomous systems -cf. RG flow space of field theories.
The analysis of higher dimensional fluids such as the Gubser flow [40,41] as well as less symmetric flow profiles in 2 + 1D [42] and 3 + 1D [43] indicate that there are attractors for these PDEs that might be subset of complete flow lines. The analysis of attractors for higher dimensional systems is more involved than the 1D longitudinal boost-invariant case so the original approach of the authors of Ref. [26] falls short to describe them. For instance, for the Gubser flow, two of the current authors showed by studying that flow as a dynamical system and performing stability analysis of IS theory [40] that the form of attractor is not described by the asymptotic perturbative gradient expansion in terms of the Knudsen number. Although this example was just a simple exercise to show the power of dynamical systems in lower spatial dimensions, this approach in essence is way more generic and could apply to fluids expanding in higher dimensions for probing not only the attractors, but also interesting phenomena such as the existence of saddle-type steady states in the mid-range, fluid turbulence, chaos, etc.
In this work, we investigate and generalize our previous findings [40] and focus on constructing all the possible links between UV and IR points using topological properties of the Gubser flow as a nonautonomous dynamical system based on the Morse-Smale decomposition and the transseries analysis. In practice, it is possible to describe the story using just the Morse decomposition, in which case we should focus on the Morse-Lyapunov functions for describing the flows. So for us, the two decompositions are interchangeable even though the demonstrations will be done without referring to a specific Morse function.
To start things off, we consider the global structure by giving an account of the stability analysis of the fixed points [44,45] and define rigorously what exactly an attractor is in a nonautonomous system. This comes handy when we discuss different types of attractors in this system. We consider the time evolution of flow lines (configuration-space trajectories) attached to individual fixed points portrayed in lower dimensions in the space of moments of distribution function and temperature, and reveal the global structure from the topology of its invariant set and constituent invariant subsets. It is shown that there exists two bounded invariant submanifolds in the configuration space and one of them has a nontrivial flow structure in the sense that it contains a local future attractor, but the system lacks a local past attractor. As for the analytic solutions to the evolution equations, we construct a formal transseries solution by expanding around the stable IR fixed point and show the asymptotic behavior of hydrodynamical quantities and the Boltzmann distribution function. Furthermore, we prove that the formal transseries is a convergent series with a finite radius of convergence 7 , and hence why the distribution does not converge to a thermal equilibrium. Indeed, one cannot find a well-defined hydrodynamic description of Gubser flow due to the fact that the Knudsen number is large and a slow manifold does not exist. This is backed by the failure of Chapman-Enskog expansion, which stems from the nonvanishing derivative of distribution at late times, implying that the system still undergoes expansion near the local equilibrium. We conclude by giving final remarks.

Stability analysis of conformal Israel-Stewart theory
We start our discussion of the global flow structure and stability analysis for the Gubser flow by considering the Israel-Stewart (IS) theory [46]. As a warm-up exercise, we use this hydrodynamic truncation scheme to get things started and introduce a set of mathematical tools to analyze the structure of the solutions near its asymptotic early-time (UV) and latetime (IR) fixed points. Some of the results presented in this section were partially derived in a previous publication [40]. We are going to build some aspects of our original approach into this paper and extend our previous findings. The generalization of the methods described in this section to the case of the relativistic Boltzmann equation within the relaxation time approximation (RTA) is presented in forthcoming sections.
The Gubser flow describes a conformal fluid undergoing longitudinally boost-invariant and azimuthally symmetric expansion on the transverse plane in Minkowski space [47,48]. The symmetry group of this flow, SO(3) q × SO(1, 1) × Z 2 is then manifestly inherited from the symmetries of dS 3 × R spacetime and thus, the analysis of dissipative equations is more conveniently described on this manifold. The coordinates, the metric g µν , and the flow velocity profile in this curved space are respectively given by From the conservation law and the IS theory, the equations of motion for the effective temperature T and the inverse Reynolds number π := π ςς /( + p) 8 are derived to be [50] dT where c > 0 is an arbitrary dimensionless real scale that relates the shear relaxation time scale τ π and the ratio of shear viscosity over entropy η/s, which in the conformal case is simply τ π T = c η/s. Eqs. (2.4) constitute what is called in the math literature a nonautonomous dynamical system since the de Sitter time variable ρ appears explicitly on the r.h.s. through the expansion rate D µ u µ = tanh ρ. In general, a nonautonomous system has some subtleties when it comes to defining and analyzing its dynamical equilibrium points, i.e. fixed points, and their stability. For instance, Eqs. In nonautonomous systems, a smooth change in the flow parameter can induce topological changes such as bifurcation that are probed in time-independent systems by tuning a control parameter. There are other differences that we will touch on as we continue.

Flow structure
In this subsection, we classify all the fixed points of the IS theory in both UV and IR limits by calculating their instability index in the extended configuration space to lay out the basic local structure of the flow. We use the skew product flow formalism of this nonautonomous system to portray both the local structure of the configuration space (T (ρ),π(ρ)) and its time-dependent change in the global picture. Following this construction and Morse decomposition, we show all the possible flow lines connecting any pair of UV/IR fixed points and discuss the existence or nonexistence of pullback attractors. 8 Here, π ςς is the only independent component of the shear viscous tensor, is the energy density, and p is the equilibrium pressure. In our previous works [36,40,49], we called π the effective normalized shear viscous tensor. We will stick to this name throughout the present paper as well. 9 For this reason, if the time direction is non-compact, it is more appropriate to use the term "fixed points at infinity" when referring to the fixed points of nonautonomous systems; see for instance [51]. It should be reminded that throughout this paper, we will assume the existence of fixed points for a fixed time interval [52].

Fixed point and stability
This subsection overlaps with previous findings [40] and recaps some of those results. We consider the analysis of fixed points in the UV and IR for Eqs. (2.4) determined by finding the zeros of the expressions on the right-hand sides. This basically amounts to solving for points at which d • /dρ = 0 and taking the limit ρ → ∓∞ in the solutions. The original nonautonomous dynamical system (2.4) can be made autonomous by defining a flow time variable τ = tanh ρ. The resulting autonomous system reads 10 where we have defined X := (T,π, τ ). One of the major advantages of performing this transformation is that the fixed points are no longer at infinity since τ = ∓1 [40]. The fixed points X c defined by F i (X) = 0 for any index i are given by The stability for each fixed point can be determined by slightly perturbing the solution in the vicinity of it. This leads to a linearized form of the ODEs which captures whether a small deviation around the fixed point X → X c + δX is attractive (a sink), repulsive (a source), a combination of both (a saddle), etc. The importance of these qualities from a physical standpoint is that we can label the points at which the flow lines are starting to diverge from or approach ultimately in the UV and IR limit by identifying these points as either stable, marginally stable, unstable, etc. For example, hydrodynamic equilibrium is always a stable point, which means that every flow line in the vicinity of this point will be attracted to it. The generic form of a linearized ODE reads where J a is called Jacobian matrix. For the autonomous system (2.5)-(2.7), J a is defined by (2.10) One can easily obtain the eigenvalue of J a as for Σ U,I(C) . • One source if all the eigenvalues are positive; • One sink if none of the eigenvalues are positive; • One index-1 saddle point if all eigenvalues are real with only one being positive; • One index-2 saddle point if all eigenvalues are real with only two being negative; • One spiral attractor if there is one real eigenvalue and one complex-conjugate pair with all having negative real parts; • One spiral repeller if there is one real eigenvalue and one complex-conjugate pair which have positive real parts; • One index-1 spiral saddle point with one positive real eigenvalue and two complex conjugate eigenvalues with negative real parts; • One index-2 spiral saddle point with one negative real eigenvalue and two complex conjugate eigenvalues with positive real parts.
The value of the eigenvalues is directly linked with the behavior of the solutions of the linearized equation around the fixed point. For example, for Σ U,I(A) one finds that where z ± = e ±2ρ , σ i is an integration constant and the Lyapunov exponents correspond to the eigenvalues governing the damping and/or growth of the perturbation close to the fixed points. It is important to emphasize that on the one hand, the solutions of the linearized ODEs give the transmonomial of the formal transseries, which are usually identified as a non-perturbative term in the sense that exponentially suppressed terms do not admit a perturbative expansion [36,49,55,56]. Since the linearized solutions (2.16) do not have the typical exponential transmonomials when written in terms of z ∓ , the solutions are deemed merely perturbative even though there might be power-law corrections. On the other hand, the asymptotics of IS theory describes a far-from-equilibrium state due to the curvature of dS 3 × R when ρ → ±∞, and it thus cannot have a hydrodynamic gradient expansion. Another way to see this is that the solutions of IS theory for the Gubser flow do not scale asymptotically with the Knudsen number Kn ∼ tanh ρ/T (ρ). In the next subsection, we will explain why the expansion parameter for the formal transseries must be z ∓ = e ∓2ρ 12 .

Global structure
The nonautonomous system is a group of ODEs which depend explicitly on time. The most distinguishing aspect of this type of system compared to an autonomous system is that in the latter, a flow is invariant under time translations, whereas in the former, it does care about the initial value of time. In other words, an explicitly time-dependent dynamical system does entail a time memory for every flow line on the contrary to an autonomous system in which the flow dynamics does not keep track of the initial time at which we considered the flow. Therefore, a flow which is sufficiently well-characterized on a time-independent smooth manifold of configuration-space variables, needs to be generalized to a smooth fiber bundle. This basically makes the space of variables in the configuration space a spacelike hypersurface, e.g., fiber or transverse manifold, that evolves over the time direction -also called base manifold 13 . Additionally, every spacelike hypersurface is labeled by the time step it is retrieved at. This is reminiscent of a spacelike hypesurface in relativity where the normal vector determines the time state throughout the hypersurface. Hence, time direction indeed becomes part of the manifold over which we let the flow lines be defined in a nonautonomous dynamical system such as Gubser flow. For the purpose of solving the whole system, it is then of utmost importance to try to get a handle on both global and local structure of the flow lines over this fibrated manifold. In order to understand the global structure, we employ skew product flow formalism (SPFF) [45,53]. This approach will require us to first construct a vector bundle with the fiber being the configuration space itself and a base space of the flow time as stated earlier, then introduce the flow lines on both spaces independently. The original nonautonomous dynamical system is realized using flows on the total space by taking an appropriate flow equation on the base space. If we fix the time, the resulting hypersurface (fiber) gives a dynamical system on the fiber space, which is a projection of the flows defined on the total system onto the fiber space. This means that on every fiber, time acts as a fixed parameter so we can consider a parameter-dependent flow structure on the fiber space and study how it transforms under any variation of the time parameter -that is not a variable anymore. One interesting phenomenon that could possibly happen is bifurcation, namely whether or not the flow hits a "turn-off" where it could go either way. This provides a powerful method for understanding the entire structure of original nonautonomous system through the lens of an autonomous system [40].
Let us define a trivial vector bundle Θ = (M, π, B, F ) for the nonautonomous dynamical system at hand. We define 14 (Total space is sometimes called the extended configuration space or extended phase space.) From a technical point of view, a fiber of a map f : X → Y is the preimage of an element y ∈ Y . Similarly, to define a τ -fiber of the projection map π defined above, we view M τ as the preimage of τ ∈ B, that is π −1 (τ ) =: {τ } × M τ . Then, we define the nonautonomous dynamical system by the SPFF [45,53].
Definition 1 A nonautonomous dynamical system (θ, φ) is defined by a continuous cocycole mapping φ on the fiber space F which is driven by an autonomous dynamical system φ acting on the base space B and time set T = R λ. The dynamical system θ on B is a group of homeomorphisms θ λ∈T with the properties that 24) and the cocycle mapping φ : In this definition, the cocycle mapping 15 has a semi-group structure, but replacing T + 0 → T yields a group instead since the system is invertible along the time direction. Namely, under time translations τ → τ + τ 0 , every positive time τ > 0 is mapped to an inverse time represented by −τ < 0 and vice versa. M is called a nonautonomous set. If there exists a submanifold S ⊂ M such that (θ λ (τ ), φ(λ, S)) = S for any λ ∈ T, we call S an invariant (sub)manifold of the flow. Strictly speaking, the IS theory is only invertible in the invariant manifold because it does not intersect any other submanifold.
As explained in the introduction, the explicit time dependence does not allow for the dynamical equilibrium state to form unless the flow continues as time tends to infinity. This means that we should distinguish between the future and past domain of M. Let us define a subspace given by M −(+) := τ ≤0(≥0) {τ } × M τ as the past(future) domain of the dynamical system.
The dynamical system (θ, φ) can be identical to the solution of the ODEs (2.5)-(2.7) for a given initial condition (τ 0 , X 0 ) at ρ = ρ 0 as θ 0 (τ 0 ) := τ (ρ 0 ) = tanh ρ 0 , φ(0, τ 0 , X 0 ) := X(ρ 0 ). Moreover, the flow parameter λ and the flow time ρ can be related to each others as θ λ (τ 0 ) = tanh(λ + arctanh τ 0 ), hence, ρ = λ + arctanh τ 0 . By setting an initial condition (τ 0 , X 0 ) ∈ M, the dynamical system (θ, φ) gives a one-dimensional orbit on the vector bundle γ τ 0 (X 0 ) ⊂ M, which is a section determined by the ODEs, as It is notable that due to the uniqueness theorem, all the flows are distinguished only by the initial condition X 0 once τ 0 is fixed. This property is not dependent on τ 0 in the autonomous flows, tying the fate of flow lines solely to the initial condition X 0 . We would like to investigate the flow structure on the τ -fiber M τ and its τ -dependence. Instead of Eq. (2.7), if one takes θ λ (τ ) to be an identity mapping, i.e. θ λ (τ ) = τ for any λ ∈ T, the resulting flows will live be on the τ -fiber M τ . What this means is that a flow on M τ determined by Eqs. (2.5) and (2.6) is now treated as if τ is a free control parameter: The solution depends on both ρ and τ , so it is written as X = X(ρ; τ ). Note that for a given initial condition X 0 (τ ) := X(τ, ρ 0 ), ρ relates to some λ ∈ T via a translation ρ = λ + ρ 0 . The fixed points on the τ -fiber are determined from Eqs. (2.5) and (2.6) with dT /dρ = dπ/dρ = 0 to be One can immediately find that the solutions (2.30) converge to the fixed points given by Eq.
. The stability around each fixed point is found again using the standard method to be given by Also, it is obvious from these formulas that we can, for example, read the instability index for each individual fixed point in the UV by calculating those in the IR and using the symmetry relation across the boundary of past and future domains, i.e., τ = 0. We find  which suggests that the Hartman-Grobman theorem [57,58] 16 does not apply at this point where the kernel of the Jacobian matrix is nontrivial. Despite the fact that this, roughly speaking, tells us that the linearization of the dynamical system at τ = 0 will fail to capture the flow structure both qualitatively and quantitatively, nonetheless the flow still exists at this point in the extended configuration space simply because the dynamical system is regular there. Hence, there should exist a bounded invariant manifold S ⊂ M for the flow such that To investigate the relation between the invariant manifold and the associated fixed points on a more microscopic level, we define the projection map π S : S → B from S onto the base space, and the τ -fiber S τ is defined by π −1 S (τ ) =: {τ } × S τ for τ ∈ B 17 . We observe that since Σ U,I (A) , Σ U,I (B) , Σ U,I (C) are isolated on M τ =∓1 , one can decompose the invariant manifold S into disconnected invariant submanifolds associated with the fixed points Σ U,I(•) and classify them. Formally, this classification might be carried out by introducing a structure group acting on the τ -fiber S τ . The procedure is explained as follows. Suppose that U α and U β are open sets on B with U α ∩ U β = ∅ and define a section as where '=' stems from the fact that the fiber bundle is globally trivial. By constructing the composite mapping on (2.40) The structure group Diff(S) can give away all there is to know about the disconnected submanifolds on the invariant manifold of the flows S, the number of which can be calculated from the topological properties of the quotient space S/Diff(S). According to the definition of S in Eq. (2.38), overall there might exist at most nine disconnected invariant submanifolds such that where a, b = A, B, C and Σ U,I(•) are given by Eq. (2.8). Since the structure of fiber space changes at τ = 0, we decompose the invariant manifold S into the past and future domains in the following way: We consider the classification of the invariant manifold. This procedure essentially requites determining a partially ordered set -so-called poset -by introducing a partial order denoted by < into an arbitrary set of fixed points.
In order to understand the flow structure of the Gubser flow, we would like to follow two crucial steps that will be fundamental to the discussions of the existence of complete UV/IR flows and attractors to come: 1. We first give a basic mathematical definition for a pullback/forward attractor or repeller on the past and future domains. A pullback/forward attractor or repeller constitutes a Morse set either on the past or future domain, which explains the (compact) invariant submanifolds of a nonautonomous flow. 17 In simpler terms, Sτ is the evolving spacelike hypersurface of the invariant manifold S, which is nothing but the material hypersurface of the flow described by the IS theory. This follows the language of fluid dynamics, in particular the theory of Lagrangian coherent structures (LCSs) which are codimension-one invariant manifolds in the extended configuration space of the time-aperiodic and finite-time nonautonomous systems [59].
2. By gluing flows on the past and future domains using the notion of connecting sets introduced in this paper, we establish the existence of flows on the invariant submanifolds satisfying Eq. (2.41).
Before doing this analysis, for later purposes we split the nonautonomous set into two subspaces One can see that the subspace with T = 0 is closed under the flow which means that there exists no flow crossing the subspace. Thus, M T ≥0 and M T ≤0 are also closed under the flows independently and one needs to look at only the one side. From now on, we would concentrate on the flow structure on M T ≥0 , but the structure of M T ≤0 can be obtained in a similar way described below. As a first step, we give a general definition for a local attractor and repeller [45,53].
Definition 2 (Local pullback and forward attractor(repeller)) Let (θ, φ) be an invertible skew product flow. A non-empty compact and invariant nonautonomous set holds for all non-positive(negative) τ ∈ B; 2. a local forward attractor(repeller) if there exists ∈ R + such that where D (A τ ) denotes an open disk of radius on the fiber space F and dist F is defined In what follows, we take an Euclid metric to define a distance on the fiber space. The pullback attractor R and forward repeller R * constitute an attractor-repeller(AR) pair (R, R * ) which gives a Morse decomposition on the past domain [53]. In the same sense, the pullback repeller A and forward attractor A * gives a pair (A, A * ) providing a Morse decomposition on the future domain. If one can find a pullback attractor (or repeller), it automatically means that there is a flow solution for which we can fix the current time and go back to history by sending the initial time τ 0 → −1 or equivalently ρ 0 → −∞.
In the case of M T ≥0 , one can immediately find from where R a ⊆ S a− and A a ⊆ S a+ satisfying Def. 2. These pairs naturally define the partial ordering among the fixed points, where < ∓ denotes the ordering on the past and future domains. It is next required to glue the flow structure on the past and future domains. To this end, we define the partial order for invariant submanifolds on the past and future domains as follows (2.50) Normally, computing the intersection in Eq. (2.50) begs to find the local attractors and repellers, but instead we take a different approach in this study due to the difficulty of solving the ODEs generally. But before divulging the technical details, at least we can trivially verify that the existence of constant solutions (T,π) In order to consider the other possibilities, we define a compact and connected subspace B ±δ on M T ≥0 ±δ with 0 < δ 1 as To incorporate the flow reaching as far as the point Σ U(C) , we choose to be large enough such that Σ U(C) ∈ Int(B −δ ). Then, we introduce a connecting setÑ bỹ −δ , it can be deduced thatÑ < R a . We should notice that any flow for −1 < τ ≤ −δ starting at a Σ U(a) can be continuously deformed into Σ a (−δ) due to the fact that the topology of N − does not change for any τ < 0, and therefore there cannot exist any singularities on the τ -fiber. In a similar way, if there exists a flow which connects Σ a (+δ) andÑ on M T ≥0 +δ , it can immediately lead to A a <Ñ . Although the connecting set N in M T ≥0 −δ connects all Σ A,B,C (−δ) through flows at τ = −δ, only Σ A (+δ) connects toÑ at τ = +δ.
Following the discussion above, the nontrivial posets are then given by A A < R B and A A < R C . Hence, we obtain the flow diagram on M T ≥0 as where the dashed and solid arrows denote < ∓ and <, respectively. In a similar way, the flow diagram on M T ≤0 is given by 18 Finally, the invariant manifold S can be decomposed into six submanifolds as For M T ≤0 , we define the compact and connected subspace B ±δ by and choose large enough such that Σ U(C) (+δ) ∈ Int(B +δ ).
One of the interesting questions that one might address is how do we find the invariant submanifold S ab having a nonempty set? This might be a very nontrivial question at first sight, but the transseries constructed in Sec. 2.2 does in fact provide a topological interpretation of the flow structure, thus indirectly answering our question. This follows next.

Transseries construction
In this subsection, we will consider the transseries around the stable IR fixed point (nonthermal equilibrium) Σ I(A) , but the transseries around Σ U(A) (free streaming point) can be constructed with a similar strategy. We refer to Ref. [55] for technical details.
Let us identify first the appropriate expansion parameter. Since the radius of convergence for the expansion rate of the Gubser flow D µ u µ = tanh ρ is infinity, an asymptotic power expansion in terms of 1/ρ would not work. Furthermore, as we have shown in the previous section, the leading-order behavior of the asymptotic solutions for the IS theory is determined by Eqs. (2.16) and (2.17). Thus, the correct expansion parameter which captures the IR limit is obtained by z → +∞. Rewriting the ODEs of (2.4a) and (2.4b) in z variable, we get For the sake of simplicity, we shift the variables T → T + T c andπ →π +π c where the constantsπ c , T c stand for the coordinates of the fixed point Σ I(A) given by Eq. (2.8). To our benefit, this transformation has brought the fixed point to the origin of configuration space (T → 0,π → 0 as z → +∞). As a result, the modified version of the ODEs (2.61) is where X = (T,π) , and the matrices are given by . (2.65) Next, we diagonalize B by way of multiplying the ODE in (2.62) by an invertible matrix U from the left: where the invertible matrix U is given by Notice that T andπ inC(X) andD(X) can be simply obtained using X = U −1X . The transseries ansatz is constructed from the solution of the linearized ODEs (2.16) and (2.17). These solutions can then be directly taken as transmonomials of the transseries, which essentially encode the information of nonperturbative contributions to the perturbative expansion in 1/z. Hence, the transseries ansatz is given bỹ where m ∈ N 2 0 , β ∈ R 2 , and σ ∈ R 2 are the integration constants. The ansatz (2.72) can be regarded as an element of a vector space V over R equipped with the basis Φ m k (k ∈ Z). Since the basis Φ m k satisfies the properties the vector space V is said to form an R-polynomial ring R[z −1 , σ 1 z β 1 , σ 2 z β 2 ], that is closed under all the operations allowed in the ODE. Notice thatC(X) andD(X) in (2.66) might be expanded in the same basis as where b i is already written in Eq. (2.68). We next normalize the integration constants as X Once again since X (0) 0 = 0, the fact that the formal transseries does not include exponential transmonomials is sufficient to prove that we have a convergent series. Now the formal transseries expanded around the UV fixed point Σ U(A) given by Eq.(2.8) is calculated using the time coordinatez = 1/z. The UV limit is attained byz → +∞, and the ODE is given by Note that the difference between Eqs. (2.61) and Eqs. (2.81b) is that the sign of last term in the second line of the former. Therefore, given the formal transseries (2.79), it would possible to directly turn it into the one around the UV fixed point Σ U(A) by flipping the sign of relaxation scale c → −c. Indeed, the transformation z →z together with c → −c gives an isomorphism between the transseries, i.e.

Double expansion of the transseries
The expansion parameter of the UV/IR transseries derived previously is determined by the asymptotic properties of the IS theory. However, it is instructive to know the transseries in Knudsen number defined as w = tanh ρ/T (ρ). In this case, we rewrite the evolution equation from the original ODEs in (2.4a) and (2.4b) in terms of w. As shown in Ref. [40], this variable does not lead to the correct asymptotic behavior. However, it is possible to construct a formal transseries by following the procedure outlined below. Instead of trying to solve Eqs. (2.4a) and (2.4b), one can first think of completely eliminating T (ρ) from these ODEs by introducing the parameter w = tanh ρ/T (ρ), that gives rise to the dynamical system [40] ∂A(w, ρ) ∂w = F A (A(w, ρ), w, ρ), (2.82) ,π(w, ρ) = 3A(w, ρ) + 2. (2.84) Trying to construct the transseries solution for the nonlinear ODE (2.82), we notice the appearance of both ρ and w in (2.82) simultaneously. Although a transseries solution depending on two variables is available as asymptotic double expansion in terms of two generators 1/z = e −2ρ and 1/w, we need to include a log-type transmonomial as well.
Let us suppose that the flows converge to the IR fixed point in the w, z → +∞ limit. Given the IR fixed point of Eq. (2.82), A c = 1/(3 √ 5) − 2/3, we can linearize the ODE of (2.82) at this value using (2.86) shows that the transseries of A(w, z) is an asymptotic expansion with leading-order term w β followed by corrections characterized by 1/w, log w and 1/z. Here, the initial condition is decoded in the integration constant σ. This result indicates that the formal transseries solution of Eq. (2.82) can be expressed in the form where A (n) (k,s,h) ∈ R are the expansion coefficient and σ ∈ R. It is important to emphasize that the Borel transform of (2.87) even in the presence of log w does not have singularities on the Borel plane. Again, we want to stress that having no exponential corrections in the expansion is a direct indicator of being perturbative in nature and therefore (2.87) should be convergent. To check this out quickly, we take z → +∞ but keep w finite, i.e. fix h = 0 in the solution (2.87). Then, the formal solution becomes (2.88) One can readily find that this series is convergent. As far as Borel summability of A(w, +∞) goes, a useful way is to take the inverse Laplace transform of the ODE (2.82). Suppose that A(w, +∞) is an asymptotic power series about the IR fixed point, i.e. A(w, +∞) =  The original asymptotic series A(w, +∞) can be reproduced by computing the Laplace integral of Y (p), 19 Notice that (c1 * c2)(p) = O(p) for constant c1 and c2, because If Y (p) is not Laplace transformable for some θ, all the terms on the r.h.s. of Eq. (2.95) must entail a singularity across the board, but here they do not for any arbitrary value of θ. Therefore, Y (p) is an entire function on the Borel plane and thus, the asymptotic series of A(w, +∞) converges 20 .

Moment expansion method for the kinetic Gubser flow
In this section, we generalize the ideas of nonautonomous dynamical systems to the case of the Boltzmann equation within the relaxation time approximation (RTA) for the Gubser flow. In what follows, the problem of solving the Boltzmann equation is basically mapped to the study of a dynamical system over the space of moments of the distribution function and de Sitter time.
The symmetry group of the Gubser flow, SO(3) q × SO(1, 1) × Z 2 , restricts the number of variables on which the distribution function depends [6,9]. It was demonstrated that the on-shell distribution function for this flow is given by f (x µ , p i ) = f (ρ, p 2 Ω , |p ς |) [6,9]. Thus, the RTA Boltzmann equation is reduced to the following relaxation equation [6,9] ∂f (ρ, is the relaxation time scale, which is taken to be a positive constant. In Eq. (3.1), the angular momentum p 2 Ω and the energy p ρ are given by in which p ς is the longitudinal momentum in the boosted frame. We assume that the equilibrium distribution function has a functional form of the Maxwell-Boltzmann type, i.e. f eq (x) = e −x . In statistical mechanics it is known that for systems invariant under a certain scaling, the zeros of the collisional kernel C[f ] = 0 admit time-independent solutions, often called steady states, which do not correspond to the global thermal equilibrium state for 20 For example, if one considers where Bf (w) = F (p). In this case, there exists a singularity at p = λ because Indeed, the Borel transformed asymptotic series of the ODE (2.96) is Borel nonsummable. This explains why the asymptotic expansion series for the Bjorken flow is Borel nonsummable since the linearized differential equation for the function A is of the form (2.96) [36,49].
which f (x, p) = f eq (p µ β µ (x)). As we shall see, the stability analysis for the moments of the distribution function indicates the dynamical equilibrium state is not a thermal one; rather it is a non-thermal steady state. This is at odds with the equilibrium state for the Bjorken flow [36]. Hence, even though the equilibrium distribution function satisfies C[f ] = 0, there might be another zero of the collisional kernel that defines a global equilibrium state [60] 21 .
We expand the distribution function in terms of orthogonal polynomials 22 where P (x) and L with θ(x) being the step function. Notice that for the Maxwell-Boltzmann type distribution, we have to set c 00 (ρ) = const. ∈ R + and c n (ρ) = 0 for ( , n) = (0, 0). Knowing the distribution function allows one to calculate the energy-momentum tensor T µν := p µ p ν f . For the ansatz (3.3), the non-zero components are found to only depend 21 In the RTA approximation, temperature is an overall function in the collisional kernel. So even if f −feq does not vanish as time tends to infinity, the vanishing of temperature guarantees that we have a global non-thermal equilibrium for the Boltzmann equation. 22 In our previous publications [36,49] (2), there should be some symmetry breaking conditions in going to the momentum space that would reduce the SO(2, 1) down to SO(1, 1), which is in part fulfilled with constraints coming from the geometry of spacetime and conservation laws. It would then be more appropriate to say that the symmetries of f are in general those of the phase space from the perspective of Hamiltonian dynamics. Yet, here the expansion in (3.3) corresponds to the genuine representation theory of the momentum space symmetries whose generators commute with a Hamiltonian that should come from the Boltzmann equation. But due to the collisional kernel being a nonlinear operator, we cannot get this Hamiltonian for the problem at hand. Thus, we are at a loss to say anything on the representation theoretical aspects of (3.3) and would rather trust the fact that it happens to solve the Boltzmann equation.
(3.16) -25 -Eqs. (3.8) as well as (3.9) make up together what we will refer to as the Gubser dynamical system. Practically speaking, it is not possible to solve the Gubser dynamical system in its full-fledged form, which is basically an infinite-dimensional system of coupled differential equations. But one can solve the system systematically with the truncation of the upper bounds and n as ≤ L ∈ N, n ≤ N ∈ N 0 in the ansatz (3.3). A surprising aspect of Gubser dynamical system is that even though it does admit a formal transseries solution, but there are no exponential transmonomials in contrast to the transseries solution of the Bjorken dynamical system. Therefore, the theory is completely perturbative. However, in both systems the evolution of the energy momentum tensor is fully determined by the moments c 0l := c l [36,49]. For this reason, we narrow down our search for the transseries structure and properties of the dynamical sub-system where the coefficients are given by , (3.20) . (3.21)

Flow structure
In this subsection, we investigate the local and global flow structure of the dynamical system (3.17) and (3.18).

Fixed point and stability
Similarly to the IS theory studied previously, we take τ to be a variable of flow time ρ in the nonautonomous system (3.17) and (3.18), and introduce an ODE for τ , which results in the autonomous dynamical system Here, we have defined the vector C whose components stand for the variables of the system, i.e., C i := (c 1 , · · · , c L , T, τ ) where 1 ≤ i ≤ L + 2. This vector essentially unveils the structure of configuration space at τ . The general solution of Eq. (3.24) is given by τ (ρ) = tanh(ρ − ρ 0 ), where ρ 0 is the integration constant. Evidently, the original nonautonomous dynamical system is reproduced by setting ρ 0 = 0. LetC i := (c 1 , · · · ,c L ,T ,τ ) be a fixed point of the above dynamical system, satisfying dC/dρ = 0, i.e.,T , for L = 21. (3.29) The stability of the dynamical system at every fixed point is demonstrated by calculating the Lyapunov exponents or simply the eigenvalues of the Jacobian matrix where a stands for the fixed point a. We then find the instability index for each fixed point to be   . . . , b L+2 ) on the complex plane. One can see that both the UV and IR fixed points in the truncated system for L ∈ 2N are unstable. Therefore, even L does not contain an invariant manifold for the flow space and would not reproduce a physical dynamical system as L → ∞. A similar conclusion was drawn for the truncation of Bjorken dynamical system in [36]. In contrast, the case of odd L stands  out in that we have a variety of different fixed points just as in the IS theory, with Σ U(A) being a spiral source, Σ I(A) being a spiral sink, and rest of fixed points being saddle points. Hence, a closed invariant manifold might be found in the flow space. For this reason, we always assume hereafter that L is odd. By calculating the expansion in (3.3) at the fixed point (3.29), one can immediately deduce that the Gubser flow does not reach a thermal equilibrium point in the IR limit, However, the distribution function correctly converges to a steady state because our collisional kernel C[f ] is proportional to T (ρ), and thus, This result will remain the same for any L and N .

Global structure
We consider the global flow structure of the system in Eqs. (3.17) and (3.18). The analysis is similar to the what was done for the IS theory in Sec. 2.1.2. Since the results will also be similar to those of IS theory, we are going to avoid divulging the technicalities and results and rather focus on the conclusions and main differences distinguishing the dynamical current system from IS theory. Moving forward, we will again assume L ∈ 2N − 1 for obtaining a source (to mimic the high energetic streaming of particles) and a sink (nonthermal equilibrium) in the UV and IR limits respectively, as explained in Sec. 3.1.1. The following is a quick sketch of the setup and what needs to be understood from the resulting flow analysis.
1. The vector bundle Θ = (M, π, B, F ) is given by Fiber space : F = R L+1 c = (c 1 , · · · , c L , T ). (3.38) The nonautonomous system (θ, φ) is given by the solution of in which we set τ (ρ) = tanh ρ for θ ρ . Furthermore, from the perspective of the Morse decomposition, the global structure of the flow space is unfolded in a seemingly similar way as in discussion of IS theory, giving rise to the following connection diagram

Transseries construction
In this subsection, we construct the transseries solution expanded around the IR fixed point Σ I(A) . The transseries analysis for the Gubser flow shows that not only it is not possible to interpret hydrodynamic or non-hydrodynamic modes, but also provides yet another evidence for the non-existence of the attractor solution. Later on, we also construct the transseries solution around the UV fixed point Σ U(A) , which can be identically done for all the other fixed points.
For L ∈ 2N − 1 and N = 0, the ODEs are given by Note that c 0 = 1 due to the Landau matching condition. The fixed points at the IR limit can be obtained by solving dc /dz = dT /dz = 0 at z → ∞. Let us denote the fixed points byc = (c ,T ). Notice thatc has an L-dependent nonzero value for any . In our analysis, we choose Σ I(A) to represent the fixed point for whichT = 0. Shifting c → c +c , we get where A T = 2 −c 1 /10 (the shift c → c +c impliesc 0 = 1 and c 0 = 0.). Next, we define the (L + 1)-component vector c = (c 1 , c 2 , · · · , c L , T ) . This will cast (3.47) in the form where A (±) is an (L + 1)-by-(L + 1) matrix, decomposed into a constant part and a cdependent part expressed by A (±) = B (±) + D (±) (c). The matrices B (±) and D (±) are respectively . . .
To simplify the matters, we are going to seek a diagonalized form of B (+) by acting on c by an invertible matrix U , and following the diagonalization procedure as In terms of the variables with a tilde, the dynamical system takes the form The formal transseries ansatz to solve these equations is given bỹ where m ∈ N L+1 0 , β ∈ C L+1 are the vectorial indices and eigenvalues characterizing the first-order corrections in the transseries, and σ ∈ C L+1 are the integration constants. Notice thatũ We note that u (m) k in Eq. (3.57) can be calculated similarly using u Moreover,D (±) (c) might be expanded in the same basis as Substituting the ansatz (3.53) into Eq. (3.52) gives the evolution equation for the coefficients as We can further fixũ ,k and β , the formal transseries solution is written as andc is the constant part appearing due to the shift implemented earlier. Comparing Eq. (3.52) against the prepared form of a generic multi-dimensional rank-1, level-1 vector differential equation around a fixed point at infinity given in [62], we are missing the term proportional toc due to the factor (z + 1) on the l.h.s. of that equation. This means that the system does not accept any further exponential corrections so the formal solution in Eq. (3.61) is merely a perturbative series, and in fact a convergent one 23 .
Following the construction of transseries solutions for the IS theory, the formal transseries around the UV fixed point Σ U(A) is generated by applying z →z := 1/z, A (+) → A (−) in the system of (3.48), making the UV limit be defined asz → +∞. Repeating the analysis leading to Eq. (3.61), one can calculate the early-time formal transseries. This might well be achieved by simply following an isomorohism between the transseries expanded around the UV and IR fixed points, of the form , and replacing θ 0 → −θ 0 .
3.3 Stability analysis when N > 0, L > 0 The Gubser dynamical system has a chain structure, i.e., the solution for the moment c nl depends on the higher-order moments c n+1,l+1 . The moments in the N = 0 sector only couple with the moments in the same sector, and one can express the energy-momentum tensor in terms of a linear combination of them. The higher moments in the N > 0 sector are known to express the energy tail for the distribution function [36,49]. Therefore, the stability analysis of the system including N > 0 sectors in both UV and IR limits gives a global picture of the nonlinear relaxation processes of the Gubser dynamical system. Suppose that we have a system l < L, n < N , the total dimension of the configuration space is I + 2, where I = (N + 1)(L + 1) − 1. Again, the stability analysis can unfold the instability index for each fixed point. Following the same procedure outlined for the IS theory and the sector N = 0, we find that Ind[J U ] around all the fixed points are These values seem to dismiss the existence of pullback(forward) attractors in the UV(IR) limit. However, by reducing the dimension of the fiber space and choosing particular initial conditions for the flows, one can find such a pullback(forward) attractor. This can also be achieved by constraining the integration constants of the transseries solutions to the dynamical system around both UV and IR fixed points.
We conclude this section by presenting in Fig. 5 a comparison between the numerical solutions of the truncated dynamical system (with N = 0, L = 3) (black lines), multiparameter transseries in the UV (green line), IR (red line) and a Taylor series expansion around ρ = 0 up to 10000th order (blue line). The discrepancy in the ρ → ±∞ regime is expected between exact and transseries since z → 1 as ρ → 0. It is also observed that there is a 23 The existence of a transseries solution of the dynamical system written fully in the variable z is enough to have a definite conclusion just on the status of convergence, no matter what the radius of convergence would be. Had we started off with ρ instead, this conclusion would have not been straightforward to draw. good agreement between the Taylor series expansion at ρ = 0 and the exact result which naturally diverges at large |ρ| due to a finite radius of convergence. In all these cases, it is always possible to analytically continue the solutions beyond the original finite domain.

Comparison with the Bjorken flow
Even though the Gubser flow manifests scale invariance and originates from the conformal subgroup of the isometry of the AdS 5 spacetime [47,48], its global flow structure and formal transseries are quite different from those of Bjorken flow within the kinetic theory framework [36,49]. In Bjorken flow, the distribution function is invariant under ISO(2) × SO(1, 1) × Z 2 symmetry and one can derive the evolution equations of the moments in the RTA approximation using the moment method with the orthogonal polynomial functions identical to the ones used in the ansatz (3.53) (see for instance Eq. (2) in Ref. [49]). For N = 0, the Bjorken dynamical system is given by [36,49] where the coefficients A ,B , C are given by Eqs. (3.19)-(3.21) and τ ∈ R + is the Milne time. The UV and IR limits are described by τ → 0 and τ → +∞, respectively. The global structure of the Bjorken dynamical system can be investigated by following the same procedure described in the previous sections. For example, if one takes L = 21, the UV/IR fixed pointsc = (c 1 , · · ·c L ,T ) arē .
In this nonautonomous system, we can obtain two disconnected regions of the flow, and those are characterized by either (i) T (τ ) = 0 for all τ ; or (ii) T (τ ) = 0. If one takes T (τ ) = 0 for all τ , the ODEs of (4.1) reduce to and thus, which is essentially an autonomous system. In generic autonomous systems, the flow structure such as the location of fixed points and their stability do not depend on the flow time. In the case of Eq. When turning on the temperature, however, the ODEs (4.1) and (4.2) describe a nonautonomous system. In this case, there is only one IR fixed point Σ I(C) for the flows, while there exists no UV fixed point because the temperature increases as we approach the UV limit (its time derivative never settles at 0). Since all the components ofc vanish at Σ I(C) , 24 In an autonomous system, the smallest ingredient of the Morse decomposition must be a fixed point, unlike the nonautonomous systems the distribution converges to the thermal equilibrium and all the flows reach Σ I(C) . If one takes w := τ T (τ ) ∈ R + 0 for the flow time, the following dynamical system is derived [36,49] This has two UV fixed points Σ U(A,B) and one IR fixed point Σ I(C) , which are exactly the same as those of the original ODEs (4.1) and (4.2). Because Σ U(A) , Σ U(B) and Σ I(C) are a source, a saddle, and a sink, respectively, the maximal invariant submanifold S might be decomposed into Indeed, Σ U(A) and Σ U(B) can be regarded as an end point of a local forward and pullback attractors, respectively, and the Morse decomposition is given by The AR pair on the past domain naturally determines a flow line, which is given the name "attractor solution" [26]. Next, we would discuss the formal transseries of the Bjorken flow. Assuming the flow time to be τ , one can obtain the formal transseries around both Σ I(B) and Σ I(C) . For the transseries around Σ I(C) , we cannot take 1/τ as the expansion parameter. This fact can be easily seen from Eq. (4.6) by obtaining the leading order behavior of the temperature, T (τ ) ∼ T 0 /τ 1/3 with T 0 ∈ R + . Thus, in order to construct a transseries that is closed under all operations in the ODEs, one has to take τ 2/3 as the expansion parameter. Using z := τ 2/3 , the ODEs are modified as follows [36] where the temperature is retrieved by T (z) =T (z)T 0 /z 1/2 . The formal transseries for the solution of the Bjorken dynamical system is defined over the C-polynomial ring C[z −1 , σ 1 e −Sz z β 1 , . . . , σ L e −Sz z β L ] with S = 3T 0 /(2θ 0 ) and it turns out to be a divergent series. The so-called non-hydrodynamic modes are associated with the exponential part of the polynomial ring, that clearly shows the nonperturbative nature of the solution. Calculating the Borel transform of the asymptotic power expansion, we can find singular points on the real axis of the Borel plane, verifying the Borel non-summability of the series. On the hand, the formal transseries (3.61) solution of the Gubser dynamical system is a convergent series with a finite radius of convergence and there are no exponential terms in the polynomial ring as proven in this paper.
In the scenario where T (τ ) is set to zero for every τ , the formal transseries given by Eq. (4.4) over the polynomial ring C[ σ 1 z β 1 , . . . , σ L z β L ] with z := log τ . Notice that this transseries looks very much like the Gubser solution in that it lacks the exponential profile that existed in the general solution for the case T = 0. Hence, understandably the conclusion on the convergence of the series is again the same as before. The polynomial structure C[ σ 1 z β 1 , . . . , σ L z β L ] entering a transseries universally appears in a generic autonomous system such as RG flow equation of quantum field theories. We invite the interested readers to check Ref. [36] for technical details on the formal transseries solution and flow structure of the Bjorken dynamical system. We summarize the comparison between the Gubser and Bjorken dynamical systems and their transseries solutions in Tab. 1.

Gubser Bjorken
Spacetime Table 1: Comparison between the Gubser and Bjorken flows. In the first column from the right, τ 2/3 and w each stand for a separate candidate for the variable z.

Conclusions and final remarks
As we have shown in this work, the asymptotic series around the UV/IR fixed points for the Gubser flow is a purely perturbative (or power-law asymptotic) with a finite radius of convergence. This literally means that the asymptotic series converges to the exact solution by adding higher order terms within a suitable truncation scheme.
In the case of a divergent series, there exists a truncation order that minimizes the error related to the deviance from the exact solution. This error is controlled by the value of the expansion parameter in the transseries solution. Therefore, if one keeps adding higherorder terms beyond a certain order, the formal transseries starts to veer off from the exact solution, a fact which was observed in the IR asymptotic expansion of the Bjorken flow [36]. It should be noted that calculating the error systematically is a well-established subject in the superasymptotic and hyperasymptotic approximation methods [63].
It is also possible to obtain the asymptotic series with the expansion parameter being τ := tanh ρ ∈ (−1, +1) and setting the initial condition c(τ = 0) =c wherec has to lie in the basin of attraction of the IR fixed point (non-thermal equilibrium) at τ → 1, even thoughc is not a fixed point. We then find the evolution equations in terms of τ from Eqs. (3.17) and (3.18) as and show that the formal transseries is has a power series form, i.e.
where u k ∈ R L+1 . The previous expression is, of course, a convergent series. We should notice that this power series can not directly encode the information of the UV/IR fixed points because c(τ ) fails to be analytic at τ = ±1 (Eqs. (5.1)-(5.2) are not obviously regular at these limits). However, one can in principle translate the power series expansion (5.3) into a transseries solution by expanding around either the UV or IR fixed point and perform an analytic continuation. For the Bjorken flow this procedure was previously proposed in Refs. [36,64] and was analyzed more recently in [43]. It must be emphasized that here we did not consider analytic continuation of the series (5.3) albeit being a very straightforward calculation. We should mention that finding the expansion parameter for a formal transseries solution of a nonlinear ODE is not an easy task. The same goes for the Gubser (or Bjorken) dynamical system, in which the expansion parameter is not an arbitrary variable and should rather be derived from the ODEs themselves. This is not necessarily limited to the transseries and must in fact be taken into account even for the naive asymptotic expansion. For instance, in the Bjorken flow [36], a naive ansatz defined as a power series in terms of τ −1 would not work for solving the dynamical system. Instead, one should use τ −2/3 as the expansion parameter. In general, for a given ODE the appropriate ansatz for the expansion parameter is constructed by ensuring that all the properties of the original theory are kept intact. So it is very risky to introduce an arbitrary expansion parameter, as in bookkeeping parameter methods, due to the basic fact the modified ODE might break the well-defined asymptotic properties of the original theory. When one assumes the ansatz for the distribution to be of the form f = f eq + αf (1) + α 2 f (2) + · · · with a bookkeeping parameter α ∈ R + , it will eventually be necessary to validate whether or not the asymptotics satisfies |f eq | |f (1) | |f (2) | · · · in the IR limit. This validation test turns out to only pass for the Bjorken flow.
As the last comment, we would like to say a few words on the so-called "attractor solution". This term has been frequently used in many recent papers on hydrodynamic attractors even though there is no clear definition of the term nor are there any established conditions for its existence. At least we can argue that a couple of implicit properties in this regard might come to light by studying these papers which we can relate to: (i) An "attractor solution" is a one-dimensional flow line approaching the local equilibrium in the late-time limit; and (ii) (Most of if not) all the flows get attracted to this flow line at large finite flow time irregardless of what the initial condition is. The techniques used in this paper and Refs. [36,49] based on the Morse decomposition and the transseries construction gives the correct mathematical interpretation for these properties and beyond. In order to obtain the "attractor solution", (i) There must exist a pullback attractor R with Ind(Σ R ) = 1 a and a forward attractor A with Ind(Σ A ) = 0 on the same invariant submanifold (i.e. A < R) where Σ R(A) is the associated fixed point for R(A), (ii) The lowest level transmonomial in the transseries centered at Σ A must be decoupled from the integration constants (e.g. the transseries be of the Gevrey-1 class).
a This condition for the pullback attractor can be generally relaxed if the transseries around the fixed point always couples with the integration constants. In such cases, one can constrain the fiber space by eliminating the extra dimensions associated to the repelling directions.
Only the Bjorken flow (in the w-coordinate) satisfies the conditions (i),(ii) all at once. However, the unfortunate fact is that under this time coordinate transformation, the global topology of the flow for the Bjorken dynamical system is changed dramatically in three ways: (1) The most important thing is that the dimension of the fiber space is lowered by one; (2) There is a singularity that appears at c 1 = 20 as seen in Eq. (4.6); and (3) The flow at T = 0 originally defines an autonomous system in (4.4), suggesting that there is a disconnected sector in the configuration space that disappears once the w-coordinate is used [36]. Technically speaking, the mapping τ → w for the Bjorken flow cannot continuously deform the topology of the w-fiber bundle to that of the τ -fiber bundle. Therefore, when T = 0, namely w = 0 for any τ > 0, the independent autonomous dynamical structure does not appear when writing the dynamical system in terms of the variable w. For T > 0, focusing only on the subfiber space over the mode configuration c i by excluding the Tdirection, does mean that c i can stay bounded despite T being divergent in the UV limit. This aspect of the original theory is similar to the flow structure of the transformed system. It is extremely important to note that if a global attractor is meant by the "attractor solution", however, it would basically ensure the existence of a pullback attractor in which case there is a universal attraction by going to the history of the flow; see, for example, [65,66]. We should notice that in the standard textbook definition (e.g. Chap. 3.2 of [53]), the global attractor possibly does not satisfy both (i) and (ii) above. In such a case, formal transseries around Σ U(A) and Σ I(A) might have, for example, a form similar to transseries of Gubser flow in the UV and IR limits, respectively. For these reasons, we have rather used the term critical line in [36] for putting the emphasis on the fact that this line is just a UV-IR complete flow line and that it lies on the boundary of an invariant submanifold.