Asymptotic Behavior of a Viscoelastic Fluid in a Closed Loop Thermosyphon: Physical Derivation, Asymptotic Analysis, and Numerical Experiments

Abstract

Fluids subject to thermal gradients produce complex behaviors that arise from the competition with gravitational effects. Although such sort of systems have been widely studied in the literature for simple (Newtonian) fluids, the behavior of viscoelastic fluids has not been explored thus far. We present a theoretical study of the dynamics of a Maxwell viscoelastic fluid in a closed-loop thermosyphon. This sort of fluid presents elastic-like behavior and memory effects. We study the asymptotic properties of the fluid inside the thermosyphon and the exact equations of motion in the inertial manifold that characterizes the asymptotic behavior. We derive, for the first time, the mathematical derivations of the motion of a viscoelastic fluid in the interior of a closed-loop thermosyphon under the effects of natural convection and a given external temperature gradient.

(ProQuest: ... denotes non-US-ASCII text omitted.)

Justine Yasappan 1 and Ángela Jiménez-Casas 1 and Mario Castro 2

Recommended by Douglas Anderson

1, Grupo de Dinámica No Lineal (DNL), Escuela Técnica Superior de Ingeniería (ICAI), Universidad Pontificia Comillas, 28015 Madrid, Spain
2, ICAI, Grupo Interdisciplinar de Sistemas Complejos (GISC) and DNL, Universidad Pontificia Comillas, 28015 Madrid, Spain

Received 11 February 2013; Accepted 19 April 2013

This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

1. Introduction

Instabilities and chaos in fluids subject to temperature gradients have been the subject of intense work for its applications in engineering and in atmospheric sciences. In such sort of systems, the fluid displays nontrivial behaviors (as turbulence or the formation of convective rolls) when subject to a heating that competes with buoyancy effects. A traditional approach that goes back to the pioneering work by Lorenz consists of the study of the system under some simplifications. Another approach is to study the controlled setups that capture the underlying complexity of the full system, being a thermosyphon one of those simpler cases [1].

In the engineering literature, a thermosyphon is a device composed of a closed-loop pipe containing a fluid where some soluble substance has been dissolved [2, 3]. The motion of the fluid is driven by the action of several forces such as gravity and natural convection. The flow inside the loop is driven by an energetic balance between thermal energy and mechanical energy. The interest on this system comes both from engineering and as a toy model of natural convection (for instance, to understand the origin of chaos in atmospheric systems).

Thus far, all the works on thermosyphons analyze the behavior of a Newtonian fluid inside the loop, consequently neglecting elastic effects in the system coming from either the fluid itself or the elastic walls of the loop. However, many interesting fluids are known to behave slightly different from the common (Newtonian) fluids in terms of their response to an applied stress and are commonly referred to as viscoelastic. Among them, it is worth emphasizing volcanic lavas, snow avalanches, flowing paint, or biological mucosas membranes.

Here, we consider a thermosyphon model in which the confined fluid is viscoelastic. This has some a priori interesting peculiarities that could affect the dynamics with respect to the case of a Newtonian fluid. On the one hand, the dynamics has memory and so its behavior depends on the whole past history, and, on the other hand, at small perturbations the fluid behaves like an elastic solid and a characteristic resonance frequency could, eventually, be relevant (e.g., consider the behavior of jelly or toothpaste).

The simplest approach to viscoelasticity comes from the so-called Maxwell constitutive equation [4]. In this model, both Newton's law of viscosity and Hooke's law of elasticity are generalized and complemented through an evolution equation for the stress tensor, σ .

The stress tensor comes into play in the equation for the conservation of momentum: [figure omitted; refer to PDF]

For a Maxwellian fluid, the stress tensor takes the form [figure omitted; refer to PDF] where μ is the fluid viscosity, E the Young's modulus, and γ the shear strain rate (or rate at which the fluid deforms). Under stationary flow, (2) reduces to Newton's law, and consequently, (1) reduces to the celebrated Navier-Stokes equation. On the contrary, for short times (where impulsive behavior from rest can be expected) equation (2) reduces to Hooke's law of elasticity.

Memory effects can be well understood from (2) after performing a separation of variables and integrating. Thus, we can rewrite [figure omitted; refer to PDF] where it is clear that the local state of stress σ ( t ) is calculated from the present and values of γ ( t ) with a memory time scale of order μ / E .

In a thermosyphon, the equations of motion can be greatly simplified because of the quasi-one-dimensional geometry of the loop. Thus, we assume that the section of the loop is constant and small compared with the dimensions of the physical device, so that the arc length coordinate along the loop ( x ) gives the position in the circuit. The velocity of the fluid is assumed to be independent of the position in the circuit; that is, it is assumed to be a scalar quantity depending only on time. This approximation comes from the fact that the fluid is assumed to be incompressible, and so [figure omitted; refer to PDF] besides the quasi-one-dimensional assumption. On the contrary, temperature is assumed to depend both on time and position along the loop.

The derivation of the thermosyphon equations of motion is similar to that in [2, 3]. The simplest way to incorporate (2) is by differentiating (1) with respect to time and replacing the resulting time derivative of σ with (2). This way to incorporate the constitutive equation allows to reduce the number of unknowns (we remove σ from the system of equations) at the cost of increasing the order of the time derivatives to second order.

The resulting second-order equation is then averaged along the loop section (as in [2]). Hence, after nondimensionalizing the variables (to reduce the number of free parameters) we arrive at our main system of equations [figure omitted; refer to PDF] where h ( x , v , T ) = l ( v ) ( _{T a} - T ) ; that is, we consider Newton's linear cooling law as in [3, 5-10] also in this paper we consider the diffusion of temperature given by the term ν ( ^{∂ 2} T / ∂ ^{x 2} ) . The parameter [varepsilon] in (5) is the nondimensional version of μ / E which has dimensions of time. Roughly speaking, it gives the (nondimensional) time scale in which the transition from elastic to fluid-like occurs in the fluid.

The model we will take into consideration forms an ODE/PDE system for the velocity v ( t ) , the distribution of the temperature T ( t , x ) of the fluid into the loop, (5) with ν ...5; 0 , where [contour integral operator] = ^{∫ 0 1} ... d x denotes integration along the closed path of the circuit. The function f describes the geometry of the loop and the distribution of gravitational forces [2, 3]. Note that [contour integral operator] f = 0 .

The system of (5) is not a trivial extension of the Newtonian model in [11] due to the first term in the differential equation for the velocity. Specifically, the addition of a term proportional to the second derivative of v is singular , in the sense that it changes qualitatively the character of the equations. The implications of this singular perturbation cannot be ascertained using a standard boundary layer analysis (see the Appendix for details). So a complementary theoretical and numerical approach is mandatory.

We assume that G ( v ) which specifies the friction law at the inner wall of the loop is positive and bounded away from zero. This function has been usually taken to be G ( v ) = G , a positive constant for the linear friction case [2] (Stokes flow), or G ( v ) = | v | for the quadratic law [12, 13], or even a rather general function given by G ( v ) = g ~ ( Re ) | v | , where Re is the Reynolds number, Re = ρ v L / μ . Here we will consider a general function of the velocity assumed to be large [9, 14]. The functions G , f , l , and h incorporate relevant physical constants of the model, such as the cross-sectional area, D , the length of the loop, L , and the Prandtl, Rayleigh or Reynolds number; see [9]. Here, we consider Newton's linear cooling law l ( v ) ( _{T a} - T ) where l represents the heat transfer law across the loop wall, and is a positive quantity depending on the velocity, and _{T a} is the (given) ambient temperature distribution, see [3, 5, 9, 10].

Hereafter, we consider G and l as continuous functions, such that G ( v ) ...5; _{G 0} > 0 , and l ( v ) ...5; _{l 0} > 0 , for _{G 0} and _{l 0} positive constants.

Our contributions in this paper are the following.

(i) To obtain the system of (5) governing a closed loop thermosyphon model with a viscoelastic fluid and to study the asymptotic behavior of this system which is a generalization of the previous model [11]. Thus, although the dynamics of viscoelastic fluids has been extensively studied in the engineering literature, it has scarcely been considered in the applied mathematics community.

(ii) To present an analysis beginning with the well-posedness and boundedness of solution. The existence of an attractor and an inertial manifold is shown and an explicit reduction to low-dimensional systems is obtained. It is noteworthy that we are able to obtain an exact finite-dimensional reduction (87) that may have a much lower number of degrees of freedom.

(iii): To provide a detailed numerical analysis of the behavior of acceleration, velocity, and temperature which includes a thorough study of the various behaviors of the system for different values of viscoelastic fluid and ambient temperature distribution.

(iv) The numerical analysis will show that viscoelasticity induces a chaotic behavior that is not captured by a boundary layer analysis (that would predict the same qualitative behaviors as in the original model in [11] see the Appendix for details) being the new (nontrivial) emergent behaviors induced by the viscoelasticity worth characterizing.

The structure of the paper is as follows. The first section provides a general introduction to the system, explaining briefly the dynamics of a thermosyphon, viscoelastic fluids, and the objectives of this work. In Section 2, we prove the existence, uniqueness, and boundedness of the solutions. In Section 3, we provide a detailed derivation of the dynamics of the system in the inertial manifold as a reduced dimensionality version of the full system of (5). In Section 4, we present the numerical integration of the reduced system of equations valid in the manifold to understand the role of the main parameters of the physical system. Finally, in Section 5, we conclude with a proposal for future works.

2. Well-Posedness and Boundedness: Global Attractor

2.1. Existence and Uniqueness of Solutions

In this section we prove the existence and uniqueness of solutions of the thermosyphon model (5).

First, we observe that for ν ...5; 0 , if we integrate the equation for the temperature along the loop taking into account the periodicity of T , that is, [contour integral operator] ( ∂ T / ∂ x ) = [contour integral operator] ( ^{∂ 2} T / ∂ ^{x 2} ) = 0 , we have ( d / d t ) ( [contour integral operator] T ) = l ( v ) ( [contour integral operator] _{T a} - [contour integral operator] T ) . Therefore, [contour integral operator] T [arrow right] [contour integral operator] _{T a} exponentially as time goes to infinity for every [contour integral operator] _{T 0} .

Moreover, if we consider τ = T - [contour integral operator] T , then from the second equation of system (5), we obtain that τ verifies the equation [figure omitted; refer to PDF] where _{τ a} = _{T a} - [contour integral operator] _{T a} .

Finally, since [contour integral operator] f = 0 , we have [contour integral operator] T f = [contour integral operator] τ f and the equation for v reads [figure omitted; refer to PDF]

