Viscoelastic rate type fluids with temperature dependent material parameters -- stability of the rest state
Judith Stein, V\'it Pr\r{u}\v{s}a

TL;DR
This paper investigates the stability of the rest state in viscoelastic fluids with temperature-dependent properties, emphasizing the role of thermodynamic analysis in model development.
Contribution
It demonstrates conditions under which the rest state remains stable, highlighting the importance of thermodynamic considerations in complex fluid models.
Findings
Stable rest state achieved with appropriate material parameters
Thermodynamic analysis is crucial for model stability
Provides conditions for stability in temperature-dependent viscoelastic fluids
Abstract
We study the dynamics of small perturbations to the rest state of a viscoelastic rate type fluid with temperature dependent material parameters. We show that if the material parameters are chosen appropriately, then the quiescent state of the fluid filling an isolated (mechanically, thermally) vessel is a stable state. The outlined analysis explicitly documents the importance of thermodynamic analysis in the development of advanced models for complex fluids.
Peer Reviews
No public reviews on file for this paper yet. If you reviewed it on a platform where reviews are public (OpenReview, ICLR, NeurIPS, ICML), you can paste yours below so the community can read it here.
Videos
No videos yet. Explain this paper in a talk, walkthrough, or lecture? Add one.
Viscoelastic Rate Type Fluids with Temperature Dependent Material Parameters – Stability of the Rest State
Judith Stein1, a)
Vít Průša2, b)
Abstract
We study the dynamics of small perturbations to the rest state of a viscoelastic rate type fluid with temperature dependent material parameters. We show that if the material parameters are chosen appropriately, then the quiescent state of the fluid filling an isolated (mechanically, thermally) vessel is a stable state. The outlined analysis explicitly documents the importance of thermodynamic analysis in the development of advanced models for complex fluids.
Keywords:
viscoelastic fluids, temperature, linearised stability
:
83.10.Gr, 47.50.Cd, 65.20.De
1 Introduction
The non-isothermal processes in viscoelastic fluids are of interest in engineering applications. A prominent example thereof is the injection molding process. Consequently, the question “how to account for the temperature changes” is perceived to be an important issue in rheology, see Tanner (2009). Regarding the viscoelastic fluids, the answer to the question basically requires one to formulate an evolution equation for the temperature field.
Various thermodynamical approaches/frameworks have been developed in order to address this issue, see in particular Leonov (1976, 1987), Wapperom and Hulsen (1998), Mattos (1998), Dressler et al. (1999) or Guaily (2015). (The treatise by Dressler et al. (1999) contains a rich bibliography on subject matter.) Recently, Hron et al. (2016) have provided a derivation of variants of incompressible Maxwell and Oldroyd-B rate type models with temperature dependent material parameters. In what follows, we elucidate the features of the model by investigating the problem of the stability of the rest state.
This is an important issue, since the requirement on the stability of the rest state gives one additional restrictions on the energetic/entropic equation of state. Such restriction are well-known in the case of compressible Navier–Stokes fluid, see for example Callen (1985) or Müller (1985), and it is worthwhile to explicitly see how a similar concept works in the case of complex viscoelastic fluids.
2 Model
The model that is used in the subsequent considerations on the stability of the rate state is the model presented by Hron et al. (2016). The derivation of incompressible Maxwell and Oldroyd-B rate type models with temperature dependent material parameters as presented by Hron et al. (2016) is a very simple one, and it is based on purely macroscopic/phenomenological arguments.
In particular, there is no need to build the extensive apparatus of GENERIC framework, see Grmela (2010), Grmela and Öttinger (1997) and Öttinger and Grmela (1997), as in the approach by Dressler et al. (1999). Further, there is no need to deal with the concept of conformation tensor, see Grmela and Carreau (1987), and the related microstructural interpretation. The extra tensorial variable is interpreted, following Leonov (1976, 1987) and Rajagopal and Srinivasa (2000) as a strain measure characterising the response of the material to the instantaneous unloading. This allows one, amongst others, to directly interpret the outcomes of the creep and stress relaxation experiments, see Průša et al. (2017) and Řehoř et al. (2016).
Note also that the concept of polymer “chains” that is an underlying concept leading to the conformation tensor approach is hardly applicable to complex structures. For example aphaltenes found in crude oil products such as asphalt binders are far more complicated than “chains”, yet viscoelastic rate type models are still applicable in the description of such materials, see Narayan et al. (2015), Málek et al. (2015a); Málek et al. (2016), Narayan et al. (2016) or Nivedya et al. (2016). Consequently, a purely phenomenological approach that does not refer, even implicitly, to the underlying microstructure of the material is most desirable.
Appealing to nearly the same concepts as Hron et al. (2016), Rao and Rajagopal (2002) have earlier developed models for non-isothermal flows of incompressible viscoelastic fluids in their effort to model complex phenomena such as crystallisation in polymers. (See also the follow-up works by Kannan et al. (2002), Kannan and Rajagopal (2004, 2005) and Kannan et al. (2006).) However, these works are restricted to a class of materials where the entropic equation of state has a special additive structure.
The assumed additive structure effectively lead to the internal energy that is a function of only the temperature. Consequently, the temperature evolution equation has a simple structure and the shear modulus is allowed to be only a linear function of temperature. Further, the governing equations for the mechanical quantities do not match the standard governing equations used in Maxwell/Oldroyd-B type models. These are substantial limitations and they have been relaxed by Hron et al. (2016).
Further, Hron et al. (2016) have also studied dynamical consequences of the derived constitutive relations in simple initial/boundary value problems. One of the interesting outcomes is the observation that the model is capable of describing the behaviour known from solid materials. Namely the phenomena that the tension can induce either cooling or heating of the material, see Joule (1859) and the related modern treatises on solid materials. The reader interested in details is kindly referred to Hron et al. (2016) for the in-depth discussion thereof.
2.1 Free energy and entropy production – general model
The derivation of the model presented in Hron et al. (2016) basically follows the approach by Rajagopal and Srinivasa (2000) and Málek et al. (2015b), but the whole procedure is extended to the non-isothermal setting. (See also Málek and Průša (2016) for an outline of the current procedure.) In principle, the governing equations are seen as consequences of the specification of material energy storage ability and entropy production ability. The specification of energy storage ability is tantamount to the specification of the internal energy or another equivalent themodynamical potential. In particular, Hron et al. (2016) have used the specific Helmholtz free energy in the form111The Helmholtz free energy is defined as a specific quantity. In other words it is the Helmholtz free energy per unit mass, and its physical dimension is .
[TABLE]
where denotes the temperature and is a quantity characterising the elastic response of the fluid, see Hron et al. (2016) or Málek et al. (2015b) for details. Symbol denotes the purely thermal contribution to the free energy and is a material parameter that can depend on temperature.
The entropy production , that is the right hand side in the evolution equation for the entropy , is assumed to take the form
[TABLE]
where , and are material parameters that depend on temperature, and, if necessary, also on the other quantities such as the pressure , symmetric part of the velocity gradient and so forth. Symbol denotes the traceless part (deviatoric part) of the corresponding tensor.
Note that if the functions , and are positive, then the entropy production is positive, and the second law of thermodynamics automatically holds for the corresponding system of governing equations.
2.2 Governing equations
Helmholtz free energy and entropy production in the form (1) and (2) accompanied by an insight into the kinematics of the tensorial quantity lead, see Hron et al. (2016) for details, to the following system of governing equations.
[TABLE]
Symbol denotes the upper convected derivative, see Oldroyd (1950),
[TABLE]
where denotes the material time derivative. Further, symbol denotes the scalar product on the space of matrices, and this scalar product also defines the norm of a matrix, .
System of governing equations (3) is a system of equations for unknown quantities , and it provides a generalisation of Maxwell/Oldroyd-B model to the non-isothermal setting. Once one provides boundary conditions for the velocity and temperature field, and the initial conditions for the velocity, temperature and field, and fixes a reference pressure value, system (3) provides a complete description of the dynamics of the quantities of interest.
Note that is in the approach byRajagopal and Srinivasa (2000) interpreted as the left Cauchy–Green tensor associated to the instantaneous elastic response, hence it is a symmetric and positive definite matrix. These properties hold also in the viscoelastic rate type models derived by appealing to the concept of conformation tensor.
2.3 Free energy and entropy production – fluid with constant specific heat capacity
Note that if the specific heat capacity at constant volume is required to be a constant, then the undetermined function in the Helmholtz free energy ansazt (1) must take the form
[TABLE]
where is the desired constant value of the specific heat capacity at constant volume , and is a reference temperature. The complete formula for the Helmholtz free energy in this case reads
[TABLE]
where is a function of temperature. Note that once the free energy as a function of and is known, one can obtain explicit formulae for the internal energy and the entropy ,
[TABLE]
as functions of the temperature and the left Cauchy–Green field . We shall occasionally use these specific explicit formulae in the ongoing analysis.
3 Stability of the rest state
3.1 General remarks
The experience shows that the isolated systems, that is the systems that do not interact with its surrounding, inherently tend to an equilibrium with homogeneous distribution of the quantities characterising the system. This observation is used when one wants to obtain further restrictions on the form of entropic/energetic equation of state.
In particular, if one investigates the compressible Navier–Stokes fluid, then the second law of thermodynamics holds provided that the viscosities and in the formula for the Cauchy stress tensor
[TABLE]
are chosen accordingly, that is and . However, even if the entropy growth is granted, additional requirements on the constitutive relations must be fulfilled provided that the preferred terminal state in the isolated system is the equilibrium with homogeneous temperature and density distribution and zero velocity. The classical construction, see for example Callen (1985), Müller (1985) or Glansdorff and Prigogine (1971), reveals that the inequalities
[TABLE]
are necessary in order to enforce the tendency to recover the homogeneous equilibrium state after the introduction of a small perturbation. Since the specific heat at constant volume and the thermodynamic pressure are obtained by manipulating the Helmholtz free energy, inequalities (9) are in fact restrictions on the form of Helmholtz free energy.
In the present case of incompressible Maxwell/Oldroyd-B fluid, the rest state is represented by the quadruple , where is the reference pressure and is the reference temperature. Indeed, is a solution to (3) with the boundary conditions and , where is the domain of interest and is the unit outward normal. These conditions mean that the system is mechanically and thermally isolated.
The second law of thermodynamics is in the present case satisfied provided that the functions , and in the entropy production ansatz (2) are positive. The question on the sign of function is however undecided by the requirement on the positivity of the entropy production. Let us now investigate whether some restrictions can be obtained by the requirement on the stability of the rest state of the isolated system.
3.2 Stability of the rest state via entropy maximisation at fixed energy
Let us first present a heuristic argument that does not deal with spatial distribution of the quantities of interest. We will focus, for the sake of simplicity, to the model with constant specific heat capacity . This means that we will deal with the Helmholtz free energy in the form (6).
We shall exploit the famous formulation of the first and second law of thermodynamics by Clausius (1865), namely the statement: “The energy of the world is constant. The entropy of the world strives to a maximum.” In the present case, the world is the isolated system kept at the fixed energy denoted as ,
[TABLE]
Apparently, the given value of the internal energy can be attained by different combinations of the quantities and , the rest state , is only one of them. If we want the rest state to be the terminal state in the isolated system, we need to show that the entropy attains its maximum at the rest state. Fortunately, we have an explicit formulae both for the internal energy and the entropy, see (7), hence we can explicitly verify whether this is true or false. Or, more precisely, we can ask what are the conditions that guarantee that the maximal entropy is attained at the rest state.
If the internal energy is kept constant , then the temperature is an implicit function of and . Formula reads
[TABLE]
We can use the implicit function theorem and differentiate (11) with respect to keeping in mind that , where plays the role of a parameter. The differentiation of (11) in the direction yields
[TABLE]
Evaluating explicitly the derivative with respect to in the last term, see for example Šilhavý (1997) for the corresponding formulae and notation, leads to
[TABLE]
Consequently, we see that
[TABLE]
In particular, we can observe that the derivative of the termperature with respect to at the rest state reads
[TABLE]
Having determined the derivative of the temperature with respect to the field , we can proceed with the investigation of the entropy. The entropy is given by (7a), where is in virtue of (11) a known function of and , hence we can express the entropy as a function of only,
[TABLE]
This function should attain maximum at the rest state .
The first derivative of the entropy with respect to in the direction can be found using the chain rule,
[TABLE]
Since the partial derivative of the entropy as a function of and with respect to reads
[TABLE]
we see that . Consequently, evaluating the right hand side of (17) yields, in virtue of (15), the equality
[TABLE]
This means that the entropy attains extreme value at the rest state . No restriction on the form of material function has arisen so far. It remains to check that the extremum is the maximum, which requires one to find the second derivative of the entropy as a function of .
The chain rule implies that the second derivative of the entropy in the directions and reads222Later, we will need only the derivative , but using different directions in the calculation helps to keep the calculation more transparent.
[TABLE]
We are interested in the sign of the second derivative at , which means that many terms in (20) in fact vanish in virtue of (15). Indeed
[TABLE]
The formulae for the derivatives of the entropy in (21) read
[TABLE]
Fortunately, we need only the values of the derivatives at the rest state which in virtue of (15) yields
[TABLE]
Consequently, we see that (23) and (21) imply
[TABLE]
We can therefore conclude that if , then the second derivative of the entropy is negative for arbitrary . In such a case the extremum is indeed a maximum, and the growth of entropy drives the system to the rest state as desired.
Note that if the term in the Helmholtz free energy ansatz were replaced by a function that vanishes at and that for all satisfies
[TABLE]
then the outcome of the presented analysis would be the same.
3.3 Stability of the rest state via direct application of the governing equations
Let us now solve a more demanding problem. We are given an isolated system, that is a vessel filled with the fluid of interest. There is no heat flux through the vessel wall, and the fluid sticks to the vessel wall, that is we have the no-slip and no-penetration boundary condition for the velocity. The rest state , where is the reference pressure and is the reference temperature is a solution to the governing equations.
The question is what happens if we slightly perturb the rest state. Again, we would like the governing equations to predict the decay of perturbations and the recovery of the rest state in the long run. This behaviour is granted provided the material coefficients/functions are chosen appropriately.
Since we are interested in small perturbations, we will investigate the dynamics dictated by the linearised governing equations in an isolated fixed domain . Note that once the stability of the linearised system is proved, then it is in principle possible to use the conditional stability theorem. (See for example Iooss and Joseph (1990) or Sattinger (1970) for application thereof in hydrodynamics.) The theorem basically claims that if a system described by the linearised governing equations is stable, then the stability is granted also for the system governed by the full system of nonlinear equations, provided that the initial perturbation is sufficiently small.
Concerning the restrictions on the material functions/coefficients, the ongoing analysis does not bring anything new compared to the simple analysis outlined above. The added value of the analysis of the linearised system is the new piece of information regarding the rates of the decay to the rest state. Further, the analysis will also confirm that the initial inhomogeneity in the spatial distribution of the quantities of interest is gradually smeared out.
Finally, let us remark that the analysis outlined below is purely formal, no effort is invested in the investigation of the actual smoothness of the involved functions. (It is assumed that the involved fields are as smooth as necessary.) Once one actually wants to work only with the smoothness dictated by the equations themselves and the concept of the weak solution, the issue gets more complicated. See Bulíček et al. (2009) for the discussion on the formulation of the evolution equation for the energy in the case of the standard Navier–Stokes–Fourier system.
Let us now proceed with the analysis. Let the quadruple denotes a reference state with given velocity, pressure, left Cauchy–Green tensor associated to the elastic part of the response and temperature field. Further, let the quadruple denotes the perturbation of the reference state. The total velocity, pressure, left Cauchy–Green tensor and temperature field are then given by the formulae
[TABLE]
In the present case, the reference state is the rest state, that is
[TABLE]
where and denote constant spatially uniform pressure and temperature field respectively. We again recall that if we substitute into (3), then all the equations are automatically satisfied.
Let us now find the evolution equations for the perturbation of the rest state (27). The perturbation is considered to be small, hence we are interested in the linearisation of the corresponding governing equations. If is a small quantity, then one gets, up to the first/second order,
[TABLE]
Using the expansions (28) we are now in the position to formulate the linearised governing equations for the perturbation.
3.3.1 Linearised governing equations and initial-boundary value problem for the perturbation
The linearised governing equations for the perturbation read
[TABLE]
where we have denoted , , , and . The boundary and initial conditions read
[TABLE]
where denotes the boundary of the domain and is the unit outward normal to the boundary of .
3.3.2 Time evolution of mechanical quantities and the decay of the mechanical energy to zero
Let us now multiply the second equation by and integrate over the domain of interest. (We are essentially following, in the linearised setting, the procedure popularised and elaborated by Serrin (1959) for the standard incompressible Navier–Stokes fluid. See also Straughan (2004).) The suggested manipulation yields
[TABLE]
which upon integrating by parts yields
[TABLE]
This follows from the fact that vanishes on the boundary of and that has zero divergence. Further is a symmetric tensor, hence .
The evolution equation for can be treated in the same way. The multiplication by yields after some manipulation
[TABLE]
Equation (32) is in fact the evolution equation for the net kinetic energy of the fluid in domain . Since is positive, the first term on the right hand side of (32) is nonpositive and it represents the decay of the kinetic energy due to viscosity. The second term may be either positive or negative and it represents the energy transfer between the net kinetic energy and the net elastic stored energy in the fluid.
Taking the sum of (32) and (33) then yields
[TABLE]
This equality can be further manipulated using the Korn equality
[TABLE]
and Poincaré inequality
[TABLE]
see for example Rektorys (1980), Gilbarg and Trudinger (2001) or Evans (1998), where is a domain dependent Poincaré constant. (Note that the (in)equalities hold provided that the boundary condition for the velocity perturbation is the zero Dirichlet boundary condition, see (30a).) Korn and Poincaré (in)equality allow one to rewrite (34) as
[TABLE]
The advantage of this form is the fact that the integrals on the left hand side and on the right hand side contain the same quantities and not their gradients. Performing just one manipulation and assuming that is a positive quantity finally yields the inequality
[TABLE]
This inequality has the structure
[TABLE]
where is a constant and
[TABLE]
denotes the net mechanical energy of the fluid. Consequently, the net mechanical energy decays exponentially in time and it tends to zero,
[TABLE]
where denotes the net mechanical energy at the initial time . This means that the mechanical quantities and indeed tend to the reference rest state (27) as desired.
3.3.3 Time evolution of the temperature field and the decay to a spatially homogeneous state
Concerning the temperature field, the linearised governing equation, see (29d), reads
[TABLE]
As before, the aim is to multiply the evolution equation by the perturbation and get an estimate on the gradient term in terms of the perturbation itself. However, the boundary condition for the temperature perturbation is not the zero Dirichlet condition as in the case of the velocity perturbation. Instead of (30a) we have zero Neumann boundary condition (30b). This prohibits the application of Poincaré inequality in the form (36). But we can still exploit another variant of Poincareé inequality, namely
[TABLE]
where
[TABLE]
and denotes the volume of the domain and is a domain dependent constant. (See for example Evans (1998), Gilbarg and Trudinger (2001) or Rektorys (1980).) The inequality can be exploited as follows.
First, we rewrite (42) as
[TABLE]
Then we multiply (45) by , integrate over the domain and use the integration by part on the right hand side. Since vanishes on the boundary of , the gradient of the spatially averaged quantity is zero, , and we get
[TABLE]
This equality tells us that the difference between the perturbation and its spatial average decays. It remains to show that the difference decays to zero. As before, Poincaré inequality, this time in the form (43), implies that
[TABLE]
This inequality has again the structure
[TABLE]
where is a constant and
[TABLE]
denotes the norm of the difference between the actual temperature perturbation field and its spatial average .
Since the spatial averaged temperature field is by definition spatially homogeneous, we see that the perturbed temperature field tends to a spatially homogeneous temperature field. Moreover, we have an upper bound on the rate of convergence, namely
[TABLE]
where denotes the norm of the difference at the initial time . We can observe, that the positivity of has been essential in all the manipulations. We can therefore conclude that the stability of the rest state is granted provided that . It remains to investigate what is the temperature value once the (new) steady state is reached.
This question can not be answered by appealing exclusively to the linearised system (29). However, we know what happens in the system, provided that we are willing to appeal to the full evolution equation for the energy. The evolution equation for the total energy reads
[TABLE]
where denotes the heat flux. Integrating (51) over the domain , using the Stokes theorem, and the fact that the velocity vanishes on the boundary and that there is no heat flux through the boundary, , we see that
[TABLE]
This means that the total net energy
[TABLE]
is conserved in the isolated system of interest.
Consequently, it must hold
[TABLE]
However, as the velocity and the left Cauchy–Green tensor associated to the elastic part of the response decay to rest state values, hence
[TABLE]
The terminal perturbation temperature value can be therefore obtain by solving the equation
[TABLE]
This means that the terminal perturbation temperature value is obtained by converting all the initial energy exclusively into the thermal energy.
Moreover, if the initial perturbation is chosen in such a way that it does not alter the energy of the rest state, that is if
[TABLE]
then the terminal perturbation temperature value is equal to zero and the original rest state is exactly recovered as .
4 Conclusion
We have considered a model for a viscoelastic rate type fluid with temperature dependent material parameters developed by Hron et al. (2016). The second law of thermodynamics places restrictions on the sign of material coefficients/functions characterising the dissipative phenomena (viscosity, thermal conductivity). On the other hand, the requirement on the stability of the rest state places restrictions on the material coefficients/functions characterising energy storage mechanisms (heat capacity, shear modulus).
Having explicitly investigated the implications of the requirement on the stability of the rest state, we have shown that the positivity of the coefficient/material function in the Helmholtz free energy (1) and the positivity of the specific heat capacity defined through the Helmholtz free energy as (3g) grants the stability of the rest state. Further, using the linearised governing equations, we have also provided an estimate on the decay rate of the (sufficiently small) perturbations to the rest state.
The outlined analysis explicitly demonstrates the consistency of the model by Hron et al. (2016), and the outlined methodology can be applied to any class of viscoelastic rate type phenomenological models, for example those obtained by Dressler et al. (1999). Moreover, since the model by Hron et al. (2016) provides explicit formulae for the energy and the entropy, it in principle allows one to easily construct functionals of the type
[TABLE]
that are of potential interest as Lyapunov functionals in the analysis of the evolution of finite perturbations, see for example Šilhavý (1997).
Finally, let us note that despite the practical importance of viscoelastic rate type models with temperature dependent material parameters, very few computational results for full models are available. (At least in comparison to the extensive literature dealing with purely mechanical viscoelastic rate type models.) Moreover, the existing implementations are either quite simple with respect to the used numerical methods, see for example Ionescu (2006) or Ireka (2015), or they deal only with an approximation of the true evolution equation for the temperature, see for example Damanik (2011). Given the available sophisticated numerical methods for purely mechanical viscoelastic models, see for example Kwack et al. (2016) or Sousa et al. (2016) and references therein, or incompressible Navier–Stokes fluid, see Elman et al. (2014), it is clear that the development of effective algorithms for full system of governing equations presents an important challenge for the numerical analysis.
The authors thank to the Ministry of Education, Youth and Sports of the Czech Republic (project LL1202 in the programme ERC-CZ) for its support.
The reference list from the paper itself. Each links out to its DOI / PubMed record.
- 1(1)
- 2Tanner (2009) R. I. Tanner, J. Non-Newton. Fluid Mech. 157 , 141–144 (2009).
- 3Leonov (1976) A. I. Leonov, Rheol. Acta 15 , 85–98 (1976).
- 4Leonov (1987) A. I. Leonov, J. Non-Newton. Fluid Mech. 25 , 1–59 (1987).
- 5Wapperom and Hulsen (1998) P. Wapperom, and M. A. Hulsen, J. Rheol. 42 , 999–1019 (1998).
- 6Mattos (1998) H. S. C. Mattos, Int. J. Non-Linear Mech. 33 , 97–110 (1998).
- 7Dressler et al. (1999) M. Dressler, B. J. Edwards, and H. C. Öttinger, Rheol. Acta 38 , 117–136 (1999).
- 8Guaily (2015) A. Guaily, Polymer 77 , 305–311 (2015).