Thus, with w = d v / d t , we get ( w , v , τ ) verifying the system (5) with _{τ a} , _{τ 0} replacing _{T a} , _{T 0} , respectively, and now [contour integral operator] τ = [contour integral operator] _{τ a} = [contour integral operator] _{τ 0} = [contour integral operator] f = 0 . Therefore, hereafter we consider the system (5) where all functions have zero average.

Also, if ν > 0 , the operator ν A = - ν ( ^{∂ 2} / ∂ ^{x 2} ) , together with periodic boundary conditions, is an unbounded, self-adjoint operator with compact resolvent in ^{L per 2} ( 0,1 ) that is positive when restricted to the space of zero average functions ^{L per 2} ( 0,1 ) . Hence, the equation for the temperature T in (5) is of parabolic type for ν > 0 .

2.1.1. The Case with Diffusion: ν > 0

We consider the acceleration w = d v / d t and write the system (5) as the following evolution system for the acceleration, velocity, and temperature: [figure omitted; refer to PDF] that is, [figure omitted; refer to PDF] with _{F 1} ( w , v , T ) = - ( 1 / [varepsilon] ) G ( v ) v + ( 1 / [varepsilon] ) [contour integral operator] T f , _{F 2} ( w , v , T ) = w and _{F 3} ( w , v , T ) = - v ( ∂ T / ∂ x ) + l ( v ) ( _{T a} - T ) and the initial data ( w v T ) ( 0 ) = ( _{w 0 v 0 T 0} ) .

The operator B = ( 1 / [varepsilon] 0 0 0 0 0 0 0 - ν ( ^{∂ 2} / ∂ ^{x 2} ) ) is a sectorial operator in ...B4; = ^{... 2} × ^{H per 1} ( 0,1 ) with domain D ( B ) = ^{... 2} × ^{H per 3} ( 0,1 ) and has compact resolvent, where [figure omitted; refer to PDF]

Using the results and techniques of sectorial operator of [15], we obtain Theorem 1.

Theorem 1.

We assume that H ( r ) = r G ( r ) and l ( v ) are locally Lipschitz, f ∈ ^{L per 2} ( 0,1 ) , _{T a} ∈ ^{H per 1} ( 0,1 ) , and l ( v ) ...5; _{l 0} > 0 . Then, given ( _{w 0} , _{v 0} , _{T 0} ) ∈ ...B4; = ^{... 2} × ^{H per 1} ( 0,1 ) , there exists a unique solution of (5) satisfying [figure omitted; refer to PDF] where w = v = d v / d t and w = ^{d 2} v / d ^{t 2} for every δ > 0 . In particular, (5) defines a nonlinear semigroup, S ( t ) in ...B4; = ^{... 2} × ^{H per 1} ( 0,1 ) , with S ( t ) ( _{w 0} , _{v 0} , _{T 0} ) = ( w ( t ) , v ( t ) , T ( t ) ) .

Proof.

We cover several steps.

Step 1. We prove the local existence and regularity. This follows easily from the variation of constants formula of [15]. In order to prove this we write the system as (9), and we have [figure omitted; refer to PDF] where the operator B is a sectorial operator in ...B4; = ^{... 2} × ^{H per 1} ( 0,1 ) with domain D ( B ) = ^{... 2} × ^{H per 3} ( 0,1 ) and has compact resolvent. In this context, the operator A = - ^{∂ 2} / ∂ ^{x 2} must be understood in the variational sense; that is, for every T , [straight phi] ∈ ^{H per 1} ( 0,1 ) , [figure omitted; refer to PDF] and ^{L per 2} ( 0,1 ) coincides with the fractional space of exponent 1 / 2 [15]. Hereafter we denote by || · || the norm on the space ^{L per 2} ( 0,1 ) . Now, if we prove that the nonlinearity F : ...B4; = ^{... 2} × ^{H per 1} ( 0,1 ) ... ^{...B4; - 1 / 2} = ^{... 2} × ^{L per 2} ( 0,1 ) is well defined and is Lipschitz and bounded on bounded sets, we obtain the local existence for the initial data in ...B4; = ^{... 2} × ^{H per 1} ( 0,1 ) .

Using H ( v ) = G ( v ) v and l ( v ) being locally Lipschitz together with f ∈ ^{L per 2} ( 0,1 ) and _{T a} ∈ ^{H per 1} ( 0,1 ) , we will prove the nonlinear terms, _{F 1} ( w , v , T ) = - ( 1 / [varepsilon] ) G ( v ) v + ( 1 / [varepsilon] ) [contour integral operator] T f , _{F 2} ( w , v , T ) = w , and _{F 3} ( w , v , T ) = - v ( ∂ T / ∂ x ) + l ( v ) ( _{T a} - T ) and satisfy _{F 1} : ^{... 2} × ^{H per 1} ( 0,1 ) ... ... , _{F 2} : ^{... 2} × ^{H per 1} ( 0,1 ) ... ... , and _{F 3} : ^{... 2} × ^{H per 1} ( 0,1 ) ... ^{L per 2} ( 0,1 ) ; that is, F : ...B4; ... ^{...B4; - 1 / 2} is well defined, Lipschitz, and bounded on bounded sets. It is possible to prove this by considering _{T a} ∈ ^{L per 2} ( 0,1 ) . In order to prove these properties of the nonlinearity F , let _{U i} = ( _{w i} , _{v i} , _{T i} ^{) t} , and we note that [figure omitted; refer to PDF] where [figure omitted; refer to PDF] and adding ± _{v 1} ( ∂ _{T 2} / ∂ x ), - or + _{v 2} ( ∂ _{T 1} / ∂ x ) and - or + _{v 1} ( ∂ _{T 1} / ∂ x ) in ( 1 ) , we have [figure omitted; refer to PDF] and adding ± l ( _{v 2} ) _{T 1} in ( 2 ) , we get [figure omitted; refer to PDF] and from the previous hypothesis on function l ( v ) , there exists M > 0 such that || _{F 3} ( _{U 1} ) - _{F 3} ( _{U 2} ) || ...4; M ( | _{v 1} - _{v 2} | + _{|| T 1 - T 2 || ^{H per 1} ( 0,1 )} ) ...4; C _{|| U 1 - U 2 || ^{... 2} × ^{H per 1}} and the rest is obvious.

Therefore, using the techniques of variations of constants formula of [15], we obtain the unique local solutions ( w , v , T ) ∈ C ( [ 0 , τ ] , ...B4; ) of (8) which are given by [figure omitted; refer to PDF]

with H ( r ) = H ( v ( r ) ) , [figure omitted; refer to PDF] and using again the results of [15], we get the regularity of solutions. In fact, from the smoothing effect of the equations, we have ( w , v , T ) ∈ C ( [ 0 , τ ] , ...B4; = ^{... 2} × ^{H per 1} ( 0,1 ) ) ∩ C ( ( 0 , τ ) , ^{... 2} × ^{H per 2} ( 0,1 ) ) and ( w , w , ( ∂ T / ∂ t ) ) ∈ C ( ( 0 , τ ) , ^{... 2} × ^{H per 2 - δ} ( 0,1 ) ) , for some positive τ and any δ > 0 . Now, for ... > 0 we have ( w ( ... ) , v ( ... ) , T ( ... ) ) ∈ ^{... 2} × ^{H per 2} ( 0,1 ) , and since F : ^{... 2} × ^{H per 2} ( 0,1 ) ... ^{... 2} × ^{H per 1} ( 0,1 ) is well defined, Lipschitz, and bounded on bounded sets, we have ( w , v , T ) ∈ C ( [ ... , τ ] , ^{... 2} × ^{H per 2} ( 0,1 ) ) ∩ C ( ( ... , τ ) , ^{... 2} × ^{H per 3} ( 0,1 ) ) and ( w , w , ( ∂ T / ∂ t ) ) ∈ C ( ( ... , τ ) , ^{... 2} × ^{H per 3 - δ} ( 0,1 ) ) . Since ... is arbitrary, we obtain the regularity of the local solution.

Step 2. Now, we prove the solutions of (8) for every time t ...5; 0 .

To prove the global existence, we must show that the solutions are bounded in ...B4; = ^{... 2} × ^{H per 1} ( 0,1 ) norm on finite time intervals. First, to obtain that the norm of T is bounded in finite time, we note that multiplying the equations for the temperature by T in ^{L per 2} ( 0,1 ) and integrating by parts, we have: [figure omitted; refer to PDF] since [contour integral operator] T ( ∂ T / ∂ x ) = ( 1 / 2 ) [contour integral operator] ( ∂ / ∂ x ) ( ^{T 2} ) = 0 .

Using Cauchy-Schwartz and the Young inequality and then the Poincaré inequality, since [contour integral operator] T = 0 together with ^{π 2} is the first nonzero eigenvalue of A = - ^{∂ 2} / ∂ ^{x 2} in ^{L per 2} ( 0,1 ) , we obtain [figure omitted; refer to PDF] and using l ( v ) ...5; _{l 0} > 0 , we get [figure omitted; refer to PDF] and we conclude that the norm of T in ^{L per 2} ( 0,1 ) remains bounded in finite time.

Now, we note that differentiating the second equation of (5) with respect to x , we obtain the same equations for || ∂ T / ∂ x || , and considering now || ∂ _{T a} / ∂ x || , we obtain [figure omitted; refer to PDF]

Thus, we show that the norm of T in ^{H per 1} ( 0,1 ) remains bounded in finite time. Then, using || T || bounded for finite time, we prove that | w ( t ) | and | v ( t ) | remain bounded in finite time and we conclude.

2.1.2. The Case with No Diffusion: ν = 0

The system now reads [figure omitted; refer to PDF] where h ( x , v , T ) = l ( v ) ( _{T a} - T ) ; that is, we consider Newton's linear cooling law as in [6-8, 10], and it is no longer of a parabolic type system and is given by [figure omitted; refer to PDF]

To prove the system is well posed, we use the techniques from [16] considering the same transport equation for temperature in different thermosyphon models as [6-8, 10, 16].

We note that if v ( t ) is a given continuous function, then the equation for the temperature can be integrated along characteristics to obtain [figure omitted; refer to PDF] and plugging this into the nonlocal differential equation for the acceleration and into the equation for the velocity yields [figure omitted; refer to PDF]

We note that for _{T 0} , _{T a} ∈ ^{L per 2} ( 0,1 ) and since in this space the translations are continuous isometries, (27) defines a continuous function of time with values in this space. Although we restrict ourselves to ^{L per 2} ( 0,1 ) , many other choices of space are possible for solving problem (26). In fact any Banach space of 1-periodic functions of x having zero mean and in which translations are continuous isometries can be used as an "admissible space" X ; see [16]. In particular ^{W per m , p} ( 0,1 ) , ^{C per k} ( 0,1 ) are admissible spaces between others. Then, we can prove Lemma 2.

Lemma 2.

Let τ > 0 , fix v ∈ C [ 0 , τ ] , and assume that _{T 0} , _{T a} ∈ X where X is an "admissible space"; see [16], in particular _{T 0} , _{T a} ∈ ^{H per 1} ( 0,1 ) . Then, the function given in (27), ^{T v} ∈ C ( [ 0 , τ ] , X ) , is an integral solution of the PDE which is satisfied only if _{T 0} and _{T a} are differentiable. In particular, if _{T 0} , _{T a} ∈ ^{H per 1} ( 0,1 ) , then ^{T v} is continuous with values in ^{H per 1} ( 0,1 ) and satisfies the PDE as an equality in ^{L per 2} ( 0,1 ) , a.e. in time. Moreover, (27) satisfies the following properties:

(i) [figure omitted; refer to PDF]

(ii) if there exist positive constants _{c d} , d = a and d = 0 such that _{T 0} , _{T a} satisfy _{|| T d ( · + h ) - T d ( · ) || X} ...4; _{c d} | k | for all k , then ^{T v} satisfies [figure omitted; refer to PDF]

: positive constant independent on time;

(iii): we assume that X ⊂ ^{L per 2} ( 0,1 ) and there exist positive constants _{c d} , d = a and d = 0 such that _{T 0} , _{T a} satisfy _{|| T d ( · + h ) - T d ( · ) || ^{L per 2}} ...4; _{c d} | k | for all k . If we also assume that _{v i} , i = 1,2 are continuous in t ∈ [ 0 , τ ] , then [figure omitted; refer to PDF]

: K is a positive constant, and _{|| v 1 - v 2 || ∞} = _{sup r ∈ [ 0 , τ ]} | _{v 1} ( r ) - _{v 2} ( r ) | .

Proof.

See [6-8, 10, 16]. Then, we have Theorem 3.

Theorem 3.

Assume that G ( v ) v is locally Lipschitz, f ∈ ^{L per 2} ( 0,1 ) , _{T 0} , _{T a} ∈ ^{H per 1} ( 0,1 ) , and _{w 0} , _{v 0} ∈ ^{... 2} . Then there exists a unique solution of (26) satisfying [figure omitted; refer to PDF] and T satisfies the PDE in the sense of (27). Moreover ( w , v , ( ∂ T / ∂ t ) ) ∈ C ( [ 0 , ∞ ) , ^{... 2} × ^{L per 2} ( 0,1 ) ) .

Proof.

As noted earlier, we need to solve the fixed point problem [figure omitted; refer to PDF] on a space of continuous functions. More precisely, we take W = { v ∈ C [ 0 , L ] , v ( 0 ) = _{v 0} , | v ( t ) - _{v 0} | ...4; M } , endowed with the sup norm, with L and M to be chosen and prove that ... is a contraction on W .

From (29) in Lemma 2 for ^{T v} , we have || ^{T v} || ...4; max { || _{T 0} || || _{T a} || } , and this, together with the local Lipschitz property of G ( v ) v shows that for fixed M , ... ( [physics M-matrix] ) ⊂ [physics M-matrix] if L is sufficiently small.

To show that ... is a contraction, it is clear that we must prove some Lipschitz dependence on [contour integral operator] ^{T v} f with respect to v ∈ W .

First, we note that from _{T 0} , _{T a} ∈ ^{H per 1} ( 0,1 ) , then verify (30), that is, _{|| T d ( · + h ) - T d ( · ) || ^{L per 2}} ...4; _{c d} | k | for all k with d = 0 , a , and given _{v i} ∈ W , again from (31) in Lemma 2, we have [figure omitted; refer to PDF] With these we find that ... is Lipschitz on W with a Lipschitz constant depending on L and M that tends to zero as L [arrow right] 0 and then ... is a contraction for small enough L . Therefore, local well-posedness follows.

To prove the global existence, it is sufficient to prove that ( w ( t ) , v ( t ) ) is bounded on finite time intervals, since from [figure omitted; refer to PDF] we find that ( w ( t ) , v ( t ) ) is of Cauchy type as t [arrow right] _{t 0} for finite _{t 0} . Consequently, the limit of ( w ( t ) , v ( t ) , T ( t ) ) exists in ^{... 2} × ^{L per 2} ( 0,1 ) and the solution can be prolonged.

But again from (29) together with (18) and (19), we obtain boundedness on finite time intervals and global existence follows.

As noted earlier, if _{T 0} , _{T a} ∈ ^{H per 1} ( 0,1 ) , then T satisfies the PDE equation as an equality in ^{L per 2} ( 0,1 ) a.e. in time. In particular, we have ( ∂ T / ∂ t ) ∈ C ( ( 0 , ∞ ) , ^{L per 2} ( 0,1 ) ) .

2.2. Boundedness of the Solutions and Global Attractor

In order to obtain asymptotic bounds on the solutions as t [arrow right] ∞ , we consider the friction function G satisfying the hypotheses from the previous section and we also assume that there exists a constant _{h 0} ...5; 0 such that [figure omitted; refer to PDF] We make use of L'Hopital's lemma proved in [11] to prove several results in this section.

Lemma 4 (L'Hopital's lemma).

Assume that f and g are real differentiable functions on ( a , b ) , b ...4; ∞ , ^{g [variant prime]} ( x ) ...0; 0 on ( a , b ) and _{lim x [arrow right] b} g ( x ) = ∞ .

(i) If _{limsup x [arrow right] b} ( ^{f [variant prime]} ( x ) / ^{g [variant prime]} ( x ) ) = L , then _{limsup x [arrow right] b} ( f ( x ) / g ( x ) ) ...4; L .

(ii) If _{liminf x [arrow right] b} ( ^{f [variant prime]} ( x ) / ^{g [variant prime]} ( x ) ) = L , then _{liminf x [arrow right] b} ( f ( x ) / g ( x ) ) ...5; L .

With this, we have Lemma 5.

Lemma 5.

If one assumes that G ( r ) and H ( r ) = r G ( r ) satisfy the hypothesis from Theorem 1 or Theorem 3 together with (36), then [figure omitted; refer to PDF] with _{H 0} = ( 1 + _{h 0} ) [varepsilon] a positive constant such that _{H 0} [arrow right] 0 if [varepsilon] [arrow right] 0 .

Proof.

Integrating by parts we have [figure omitted; refer to PDF] and using Lemma 4 we obtain [figure omitted; refer to PDF] and from (36) we conclude.

Remark 6.

We note that the conditions (36) are satisfied for all friction functions G considered in the previous works, that is, the thermosyphon models, where G is constant or linear or quadratic law. Moreover, the conditions (36) are also true for G ( s ) [approximate] A | s ^{| n} , as s [arrow right] ∞ .

Now, we use the asymptotic bounded for temperature to obtain the asymptotic bounded for the velocity and the acceleration functions.

Theorem 7.

Under the previous notations and hypothesis of Theorem 1 or Theorem 3, if one assumes also that G satisfies (37) for some constant _{H 0} ...5; 0 , then one has the following.

Part (I). General case:

(i) [figure omitted; refer to PDF]

: In particular: If _{limsup t [arrow right] ∞} || T || ∈ ... , then [figure omitted; refer to PDF]

(ii) If _{limsup t [arrow right] ∞} || T || ∈ ... and ^{G 0 *} = _{limsup t [arrow right] ∞} G ( v ( t ) ) , then [figure omitted; refer to PDF]

Part (II). If ν ...0; 0 and also assume that there exists _{L 0} a positive constant such that _{L 0} ...5; l ( v ) ...5; _{l 0} , then for any solution of (5) in the space ...B4; = ^{... 2} × ^{H per 1} ( 0,1 ) , one has

(i) [figure omitted; refer to PDF]

(ii) [figure omitted; refer to PDF]

(iii): if ^{G 0 *} = _{limsup t [arrow right] ∞} G ( v ( t ) ) , [figure omitted; refer to PDF]

In particular, (5) has a global compact and connected attractor, ...9C; , in ...B4; = ^{... 2} × ^{H per 1} ( 0,1 ) .

Proof.

Part (I). General case.

(i) From (8) we have that [figure omitted; refer to PDF] and w ( t ) = d v / d t satisfies [figure omitted; refer to PDF] where H ( r ) = H ( v ( r ) ) = v ( r ) G ( v ( r ) ) . First, we rewrite (48) as [figure omitted; refer to PDF] with [figure omitted; refer to PDF] Next, for any δ > 0 there exists _{t 0} > 0 such that δ ( s ) = w ( 0 ) ^{e - 1 / [varepsilon]} < δ for any s ...5; _{t 0} and integrating with t ...5; _{t 0} , we obtain [figure omitted; refer to PDF] Using L'Hopital's Lemma 4 proved in [11], we get [figure omitted; refer to PDF] Moreover, using again the L'Hopital's Lemma 4 proved in [11], we get [figure omitted; refer to PDF] and from (51) together with (37), we conclude that [figure omitted; refer to PDF] for any δ .

(ii) From (47) together with Gronwall's lemma, we get [figure omitted; refer to PDF] where G ( r ) = G ( v ( r ) ) . Consequently, for any δ > 0 , there exists _{t 0} such that for any t ...5; _{t 0} , [figure omitted; refer to PDF] that is, [figure omitted; refer to PDF] for any δ > 0 , and using the previous results (i) we get (42).

Part (II). (i) From (23) together with (24), we get [figure omitted; refer to PDF] and by elementary integration we obtain (44). Using Part (I), the rest (ii) and (iii) are obvious. Since the sectorial operator B defined in Section 2.1.1 has compact resolvent, the rest follows from [17, Theorems 4.2 . 2 and 3.4 . 8 ].

Remark 8.

First, we note that the hypothesis about the function l ( v ) in Theorem 7, l ( v ) ...4; _{L 0} is satisfied when we consider Newton's linear cooling law h = k ( _{T a} - T ) , where k is a positive quantity; that is, l ( v ) = k = _{L 0} , as [18]. Moreover, this condition is also satisfied if we consider h = l ( v ) ( _{T a} - T ) where l ( v ) is a positive upper bounded function.

Second, it is important to note that we prove in the next section the existence of the global compact and connected attractor and the inertial manifold for the system (8), when we consider the general Newton's linear cooling law without the additional previous hypothesis on l ( v ) ; but we assume that the friction function G ( v ) always satisfies (36).

In order to get this, we consider the Fourier expansions and observe the dynamics of each coefficient of Fourier expansions to improve the asymptotic bounded of temperature. In particular, we will prove _{limsup t [arrow right] ∞} || T ( t ) || ...4; || _{T a} || for every locally Lipschitz and positive function l ( v ) and also for every ν ...5; 0 (see (72) in Proposition 9) and for every friction function G ( v ) satisfying (36).

3. Asymptotic Behavior: Reduction to Finite-Dimensional Systems

We take a close look at the dynamics of (5) by considering the Fourier expansions of each function and observing the dynamics of each Fourier mode. Assume that _{T a} ∈ ^{H per 1} ( 0,1 ) and f ∈ ^{L per 2} ( 0,1 ) are given by the following Fourier expansions: [figure omitted; refer to PDF] while the initial data _{T 0} ∈ ^{H per 1} ( 0,1 ) is given by _{T 0} ( x ) = _{∑ k ∈ ^{... *}} ... _{a k 0} ^{e 2 π k i x} .

Assume that T ( t , x ) ∈ ^{H per 1} ( 0,1 ) is given by [figure omitted; refer to PDF] Then, we find that the coefficient _{a k} ( t ) in (60) is a solution of [figure omitted; refer to PDF] Since all the functions involved are real, we have _{a - k} = _{a - k} , _{b - k} = _{b - k} , and _{c - k} = _{c - k} . Therefore, (5) is equivalent to the infinite system of ODEs consisting of (61) coupled with [figure omitted; refer to PDF]

The system of (5) reflects two of the main features: (i) the coupling between the modes enter only through the velocity, while diffusion acts as a linear damping term, and (ii) it is important to note in this model that we have also the nonlinear term given by Newton's linear cooling law.

We note that the system (5) is equivalent to the system (8) for acceleration, velocity and temperature. It is equivalent to the following infinite system of ODEs: [figure omitted; refer to PDF]

In what follows, we will exploit this explicit equation for the temperature modes to analyze the asymptotic behavior of the system and to obtain the explicit low-dimensional models.

3.1. Attractors and Inertial Manifolds

The existence of an inertial manifold does not rely, in this case, on the existence of large gaps in the spectrum of the elliptic operator but on the invariance of certain sets of Fourier modes.

A similar explicit construction was given by Bloch and Titi in [19] for a nonlinear beam equation where the nonlinearity occurs only through the appearance of the ^{L 2} norm of the unknown. A related construction was given by Stuart in [20] for a nonlocal reaction-diffusion equation.

In order to get the inertial manifold for this system, we first improve the bounds on acceleration, velocity, and temperature of the previous section for all situations with ν ...5; 0 and l ( v ) ...5; _{l 0} > 0 general locally Lipschitz function under the hypotheses of Lemma 5 for some _{H 0} ...5; 0 , in particular G ( v ) satisfying (36).

We will prove in Proposition 9 that we have always an upper bounded for the temperature in ^{L 2} ( 0,1 ) independent of the velocity and the function l ( v ) , considered in Newton's linear cooling law and also independent of the diffusion coefficient. That is, [figure omitted; refer to PDF]

Proposition 9.

Under the previous notations, for every solution of the system (5), ( w , v , T ) , and for every k ∈ ^{... *} , one has

(i) [figure omitted; refer to PDF]

(ii) [figure omitted; refer to PDF]

: and _{G 0} a positive constant such that G ( v ) ...5; _{G 0} ;

(iii): [figure omitted; refer to PDF]

In particular, one has a global compact and connected attractor ...9C; ⊂ [ - M , M ] × [ - N , N ] × ...9E; where M , N are the upper bounds for acceleration and velocity as given in (67) and (66) and _{T 0} ∈ ...9E; = { R ( x ) = _{∑ k ∈ ^{... *}} ... _{r k} ^{e 2 π k i x} , | _{r k} | ...4; | _{b k} | } .

Proof.

From (61), we have [figure omitted; refer to PDF] with [figure omitted; refer to PDF] Thus, we obtain [figure omitted; refer to PDF] and we get _{limsup t [arrow right] ∞} | _{a k} ( t ) | ...4; | _{b k} | . Using Theorem 7 together with [contour integral operator] T f = _{∑ k ∈ ^{... *}} ... _{a k} ( t ) _{c - k} , we get [figure omitted; refer to PDF] From this upper bounded for the velocity, we also have _{L 0} the upper bound for the continuous positive function l ( v ) , and using Part (II) from Theorem 7 we conclude.

We note from the previous result that we have always the upper bound for || T || and from Theorem 7 for the velocity. Therefore, we can consider _{L 0} the upper bound for the continuous positive function l ( v ) ; we note that _{L 0} = _{limsup t [arrow right] ∞} l ( v ) and we prove in Proposition 10 the bound of solutions to show the influence of diffusion coefficient ν .

Proposition 10.

Under the previous notations, for every solution of the system (5), ( w , v , T ) , and for every k ∈ ^{... *} , one has [figure omitted; refer to PDF] Moreover, one has that [figure omitted; refer to PDF] and _{G 0} positive constant such that G ( v ) ...5; _{G 0} ; [figure omitted; refer to PDF]

Proof.

Using again [figure omitted; refer to PDF] if we assume that 0 < _{l 0} ...4; l ( v ) ...4; _{L 0} , then we obtain [figure omitted; refer to PDF] and we get [figure omitted; refer to PDF] that is, if | _{a k 0} | ...4; ( _{L 0} / ( 4 ν ^{π 2 k 2} + _{l 0} ) ) | _{b k} | , then | _{a k} ( t ) | ...4; ( _{L 0} / ( 4 ν ^{π 2 k 2} + _{l 0} ) ) | _{b k} | .

In particular | _{a k} ( t ) | ...4; ( _{L 0} / ( 4 ν ^{π 2 k 2} + _{l 0} ) ) | _{b k} | ...4; ( _{L 0} / ( 4 ν ^{π 2} + _{l 0} ) ) | _{b k} | for every k ∈ ^{... *} and also | _{a k} ( t ) | ...4; ( _{L 0} / ( 4 ν ^{π 2 k 2} + _{l 0} ) ) | _{b k} | ...4; ( _{L 0} / ( 4 ν ^{π 2 k 2} ) ) | _{b k} | ; that is, | _{a k} ( t ) | ...4; ( _{L 0} / 4 ν ^{π 2} ) | _{b k} | .

Then, we also have that [figure omitted; refer to PDF]

Finally, we note that if _{L 0} = _{limsup t [arrow right] ∞} l ( v ( t ) ) , then given δ > 0 , there exists _{t 0} such that L ( v ( t ) ) ...4; _{L 0} + δ for every t ...5; _{t 0} , and integrating in t ...5; _{t 0} , we obtain _{limsup t [arrow right] ∞} | _{a k} ( t ) | ...4; ( _{L 0} / ( 4 ν ^{π 2 k 2} + _{l 0} ) ) | _{b k} | + δ , for any δ > 0 , and [figure omitted; refer to PDF] Using again Theorem 7, we conclude.

We note that if ν = 0 , then | _{a k} | ...4; | _{b k} | ...4; ( _{L 0} / _{l 0} ) | _{b k} | . If ν > 0 we observe in the numerical experiments, as ν is bigger, that the solution becomes stable or periodic.

As a consequence, we have the following result on the smoothness of the attractor of (5).

Corollary 11.

(i) If | _{a k 0} | ...4; ( _{L 0} / ( 4 ν ^{π 2 k 2} + _{l 0} ) ) | _{b k} | , then | _{a k} ( t ) | ...4; ( _{L 0} / ( 4 ν ^{π 2 k 2} + _{l 0} ) ) | _{b k} | for every t ...5; 0 .

(ii) If ...9C; is the global attractor in the space ...B4; = ^{... 2} × ^{H per 1} ( 0,1 ) , then for every ( _{w 0} , _{v 0} , _{T 0} ) ∈ ...9C; , with _{T 0} ( x ) = _{∑ k ∈ ^{... *}} ... _{a k} ^{e 2 π k i x} , one gets [figure omitted; refer to PDF] In particular, if _{T a} ∈ ^{H per m} ( 0,1 ) with m ...5; 1 , the global attractor ...9C; ... ^{... 2} × ^{H per m + 2} ( 0,1 ) and is compact in this space.

Proof.

(i) From (77) we have | _{a k} ( t ) | ...4; ( _{L 0} / ( 4 ν ^{π 2 k 2} + _{l 0} ) ) | _{b k} | + ( | _{a k 0} | - ( _{L 0} / ( 4 ν ^{π 2 k 2} + _{l 0} ) ) | _{b k} | _{) +} ^{e - ∫ 0 t ... l ( v )} , therefore, if | _{a k 0} | ...4; ( _{L 0} / ( 4 ν ^{π 2 k 2} + _{l 0} ) ) | _{b k} | , then | _{a k} ( t ) | ...4; ( _{L 0} / ( 4 ν ^{π 2 k 2} + _{l 0} ) ) | _{b k} | for every t ...5; 0 and k ∈ ^{... *} .

(ii) We note that from (i), if _{T a} ( x ) = _{∑ k ∈ ^{... *}} ... _{b k} ^{e 2 π k i x} ∈ ^{H per m} , then _{∑ k ∈ ^{... *}} ... ^{k 2 m} | _{b k} ^{| 2} < ∞ , and therefore _{T 0} ∈ ...9E; = { R ( x ) = _{∑ k ∈ ^{... *}} ... _{r k} ^{e 2 π k i x} ∈ ^{H per m + 2} , 4 ν ^{π 2 k 2} | _{r k} | ...4; | _{b k} | } ; that is, ...9C; ⊂ [ - M , M ] × [ - N , N ] × ...9E; where M , N are the upper bounds for acceleration and velocity as given in (74) and (73). But, the set ...9E; is compact in ^{H per m + 2} since for any sequence { ^{T n} } in ...9E; we can extract a subsequence that we still denote { ^{T n} } such that it converges weakly to a function T and such that for any k ∈ ^{... *} , the Fourier coefficients verify ^{a k n} [arrow right] _{a k} as n [arrow right] ∞ , where _{a k} is the k th Fourier coefficient of T . Therefore, 4 ν ^{π 2 k 2} | _{a k} | ...4; | _{b k} | and for every integer E , [figure omitted; refer to PDF] where _{|| · || m + 2} denotes the norm in ^{H per m + 2} . Hence, the first term goes to zero as n [arrow right] ∞ and the second term can be made arbitrarily small as E [arrow right] ∞ . Consequently, T ∈ ...9E; and ^{T n} [arrow right] T in ^{H per m + 2} , and the result follows.

Note that this result reveals in particular the asymptotic smoothing of (5). In the next result we will prove that the dynamical system induced by (8) in the phase space ...B4; = ^{... 2} × ^{H per m} ( 0,1 ) , m ...5; 1 , has an inertial manifold. According to [21], we have the following definition.

Definition 12.

Let S ( t ) , t ...5; 0 , be a nonlinear semigroup in a Banach space in ...B4; that has a global attractor ...9C; . Then, a smooth manifold [physics M-matrix] ⊂ ...B4; is called an inertial manifold if

(i) [physics M-matrix] is positively invariant, that is, S ( t ) [physics M-matrix] ⊂ [physics M-matrix] for every t ...5; 0 ;

(ii) [physics M-matrix] contains the attractor, that is, ...9C; ⊂ [physics M-matrix] ;

(iii): [physics M-matrix] is exponentially attracting in the sense that there exists a constant δ > 0 such that for every bounded set B ⊂ ...B4; there exists C = C ( B ) ...5; 0 such that [figure omitted; refer to PDF]

for every t ...5; 0 .

See, for example, [21, 22].

Assume the ambient temperature given by [figure omitted; refer to PDF] where K ⊂ ... , that is, with _{b k} ...0; 0 for every k ∈ K ⊂ ... with 0 ∉ K , since [contour integral operator] ... _{T a} = 0 .

Then, we denote by _{V m} the closed linear subspace of ^{H per m} ( 0,1 ) spanned by { ^{e 2 π k i x} , k ∈ K } and consider the following spectral decomposition in ^{H per m} ( 0,1 ) : T = ^{T 1} + ^{T 2} , where ^{T 1} denotes the projection of T onto V and ^{T 2} the projection onto the space generated by { ^{e 2 π k i x} , k ∉ K } . Note that (8) is equivalent to [figure omitted; refer to PDF]

Note that from (61), if _{b k} = 0 , then the k th mode for the temperature is damped out exponentially, and therefore the space V attracts the dynamics for the temperature. This is precisely stated in the following result.

Theorem 13.

Assume that _{T a} ∈ ^{H per m} ( 0,1 ) and f ∈ ^{L per 2} ( 0,1 ) . Then, the set [physics M-matrix] = ^{... 2} × _{V m} is an inertial manifold for the flow of S ( t ) ( _{w 0} , _{v 0} , _{T 0} ) = ( w ( t ) , v ( t ) , T ( t ) ) in the space ...B4; = ^{... 2} × ^{H per m} ( 0,1 ) . Moreover, if f ∈ _{V m} the inertial manifold [physics M-matrix] has the exponential tracking property; that is, for every ( _{w 0} , _{v 0} , _{T 0} ) ∈ ^{... 2} × ^{H per m} ( 0,1 ) , there exists ( _{w 1} , _{v 1} , _{T 1} ) ∈ [physics M-matrix] such that if ( _{w i} ( t ) , _{v i} ( t ) , _{T i} ( t ) ) , i = 0,1 , are the corresponding solutions of (8), then ( _{w 0} ( t ) , _{v 0} ( t ) , _{T 0} ( t ) ) - ( _{w 1} ( t ) , _{v 1} ( t ) , _{T 1} ( t ) ) [arrow right] 0 in ^{... 2} × ^{H per m} ( 0,1 ) . In particular if K is a finite set, the dimension of [physics M-matrix] is | K | + 2 , where | K | is the number of elements in K .

Proof.

In order to prove that [physics M-matrix] is invariant, we note if k ∉ K , then _{b k} = 0 , and therefore _{a k 0} = 0 , ( ^{T 0 2} = 0 ) , from (70), we get that _{a k} ( t ) = 0 , for every t ; that is, T ( t , x ) = _{∑ k ∈ K} ... _{a k} ( t ) ^{e 2 π k i x} = ^{T 1} . Thus, if ( _{w 0} , _{v 0} , _{T 0} ) ∈ [physics M-matrix] , then ( w ( t ) , v ( t ) , T ( t ) ) ∈ [physics M-matrix] for every t , that is, invariant manifold.

We consider the decomposition in ^{H per m} , T = ^{T 1} + ^{T 2} , where ^{T 1} is the projection of T on _{V m} and ^{T 2} is the projection of T on the subspace generated by { ^{e 2 π k i x} , k ∈ ^{... *} \ K } ; that is, ^{T 1} = _{∑ k ∈ K} ... _{a k} ^{e 2 π k i x} and ^{T 2} = _{∑ k ∈ ^{... *} \ K} ... _{a k} ^{e 2 π k i x} = T - ^{T 1} .

From (77) taking into account that b k = 0 for k ∈ ... * \ K , we have that | a k ( t ) | ...4; | a k 0 | e - ( 4 ν π 2 k 2 + l 0 ) t for every k ∈ ... * implies that there exist positive constants C i , such that || T 2 ( t ) || H per 2 m ...4; C 1 || T 2 ( t ) || H per m ...4; C 2 || T 0 2 || H per m e - ( 4 ν π 2 + l 0 ) t , that is, T 2 ( t ) [arrow right] 0 in H per 2 m ( 0,1 ) if t [arrow right] ∞ .

Therefore, we have in particular that || T 2 ( t ) || H per m [arrow right] 0 as t [arrow right] ∞ with exponential decay rate e - ( 4 ν π 2 + l 0 ) t . Thus, [physics M-matrix] also attracts ( w ( t ) , v ( t ) , T ( t ) ) with exponential rate e - 4 ν π 2 t , since dist ...B4; ( ( w ( t ) , v ( t ) , T ( t ) ) , [physics M-matrix] ) = dist H per m ( T ( t ) , V ) = || T 2 ( t ) || H per m ...4; C 2 || T 0 2 || H per m e - ( 4 ν π 2 + l 0 ) t .

To prove the exponential tracking property just note that the flow inside [physics M-matrix] is given by setting ^{T 2} = 0 where [contour integral operator] T ( x ) · f = _{∑ k ∈ K} ... _{a k} ( t ) · _{c - k} ; that is, [figure omitted; refer to PDF] Therefore, if f ∈ V , then [contour integral operator] T f = [contour integral operator] _{T 1} f , and given ( _{w 0} , _{v 0} , _{T 0} ) ∈ ...B4; and ( w ( t ) , v ( t ) , T ( t ) ) , the solution of (84), we decompose _{T 0} = ^{T 0 1} + ^{T 0 2} and T ( t ) = ^{T 1} ( t ) + ^{T 2} ( t ) . Then, we consider ( w ( t ) , v ( t ) , ^{T 1} ( t ) ) ∈ [physics M-matrix] and it is still a solution of (85). Hence, ( w ( t ) , v ( t ) , T ( t ) ) - ( w ( t ) , v ( t ) , ^{T 1} ( t ) ) = ( 0,0 , ^{T 2} ( t ) ) and the right-hand side is of order ^{e - ( 4 ν π 2 + _{l 0} ) t} . In particular, if the set K is finite, then the inertial manifold [physics M-matrix] is of finite dimension and the flow inside is equivalent to the finite system of ODEs given by (85). Thus, the theorem is proved.

3.2. The Reduced Subsystem

Under the hypotheses and notations of Theorem 13, we suppose that [figure omitted; refer to PDF] with _{c k} ...0; 0 for every k ∈ J ⊂ ... . Then, [contour integral operator] ( T · f ) = _{∑ k ∈ K ∩ J} ... _{a k} ( t ) _{c - k} . So, the evolution of velocity v and acceleration w depend only on the coefficients of T which belong to the set K ∩ J . Note that in (85) the set of equations for _{a k} with k ∈ K ∩ J together with the equation for v and w is a subsystem of coupled equations denoted by the reduced subsystem.

Thus, we will reduce the asymptotic behavior of the initial system (5) to the dynamics of the reduced explicit system (87) when we consider the relevant modes of temperature _{a k} , k ∈ K ∩ J . Consider [figure omitted; refer to PDF] where _{a - k} = _{a - k} , _{b - k} = _{b - k} , and _{c - k} = _{c - k} as we consider only the real functions. After solving this, we must solve the equations for k ∉ K ∩ J which are linear autonomous equations.

Now, we will show the modes in k ∈ K ∩ J that will play an essential role in the dynamics. With the previous notations, we further decompose _{T 1} as follows: [figure omitted; refer to PDF] where τ is the projection onto the space generated by { ^{e 2 π k i x} , k ∈ K ∩ J } and θ is the projection onto the space generated by { ^{e 2 π k x} , k ∈ K \ J } . We denote by P the projection P ( w , v , T ) = ( w , v , τ ) and Q = I - P . With these notations and decomposing _{T a} as _{T a} = _{τ T a} + _{θ T a} , (8) and (84) can be decomposed as a system of the form [figure omitted; refer to PDF] Since [contour integral operator] θ f = 0 and setting _{T 2} = 0 , the first four equations give the flow inside the inertial manifold [physics M-matrix] ; that is, they are equivalent to (84) while the first three are the only nonlinearity coupled equations. Therefore, once this subsystem is solved, the other unknowns are determined through linear nonhomogeneous equations.

To make this idea more precise in terms of semigroup and attractors, we proceed as in [14]. We denote by S ( t ) the semigroup generated by (8) on ...B4; : = ^{... 2} × ^{H per 1} ( 0,1 ) and by _{S M} ( t ) its restriction to the inertial manifold [physics M-matrix] , that is, the semigroup generated by (90). Consider [figure omitted; refer to PDF]

We will find a reduced semigroup on the reduced space _{...B4; ...} : = P ( ...B4; ) , denoted as _{S R} ( t ) , that, in a sense, determines the asymptotic behavior of _{S M} ( t ) and therefore that of S ( t ) . Note that S ( t ) and _{S M} ( t ) have the same attractor, while the dimension of the space _{...B4; ...} might be much smaller than that of [physics M-matrix] . The next result states, in particular, that the attractor of the full system can be reconstructed from the attractor of the reduced one.

Proposition 14.

With the previous notation, one has the following conditions.

(i) The system of equations [figure omitted; refer to PDF]

: defines a nonlinear semigroup, denoted as _{S R} ( t ) , on _{...B4; ...} : = P ( ...B4; ) that can be identified with P _{S M} ( t ) P = P _{S M} ( t ) restricted to _{...B4; ...} .

(ii) : If ...9C; denotes the maximal attractor of (8), then _{...9C; R} = P ( ...9C; ) is the maximal attractor of (91). Moreover, [figure omitted; refer to PDF]

: where ...A2; : _{...9C; ...} ... ...9C; is continuous.

(iii) : If the set K ∩ J is finite, (91) is equivalent to a system of complex ODEs of the form (87). Consequently, the asymptotic behavior of (8) is described by an explicit system of ODEs in ^{... N} with N = | K ∩ J | + 2 an even number. In particular, if K ∩ J = ∅ , l ( v ) = _{l 0} , and G ( v ) = _{G 0} for every ( _{w 0} , _{v 0} , _{T 0} ) ∈ ^{... 2} × ^{H p e r 1} ( 0,1 ) , one has that the associated solution verifies v ( t ) [arrow right] 0 and T ( t ) [arrow right] _{θ ∞} in ^{H p e r 1} ( 0,1 ) , where _{θ ∞} ( x ) is the unique solution in ^{H p e r 2} ( 0,1 ) of the equation [figure omitted; refer to PDF]

Moreover, if ν = 0 , one gets v ( t ) [arrow right] 0 and T ( t ) [arrow right] _{T a} .

Proof.

(i) Working as aforementioned, we prove that the semigroups _{S M} ( t ) and _{S R} ( t ) are well defined and prove the existence of attractor _{...9C; R} . Using the techniques of [16] we can work as in [14] to prove (ii). Then (iii) we note that 0 ∉ K ∩ J , and since K = - K and J = - J , then the set K ∩ J is a symmetric set and has an even number of elements that we denote by 2 _{n 0} . Therefore, the number of the positive elements of K ∩ J , ( K ∩ J _{) +} , is _{n 0} .

Note that [contour integral operator] T · f = _{∑ k ∈ ^{... *}} ... _{a k} ( t ) _{c - k} = _{∑ k ∈ K ∩ J} ... _{a k} ( t ) · _{c - k} . Thus, the dynamics of the system depends only on the coefficients in K ∩ J . Moreover, the equations for _{a - k} are conjugates of the equations for _{a k} , and therefore we have that [figure omitted; refer to PDF] Taking real and imaginary parts of _{a k} , k ∈ ( K ∩ J _{) +} , that is, employing real variables, _{a k} = _{x k} + i _{y k} , we have a system in ^{... N} with N = 2 _{n 0} + 2 .

If K ∩ J = ∅ , then from the equation for the velocity [figure omitted; refer to PDF] we have [figure omitted; refer to PDF]

Moreover, from the equation for the temperature in (5), we have that the function θ = T - _{θ ∞} satisfies the equation [figure omitted; refer to PDF] We can multiply by θ in ^{L per 2} ( 0,1 ) , and taking into account that [contour integral operator] ( ∂ θ / ∂ x ) θ = ( 1 / 2 ) [contour integral operator] ( ∂ ( ^{θ 2} ) / ∂ x ) = 0 , since θ is periodic, we have [figure omitted; refer to PDF] and using Cauchy-Schwartz and Young inequality with δ , _{C δ} = 1 / 4 δ and then Poincaré inequality, since [contour integral operator] θ = 0 , together with - _{l 0} [contour integral operator] ^{θ 2} ...4; 0 , we have that [figure omitted; refer to PDF] Next, using v ( t ) [arrow right] 0 , we prove that θ ( t ) [arrow right] 0 in ^{L per 2} ( 0,1 ) .

Now, we multiply (97) by - ^{∂ 2} θ / ∂ ^{x 2} in ^{L per 2} . Integrating by parts, applying Young inequality and taking into account again that [contour integral operator] ( ∂ θ / ∂ x ) ( ^{∂ 2} θ / ∂ ^{x 2} ) = 0 since ∂ θ / ∂ x is periodic, together with - _{l 0} [contour integral operator] θ ( - ^{∂ 2} θ / ∂ ^{x 2} ) = - _{l 0} [contour integral operator] ( ∂ θ / ∂ x ^{) 2} ...4; 0 , we obtain that [figure omitted; refer to PDF] for every δ > 0 with _{C δ} = 1 / 4 δ . Thus, working as aforementiond and taking into account that | v ( t ) | [arrow right] 0 , we get ( ∂ θ / ∂ x ) ( t ) [arrow right] 0 in ^{L per 2} ( 0,1 ) ; that is, θ [arrow right] 0 ∈ ^{H per 1} ( 0,1 ) .

Next, we pay attention to the other modes for the temperature _{a k} where k ∉ ( K ∩ J ) . Also, note that these modes are determined as solution of the linear nonhomogeneous equations [figure omitted; refer to PDF] with initial data _{a k} ( 0 ) ∈ C . Therefore, we call these the slave modes.

We will show in Proposition 15 that the dynamics of these modes are completely determined by the solution of (87), in the sense that the solution will have only one asymptotic behavior as time goes to infinity.

Proposition 15.

Assume that { ( w ( t ) , v ( t ) , _{a k} ( t ) ) , k ∈ K ∩ J } is a solution of (87). Then, for any k ∉ ( K ∩ J ) , there exists a solution of (101), denoted as ^{a k *} ( t ) , such that | ^{a k *} ( t ) | ...4; | _{b k} | for every t ...5; 0 and for any other solution of (101), [figure omitted; refer to PDF] at an exponential rate independent of k , as t [arrow right] ∞ . Moreover, if k ∉ K , that is, if _{b k} = 0 , then ^{a k *} ( t ) = 0 ; that is, this subset of the slave mode is damped out exponentially. In particular, if { ( w ( t ) , v ( t ) , _{a k} ( t ) ) , k ∈ K ∩ J } is a stationary or periodic (resp., quasiperiodic, almost periodic) solution, then ^{a k *} ( t ) can be chosen such that it is stationary or periodic with the same period (resp., quasiperiodic, almost periodic with a set of frequencies contained in those of v ( t ) ) .

Proof.

Define | ^{a k *} ( t ) | as a solution of (101) with an initial condition satisfying | ^{a k *} ( t ) | ...4; | _{b k} | . In particular, if k ∉ K , that is, if _{b k} = 0 , then | ^{a k *} ( t ) | = 0 . Then for any other solution of (101), _{z k} = _{a k} ( t ) - ^{a k *} ( t ) satisfies the homogeneous equation [figure omitted; refer to PDF] and | _{z k} ( t ) | ...4; | _{z k} ( 0 ) | ^{e - ∫ 0 t ... [ 2 π k i v + 4 ν π 2 k 2 + l ( v ) ]} which proves the statement.

If the solution of (87) is stationary, that is, independent of time, then we choose ^{a k *} ( t ) to be the solution of ( 2 π k v i + 4 ν ^{π 2 k 2} + l ( v ) ) _{a k} ( t ) = l ( v ) _{b k} and the result follows.

If { ( w ( t ) , v ( t ) , _{a k} ( t ) ) , k ∈ K ∩ J } is a periodic solution of (87), then since the ^{a k *} ( t ) are the solutions of linear scalar differential equations of the form x ( t ) + A ( t ) x ( t ) = f ( t ) , with A ( t ) and f ( t ) period of the same period and the homogeneous equation is stable, the result follows from Fredholm's alternative. For the quasiperiodic or almost periodic case, the result follows from Theorem 6.6 in [23].

Remark 16.

Taking real and imaginary parts of _{a k} , _{b k} , and _{c k} as [figure omitted; refer to PDF] the asymptotic behavior of the system (5) is given by a reduced explicit system in ^{... N} , where N = 2 _{n 0} + 2 , given by [figure omitted; refer to PDF] where _{a - k} = _{a - k} , _{b - k} = _{b - k} , and _{c - k} = _{c - k} .

Observe that from the previous analysis, it is possible to design the geometry of circuit and/or the external heating by properly choosing the functions f and/or the heat flux l and the ambient temperature _{T a} so that the resulting system has an arbitrary number of equations of the form N = 2 n + 2 .

Note that the set K ∩ J can be much smaller than the set K , and therefore the reduced subsystem may possess far fewer degrees of freedom than the system on the inertial manifold. Also note that it may be the case that K and J are infinite sets, but their intersection is finite. Also, for a circular circuit, we have f ( x ) ~ a sin ( x ) + b cos ( x ) , that is, J = { ± 1 } and then K ∩ J is either { ± 1 } , or the empty set. Also, if in the original variable for (5) _{T a} is constant, we get K ∩ J = ∅ for any choice of f .

The physical and mathematical implications of the resulting system of ODEs which describe the dynamics at the inertial manifold need to be analyzed numerically. The role of the parameter [varepsilon] which contains the viscoelastic information of the fluid deserves special attention and will be the aim of the next section.

4. Numerical Experiments

In this section, we describe the results of the numerical experiments obtained using the MATHEMATICA package [24] for the resolution of the differential equations, using a fourth-order explicit Runge-Kutta method for stiffness equations following the method used in previous works [5, 15]. We will solve a system of ordinary differential equations which are the projection of the partial differential equations (5) on the inertial manifold derived in the preceding sections. All the variables and equations that we deal with are nondimensional. As the system is multidimensional, we present the results in temporal graphs (a given variable versus time) and phase-space graphs (two physical variables plot against each other).

Specifically, we are integrating the system of (87) where we consider only the coefficients of temperature _{a k} ( t ) with k ∈ K ∩ J (relevant modes). Then, [figure omitted; refer to PDF] where _{a - k} = _{a - k} , _{b - k} = _{b - k} , and _{c - k} = _{c - k} since all the physical observable are real functions. In particular, we will consider a thermosyphon with a circular geometry; so J = { ± 1 } and K ∩ J = { ± 1 } . Consequently, we take k = 1 and omit the equation for k = - 1 . Hence, [figure omitted; refer to PDF] where the unknowns are w ( t ) (the acceleration of the fluid), v ( t ) (velocity of the fluid), and _{a 1} ( t ) (the Fourier mode of the temperature). More complex geometries will result in higher dimensional dynamics on the inertial manifold. It is beyond the scope of the present work.

In order to reduce the number of parameters we make the change of variables _{a 1 c - 1} [arrow right] _{a 1} and we define the real and imaginary parts of the equations in the following way: [figure omitted; refer to PDF] with A ∈ ... , B ∈ ... . Therefore, our central results correspond to the system of equations [figure omitted; refer to PDF] Note that it is a system of four equations with four unknowns where we need to make explicit choices for the constitutive laws for both the fluid-mechanical and thermal properties. Thus, for the friction law G ( v ) and heat flux l ( v ) we will take the ones used in [5, 15]. For the numerical experiments, which are of a similar model of thermosyphon for a fluid with one component they use the functions G ( v ) = ( | v | + 1 ^{0 - 4} ) and l ( v ) = ( 1 ^{0 - 2} | v | + 1 ) . The function G ( v ) has a clear physical meaning; it interpolates between a low Reynolds number friction law (in which the overall friction G ( v ) v is linear (Stokes friction law) and high Reynolds number (in which the friction is a quadratic law).

A and B refer in this model to the position-dependant ( x ) heat flux inside the loop and will be used as tuning parameters. Without loss of generality, we will assume that A = 0 in order to simplify in analogy with Lorenz's model, as it is shown in [5, 15] (changing A and B simultaneously only results in a change in the phase of initial temperature profile).

We have carried out two different sets of numerical experiments with regard to heat diffusion. The first numerical experiments are carried out keeping the heat diffusion to zero as it was done in [11]. And the second numerical experiments are performed with heat diffusion. The initial conditions are fixed as w ( 0 ) = 0 , v ( 0 ) = 0 , _{a 1} ( 0 ) = 1 , _{a 2} ( 0 ) = 1 . This split would appear naive as diffusion tends to smooth the solution however, as the order of the equations changes in the presence of diffusion (from first to second order due to the Laplacian), it is worth studying both cases separately.

Numerical analysis has been carried out keeping [varepsilon] the viscoelastic coefficient as the tuning parameter ranging from 100 to 0.0001 and B the heat flux also as a tuning parameter ranging from 1 to 10000. The impact of [varepsilon] on the system has been keenly observed for various parameters of time t , as short as 50 time units and as long as 5000 time units. We will show that in analogy with the classical Lorenz system, as [varepsilon] varies, the dynamics of the model undergoes various transformations including steady asymptotic behavior, meta-stable chaos, that is, transient irregular behavior followed by convergence to equilibria, periodic behaviors, and chaotic progressions.

4.1. Dynamics of the Thermosyphon without Heat Diffusion ( ν = 0 )

In this section, we summarize some of the outcomes of the model equations. As we have mentioned earlier, this behavior is highly sensitive to the choice of parameters. Thus, we present those results in different subsections accounting for the most relevant signature for each set of numerical experiments.

4.1.1. Chaotic Behavior of the Model for Large Values of [varepsilon]

The simulations of the numerical experiments done for large values of [varepsilon] , for instance, [varepsilon] ranging from 2 to 1000, show that the system exhibits a pattern of chaotic behavior. For relatively large values of viscoelastic components, chaotic behaviors are observed. For all the values of heat flux B , starting from 1 to 10000, this chaotic behavior is observed (see Table 2). To illustrate the previously stated chaotic behavior, we observe the plots obtained for the values of [varepsilon] = 10 and B = 50 .

In Figure 1, we show a time graph of acceleration for a large value of the viscoelastic parameter [varepsilon] . The acceleration ranges from -3.5 to 3.5. Since the very beginning, it displays a chaotic behavior. The curve is not very erratic although it does not show any sort of periodicity. As velocity is the time integral of acceleration, in Figure 2, the curve does not present abrupt changes close to some maxima and minima, but the nonperiodic features are also captured by this observable.

Figure 1: The chaotic progress of acceleration for [varepsilon] = 10 , B = 50 , ν = 0 .

[figure omitted; refer to PDF]

Figure 2: The inconsistent behavior of velocity for [varepsilon] = 10 , B = 50 , ν = 0 .

[figure omitted; refer to PDF]

In Figure 3, we show a phase-diagram plot for the real and imaginary parts of the temperature. As expected, it also exhibits a nonperiodic pattern in which the trajectory in this phase-plane moves inwards and outwards the graph. This graph illustrates the complex underlying dynamics of the attractor (of which Figure 3 is a two-dimensional projection).

Figure 3: A chaotic global attractor of real and complex temperature for [varepsilon] = 10 , B = 50 , ν = 0 .

[figure omitted; refer to PDF]

This sort of behavior remains similar for other values of B as long as the viscoelastic parameter takes larger than unity values as summarized in Table 2. In other words, the elastic effects introduce a memory effect in the dynamics which avoids the system to fully stabilize but, rather, viscoelasticity sustains the chaotic pattern. This memory effect can be understood from (3).

To sum up this section, large values of the viscoelastic parameter [varepsilon] result on sustained chaotic behaviors. The dynamics becomes more complex, characterized in all the cases by periods of chaos and of violent oscillations, giving an idea of the complexity of the solutions of the system under these variables. In the detailed analysis of the evolution of the acceleration, velocity, and temperature, we say that the chaotic behavior of the system reveals the chaotic nature of the viscoelastic fluids.

4.1.2. Transient Irregular Behavior Followed by Stable Behavior for [varepsilon] = 0.1 to 1

For values of [varepsilon] ranging from 0.1 to 1 , the system tends towards a stable fixed point, although it is still chaotic in the initial stages. This transient irregular behavior followed by equilibrium is shown in Figures 4, 5, and 6.

Figure 4: The stabilizing progress of acceleration for [varepsilon] = 0.1 , B = 100 , ν = 0 .

[figure omitted; refer to PDF]

Figure 5: Velocity stabilizes at 3.19981 for [varepsilon] = 0.1 , B = 100 , ν = 0 .

[figure omitted; refer to PDF]

Figure 6: Equilibrium velocity scale.

[figure omitted; refer to PDF]

The general behavior of acceleration is that it has a chaotic outburst in the initial stages. But as time progresses it tends to stabilize, attaining equilibria. The velocity too, in the initial stages, when the time period is less than 20 units, is very inconsistent and at times unpredictable. But as the time progresses, velocity converges to a stable fixed point. Interestingly, this fixed point for the velocity is not trivial ( v ...0; 0 ) but, on the contrary, it depends strongly on the choice of the parameters. Specifically, in Table 1, we summarize this asymptotic value. Physically, this means that the fluid inside the thermosyphon performs sustained flow in the same direction over time. This stage could be identified by the existence of convective rolls in a fully spatially extended system.

Table 1: Equilibrium values of velocity for different values of heat flux B .

B	[straight epsilon] = 1	[straight epsilon] = 0.1
1	0.6718	0.6718
10	1.4723	1.4723
20	1.8601	1.8601
30	2.1325	2.1325
40	2.3495	2.3495
50	2.5329	2.5329
100	3.1998	3.1998
1000	6.9800	6.9800
10000	15.430	15.430

Table 2: Behavior of solutions without diffusion ( ν = 0) for different values of the viscoelastic characteristic time, [straight epsilon] (rows), and the temperature gradient, B (columns). We introduce the following notation to account for the obtained numerical results: "C" denotes a fully chaotic behavior, "CS" a transition from chaotic outburst to stable equilibria, "P" a stable periodic orbit, and "CP" a transitional behavior from chaotic to periodic.

B / [straight epsilon]	10 - 4	10 - 3	10 - 2	10 - 1	10 0	10 1	10 2
1	CS	CS	CS	CS	CS	C	C
10	CS	CS	CS	CS	CS	C	C
20	CS	CS	CS	CS	CS	C	C
30	CP	CP	CP	CS	CS	C	C
40	CP	CP	CP	CS	CS	C	C
50	CP	CP	CP	CS	CS	C	C
100	P	P	CP	CS	CS	C	C
1000	CP	CP	CP	CS	CS	C	C
10000	CP	CP	CP	CS	CS	C	C

It is worth noting in Table 1 that both columns have the same value for the asymptotic velocity. This is a signature that viscoelastic effects do not play any role in this case (so memory effects are damped out after the early chaotic transient).

Another comment regarding Table 1 concerns the role of the parameter B . Roughly, B accounts for the scale of temperature gradients inside the thermosyphon. These results suggest that higher temperature gradients produce higher values of the sustained stationary velocity. Actually, as shown in Figure 6, the equilibrium velocity scales with B nontrivially (as a power law with exponent 1 / 3 approximately).

To sum up this section, we conclude that although the system has a chaotic initial transient, it tends to stabilize at longer times reaching a temperature gradient dependent equilibrium velocity. Notwithstanding, this asymptotic velocity depends nontrivially on temperature as a power law, being this a signature of the underlying nonlinearity of the equations.

4.1.3. Transition to Periodic Pattern of Behaviors for Small Values of [varepsilon]

When the viscoelastic effects are gradually less important (values of [varepsilon] between 0.01 and 0.0001), the system exhibits different behaviors as a function of the temperature gradient (see Table 2). For lower values of B , the system behaves in a similar fashion as for higher values of [varepsilon] . On the contrary, for larger values of B , the system displays a periodic pattern. As in the previous case, both long-term behaviors may be preceded by an initial chaotic transient (probably caused by viscoelasticity).

To illustrate this, in Figure 7, we show a periodic (nontrivial) behavior for [varepsilon] = 0.001 . In some cases, the initial transient cannot be distinguished from the periodic one and we refer it simply to periodic type in Table 2.

Figure 7: The periodic progress of velocity for [varepsilon] = 0.001 , B = 30 , ν = 0 .

[figure omitted; refer to PDF]

In Figure 8, we show the phase-space of the complex components of the temperature. Although at first sight it resembles a typical chaotic attractor motif, after the mentioned initial transient, the trajectories are overlapped for longer times.

Figure 8: The chaotic but periodic plot of real and complex temperature for [varepsilon] = 0.0001 , B = 40 , ν = 0 .

[figure omitted; refer to PDF]

Thus, in order to summarize the information covered in the last three subsections, we collect all the outcomes of the model in Table 2.

4.2. Dynamics of the Thermosyphon with Heat Diffusion ( ν ...0; 0 )

In this case, to avoid unnecessary repetitions in the text, we focus on the main differences between this case and that in Section 4.1.

Thus, in this second set of numerical experiments we introduce a nonzero value for the thermal diffusivity, ν . The role of diffusion is to reduce temperature gradients. This can be seen in a hand-waving way by realizing that [figure omitted; refer to PDF] namely, the Laplazian can be understood as the flux of temperature created by a temperature current , J = - ∇ T . This current is larger in those regions where the temperature variations are also larger. Thus, the system tends to reduce those differences. As shown in Table 3, the variety of the behaviors is clearly less rich than in the case when ν = 0 .

Table 3: Behavior of solutions with diffusion ( ν ...0; 0 ) for different values of the viscoelastic characteristic time, [straight epsilon] (rows), and the temperature gradient, B (columns). We introduce the following notation to account for the obtained numerical results: "C" denotes a fully chaotic behavior, "CS" a transition from chaotic outburst to stable equilibria, and "P" a periodic orbit.

B / [straight epsilon]	10 - 4	10 - 3	10 - 2	10 - 1	10 0	10 1	10 2
1	CS	CS	CS	CS	CS	CS	CS
10	CS	CS	CS	CS	CS	CS	CS
20	CS	CS	CS	CS	CS	CS	CS
30	CS	CS	CS	CS	CS	CS	CS
40	P	P	P	CS	CS	CS	CS
50	P	P	P	CS	CS	CS	CS
100	P	P	P	CS	CS	CS	CS
1000	P	P	P	CS	CS	CS	CS
10000	P	P	P	CS	CS	CS	CS

So, here we will only illustrate the most interesting behaviors with two examples. The first one takes the values of heat diffusion ν = 1 , [varepsilon] = 5 , and B = 1000 . As shown in Figure 9, the acceleration performs a series of damped oscillations that eventually stabilizes.

Figure 9: The stabilizing process of acceleration for ν = 1 , [varepsilon] = 5 , B = 1000 , ν ...0; 0 .

[figure omitted; refer to PDF]

Similarly, the second example (Figure 10) takes the values ν = 2 , [varepsilon] = 5 , and B = 1000 . The behavior is qualitatively equal but the period of the oscillations is enlarged. In Table 3, we summarize the interaction between the tuning parameters.

Figure 10: The fast stabilization of acceleration for ν = 2 , [varepsilon] = 5 , B = 1000 , ν ...0; 0 .

[figure omitted; refer to PDF]

To sum up this section, we have found that greater values of the heat diffusion, ν , smoothen the dynamics of the system which, invariably, tends to stabilize, either reaching an equilibrium steady state or stable periodic orbits. These two behaviors are governed by the value of the temperature gradient which in a similar fashion as in Section 4.1 tends to produce richer behaviors for large values.

5. Conclusion

In this work, we have derived a novel system of equations to study the behavior of a viscoelastic material inside a thermosyphon. This model serves as a preliminary simplification of a more complex fully spatially extend system. This model illustrates the presence/absence of complex chaotic behaviors and also enables to relate them with the underlying viscoelastic (memory effects).

The main result is that we are able to prove that the original system (which involves both ordinary and partial different equations) possesses an inertial manifold in which the dynamics can be accurately described by a system of ODEs. By numerical integration of the reduced equations, we have been able to better understand the role of viscoelasticity (as opposed to a simpler Newtonian fluid) through the parameter [varepsilon] . This parameter is a nondimensional version of the so-called Maxwellian viscoelastic time [4] which accounts for memory effects.

Our results suggest that in the absence of heat diffusion ( ν = 0 ), a large value of [varepsilon] drives the dynamics to chaotic behaviors for all the physical observable (acceleration, velocity, and temperature). As the value of [varepsilon] gradually decreases, the system is no longer chaotic but stable or periodic. Notably, these results cannot be understood in terms of a boundary layer theory as the attractor of the dynamics changes dramatically when the second-order derivative term [varepsilon] ^{d 2} v / d ^{t 2} is introduced in (5).

Physically, this induction of chaotic behaviors can be rationalized in terms of the memory effects inherent to viscoelastic models explicitly shown in (3). Thus, in the same way as delayed equations are known to produce chaos, even in the simplest situations, viscoelasticity apparently produces the same kind of transition (see, e.g., [25]).

Other interesting results are related to the effect of heat diffusion. We have found that as the heat diffusion is not zero, the system tends to stabilize either to a fixed equilibrium point or to a ( ν -dependent periodicity) periodic orbit. This abrupt change in the qualitative behavior of the system is, again, due to the fact that the term ν ^{d 2} T / d ^{x 2} is a singular perturbation that changes the type of equation for the temperature from hyperbolic to parabolic, which justifies the explicit analysis of each situation separately in Section 2.

Our work can be generalized in many different ways, from changing the constitutive equation (from Maxwellian to other more complex situations) or to include shear-thinning effects [4] common to many non-Newtonian materials. Shear thinning is the manifestation of a shear-rate-dependent viscosity. Thus, it is commonly observed that many fluids reduce their resistance to flow for large enough imposed stresses (in our case, temperature gradients), for instance, tooth paste, paint, or lava. In principle, shear thinning means that the viscosity is smaller for higher velocities what could help to sustain oscillations once they have been formed inside the loop. However, as shown in the present work, this sort of hand waving analysis might lead to erroneous conclusions, and a detailed consideration of those effects deserves further (and deeper) investigation. These generalizations will be the aim of future works.

Acknowledgments

This work is partially supported by projects MTM2009-07540, GR58/08 Grupo 920894 BSCH-UCM, Grupo de Investigación CADEDIF, and FIS2009-12964-C05-03, Spain, and partially supported by Grant MTM2012-31298 from Ministerio de Economia y Competitividad, Spain.

References

[1] M. C. Cross, H. Greenside Pattern Formation and Dynamics in Nonequilibrium Systems , Cambridge University Press, New York, NY, USA, 2009.

[2] J. B. Keller, "Periodic oscillations in a model of thermal convection," Journal of Fluid Mechanics , vol. 26, no. 3, pp. 599-606, 1966.

[3] P. Welander, "On the oscillatory instability of a differentially heated fluid loop," Journal of Fluid Mechanics , vol. 29, no. 1, pp. 17-30, 1967.

[4] F. Morrison Understanding Rheology , Oxford University Press, New York, NY, USA, 2001.

[5] R. Greif, Y. Zvirin, A. Mertol, "The transient and stability behavior of a natural convection loop," Journal of Heat Transfer , vol. 107, pp. 684-688, 1987.

[6] A. Jiménez Casas, A. Rodríguez-Bernal Dinámica no Lineal: Modelos de Campo de Fase y un Termosifon Cerrado , Editorial Acadmica Espaola, Lap Lambert Academic Publishing, 2012.

[7] A. Jiménez-Casas, A. Rodríguez-Bernal, "Finite-dimensional asymptotic behavior in a thermosyphon including the Soret effect," Mathematical Methods in the Applied Sciences , vol. 22, no. 2, pp. 117-137, 1999.

[8] A. Jiménez-Casas, "A coupled ODE/PDE system governing a thermosyphon model," Nonlinear Analysis , vol. 47, pp. 687-692, 2001.

[9] J. J. L. Velázquez, "On the dynamics of a closed thermosyphon," SIAM Journal on Applied Mathematics , vol. 54, no. 6, pp. 1561-1593, 1994.

[10] A. Jiménez Casas, A. M. L. Ovejero, "Numerical analysis of a closed-loop thermosyphon including the Soret effect," Applied Mathematics and Computation , vol. 124, no. 3, pp. 289-318, 2001.

[11] A. Rodríguez-Bernal, E. S. van Vleck, "Complex oscillations in a closed thermosyphon," International Journal of Bifurcation and Chaos in Applied Sciences and Engineering , vol. 8, no. 1, pp. 41-56, 1998.

[12] M. A. Herrero, J. J. L. Velázquez, "Stability analysis of a closed thermosyphon," European Journal of Applied Mathematics , vol. 1, no. 1, pp. 1-23, 1990.

[13] A. Linan, M. G. Velarde, C. I. Christov, "Analytical description of chaotic oscillations in a toroidal thermosyphon," Fluid Physics , of Lecture Notes of Summer Schools, pp. 507-523, World Scientific, River Edge, NJ, USA, 1994.

[14] A. Rodríguez-Bernal, E. S. Van Vleck, "Diffusion induced chaos in a closed loop thermosyphon," SIAM Journal on Applied Mathematics , vol. 58, no. 4, pp. 1072-1093, 1998.

[15] D. Henry Geometric Theory of Semilinear Parabolic Equations , vol. 840, of Lecture Notes in Mathematics, Springer, Berlin, Germany, 1981.

[16] A. Rodríguez-Bernal, "Attractors and inertial manifolds for the dynamics of a closed thermosyphon," Journal of Mathematical Analysis and Applications , vol. 193, no. 3, pp. 942-965, 1995.

[17] J. K. Hale Asymptotic Behavior of Dissipative Systems , AMS, Providence, RI, USA, 1988.

[18] F. P. Incropera, T. L. Bergman, A. S. Lavine, D. P. DeWitt Fundamentals of Heat and Mass Transfer , John Wiley & Sons, 2011.

[19] A. M. Bloch, E. S. Titi, "On the dynamics of rotating elastic beams," New Trends in Systems Theory , vol. 7, of Progress in Systems and Control Theory, pp. 128-135, Birkhäauser, Boston, Mass, USA, 1991.

[20] A. M. Stuart, M. Ainsworth, J. Levesley, W. A. Light, M. Marletta, "Pertubration theory of infinite-dimensional dynamical systems," Theory and Numerics of OrdInary and Partial Differential Equations , Oxford University Press, Oxford, UK, 1994.

[21] C. Foias, G. R. Sell, R. Temam, "Inertial manifolds for nonlinear evolutionary equations," Journal of Differential Equations , vol. 73, no. 2, pp. 309-353, 1988.

[22] A. Rodríguez-Bernal, "Inertial manifolds for dissipative semiflows in Banach spaces," Applicable Analysis , vol. 37, no. 1-4, pp. 95-141, 1990.

[23] A. M. Fink Almost Periodic Differential Equations , vol. 377, of Lecture Notes in Mathematics, Springer, Berlin, Germany, 1974.

[24] S. Wolfram The Mathematica book , Cambridge University Press, 1999.

[25] Y. Kuang Delay Differential Equations with Applications in Population Dynamics , Academic Press, 1993.

[26] F. Verhulst Methods and Applications of Singular Perturbations , Springer, New York, NY, USA, 2005.

Appendix

Appendix

Boundary Layer Theory

The system of (5) can be seen as a singular perturbation problem provided that when [varepsilon] = 0 , the order of the differential equation reduces to one. In order to provide some insight about the solutions of those equations, a standard boundary layer theory is customary (see, e.g., [26]).

In particular, one splits the problem into two parts: the inner problem and the outer problem. The inner problem is defined as the dynamics of the system for times up to 1 / [varepsilon] in which the term [varepsilon] ^{d 2} v / d ^{t 2} is dominant. This regime is strongly dominated by transient impulsive changes in the physical variables. Following [26], we define a new time scale τ = t / [varepsilon] . So, the first equation in (5), up to order [varepsilon] , is given by [figure omitted; refer to PDF] whose solution is [figure omitted; refer to PDF] with β a constant that can be determined by matching with the outer solution. For instance, if one assumes, the velocity has initial condition v ( 0 ) and the matching time is τ = O ( 1 ) ≈ 1 , then it is straightforward to see that the velocity at a matching time τ = 1 converges exponentially to _{v matching} = v ( 0 ) + β ( ^{e - 1} - 1 ) .

Besides, the outer problem is defined as the naive approximation [varepsilon] [arrow right] 0 . Thus, the system of (5) reduces to that in [11] with an effective initial condition given by _{v matching} . So one would expect the same qualitative behaviors as in that work.

For instance, if the parameter B = 10 , the model in [11] predicts a stable behavior. However, as shown in Table 2 (second row), increasing the value of [varepsilon] induces a chaotic behavior. Naturally, different values of the parameters (including B ) would provide different values of [varepsilon] for which the trivial case [varepsilon] [arrow right] 0 cannot account for the observed dynamics (even for small enough values of [varepsilon] ).

This simple analysis shows the intrinsic complexity of the physical problem when viscoelasticity is considered and, more importantly, the need to study every parameter set in detail to provide an accurate description of the type of dynamics in which the system evolves.